求解一维对流扩散方程的一种新方法

第43卷第1期2010年2月武汉大学学报(工学版)

EngineeringJournalofWuhanUniversityVol.43No.1Feb.2010

文章编号:1671-8844(2010)01-0010-04

求解一维对流扩散方程的一种新方法

陈翠霞,张小峰

(武汉大学水资源与水电工程科学国家重点实验室,湖北武汉 430072)

摘要:针对常系数对流扩散方程,基于微分算子分裂算法思想,分别对对流步与扩散步运用待定系数法,以格式的

数值振荡和数值扩散最小为目标,得出各节点的权重系数,并在格式中引入无因次系数.用对流步进行计算,并将其结果作为已知值运用到扩散步的求解中,构造出一种新的一维对流扩散方程的数值求解格式.数值试验表明,相比其他已有格式,该格式可有效控制格式的数值振荡和数值扩散问题,易于编程,精度高,数值结果令人满意,较好地实现物质输移扩散的真实物理过程.

关键词:对流扩散方程;微分算子分裂算法;待定系数法;数值格式中图分类号:TV131 文献标志码:A

Anewsolutiontoone-dimensionalconvection-diffusionequation

CHENCuixia,ZHANGXiaofeng

(StateKeyLaboratoryofWaterResourcesandHydropowerEngineeringScience,

WuhanUniversity,Wuhan430072,China)

Abstract:Anewmethodisdevelopedwiththedifferentialoperatorsplittingalgorithmsandtheundeter-minedcoefficientmethodtosolveone-dimensionalconvection-diffusionequation,whichhasattractedatten-tionbecauseofitsimportance.Withthetargetofnumericaloscillationoftheleastnumericaloscillationandnumericalproliferation,themethodgetsweightindexofeachnodeandintroducesadimensionlessco-efficient.Calculatingbyconvectionandapplyingitsresultsasknownvaluetothesolutionofstep-prolifer-ation.Theresultsofnumericalexperimentsshowthatcomparedwithotherexistingschemes,theschemecancontrolthenumericaloscillationandproliferationmoreeffectivelyandhashigheraccuracy.Theschemeproposedcansimulatetherealphysicalprocessofmaterialtransportanddiffusion.Keywords:convection-diffusionequation;differentialoperatorsplittingalgorithms;undeterminedcoeff-icientmethod;numericalscheme

对流扩散方程一般用来描述温度场、浓度场或大气、海洋、江河及地下水中的污染物运输转移扩散过程,因其重要性而备受关注.解决此类问题,常借助于数值求解.目前常用的差分格式有交替方向特征有限元法、Crank-Nicolson格式

[1]

[2]

计算中取得了较好的效果,并通过该方法,对一

些早期经典格式进行了改进,构造出HUAC2格式,推广应用到一维对流扩散方程的数值模拟中,取得了较好的结果.

本文在以往格式计算的基础上,基于微分算子分裂算法,运用待定系数法,构造一种新的一维对流扩散方程数值求解格式,并用算例验证该格式的优越性.

[10]

和Lax-Wen-

droff格式[3].此外钟万勰等[4-5]应用子域精细积分方法对一些偏微分方程进行了数值求解;忻孝康等采用微分算子理论并运用欧拉-拉格朗日分裂格式对对流扩散型方程进行了求解;张小峰等

[7-9]

[6]

提

1 一维对流扩散方程求解

一维对流扩散基本方程为

出了一种构造高阶精度的待定系数法,在一维对流

收稿日期:2009-09-09

作者简介:陈翠霞(1987-),女,硕士研究生,主要从事河流动力学研究,E-mail:whu-wzmqccx@http://www.77cn.com.cn.基金项目:国家杰出青年科学基金项目(编号:50725930).

第1期陈翠霞,等:求解一维对流扩散方程的一种新方法

11

+k=E22

其中,E>0,是一个常数.

采用微分算子分裂法,令

L1uS+k=0

2L2

2uS-E=0

2(1)

a1-a34k3$t3$x/3+(a1+a3-4

a4-a6)$x4/24]+,

4

(7)

式(7)为式(6)模拟对流方程式(2)的等价微分

(2)(3)

方程.为使该格式的数值振荡和数值扩散同时达到最小,得到下面6个相互独立的代数方程:

a1+a2+a3=1-a1+a3+a4-a6=Cr

-a1-a2-a3+a4+a5+a6=0

Cr2+2a1-a3Cr+a1+a3-a4-a6=0Cr3-3a3-a1Cr2+3a1+a3Cr+ a1-a3-a4+a6=0

Cr4+4a1-a3Cr3+6a1+a3Cr2+ 4a1-a3Cr+a1+a3-a4-a6=m

(8)

式中:m为无因次参数.这是根据数值格式稳定性分析的Fourier谱分析理论得知的.当等价微分方程右边第一项为奇数阶导数项时,格式以数值色散特性为主;当为偶数阶导数项时,格式具有数值扩散特性

