ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

一维常物性无源对流导热微分方程数值解(大概)

2020-03-04 12:07:05  阅读:316  来源: 互联网

标签:phi frac 无源 tag dx delta alpha 导热 物性


计算传热学第一次大作业

1 Taylor级数展开法

1.1 网格划分

Taylor级数展开发的网格划分使用外点法
在这里插入图片描述

1.2 离散方程表达式

ρudϕdx=Γd2ϕdx2(1)\rho u\frac{d\phi}{d x}=\Gamma \frac{d^{2}\phi}{dx^{2}}\tag{1}ρudxdϕ​=Γdx2d2ϕ​(1)
ϕ\phiϕ分别在iii+1点和iii-1点在iii点展开
ϕ(i+1)=ϕ(i)+dϕdxiδx+d2ϕdx2iδx22!+(2)\phi(i+1) = \phi(i)+\left.\frac{d\phi}{dx}\right|_{i}\delta x+\left.\frac{d^{2}\phi}{dx^{2}}\right|_{i}\frac{\delta x^{2}}{2!}+…… \tag{2}ϕ(i+1)=ϕ(i)+dxdϕ​∣∣∣∣​i​δx+dx2d2ϕ​∣∣∣∣​i​2!δx2​+……(2)
ϕ(i1)=ϕ(i)dϕdxiδx+d2ϕdx2iδx22!+(3)\phi(i-1) = \phi(i)-\left.\frac{d\phi}{dx}\right|_{i}\delta x+\left.\frac{d^{2}\phi}{dx^{2}}\right|_{i}\frac{\delta x^{2}}{2!}+…… \tag{3}ϕ(i−1)=ϕ(i)−dxdϕ​∣∣∣∣​i​δx+dx2d2ϕ​∣∣∣∣​i​2!δx2​+……(3)
联立(2)、(3)式得
dϕdxi=ϕi+1ϕi12δx,o(δx2)(4)\left.\frac{d\phi}{dx}\right|_{i}=\frac{\phi_{i+1}-\phi_{i-1}}{2\delta x},o(\delta x^{2})\tag{4}dxdϕ​∣∣∣∣​i​=2δxϕi+1​−ϕi−1​​,o(δx2)(4)
同样的,联立(2)、(3)式得
d2ϕdx2i=ϕi+12ϕi+ϕi1δx2,o(δx2)(5)\left.\frac{d^{2}\phi}{dx^{2}}\right|_{i}=\frac{\phi_{i+1}-2\phi_{i}+\phi_{i-1}}{\delta x^{2}},o(\delta x^{2})\tag{5}dx2d2ϕ​∣∣∣∣​i​=δx2ϕi+1​−2ϕi​+ϕi−1​​,o(δx2)(5)
将(4)、(5)式带入(1)式
ρuϕi+1ϕi12δx=Γϕi+12ϕi+ϕi1δx2(6)\rho u \frac{\phi_{i+1}-\phi_{i-1}}{2\delta x}=\Gamma \frac{\phi_{i+1}-2\phi_{i}+\phi_{i-1}}{\delta x^{2}}\tag{6}ρu2δxϕi+1​−ϕi−1​​=Γδx2ϕi+1​−2ϕi​+ϕi−1​​(6)
化简得
4Γϕi=(ρuδx+2Γ)ϕi1+(2Γρuδx)ϕi+1,o(δx2)(7)4\Gamma \phi_{i}=(\rho u \delta x+2\Gamma)\phi_{i-1}+(2\Gamma-\rho u \delta x)\phi_{i+1},o(\delta x^{2})\tag{7}4Γϕi​=(ρuδx+2Γ)ϕi−1​+(2Γ−ρuδx)ϕi+1​,o(δx2)(7)

2 控制容积积分法

2.1 网格划分

控制容积积分法使用内节点法划分网格

在这里插入图片描述

2.2 离散方程表达式

