对流扩散方程的离散格式
对流扩散问题有限体积法
流体仿真与应用第八讲二、对流-扩散问题的有限体积法◆中心差分格式（例子）节点增加到20个结果◆离散格式的性质在数学上，一个离散格式必须要引起很小的误差（包括离散误差和舍入误差）才能收敛于精确解，即要求离散格式必须要稳定或网格必须满足稳定性条件。
在物理上，离散格式所计算出的解必须要有物理意义，对于得到物理上不真实的解的离散方程，其数学上精度再高也没有价值。
通常，离散方程的误差都是因离散而引起，当网格步长无限小时，各种误差都会消失。
然而，在实际计算中，考虑到经济性（计算时间和所占的内存）都只能用有限个控制容积进行离散。
因此，格式需要满足一定的物理性质，计算结果才能令人满意。
主要的物理性质包括：守恒性、有界性和迁移性。
◆离散格式的性质——守恒性满足守恒性的离散方程不仅使计算结果与原问题在物理上保持一致，而且还可以使对任意体积（由许多个控制容积构成的计算区域）的计算结果具有对计算区域取单个控制容积上的格式所估计的误差。
◆离散格式的性质——迁移性③当Pe 为有限大小时，对流和扩散同时影响一个节点的上、下游相邻节点。
随着Pe 的增加，下游受的影响逐渐增大，而上游受的影响逐渐变小。
①，即纯扩散，无对流。
②，即纯对流，无扩散。
0=Pe ∞=Pe◆迎风格式迎风格式（Upwind Differencing Scheme ）在确定控制容积界面上的值时就考虑了流动的方向性，其思想为：在控制容积界面上对流项的取上游节点处的值，称之为第二类迎风格式。
中心差分格式的缺点是，它不能识别流动的方向，控制容积界面上的值取相邻上、下游节点的平均值。
当对流作用较强时，这样的处理就与其物理特征（某点的值受上游的影响，而不受下游的影响）不一致了。
φφφ◆迎风格式◆迎风格式在控制容积界面上对流项的取其上游节点处的值EW →φWw φφ=Pe φφ=()()W P w P E e W w P e D D F F φφφφφφ−−−=−()()[]()Ee W w w P w e e w w D F D F F D F D φφφ++=−+++WE →Pw φφ=Ee φφ=()()[]()Ee e W w Pw e e e w F D D F F F D D φφφ−+=−+−+◆迎风格式通用形式WW E E P P a a a φφφ+=()w e E W P F F a a a −++=EW →ww W F D a +=eE D a =W E →w W D a =ee E F D a −=◆迎风格式的特点迎风格式满足守恒性。
第五章对流扩散方程的离散格式
aP = aE + aW
aE = De – Fe / 2 aW = Dw + Fw / 2
在流场的实际求解过程中， 每一个迭代层次上，即使速度 场尚未收敛，也要保证连续方 程是满足的。
3. 对流项的中心差分与迎风差分
3.2 对流项的中心差分
三点说明：
系数 aE , aW 包含了对流 F 与扩散 D的作用的影响；
对均分网格：
2. 对流项离散格式的重要性及两种离散方式
2.2 构造对流项离散格式的两种方式 （2）控制容积积分
给出界面上被求函数的插值方式
2. 对流项离散格式的重要性及两种离散方式
2.2 构造对流项离散格式的两种方式 （3）两种定义之间的关系
对某种对流项离散格式，都可以用两种方法给出其相应 的定义；
3.1 一维对流-扩散问题模型方程的精确解
上游优势
3. 对流项的中心差分与迎风差分
3.1 一维对流-扩散问题模型方程的精确解
希望所构建的离散方程形式也具有这样的物理特性
3. 对流项的中心差分与迎风差分
3.2 对流项的中心差分 Central Scheme （CS）
分段线性
均分网格
令
对流项
----界面上的流量
1. 简 介
对流与扩散作用在物理本质上的区别
从物理过程来看，扩散作用与对流作用在传递信息或扰动方面 的特点有很大区别：
扩散是由分子的不规则热运动所致，分子不规则热运动对空间不同方向
的几率都是一样的，因而扩散过程可以把发生在某一地点上的扰动的影响 向各个方向传递。
对流是流体微团宏观的定向运动，带有强烈的方向性。在对流的作用下，
两种定义方式的截断误差阶数是一致的，均为二阶截差 （中心差分，分段线性）；
CFD大作业解析(对流扩散方程)
求解迁移方程
求解域为 初始条件为
求解当 t = 0.7 时，u 的值。
∂u ∂u ∂2u ∂t + ∂x = α ∂x2
x ∈ [0.0, 1.0]
{
u=
2.0, 0.1 < x < 0.2; 1.0, others.
2 迁移方程 α = 0 差分格式的构造
2.1 显式 Euler 法 + 空间后差
空间采取 Euler 后差。构造的差分格式如下
uni +1 − uni = − uni − uni−1
∆t
∆x
将其化为可求解的格式
uni +1
=
(1
−
∆t ∆x
)uni
+
∆t ∆x
uni−1
2.2 显式 Euler 法 + 空间前差
空间采取 Euler 前差。构造的差分格式如下
uni +1 − uni = − uni+1 − uni
∂u + ∂u = 2π cos(2πx) ∂t ∂x
边界条件为 x = 0, u = 0; x = 1, u = 0。 粘性项的验证，令 a = 0α = −1，令理论解为 u = x2 − x，则控制方程为
∂u ∂2u ∂t + ∂x2 = 2x
边界条件为 x = 0, u = 0; x = 1, u = 0。
1 N
∑精Ni=度1(测ui试−的Ui方)，法以为l空n N间, l离n e散为点坐数标N作图取，为计5算0、两1点00间、的20斜0，率进。行时间推进，直到达到理论解，计算误差
e
=
6 作业要求
• 编写程序，实现 2 中列出的四种方法，完成空间二阶导数的离散
哈尔滨工业大学 计算传热学 第五章 对流-扩散方程的离散格式-2013
aPP aEE aWW
Fe Fw exp( Pw ) aE , aW exp( Pe ) 1 exp( Pw ) 1
(D)
aP aE aW (Fe Fw )
区别就在函数 aE和aW
aE De
Pe aE De exp( Pe ) 1
aE Pe De
该格式计算量比指数小，且指数格式的解差别很小。
§ 5-3
为了在讨论中引入 PE 记
通用表达式
x
i
J*
i+1 i+1/2
x
1 界面i+ 上的值可以用界面两侧节点值表示 2
J * Bi Ai 1 （y）
系数A和B的性质的讨论 (1)当 i i 1 时，扩散量=0， J *完全由对流造成，即
即
aPP aEE aW W
显然不论那种格式，仅仅是 A(| P |) 表达式的区别。
A( P )
A(|P |)
中心 1 0.5 | P | 迎风 1 混合 [| 0,1 0.5 | P | |] 指数 | P | [exp(| P |) 1]
1.0
迎风
指数 乘方
乘方 | 0, (1 0.1| P |)5 |
中心
混合
P
§ 5-4
原始的假扩散概念
关于假扩散的讨论
一维非稳态对流方程（纯对流,没有扩散）
u t x
显示迎风差分格式
in1 in
t
u
in in 1
x
, o(x, t )
将上式在（i，n）点做Taylar级数展开，保留二阶。
上述若对任何成立，必得
B( P ) A( P ) A( P ) B( P )
数值传热第六章作业
6-3 试在直角坐标系的交错网格上，写出动量离散方程式（6-5）、（6-6）中的系数nb a （即S N W E a a a a ,,,）,n n e e A a A a ,,,的表达式。
为简便起见，设（1）流体物性为常数；（2）在x, y 方向上网格各自均匀划分。
速度e u 的邻点可参阅图6-5， 速度n υ的邻点参见图6-32.对流、扩散项的离散可采用五种三点格式之一。
解：根据课本P145式(5-13)、(5-16)、(5-18)，对流、扩散项采用指数格式计算本题 在二维直角坐标系中，对流—扩散方程的通用形式为：()()()φφφφφρυφφρρφS y y x x y u x t +⎪⎪⎭⎫⎝⎛∂∂Γ∂∂+⎪⎭⎫ ⎝⎛∂∂Γ∂∂=∂∂+∂∂+∂∂ 对于动量方程，把压力梯度项放到源项中了。
引入在x 及y 方向的对流—扩散总通量密度，上式可改写为：()⎪⎪⎭⎫ ⎝⎛∂∂+∂∂-=⎪⎪⎭⎫⎝⎛∂∂Γ-∂∂+⎪⎭⎫ ⎝⎛∂∂Γ-∂∂+∂∂y p x p S y y x u x t φρυφφφρρφφφ 即：()⎪⎪⎭⎫ ⎝⎛∂∂+∂∂-=∂∂+∂∂+∂∂y p x p S y J x J t yx ρφ （1） 其中：yJ xu J y x ∂∂Γ-=∂∂Γ-=φρυφφφρφφ将（1）式对P 控制容积做时间与空间上的积分得：e E P p P c s n w e pp A P P V S S J J J J V t)()()()()()(0-+∆+=-+-+∆∆-φρφρφ将通用变量φ换成速度u ，相应的其控制容积变为：所以上式可改写为：e E P p p c s n w e ee A P P V S S J J J J V t u u )()()()()()(0-+∆+=-+-+∆∆-φρρ （2）式（6-5）为：()e E P nb nbe e A p p b u au a -++=∑对上式用界面总通量表达式为：ee E e e EE e u a u F a J -+=)( （3）e w W w W w u F a u a J )(--= （4）n N e n N n u a u F a J -+=)( （5）e s S s S s u F a u a J )(--= （6） 把以上方程代入方程（2）得：e E p e p c e s S s S n N en N w w e w W ee E e e EE ee A P P V u S S u F a u a u a u F a u a u F a u a u F a V tu u )()()()()()(0-+∆+=-+--++--+-++∆∆-ρρ整理得：eE p ec s S n N w W ee E ep s S n N w W e EE A P P u tVV S u a u a u a u a u V S F a F a F a F a tV)(])()()()([0-+∆∆+∆++++=∆--+++-+++∆∆ρρ当对流、扩散项的离散采用指数格式时， 则上式中的系数分别为：1)ex p()(-==∆∆e ee e EE P F P A D a 1)e x p ()e x p ()(-==∆∆∆w w w w w W P P F P B D a1)ex p()(-==∆∆n n n n N P F P A D a 1)e x p ()e x p ()(-==∆∆∆s s s s s S P P F P B D a tVa e ∆∆=ρ0V S a F F F F a a a a a p e s n w e S N W EE e ∆-+-+-++++=00e e c u a V S b +∆=y A e ∆=同理对（6－6）()n N P nb nbn n A p p b aa -++=∑υυ，类似地有：1)ex p()(-==∆∆n n n n NN P F P A D a 1)e x p ()e x p ()(-==∆∆∆s s s s s S P P F P B D a 1)ex p()(-==∆∆e e e e E P F P A D a 1)e x p ()e x p ()(-==∆∆∆ww w w w W P P F P B D atVa n ∆∆=ρ0V S a F F F F a a a a a p n s e w n S E W NN n ∆-+-+-++++=000n n c u a V S b +∆=x A n ∆=6-4 对图6-11所示的二维流动情形，已知：10,0,20,50====E N s w p p v u 流动是稳态的，且密度为常数。
输运方程对流扩散方程
输运方程对流扩散方程输运方程是描述物质传输过程的数学模型，常见的有对流扩散方程。
对流扩散方程是由对流和扩散两种机制共同产生的输运过程来描述的，它的一般形式为：∂c/∂t+∇·(v*c)=∇·(D*∇c)其中，c表示物质的浓度或者响应变量，t表示时间，v表示流体的速度场，D表示物质的扩散系数，∇表示梯度运算符。
对流项描述了物质的对流运动，即物质随着流体的移动而移动。
对于三维坐标系来说，对流项可以表示为∇·(v*c)。
具体来说，对流项的每一项分别表示了物质在x、y和z方向上的携带速度与浓度梯度的乘积。
扩散项描述了物质由浓度高处至浓度低处的扩散现象，即物质自发性地从高浓度区域向低浓度区域传播。
扩散项可以表示为∇·(D*∇c)，其中D是扩散系数，表示物质扩散的速率与浓度梯度的乘积。
对流扩散方程的物理意义是描述了物质在流体中传输的速率与物质浓度梯度之间的关系。
通过对流项，方程能够描述物质随着流体的运动快速传输的现象；而通过扩散项，方程能够描述物质由浓度高处向浓度低处传输的现象。
综合考虑对流和扩散的作用，对流扩散方程能够比较准确地描述物质在流体中的传输过程。
对流扩散方程在科学和工程领域有广泛的应用。
例如，在污染物传输和扩散模拟中，对流扩散方程可用于描述污染物由源区到周围空气或水体的传输过程。
在热传导模拟中，对流扩散方程可用于描述热量由高温区域到低温区域的传导过程。
在物质传递过程中，对流扩散方程也被广泛应用于描绘物质的传输行为。
总结起来，对流扩散方程是一种常见的输运方程，它能够描述物质由流体传输并扩散的过程。
通过对流项和扩散项的综合作用，对流扩散方程能够比较准确地描述物质在流体中的传输行为，所以在科学和工程领域有着广泛的应用。
2.4常用的离散格式
低阶格式的假扩散特性
迎风格式,指数格式,混合格式及乘方格式等 一阶格式应用于实际问题时都可能引起较严 重的假扩散,这在HVAC领域的高大空间流体 流动及传热计算中尤为明显. 因此,为了有效地克服或减轻假扩散所带来的 计算误差,空间导数应当采用二阶或更高阶的 格式(如QUICK格式,二阶迎风差分格式等).
离散格式
假设速度场已知，则为求解离散方程，需计算广义未 知量在边界e和w处的值。
为完成这一任务，必须决定界面物理量如何通过节点 物理量的插值表示。
各种不同的插值方法就构成了不同的离散格式。
中心差分格式
一阶迎风格式
混合格式
指数格式
乘方格式
1
2.4.1术语的约定
对离散格式的讨论以一维稳态对流扩散方程为例，不 涉及瞬态项。
3
Central differencing scheme 中心差分格式
(x) P
P
interpolated value
e E
eE
We determine the value of at the face by linear
interpolation between the cell centered values.就是界 面上的物理量采用线性插值公式来计算。
基于此限制，中心差分格式不能作为对于一般 流动问题的离散格式，需创建其它更合适的格 式（对纯扩散稳态，如热传导是适用的）。
5
对流扩散方程的精确解
6
精确解随Pe数的变化
（Pe＝0纯扩散，Pe增大对流增强）
7
具体算例
（不同计算工况意味着不同Pe数）
8
第一种工况Pe=0.2
尽管网格粗糙，但数值解与精确解非常接近。
对流-扩散方程的离散格式
第5章 对流-扩散方程的离散格式
2009年3月13日
1/59
传热与流体流动的数值计算
§5.1 对流项离散格式的重要性 及两种离散方式
一、对流项离散格式的重要性
1、数值解的准确性（假扩散） 2、数值解的稳定性 3、数值解的经济性
二、构造离散格式的两种方式
1、Taylor展开法 2、控制容积积分法
u e
Fee
P
max Fe,0
E
max Fe,0
w界面
P Fe ,0 E Fe ,0
uw 0 , W ; uw 0 , P
u w
Fww
W
max Fw,0
P
max Fw,0
W Fw,0 P Fw,0
9/59
传热与流体流动的数值计算
三、对流项的迎风格式（续）
3、一阶迎风格式截差阶数低，除非采用相当密的网格， 否则计算结果的误差较大。
4、一阶迎风格式的启示：应当在迎风方向取更多的信 息构造格式，更好地反映对流过程的物理本质。
5、在调试程序或计算的中间过程仍可以采用一阶迎风 格式。
11/59
传热与流体流动的数值计算
§5.3 对流-扩散方程的混合格式及乘方格式
一、通量密度及其离散表达式
d dx
u
d dx
d
dx
总通量密度J：单位时间内、单位面积上由扩散
及对流作用而引起的某一物理量的总转移量。
J
u
d
dx
x
P
d
d
x
x
J*
J D
P
d
d
x
x
18/59
传热与流体流动的数值计算
数值分析-对流项的离散
数值分析-对流项的离散对流项的数值离散扩散项的离散在先前已经介绍，对于动量⽅程，更为棘⼿的问题在于N-S⽅程中对流项和压⼒项的离散处理。
在不可压流场的数值求解过程中，最重要的问题均是这两个⼀阶偏导项所引起的。
这⾥主要讨论关于对流-扩散项的离散⽅法。
对流项的离散⽅法⾸先是关于对流项的两种离散思路的介绍Taylor展开离散Taylor展开的⽅法，就是直接套⽤泰勒公式将偏导转化为线性离散控制容积积分通过将对流项的⼀阶导数在控制容积P中进⾏积分，利⽤相邻两个结点的值来获得控制容积边界e，w处的插值两种思路所给出的截断误差的阶数是相同的⼀维稳态对流-扩散的离散由于整个N-S⽅程还是过于复杂，这⾥我们对整个物理问题可以进⾏相应的简化，从仅含有对流项和扩散的场景出发，来讨论如何解决对流项中的⼀阶导离散问题。
因此，⾸先对于⼀维稳态⽆内热源的对流扩散问题的数学描述，我们可以很轻松的给出⽅程的守恒形式：ddx(ρuϕ)=ddx(Γdϕdx)边界条件写为：ϕ(0)=ϕ0；ϕ(L)=ϕL值得注意的是，以下内容表⽰的是对对流项（convection term）采⽤不同的差分格式进⾏描述，⽽公式中的扩散项仍然保持中⼼差分格式。
中⼼差分格式采⽤控制容积法，对流项的中⼼差分格式相当于在控制容积P的界⾯上取分段的线性函数，即(ρuϕ)e−(ρuϕ)w=(ρu)e ϕE+ϕP2−(ρu)wϕP+ϕW2最终的离散结果可以写为：ϕP[(ρu)e2−(ρu)w2+Γe(δx)e+Γw(δx)w]=ϕE(Γe(δx)e−(ρu)e2)+ϕW(Γw(δx)w+(ρu)w2)这⾥可以将界⾯上的流量记为$ F = \rho u，界⾯上单位扩散阻⼒的倒数记为D = \frac{\Gamma}{\delta x}$最终式(2)可以化简为：a PϕP=a EϕE+a WϕWpython代码可以轻松实现中⼼差分算法的描述：# 中⼼差分格式（x结点数，y结点数，x⽅向扩散，y⽅向扩散，x界⾯流量，y界⾯流量，linux下debug符号）def CDiff(x_nodes,y_nodes,D_x, D_y, F_x, F_y, Boundary, debug="n"):"""Requires:x_nodes = number of nodes along xy_nodes = number of nodes along yD_x = Diffusivity along x, type = floatD_y = Diffusivity along y, type = floatF_x = Force along X, type = matrix (x by y)F_y = Force along Y, type = matrix (x by y)Optional:debug, "y" if you want debug, off by default"""## setting up stuff (simple math)dim = x_nodes * y_nodes # specified to reduce calculationsw_bound = [0] * y_nodes # creating 0-matricies to finde_bound = [0] * y_nodes # all boundary nodess_bound = [0] * x_nodesn_bound = [0] * x_nodesif debug == 'y': # Debugger notifies progressprint("finding boundaries")for i in range (x_nodes): # boundary nodes are populateds_bound[i] = i # used to populate matriciesn_bound[i] = (dim - x_nodes + i)for i in range (y_nodes):w_bound[i] = (i * x_nodes)e_bound[i] = (i+1) * x_nodes -1Processing math: 100%if debug == 'y': # Debugger prints boundary nodesxB = "n:" + str(n_bound) + "\n"+ " s:" + str(s_bound)yB = "e:" + str(e_bound) + "\n"+ " w:" + str(w_bound)print("boundaries: \n", xB, "\n", yB)# Declaring variablesa_w = [0] * dim # Creating empty matricies for all nodesa_e = [0] * dima_s = [0] * dima_n = [0] * dimsp_p = 0su_p = 0B = [0] * dimsp_W = [0] * dimsu_W = [0] * dimsp_E = [0] * dimsu_E = [0] * dimsp_S = [0] * dimsu_S = [0] * dimsp_N = [0] * dimsu_N = [0] * dim#Boundariesphi_W = Boundary[0]phi_E = Boundary[1]phi_S = Boundary[2]phi_N = Boundary[3]## answer matricies ##a = [0] * dim # make afor i in range(dim): # 0-arraya[i] = [0] * dimB = [0] * dim# POPULATING USING CENTRAL DIFFERENCEprint("populating central difference")pbar = tqdm(range(dim)) # load progress barfor i in range(dim):pbar.update(1) # progress indicatorif i not in w_bound: # populating non-boundary conditionsa[i][i - 1] = a[i][i - 1] - (D_x + (F_x[i - 1] / 2))a[i][i] = a[i][i] + (D_x + (F_x[i - 1] / 2))if i not in e_bound:a[i][i + 1] = a[i][i + 1] - (D_x - (F_x[i + 1] / 2))a[i][i] = a[i][i] + (D_x - (F_x[i + 1] / 2))if i not in s_bound:a[i][i - x_nodes] = a[i][i - x_nodes] - (D_y + (F_y[i - x_nodes] / 2))a[i][i] = a[i][i] + (D_y + (F_y[i - x_nodes] / 2))if i not in n_bound:a[i][i + x_nodes] = a[i][i + x_nodes] - (D_y - (F_y[i + x_nodes] / 2))a[i][i] = a[i][i] + (D_y - (F_y[i + x_nodes] / 2))if i in w_bound: # populating boundary conditionssp_W[i] = -(2*D_x + F_x[i])su_W[i] = (2*D_x + F_x[i])* phi_Wa[i][i] = a[i][i] - sp_W[i]B[i] = B[i] + su_W[i]if i in e_bound:sp_E[i] = -(2*D_x - F_x[i])su_E[i] = (2*D_x - F_x[i])* phi_Ea[i][i] = a[i][i] - sp_E[i]B[i] = B[i] + su_E[i]if i in s_bound:sp_S[i] = -(2*D_y + F_y[i])su_S[i] = (2*D_y + F_y[i])* phi_Sa[i][i] = a[i][i] - sp_S[i]B[i] = B[i] + su_S[i]if i in n_bound:sp_N[i] = -(2*D_y - F_y[i])su_N[i] = (2*D_y - F_y[i]) * phi_Na[i][i] = a[i][i] - sp_N[i]B[i] = B[i] + su_N[i]return a, B迎风差分格式利⽤控制容积法，对对流项的⼀阶迎风差分格式进⾏描述：在e界⾯上u e>0,ϕ=ϕP u e<0,ϕ=ϕE 在w界⾯上u w>0,ϕ=ϕW u w<0,ϕ=ϕP 因此对流迎风差分格式可以表达为：(ρuϕ)e=F eϕe=ϕP MAX[F e,0]−ϕE MAX[−F e,0]最终的离散结果可写为：ϕP[MAX[F e,0]+MAX[−F w,0]+Γe(δx)e+Γw(δx)w]=ϕE(Γe(δx)e+MAX[−F e,0])+ϕW(Γw(δx)w+MAX[F w,0])需要注意的是：在数值解不出现振荡的范围内，中⼼差分要⽐迎风差分结果的误差⼩⼀阶迎风格式在应⽤上误差较⼤，应⽤需由限制，但新发展出的⼆阶迎风、三阶迎风都从中吸取了思想python代码中同样也可以实现⼀阶迎风差分算法的描述：# ⼀阶迎风差分格式（x结点数，y结点数，x⽅向扩散，y⽅向扩散，x界⾯流量，y界⾯流量，linux下debug符号）def UDiff(x_nodes,y_nodes,D_x, D_y, F_x, F_y, Boundary, debug="n"):"""Requires:x_nodes = number of nodes along xy_nodes = number of nodes along yD_x = Diffusivity along x, type = floatD_y = Diffusivity along y, type = floatF_x = Force along X, type = matrix (x by y)F_y = Force along Y, type = matrix (x by y)Optional:debug, "y" if you want debug, off by default"""## setting up stuff (simple math)dim = x_nodes * y_nodes # specified to reduce calculationsw_bound = [0] * y_nodes # creating 0-matricies to finde_bound = [0] * y_nodes # all boundary nodess_bound = [0] * x_nodesn_bound = [0] * x_nodesif debug == 'y': # Debugger notifies progressprint("finding boundaries")for i in range (x_nodes): # boundary nodes are populateds_bound[i] = i # used to populate matriciesn_bound[i] = (dim - x_nodes + i)for i in range (y_nodes):w_bound[i] = (i * x_nodes)e_bound[i] = (i+1) * x_nodes -1if debug == 'y': # Debugger prints boundary nodesxB = "n:" + str(n_bound) + "\n"+ " s:" + str(s_bound)yB = "e:" + str(e_bound) + "\n"+ " w:" + str(w_bound)print("boundaries: \n", xB, "\n", yB)dim = x_nodes* y_nodesa_w = [0] * dim # Creating empty matricies for all nodesa_e = [0] * dima_s = [0] * dima_n = [0] * dimsp_p = 0su_p = 0B = [0] * dimsp_W = [0] * dimsu_W = [0] * dimsp_E = [0] * dimsu_E = [0] * dimsp_S = [0] * dimsu_S = [0] * dimsp_N = [0] * dimsu_N = [0] * dim#Boundariesphi_W = Boundary[0]phi_E = Boundary[1]phi_S = Boundary[2]phi_N = Boundary[3]## answer matricies ##a = [0] * dim # make afor i in range(dim): # 0-arraya[i] = [0] * dimB = [0] * dim# POPULATING USING CENTRAL DIFFERENCEprint("populating upwind difference")pbar = tqdm(range(dim))for i in range(dim):pbar.update(1)if i not in w_bound: # populating non-boundary nodesa[i][i - 1] = a[i][i - 1] - (D_x + max((F_x[i - 1]),0))a[i][i] = a[i][i] + (D_x + max((F_x[i - 1]),0))if i not in e_bound:a[i][i + 1] = a[i][i + 1] - (D_x + max(-(F_x[i + 1]),0))a[i][i] = a[i][i] + (D_x + max(-(F_x[i + 1]),0))if i not in s_bound:a[i][i - x_nodes] = a[i][i - x_nodes] - (D_y + max((F_y[i - x_nodes]),0))a[i][i] = a[i][i] + (D_y + max((F_y[i - x_nodes]),0))if i not in n_bound:a[i][i + x_nodes] = a[i][i + x_nodes] - (D_y + max(-(F_y[i + x_nodes]),0)) a[i][i] = a[i][i] + (D_y + max(-(F_y[i + x_nodes]),0))if i in w_bound: # populating boundary nodessp_W[i] = -(2*D_x + max((F_x[i]),0))su_W[i] = (2*D_x + max((F_x[i]),0))* phi_Wa[i][i] = a[i][i] - sp_W[i]B[i] = B[i] + su_W[i]if i in e_bound:sp_E[i] = -(2*D_x + max(-(F_x[i]),0))su_E[i] = (2*D_x + max(-(F_x[i]),0))* phi_Ea[i][i] = a[i][i] - sp_E[i]B[i] = B[i] + su_E[i]if i in s_bound:sp_S[i] = -(2*D_y + max((F_y[i]),0))su_S[i] = (2*D_y + max((F_y[i]),0))* phi_Sa[i][i] = a[i][i] - sp_S[i]B[i] = B[i] + su_S[i]if i in n_bound:sp_N[i] = -(2*D_y + max(-(F_y[i]),0))su_N[i] = (2*D_x + max(-(F_y[i]),0))* phi_Na[i][i] = a[i][i] - sp_N[i]B[i] = B[i] + su_N[i]return a, B混合格式乘⽅格式指数格式假扩散。
《对流扩散方程》课件
环境科学
描述污染物在大气、水体等环境 介质中的扩散、输移和归宿。
在环境科学中，对流扩散方程用 于模拟污染物在大气、水体等环 境介质中的扩散、输移和归宿过
程。
在环境保护、污染治理等领域， 对流扩散方程具有重要的应用价
值。
化学反应动力学
描述化学反应在流体或固定床 反应器中的传递和反应过程。
在化学反应动力学中，对流 扩散方程用于模拟化学反应 在流体或固定床反应器中的
初始条件
指定在求解开始时刻的解的性质，如 常数、函数等。
03 对流扩散方程的应用
流体动力学
01
描述流体在运动状态下的物质传递和扩散现象。
02
在流体动力学中，对流扩散方程用于模拟流体中的物质传递过
程，如温度、浓度、速度等。
在航空航天、船舶、汽车等领域的流体动力学分析中，对流扩
03
散方程被广泛应用。
应用于多尺度问题
研究对流扩散方程在多尺度问题中的应用，如 微纳尺度流动、大气污染扩散等。
探索新的应用领域
将该方程应用于其他领域，如生物医学、环境科学等。
与其他领域的交叉研究
与流体动力学结合
研究对流扩散方程与流体动力学之间的相互 作用和影响，探索更深入的物理机制。
与偏微分方程理论的交叉
将对流扩散方程的研究与偏微分方程理论相 结合，推动数学理论的发展。
02
03
有限体积法
将连续的求解域离散化为有限个小的 体积，在每个体积上近似函数，将微 分方程转化为代方程进行求解。
有限差分法
向前差分法
将微分方程中的导数项用前一步的函数值近似代替，得到向前差 分方程。
向后差分法
将微分方程中的导数项用后一步的函数值近似代替，得到向后差 分方程。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
§5.1 对流项离散格式的重要性 及两种离散方式 一、对流项离散格式的重要性
1、数值解的准确性（假扩散） 2、数值解的稳定性 3、数值解的经济性
二、构造离散格式的两种方式
1、Taylor展开法 2、控制容积积分法
1/59
两种定义截差阶数一致，但截差首项系数有所不同。
传热与流体流动的数值计算
2/59
传热与流体流动的数值计算
一、一维对流-扩散问题模型方程的精确解（续）
0 eux 1 e Pex L 1 uL Pe L 0 e 1 e 1
Peclet数：
Pe
uL
0
Pe表示对流与扩散作用 的相对大小。
0
3/59
传热与流体流动的数值计算
三、对流项的迎风格式（续）
e界面
ue 0 , P ; ue 0 , E
u e Fee P max Fe ,0 E max Fe ,0
w界面
P Fe ,0 E Fe ,0
uw 0 , W ; uw 0 , P
aP aE aW Fe Fw Fe Fw , aE De , aW Dw 2 2
在数值计算过程中，如果连续性方程始终得到满足， 则： aP aE aW 在求解过程中，始终保持连续性方程满足非常重要。 常物性条件下均分网格： 1 0.5P E 1 0.5P W P
aE De Fe ,0 , aW Dw Fw ,0
对流项一阶迎风：
aW i 1 aE i 1 P ,0 1 P ,0 P D D
11/59
传热与流体流动的数值计算
2
5/59
传热与流体流动的数值计算
二、对流项的中心差分（续）
例：在一维模型方程离散求解的 均分网格中，已知W =100， E =200。试对P =0，1，2及4 四种情况按中心差分格式计算 P之值。
负系数会导致物理上不真实的解。
6/59
传热与流体流动的数值计算
三、对流项的迎风格式
Taylor展开法
§5.2 对流项的中心差分与迎风格式
一、一维对流-扩散问题模型方程的精确解
d d d u dx dx dx
边界条件： x 0 , 0 ; x L , L
d
udx
C
C2
ln C1
ux
C1eux C2
二、混合格式(Spalding,1971)
0 aE 1 0.5Pe De Pe , Pe 2 , , 2 Pe 2 Pe 2
aE Pe , 1 0.5Pe , 0 De10/59来自传热与流体流动的数值计算
§5.3 对流-扩散方程的混合格式及乘方格式
一、系数aE与aW 之间的内在联系
aE(i)与aW (i+1)共享同一个界面。 对流项中心差分：
Fe Fw aE De , aW Dw 2 2 aW i 1 aE i P P 1 1 P D D 2 2
u w Fww W max Fw ,0 P max Fw ,0
W Fw ,0 P Fw ,0
8/59
传热与流体流动的数值计算
三、对流项的迎风格式（续）
迎风格式离散形式：
aPP aEE aWW
二、对流项的中心差分
d d d u 采用控制容积积分法 对方程 dx dx dx e u e w u w P 2 2 x w x e
e w u e u w E W 2 2 x e x w
d i i 1 , ui 0 dx i x
i 1 i , ui 0 x
控制容积积分法 e界面 ue 0 , P ; ue 0 , E w界面 uw 0 , W ; uw 0 , P
7/59
传热与流体流动的数值计算
记：F = u 通过界面的流量。 D= x 界面上单位面积扩散阻力的倒数(扩导)。 F u u x = P D x 4/59
传热与流体流动的数值计算
二、对流项的中心差分（续）
Fe Fw Fe Fw De Dw P De E Dw W 2 2 2 2 aPP aEE aWW
a E De Fe , 0 aW Dw Fw , 0 a P a E aW Fe Fw
9/59
传热与流体流动的数值计算
四、中心差分与一阶迎风格式的讨论
1、对流项中心差分在不发生振荡的参数范围内，比一 阶迎风格式的误差更小。 2、一阶迎风格式离散方程系数永远大于零，不会引起 解的振荡，得到物理上看似合理的解。 3、一阶迎风格式截差阶数低，除非采用相当密的网格， 否则计算结果的误差较大。 4、一阶迎风格式的启示：应当在迎风方向取更多的信 息构造格式，更好地反映对流过程的物理本质。 5、在调试程序或计算的中间过程仍可以采用一阶迎风 格式。