(5)

[11]

其中:L1u为求解对流步;L2u为求解扩散步.

通过对方程(2)和(3)的求解来实现方程(1)的求解.在对这两个简单方程求解时,时间步长涉及到(t,t

n

n+1

,t

n+2

)3个时间层,即利用方程(2)在(t,

n

tn+1)内求得un+1,然后再利用方程(3)在(tn+1,tn+2)内求得u

n+2

,最终求得结果.

1.1 对流步的待定系数法求解

利用Crank-Nicolson格式对方程(2)离散得

+1+1n+1nninin+1i-1ini+1i-1

+k+=0(4)2$x2$x

式中:$t为时间步长;$x为空间步长.

对式(4)进行变化得

+1n+1+1

-Cruni-1/4+ui+Cruni+1/4=

.上述推导过程,式(7)中方程的右边引入一个

Crui-1+ui/4-Crui+1/4

nnn

无因次参数m,用以增强格式的稳定性.

方程组(8)的唯一解为a1=Cr-1

12

Cr-+

4Cr(1+Cr)

式中:Cr=2k$t/$x.此格式对对流项的离散采用了代数平均,没有考虑各节点在格式中所起的不同作用.为了提高格式的精度,运用待定系数法对各节点赋予不同的权重系数,同时以格式的数值振荡和数值扩散为最小作为改进目标,最后定出这些待定权重系数的恰当值.

用待定系数ai(i=1,2,,,6)代替式(5)中的各系数:

+1n+1+1nnna1uni-1+a2ui+a3uni+1=a4ui-1+a5ui+a6ui+1

2a2=1-Cr+2-62Cr-(Cr+1)

a3=

Cr+12

Cr++

4Cr(Cr-1)

a4=Cr+12Cr+-4Cr(Cr-1)

(9)

(6)

对式(6)在(xi,tn)处作泰勒展开,得(a1+a2+a3)$t$x=+-a1+a3+a4-a6#

a5=1-Cr2+2+

62Cr-(Cr+1)a6=Cr-112

Cr-2-4Cr(1+Cr)

1.2 扩散步的待定系数法求解

同对流方程式的离散,用待定系数bii=1,2,,,6代替式(5)中的各系数得

+2n+2+2n+1n+1+1b1uni-1+b2ui+b3uni+1=b4ui-1+b5ui+b6uni+1

-a1-a2-a3+a4+a5+a6u-$x/2]+

5x

2

3

3

2

[2a1+a2+a3k2$t2+2a1-a3k$t$x+a1+a3-a4-a6

(10)

2

对式(10)在xi,tn+1处作泰勒展开,得b1+b2+b3$t+

=

-b1+b3+b4-b6$x2

[4a1+a2+a3k$t/3+2a1-a3k$t$x+

2

2

3

a1+a3k$t$x+

3

a1-a3-a4+a6

4

$x/6]-5x

4

-b1-b2-b3+b4+b5+b6u-

[2a1+a2+a3k$t/3+a1-a3k$t$x3/3+

a1+a3k2$t2$x2+

b1+b3-b4-b6$x2-25x[-b1+b3+b4-b6($x3/6)+

12

3

$x$tE]-3

武汉大学学报(工学版)第43卷

+2-b1+b3

a5a4

B1=

(11)

b2b1

C=

a6a5w

a6wa4

b3b2w

b3wb1

wb2b1

b3bwa5a5

a6a2

2b1+b2+b3E$t2+b1+b3$x3$tE++,

4

b1+b3-b4-b624

4

式(11)为式(10)模拟扩散方程式(3)的等价微分方程.为使该格式的数值振动和数值扩散同时达到最小,得到下面6个相互独立的代数方程:b1+b2+b3$t=1-b1+b3+b4-b6=0

-b1-b2-b3+b4+b5+b6=0-b1+b3-b4-b6$x2/4=E-b1+b3+b4-b6$x/12+

2

2b1+b2+b3$t2E+

3

b5b4

B2=

b6b5w

b6w

b4

wb5b5

b6bT

T

T

b3-b1$x$tE=0

b1+b3$x2$tE+

4

b1+b3-b4-b6$x/24=0

(12)

Un=D1=D2=

nn

u1n,u2n,un3,u4,,,uN-1

方程组(12)的唯一解为

b1=-+2;b2=12$t6$tb3=b5=

-;b4=+2

12$t$x12$t$x2

(13)

n+1n+1

un0-a1u0,0,0,0,,,unN-a3uN

+1+2+1+2un0-b1un0,0,0,0,,,unN-b3unN

2 数值试验

算例1 矩形波问题.其初边值为Ux,0=

1,当0.05<x<0.25时

0,其他ux,0=

0,0<x<1

u0,t=0,u1,t=0,0<t<T其中k、E是常数.这个问题的理论解是Ux,t=

erf+erf222(14)

x

2

其中,erfx是误差函数,erfx=e-tdt.

n+1

-2;b6=+6$

t12$t$x$x

1.3 对流扩散方程的求解

在tn,tn+1时间步长内,利用式(9)对式(6)进行求解;在t

n+1

,t

n+2

时间步长内,利用t,t

n时

间步长内的计算结果对式(10)进行求解,便构成了对流扩散方程(1)的完整的分步离散求解过程.

结合式(6)、(10)可得如下矩阵形式:

AUn+1=B1Un+D1CU

a2a3a1

其中 A=

a2w

n+2

=B2Ua3wa1

wa2a1

n+1

+D2

取$x=0.01,$t=0.001,k=0.5,E=0.02,t=1.0.当无因次参数m取不同值时,图1给出了用待

a3a定系数法计算这一对流扩散运动过程的结果及相应的理论解.图2分别给出了Crank-Nicolson格式与迎风格式在相同初始和边界条件下的计算结果.

图1 算例1计算结果与理论解的比较

第1期陈翠霞,等:求解一维对流扩散方程的一种新方法

13

由图1中可以看出,当m=0时,格式的等价微分方程右边的第一项为5阶导数项,无足够的数值粘性,计算结果有明显的数值振动;随着m的增大,数值粘性的作用随之增大,而数值振动的作用随之减小;使数值振动和数值扩散同时达到比较理想的m值约在0.02-0.06之间,以小于0.06时为宜,在实际计算中选用m=0.02.对比图2与图1(b)得,迎风差分格式在模拟锯齿波的对流扩散时,具有很强的数值扩散效应,它使初始锯齿波顶部变尖变滑,底部向两端耗散严重,整体波形趋于矮胖,计算结果与理论解相差较大.Crank-Nicolson格式模拟时有明显的数值振荡作用,与理论解有较大的差别.相比之下,本文格式引入不为零的因数m,较好实现物质输移扩散的真实物理过程,可使数值振动和数值扩散的控制同时达到理想效果

.

图3 算例2计算结果与理论解的比较

2)运用分步离散的方法,对方程的对流步与扩散步运用待定系数法,并用对流步的计算结果对扩散步进行求解,构成对流扩散方程的完整的分步离散求解过程.

3)数值试验表明,本文格式在模拟对流扩散方程的过程中,可有效地控制数值格式的振动与扩散问题,较真实地描述波的传输和扩散过程.在实际的计算中,节点简单,边界条件易于确定,精度高,便于实际运用.参考文献:

