基于Peclet数判别法的一维对流扩散方程分类研究

基于Peclet 数判别法的一维对流扩散方程分类研究摘要：采用Peclet 数的绝对值大小来判别一维对流扩散方程为对流占优型或是扩散占优型方程，运用三种隐式差分格式—中心隐式格式、对流C-N 型格式和扩散C-N 格式，对不同Peclet 数的算例进行离散和求解。

然后，将计算区域中所有节点的解析解与数值解表示成矩阵形式，并求解出它们的矩阵2范数之后作比较，两者越接近则代表差分格式精度越高。

通过比较得出了当方程Peclet 数的绝对值小于等于0.5时，方程为扩散占优型方程。

在离散方法选取方面，针对扩散项的离散可以采用更高精度的差分格式，如扩散C-N 格式；当Peclet 数的绝对值大于等于20时，方程为对流占优型方程。

此时，针对对流项可以采用更高精度的差分格式，如对流C-N 格式；当Peclet 数的绝对值介于0.5与20之间时，无法用Peclet 数判断方程类型,不过可以选择折衷的离散格式减小误差，如中心隐式格式。

关键字：一维对流扩散方程 Peclet 数判别法 有限差分方法 数值模拟MR(2010)主题分类号：39A14；65M06 中图分类号：O242.2 文献标识码： A1.引言一维对流扩散方程是描述流体流动和传热问题的一类线性化模型方程。

土壤、大气等多孔介质中水、盐分、温度以及污染物质的对流扩散问题都会遇到此类方程。

在一维对流扩散方程的求解过程中，反映流体对流和扩散两种物理作用的分别是对流项和扩散项。

所以，根据方程中对流项还是扩散项占主导作用，通常可将方程分为对流占优型和扩散占优型两类方程。

然而，要想得到精确度较高的数值结果，这两种类型方程的离散方法不能采用相同的离散格式。

因此，需要有一种判别方法来判断方程的类型，关于对流占优型和扩散占优型方程的判别方法一直是近年来研究的热点问题。

这对研究不同类型的方程使用合适的差分格式进行离散具有实际的意义。

由于Peclet 数的绝对值表示了对流作用相对扩散作用的大小，即绝 大，扩散所起的作用就可以忽略。

反之，当Peclet 数为零时，方程就为纯扩散方程。

本文选用一维定解非稳态对流扩散方程为例，通过考察Peclet 数的绝对值大小来对方程进行分类，方程一般形式如下：2(,),,0122(1)(,0)()(,)(),(,)()12(,)u u u a f x t x x x t t x xu x g x u x t t u x t t u u x t υϕφ⎧⎪⎪⎪⎨⎪⎪⎪⎩∂∂∂+=+≤≤≥∂∂∂====其中a 和υ分别代表对流项系数和扩散项系数。

假定求解区间长度为s ， Peclet 数的绝对值计算公式为：.(2)asPe υ=从公式（2）中可以看出，当计算区间长度给定，Peclet 数是由对流和扩散系数确定的[1]。

下面介绍方程（1）的离散方法。

2. 离散方法2.1 显式格式离散对于上述方程（1），需要离散非稳态项（简称U 项）、对流项（简称C 项）和扩散项（简称D 项）。

常见的离散方法有显式格式和隐式格式两种。

显式格式有：中心显式格式、修正中心显式格式、迎风差分格式等。

比如，以中心显式格式为例，即使用向前差分格式、一阶中心差分格式与二阶中心差分格式组合分别离散U 项、C 项和D 项。

其离散形式如下：1111122(3)2n n n nn n nj jj j j j j u u u u u u u ahhυτ++-+----++=其截断误差为2()o h τ+。

然而，由von Neumann 判别条件判断此种格式将受到稳定性条件的限制，即：2221,2a h υττυ≤≤。

相应其它显式格式同样有稳定性条件限制。

所以，显示格式时间步长τ和空间步长h 取值将受到限制。

因此，若采用显式格式求解一维非稳态对流扩散方程问题，得到的数值解精度将受到限制，甚至误差很大。

所以，显式格式的离散效果欠佳，为了弥补它的缺陷，尝试采用无条件稳定的隐式格式离散（1）式。

2.2隐式格式离散常见的隐式格式有三种：向后差分格式、一阶中心差分与二阶中心差分组合（简称中心隐式格式）；向后差分格式、C-N 型一阶中心差分与二阶中心差分组合（简称对流C-N 型格式）；向后差分格式、一阶中心差分与C-N 型二阶中心差分组合（简称扩散C-N 型格式）。

三种隐式格式的离散形式如下： 1）中心隐式格式：111111111122()(4)2n n n n n n n j jj j j j j u u u u u u u a hhυτ+++++++-+----++=2）对流C-N 型格式：11111111111122()(5)44n nn n n nn n n j jj j j j j j j u u u u u u u u u a hhhυτ+++++++-+-+-----+++=3）扩散C-N 型格式：1111111111112222()(6)222n n n n n n n n n nj jj j j j j j j j u u u u u u u u u u ahhhυτ+++++++-+-+----+-++=+由von Neumann 条件判断上述三种隐式格式均为无条件稳定的格式，即在网格系统较为粗糙时，也不会产生数值震荡现象。