ρudϕdx=Γd2ϕdx2(8)\rho u\frac{d\phi}{d x}=\Gamma \frac{d^{2}\phi}{dx^{2}}\tag{8}ρudxdϕ​=Γdx2d2ϕ​(8)
将方程在空间内积分
weρudϕdxdx=weΓd2ϕdx2dx(9)\int_{w}^{e}\rho u\frac{d\phi}{d x}dx=\int_{w}^{e}\Gamma \frac{d^{2}\phi}{dx^{2}}dx\tag{9}∫we​ρudxdϕ​dx=∫we​Γdx2d2ϕ​dx(9)
其中
wedϕdxdx=ϕeϕw(10)\int_{w}^{e}\frac{d\phi}{d x}dx=\left.\phi\right|_{e}-\left.\phi\right|_{w}\tag{10}∫we​dxdϕ​dx=ϕ∣e​−ϕ∣w​(10)
假设ϕ\phiϕ对空间呈分段线性变化
ϕeϕw=ϕE+ϕP2ϕP+ϕW2=ϕEϕW2(11)\left.\phi\right|_{e}-\left.\phi\right|_{w}=\frac{\phi_{E}+\phi_{P}}{2}-\frac{\phi_{P}+\phi_{W}}{2}=\frac{\phi_{E}-\phi_{W}}{2}\tag{11}ϕ∣e​−ϕ∣w​=2ϕE​+ϕP​​−2ϕP​+ϕW​​=2ϕE​−ϕW​​(11)
对于扩散项
wed2ϕdx2dx=dϕdxedϕdxw(12)\int_{w}^{e}\frac{d^{2}\phi}{dx^{2}}dx=\left.\frac{d\phi}{dx}\right|_{e}-\left.\frac{d\phi}{dx}\right|_{w}\tag{12}∫we​dx2d2ϕ​dx=dxdϕ​∣∣∣∣​e​−dxdϕ​∣∣∣∣​w​(12)
假设dϕdx\frac{d\phi}{dx}dxdϕ​在空间呈分段线性变化
dϕdxedϕdxw=ϕEϕP(δx)eϕPϕW(δx)w(13)\left.\frac{d\phi}{dx}\right|_{e}-\left.\frac{d\phi}{dx}\right|_{w}=\frac{\phi_{E}-\phi_{P}}{(\delta x)_{e}}-\frac{\phi_{P}-\phi_{W}}{(\delta x)_{w}}\tag{13}dxdϕ​∣∣∣∣​e​−dxdϕ​∣∣∣∣​w​=(δx)e​ϕE​−ϕP​​−(δx)w​ϕP​−ϕW​​(13)
若使用均分网格
                         =ϕEϕPΔxϕPϕWΔx(14)~~~~~~~~~~~~~~~~~~~~~~~~~=\frac{\phi_{E}-\phi_{P}}{\Delta x}-\frac{\phi_{P}-\phi_{W}}{\Delta x}\tag{14}                         =ΔxϕE​−ϕP​​−ΔxϕP​−ϕW​​(14)
将式(11)、(14)带入(8)
ρu(ϕE+ϕP2ϕP+ϕW2)=Γ(ϕEϕPΔxϕPϕWΔx)(15)\rho u (\frac{\phi_{E}+\phi_{P}}{2}-\frac{\phi_{P}+\phi_{W}}{2})=\Gamma (\frac{\phi_{E}-\phi_{P}}{\Delta x}-\frac{\phi_{P}-\phi_{W}}{\Delta x})\tag{15}ρu(2ϕE​+ϕP​​−2ϕP​+ϕW​​)=Γ(ΔxϕE​−ϕP​​−ΔxϕP​−ϕW​​)(15)
整理得
4ΓϕP=(2ΓρuΔx)ϕE+(2Γ+ρuΔx)ϕW(16)4\Gamma \phi_{P}=(2\Gamma-\rho u \Delta x)\phi_{E}+(2\Gamma+\rho u \Delta x)\phi_{W}\tag{16}4ΓϕP​=(2Γ−ρuΔx)ϕE​+(2Γ+ρuΔx)ϕW​(16)

3. Gauss-Seidel迭代法