[1] 谢树森.对流扩散问题的交替方向特征有限元方法

[J].高等学校计算数学学报,1996,(3):282-291.[2] 王同科.一维对流扩散方程Crank-Nicolson特征差分

格式[J].应用数学,2001,14(4):55-60.

图2 其他格式计算结果与理论解的比较

[3] 陆金甫,关冶.偏微分方程数值解法[M].北京:清华大

学出版社,1987:165-171.

[4] 钟万勰.单点子域积分与差分[J].力学学报,1996,28

(2):159-163.

[5] 钟万勰.子域精细积分及偏微分方程数值解[J].计算

结构力学及其应用,1995,12(3):253-260.

算例2 如下对流扩散方程的初边值问题:ux,,0<x<1

u0,t=0,u1,t=1,0<t<T

其中k、E是常数.问题的理论解是ux,t=

+e-1ek

n=1

E

]

n+k/(2E)

2

n2E+n

#

[6] 忻孝康,黄光伟.对流扩散型方程的一种简单、有效的

欧拉-拉格朗日分裂格式[J].空气动力学学报,1986,4(1):65-72.

[7] ZhangXiaofeng.Anewhighorderschemeforconvec-tionequationanditsapplication[J].JournalofHy-drodynamics,2002,14(1):81-86.

[8] 张小峰,张红武.Crank-Nicolson格式精度的改进[J].

水科学进展,2001,12(1):33-38.

[9] 张小峰,陆俊卿,易灵.求解一维对流扩散方程的一种

高精度数值格式[J].武汉大学学报(工学版),2005,38(2):10-14.

[10]李炜.黏性流体的混合有限分析解法[M].北京:北京

科学出版社,2001.

[11]刘导治.计算流体力学基础[M].北京:北京航空航天

大学出版社,1989:60-72.

sinnPxe-

t

取$t=0.002,$x=0.01,k=1.0,E=0.5,t=1.0.同算例1,取m=0.02.时间步长及空间步长等参数的取值需满足稳定性条件0<Cr[1.0,0<$tk/($x)[0.5.图3绘出了本文格式模拟一

维线性对流扩散情况,计算结果表明,本文格式的计算结果和理论解吻合良好,可有效真实地描述波的传输与扩散过程.

2

3 结语

1)针对常系数对流扩散方程,基于微分算子分裂算法思想,运用待定系数法,构造出一种新的一维对流扩散方程的数值求解格式.

本文来源：https://www.bwwdw.com/article/zpti.html