下面将给出上述三种差分格式的稳定性分析[2]。

3.稳定性分析将（4）~（6）式分别按网格节点排列如下：111'11(2)(24)(2)2(4)n n n nj j j ja h u u a h u u λυλυλλυλ++++--++-+=111'1111(4)(48)(4)4(5)n n n n n nj j j j j j a h u u a h u a hu u a hu λυλυλλυλλλ++++-+--++-+=-++111'1111()(22)()(22)(6)n n n n n n j j j j j j a h u u a h u u u u λυλυλλυλυλνλυλ++++-+--++-+=+-+其中网格比。

假定2hτλ=(7)ikjh n n u v ej =其中cos()sin()ikh e kh i kh =+⋅。

把（7）式代入'(4)~'(6)式并消去公因子ikjh e ，容易求出上述四式的增长因子分别为：222(,)(8)(244cos())(2)G k kh a h τυλυλλ=+-+2224(2(1))(,)(9)(244cos())(2)a h G k kh a h λθτυλυλλθ+-=+-+2224(1)4(1)cos()(,)(10)(244cos())(2)kh G k kh a h υλθυλθτυλθυλθλ--+-=+-+可以看出，(8)(10) 式的值均小于等于1。

因此，满足von Neumann 判别条件。

所以三种隐式格式均无条件稳定[34]-。

4.数值算例为了通过Peclet 数判别法讨论一维对流扩散方程的分类，运用上述（4）~（6）式的三种离散格式进行了大量实例计算。

本文列举其中部分数值算例如下。

为了讨论的必要，所有算例的计算区间长度s 均取1m ，模拟时间取1s ；时间步长取0.1s ，每个算例的空间步长分别取0.2m ，0.1m ，0.05m 进行比较计算。

按照方程（1），算例的条件依次如下： 例 122(100)800(10010)80010,0.1,01,0(,0),(0,)(1,)0(11)20(,)(),(,)04000.1x x t a x t u x eu t u t u x t e f x t t υ⎡⎤--⎢⎥⎢⎥⎣⎦⎡⎤---⎢⎥⎢⎥⎣⎦⎧==≤≤≥⎪⎪⎪===⎨⎪⎪=⋅=⎪+⎩例 221,0.05,01,0(,0)(1),(0,)(1,)0(12)(,)(1),(,)(312)t t a x t u x x x u t u t u x t x x e f x t x x e υυ--⎧==≤≤≥⎪=-==⎨⎪=-=-++⎩例 390.0990.0990.090.1,0.01,01,0(,0),(0,),(1,)(13)(,),(,)0x tt x t a x t u x e u t eu t e u x t e f x t υ---==≤≤≥⎧⎪===⎨⎪==⎩例 4255(0.010.25)0.1,0.01,01,0(,0)sin(),(0,)(1,)0(14)(,)sin(),(,)0xx t a x t u x e x u t u t u x t ex f x t πυππ⎡⎤-+⎣⎦⎧==≤≤≥⎪⎪===⎨⎪==⎪⎩例 51.7712434446770460.091.7712434446770460.091.7712434446770460.090.1,0.02,01,0(,0),(0,)(15)(1,)(,),(,)0xt t x t a x t u x eu t e u t eu x t e f x t υ---==≤≤≥⎧⎪==⎪⎨=⎪⎪==⎩例 61111111111,1,01,0(,0)1(1),(0,)(1,)0(16)(,)(1)(,)12(1)xx txta x t u x e e x e u t u t u x t e e x e e f x t e e x e e υυυυυυυυυυ--------------==≤≤≥⎧⎪⎪=-+--==⎪⎪⎡⎤⎨=+--⎢⎥⎪⎣⎦⎪⎡⎤⎪=-+--⎢⎥⎪⎣⎦⎩例 71,1,01,0(,0)sin(),(0,)sin()(17)(1,)sin(1)(,)sin(),(,)0ttt a x t u x x u t e t u t e t u x t e x t f x t υ---==≤≤≥⎧⎪==-⎪⎨=-⎪⎪=-=⎩例 82221,1,01,0(,0),(0,)(1,)0(18)(,)(),(,)(1)t t a x t u x x x u t u t u x t x x e f x t x x e υ--=-=≤≤≥⎧⎪=-==⎨⎪=-=-++⎩例 91()1,1,01,0(,0),(0,),(1,)(19)(,),(,)x t t x t x t a x t u x e u t e u t e u x t e f x t e υ-----==≤≤≥⎧⎪===⎨⎪==⎩例 1020.250.25(0.01250.2)0.1,0.2,01,0(,0)sin(),(0,)(1,)0(20)(,)sin(),(,)0xx t a x t u x x e u t u t u x t ex f x t πυππ⎡⎤-+⎣⎦⎧==≤≤≥⎪⎪===⎨⎪==⎪⎩例 1120.220.22(0.02420.5)0.22,0.5,01,0(,0)sin(),(0,)(1,)0(21)(,)sin(),(,)0xx t a x t u x x e u t u t u x t ex f x t πυππ⎡⎤-+⎣⎦⎧==≤≤≥⎪⎪===⎨⎪==⎪⎩下表为以上11个算例[510]-在使用三种格式离散和取不同空间步长的情况下，得到的解析解与数值解矩阵2范数之差的绝对值。