原问题的代数方程
4Γϕi=(ρuδx+2Γ)ϕi1+(2Γρuδx)ϕi+14\Gamma \phi_{i}=(\rho u \delta x+2\Gamma)\phi_{i-1}+(2\Gamma-\rho u \delta x)\phi_{i+1}4Γϕi​=(ρuδx+2Γ)ϕi−1​+(2Γ−ρuδx)ϕi+1​
α=ρuΓ\alpha = \frac{\rho u}{\Gamma}α=Γρu​则
4ϕi=(αδx+2)ϕi1+(2αδx)ϕi+1(17)4 \phi_{i}=(\alpha\delta x+2)\phi_{i-1}+(2-\alpha \delta x)\phi_{i+1}\tag{17}4ϕi​=(αδx+2)ϕi−1​+(2−αδx)ϕi+1​(17)
假设有NNN个格点,则该方程得Guass-Seidel迭代格式
{ϕ1k+1=14[(αδx+2)ϕ0+(2αδx)ϕ2k]ϕ2k+1=14[(αδx+2)ϕ1k+1+(2αδx)ϕ3k]ϕNk+1=14[(αδx+2)ϕN1k+1+(2αδx)ϕN+1]\left\{ \begin{aligned} \phi_{1}^{k+1} & = &\frac{1}{4}[(\alpha\delta x+2)\phi_{0}+(2-\alpha \delta x)\phi_{2}^{k}] \\ \phi_{2}^{k+1} & = & \frac{1}{4}[(\alpha\delta x+2)\phi_{1}^{k+1}+(2-\alpha \delta x)\phi_{3}^{k}] \\ \vdots\\ \phi_{N}^{k+1} & = & \frac{1}{4}[(\alpha\delta x+2)\phi_{N-1}^{k+1}+(2-\alpha \delta x)\phi_{N+1}] \end{aligned} \right. ⎩⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎧​ϕ1k+1​ϕ2k+1​⋮ϕNk+1​​===​41​[(αδx+2)ϕ0​+(2−αδx)ϕ2k​]41​[(αδx+2)ϕ1k+1​+(2−αδx)ϕ3k​]41​[(αδx+2)ϕN−1k+1​+(2−αδx)ϕN+1​]​
写为矩阵形式:
(42αδx2+αδx2αδx2+αδx4)(ϕ1ϕ2ϕ3ϕN)((2+αδx)ϕ000(2αδx)ϕN+1)(18) \begin{pmatrix} -4 & 2-\alpha \delta x & {} & {} & {} \\ 2+\alpha \delta x & {} & {} & {} & {} \\ {} & {} & \ddots & {} & {} \\ {} & {} & {} & {} & 2-\alpha \delta x \\ {} & {} & {} & 2+\alpha \delta x & -4 \\ \end{pmatrix} \begin{pmatrix} {{\phi }_{1}} \\ {{\phi }_{2}} \\ {{\phi }_{3}} \\ \vdots \\ {{\phi }_{N}} \\ \end{pmatrix} \doteq \begin{pmatrix} -(2+\alpha \delta x){{\phi }_{0}} \\ 0 \\ 0 \\ \vdots \\ -(2-\alpha \delta x){{\phi }_{N+1}} \\ \end{pmatrix}\tag{18} ⎝⎜⎜⎜⎜⎛​−42+αδx​2−αδx​⋱​2+αδx​2−αδx−4​⎠⎟⎟⎟⎟⎞​⎝⎜⎜⎜⎜⎜⎛​ϕ1​ϕ2​ϕ3​⋮ϕN​​⎠⎟⎟⎟⎟⎟⎞​≐⎝⎜⎜⎜⎜⎜⎛​−(2+αδx)ϕ0​00⋮−(2−αδx)ϕN+1​​⎠⎟⎟⎟⎟⎟⎞​(18)
解得方程在α=(5,1,0,1,5)\alpha=(-5,-1,0,1,5)α=(−5,−1,0,1,5)时的曲线
在这里插入图片描述

4. TDMA方法

对于任意一个三对角阵
AiTi=BiTi+1+CiTi1+Di(19) A_{i}T_{i}=B_{i}T_{i+1}+C_{i}T_{i-1}+D_{i}\tag{19} Ai​Ti​=Bi​Ti+1​+Ci​Ti−1​+Di​(19)
我们都可以将其转化至
Ti1=Pi1Ti+Qi1(20) T_{i-1}=P_{i-1}T_{i}+Q_{i-1}\tag{20} Ti−1​=Pi−1​Ti​+Qi−1​(20)
联立(19)、(20)可得
TiCiPiTi=BiTi+1+Di+CiQi1(21) T_{i}-C_{i}P_{i}T_{i}=B_{i}T_{i+1}+D_{i}+C_{i}Q_{i-1}\tag{21} Ti​−Ci​Pi​Ti​=Bi​Ti+1​+Di​+Ci​Qi−1​(21)
整理得
Ti=BiAiCiPi1Ti+1+Di+CiQi1AiCiPi1(22) T_{i}=\frac{B_{i}}{A_{i}-C_{i}P_{i-1}}T_{i+1}+\frac{D_{i}+C_{i}Q_{i-1}}{A_{i}-C_{i}P_{i-1}}\tag{22} Ti​=Ai​−Ci​Pi−1​Bi​​Ti+1​+Ai​−Ci​Pi−1​Di​+Ci​Qi−1​​(22)
观察(20)、(22)式,可以发现
Pi=BiAiCiPi1;Qi=Di+CiQi1AiCiPi1(23) P_{i}=\frac{B_{i}}{A_{i}-C_{i}P_{i-1}}; Q_{i}=\frac{D_{i}+C_{i}Q_{i-1}}{A_{i}-C_{i}P_{i-1}}\tag{23} Pi​=Ai​−Ci​Pi−1​Bi​​;Qi​=Ai​−Ci​Pi−1​Di​+Ci​Qi−1​​(23)
结合代数形式(16)及矩阵形式(18)、(19)与(23)可以得到有NNN个格点的各个系数向量
A=(4   4   4   4   4   4   4   4   4         4)TB=(2αΔx  2αΔx  2αΔx    0)TC=(0  2+αΔx  2αΔx  2αΔx  )TD=((2+Δx)ϕ0  0  0    (2Δx)ϕN+1)T \begin{aligned} A & = & (4~~~4~~~4~~~4~~~4~~~4~~~4~~~4~~~4~~~~~\cdots~~~~4)^{T} \\ B & = & (2-\alpha \Delta x~~2-\alpha \Delta x~~2-\alpha \Delta x~~\cdots~~0)^{T} \\ C & = & (0~~2+\alpha \Delta x~~2-\alpha \Delta x~~2-\alpha \Delta x~~\cdots)^{T} \\ D & = & ((2+\Delta x)\phi_{0}~~0~~0~~\cdots~~(2-\Delta x)\phi_{N+1})^{T} \end{aligned} ABCD​====​(4   4   4   4   4   4   4   4   4     ⋯    4)T(2−αΔx  2−αΔx  2−αΔx  ⋯  0)T(0  2+αΔx  2−αΔx  2−αΔx  ⋯)T((2+Δx)ϕ0​  0  0  ⋯  (2−Δx)ϕN+1​)T​








解得方程在α=(5,1,0,1,5)\alpha=(-5,-1,0,1,5)α=(−5,−1,0,1,5)时的曲线

在这里插入图片描述

5. 网格独立化考核

网格数量对数值计算结果有着重要的影响。当网格足够细密以至于再进一步加密网格已对数值计算结果基本上没有影响时所得到的数值解称为网格独立解。
在这里插入图片描述

为了避免浪费计算资源,在精度足够高的同时,减少网格数,这里使用TDMA方法计算在不同网格划分情况下xxx=0.7处的值,进行网格独立化考核。

显然,当格点数大于80时,再增加格点数时ϕ\phiϕ(0.7)的值其变化已经不明显了,故我们认为格点数达到80时,此时的解为网格独立解。

标签:phi,frac,无源,tag,dx,delta,alpha,导热,物性
来源: https://blog.csdn.net/qq_41640870/article/details/104650455

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有