对流占优扩散方程的差分法.doc

《对流占优扩散方程的差分法.doc》由会员分享，可在线阅读，更多相关《对流占优扩散方程的差分法.doc（38页珍藏版）》请在沃文网上搜索。

1、摘要对流占优扩散方程主要包含对流项和扩散项，其中对流项系数远远大于扩散项系数。在数值计算中，方程中的扩散项一般采用具有优良物理特性和计算精度的中心差分离散格式，而关于对流项的处理就稍显困难，若处理不当便会产生数值震荡或数值弥散，给数值计算带来困难。因此，需要对求解的方法做出改进。本文主要讨论迎风差分格式，迎风加权差分格式，以及特征有限差分格式。三种方法都能够消除数值震荡，但各种方法间又各有差异。迎风格式计算量较小，能够消除数值震荡，但是数值解的精度不高。特征有限差分格式中含有多个未知的点，计算量特别大，从误差分析中可以看出，其数值解拥有较高的精度。迎风加权差分格式，是在迎风格式的基础上改进得到

2、的，精度较高，其数值解不仅受到时间和空间步长的影响，还受到不同参数的影响。可以选取不同的参数是迎风加权格式的一个优点。关键词：对流占优扩散方程；迎风格式；迎风加权差分格式；特征有限差分法AbstractConvection-dominated diffusion problems mainly contain convection and diffusion terms, which the convection coefficient is much larger than the diffusion coefficient. In the numerical calculation, di

3、ffusion terms in the equation commonly used central difference discretization scheme with excellent physical properties and calculation accuracy. However, the method of the convective terms slightly difficult. It would produce numerical shock or numerical dispersion if not handled properly. Therefor

4、e, we need to make some improvements.This article focuses on upwind difference scheme, upstream weighted scheme, as well as characteristic finite difference method. The numerical oscillation can be eliminated by all three methods, but there are differences between each method. Upwind difference sche

5、me has smaller amount of calculation, to eliminate the numerical oscillation, but the accuracy of numerical solution is not as good as we expect. Characteristic finite difference method which contains a number of unknown point, with a large amount of calculation, and we can see from the error analys

6、is, the accuracy of numerical solution is much higher. Upstream weighted scheme, which improved based on upwind scheme, is not only influenced by the time and space step, but also affected by different parameter of . To choose a different parameter of is also an advantage of upstream weighted scheme

7、.Key Words: Convection-dominated diffusion problem; Upwind difference scheme; Upstream weighted scheme; Characteristic finite difference method目录1、绪论11.1设计（论文）的背景及目的11.2 国内外研究现状11.3 论文主要研究内容21.4 研究思路和方法32、论文的预备知识42.1 差分法简介42.2 方法52.3 差分格式的稳定性定理63、含对流项的一维抛物型方程73.1 中心差分格式的推导73.2稳定性分析93.3中心差分格式的缺陷104、迎

8、风格式114.1 对流占优扩散方程的迎风差分格式114.2迎风差分格式的稳定性分析135、迎风加权差分格式155.1加权差分格式的建立155.2稳定性分析166、特征有限差分法176.1特征差分格式的建立176.2双线性插值197、数值算例20结论27谢辞28参考文献29附录30对流占优扩散方程的差分法1、绪论1.1设计（论文）的背景及目的对流占优扩散方程是一类基本的运动方程，它可用于环境科学、能源开发、流体力学和电子科学等许多领域，对该方程数值计算方法的研究具有重要的理论和实际意义。对流扩散方程右端第一项为扩散项，左端第二项则是对流项。由于其方程本身的特点，给建立准确有效的数值求解方法带来一

9、定的困难。对流和扩散给流体中由流体携带的某种物理量的变化过程，可以通过一个无量纲的特征参数(数)来描述。如果数较小，即对流效应相对较弱，这类问题中，扩散占主导地位，方程是椭圆型或抛物线型；如果数较大，即溶质分子的扩散相对于流体速度而言是缓慢的，这类问题中，对流占优，方程具有双曲型方程的特点。该方程表征了流动系统的质量传递规律，求解此方程可得出浓度分布。此方程系通过对系统中某空间微元体进行物料衡算而得。对于双组分系统，A组分流入某微元体的量，加上在此微元体内因化学反应生成的量，减去其流出量,即为此微元体中组分A的积累量。考虑到组分A进入和离开微元体均由扩散和对流两种作用造成,而扩散通量是用斐克定

10、律表述的。对流占优扩散方程具有一个共性，即对流占优性，由于对流项的存在给数值求解带来许多困难。因此，寻找一种有效数值解法一直是计算数学中重要研究内容。本文考虑一维对流扩散方程。用通常的差分法进行求解可能会出现数值震荡，为克服数值震荡，需要对传统的差分方程改进，如迎风格式，迎风加权差分格式，特征有限差分方法等。本文将讨论上述三种解决对流占优问题的差分算法，说明每种算法的稳定性条件，收敛性等，并结合数值算例说明。1.2 国内外研究现状年代，和等提出特征修正技术求解对流占优扩散问题，与其他方法相结合，提出了特征有限差分方法、特征有限元方法、特征混合有限元方法等，并给出理论分析；提出过一种沿线方向附加

11、人工黏性的间断有限元法，称为流线扩散方程（）。有限差分法、有限元法、有限体积法是工程应用中的主要方法。国内也有许多这一方面的文献，秦新强在对流占优扩散发方程的一种特征差分算法中，提出了解对流扩散方程的特征线法, 这一方法考虑沿着特征线(流动方向) 的离散, 利用了对流扩散问题的物理力学特征, 可以有效地克服数值震荡, 保证数值解的稳定；梁栋在对流扩散方程的一类迎风格式中，对其方程进行分析, 得到了稳定性和收敛性定理。1.3 论文主要研究内容对流占优问题的求解，若采用常规的方法，很容易出现数值震荡。为了避免求解结果产生数值振荡，获得稳定解，本文讨论一些改进的差分方法，来求解常系数的对流占优对流扩

12、散问题，使其得到稳定的数值解。 （1）迎风格式基于广义差分法，数学家们提出了数值求解的一类迎风格式。从简单的一维常系数双曲方程着手，对我们构造差分格式是很有启发的。构造对流占优扩撒问题的差分逼近，为了导出稳定性条件，通常用局部固定系数法或视变系数为常系数法。最后按照气体理学的含义，系数表示气流速度，因此，人们称其为迎风格式。再对方程进行分析，得到了稳定性和收敛性定理，并对一类模型问题进行试算，结果良好。 （2）迎风加权差分格式这也是一种解非定常对流占优扩散方程的有效方法，它实际上是对迎风格式的进一步改进，即对一般的空间中心格式和迎风格式进行加权处理。此格式适合对流项占优时求解，它是一个显示格式

13、，计算量比较适中。另外，还有一个优点是可以通过选取参数而获得差分格式的适应性。经精度分析和数值验证，可得其稳定性良好，便于求解。 （3）特征有限差分方法解对流扩散方程的特征线法，即考虑沿着特征线(流动方向) 的离散, 利用了对流扩散问题的物理力学特征, 可以有效地克服数值震荡, 保证数值解的稳定, 众多的理论和应用成果均采用了这一方法。另外，如果要进一步消除因分割区间的步长不同，而引起的数值震荡现象，则可对网格运用双线性插值的方法，构造出一种新的特征差分算法，同时由于算法构造的独特性，该算法还特别适用于求解变系数的对流占优扩散方程。1.4 研究思路和方法关于对流占优扩散问题的研究，目前国内外都

14、没有的绝对良好的方法，只能尽可能的采用适合的逼近方法，来获取较稳定数值解。本文讨论前面提到的三种改进的差分格式，主要对三种方法进行阅读，理解其精髓，学会用方法判断其稳定性，后期要编程，并用实现数值算例的求解。这就要求在看懂方法的基础上，用数值算例来做检验，看其是否符合先前所做的稳定性分析和截断误差分析，从而验证该方法的正确性。并且对三种方法做出比较，能够在不同的情形下使用相应的差分方法来解决问题。2、论文的预备知识 在进行论文写作之前，我们先要对毕业设计中可能用到的各种知识或者定理，进行说明和引述，以便在之后的论文中直接运用，不再赘述。本文中主要涉及到的知识有差分法，方法，以及差分格式相关的稳

15、定性定理。而展开式等基本的方法，将不在本节详细介绍。2.1 差分法简介差分法是解微分方程的主要数值方法。由于数字电子计算机只能存储有限个数据和作有限次运算，所以任何一种计算机解题的方法，都必须把连续问题（微分方程的边值问题、初值问题等）离散化，最终化成有限形式的线性代数方程组。用差分法将连续问题离散化的步骤是，首先对求解区域做网格剖分，用有限个网格节点代替连续区域；其次将微分算子离散化，从定解问题的微分或积分形式出发，用数值微商或数值积分公式导出相应的代数方程，从而把微分方程的定解问题化为线性代数方程组的求解问题。（1）差分法的基本问题a) 对求解域做网格剖分。一维情形是把区间分成一些等距或不

16、等距的小区间，称之为单元。二维情形则把区域分割成一些均匀或不均匀的矩形，其边与坐标轴平行。也可分割成一些三角形或凸四边形等。b) 构造逼近微分方程定解问题的差分格式。有直接差分化法、有限体积法或广义差分法、以及变分差分法。c) 差分解的存在唯一性、收敛性及稳定性的研究。这些理论问题都归结到对差分解做出先验估计。d) 差分方程的解法。由于代数方程组的某些特点，容易导致数值震荡、病态，所以求解时应采取某些特殊技巧。（2）偏微分方程差分法的初值问题许多物理现象随着时间而发生变化、如热传导过程、气体扩散过程和波的传播过程都与时间有关。描述这些过程的偏微分方程具有这样的性质：若初始时刻的解已经给定，则时

17、刻的解完全取决于初始条件和某些边界条件，利用差分法解这类问题，就是从初始值出发，通过差分格式沿时间增加的方向，逐步求出微分方程的近似解。（3）偏微分方程差分法的边值问题物理上的定常问题，如弹性力学中的平衡问题、电磁场及引力场等。其定解问题为各种边值问题，即要求解在某个区域D内满足微分方程，在边界上满足给定的边界条件。差分解法可归结为选取合理的差分网格，建立差分格式，求解代数方程组以及考察差分格式的收敛性等问题。2.2 方法方法是分析差分方法稳定性的一种方法，又被称作方法，它是在第二次世界大战期间由首先提出的。该方法是目前分析线性常系数差分方程稳定性的应用最为广泛的方法。我们用方法（包括积分和级

18、数）将差分方程中的空间变量和时间变量分离，从而将差分方程的稳定性归结为有限阶矩阵族的一致有界性。考虑一维的线性常系数抛物型方程方程的初值和周期边值条件为，其中。则逼近它的二层差分方程的一般形式可以写为， 这是空间网点处的差分方程，是包含及其邻近的正负整数的有限集合，不依赖但和有关。由于是周期边值条件，故可将周期开拓，使其对一切有意义，且方程对所有整数成立。为了应用方法，我们再将。为此，取半整数点，并用如下阶梯函数逼近初试函数：，当，再将看成在任意成立，则得具连续变量的差分解。显然仍是的周期函数。这样我们就可将方法用于具连续空间变量的差分方程： 将展成级数：， 把，再比较对应项系数，得， 其中

19、我们称为增长因子，它是判断差分格式稳定与否的重要依据，也是进行稳定性分析的重要手段。2.3 差分格式的稳定性定理定理2.3.1 差分方程按谱范数稳定的充分必要条件是，对于任意的，均成立不等式， 即一致有界。其中为一个独立的正常数。推论 差分方程按谱范数稳定的条件是，对于任意的，均成立不等式。 不等式通常也被称作条件。值得注意的是，只有当建立的差分格式为一维问题时，条件才是差分格式稳定的充要条件。然而从更普遍的情况来讲，它只是差分格式稳定的必要条件。本文中主要讨论一维的对流占优扩散方程，可以将式作为充要条件。3、含对流项的一维抛物型方程当方程中同时存在对流项与扩散项的时候，方程便被称作对流扩散方

20、程。在讨论对流占优扩散方程之前，我们暂时不限定方程的系数，仅以抛物型方程的差分法为基础，从简单的模型着手，考虑一维的含对流项的抛物型方程，并建立其差分格式。对流扩散方程描写的状态不随时间改变，称为驻定问题（定常问题）；而随时间改变的，称为非注定问题（非定常问题）。定常问题可以看成是某一非定常问题当的渐近状态，当用渐近方法求解定常问题时，只关心最终状态。相反，非定常问题的瞬时状态有物理意义，需要求解。本文中我们主要讨论非定常问题。另外，若不特别说明，文中出现的方程均为常系数方程。3.1 中心差分格式的推导考虑非定常的常系数一维含对流项的抛物型方程：， 其中是常数。而初始条件和边值条件分别为，设在

21、相应的区域光滑，且在与边值相容，使得上述问题有唯一充分光滑的解。现在考虑上述方程的初边值问题的差分逼近。取空间步长和时间步长，其中都是正整数。用两族平行直线将矩形域分割成矩形网格，定义网格节点为。并且用表示定义在的函数。推导此差分格式的主要思路，就是使用已剖分好的网格内第层上的点来差分第层上的点（图）。首先，对一阶算子进行离散。为此，对充分光滑的解，在节点处沿时间方向离散，由展开式可得， 于是有， 舍去高阶项，用差商代替微商算子，即，同理可得 。 现将二阶的微商算子在节点处沿空间方向离散化，同样由展开式可得， 将代入到，并舍去截断误差，则逼近该对流扩散方程的中心差分格式为。 3.2稳定性分析对

22、提出的差分格式进行稳定性分析，使用方法，首先对差分格式进行改写：，其中，。可以知道其增长因子为， 通过计算可得， 根据条件，一维的对流扩散方程稳定的充要条件是成立。由于中，可以略去，因而只要求即 从而得稳定性条件为 但是，这只是对充分小得到的稳定性条件，而实际计算时步长总是有限的，为得出相应的稳定性条件，应该要求令，考虑自变量为的抛物线方程容易得到该函数在区间内的值非负，不妨在端点取值：接下来我们讨论该函数的最大值。假设，则抛物线下凸，所以在端点取极大，故稳定性要求。于是用端点值可以导出，连同假设条件一起，从而得稳定性条件， 另外，假设，此时抛物线上凸，极值点位于抛物线顶点极值为易知左侧，又，

23、所以可以取最大值。为使，必须，从而。这便与假设矛盾，此时格式不稳定。总而言之，该含对流项的抛物型方程的中心差分格式稳定的充要条件是成立。3.3中心差分格式的缺陷根据我们推到出的稳定性条件，进一步对两个系数做分析，考虑对流项系数和扩散项系数在不同大小的情况下差分格式的稳定性。将， 代入，得到不等式， 可以推出步长需要满足的条件，即。当时，容易满足此条件，用中心差分格式计算不会出现数值震荡。但是，若当系数时，此条件不容易得到满足，即稳定性条件不成立。因此，求解一般的含对流项的抛物型方程，只要满足就可以运用中心差分格式。而求解对流占优扩散方程时，即系数，稳定性条件难以得到满足，如果仍然使用提到的差分

24、格式，则可能出现不应有的数值震荡。从上面的推导我们可以看到，针对于解对流占优扩散方程，式是有缺陷的，往往不能精确求解。所以，才在此基础上提出了许多改进的差分格式。4、迎风格式前面我们已经讨论过，当方程中的系数时，即对流项占优时，若简单的运用式的差分格式，则该对流占优扩散方程的数值解往往会出现数值震荡，给求解带来不便。为了克服这类数值震荡，在常规差分方法的基础上，提出了一种新的差分格式。此方法是沿着特征方向，用下一层向上一层差分的方法代替微分，看起来似乎是迎风而上的差分格式。当对流项系数是变系数时，按照气体力学的含义，此对流项系数可以表示气流速度，所以被称之为迎风差分格式。4.1 对流占优扩散方

25、程的迎风差分格式考虑非定常的一维对流占优扩散方程：， 其中是常数，。初始条件和边值条件分别为，此时，上述方程虽是含对流项的抛物型方程，但是其解却具有双曲方程的性质，我们在构造差分格式时需要考虑解的这一性质。为此，首先来处理方程的左端。我们令与相伴的特征方向为，沿的方向导数，所以可以将改写成：， 取时间步长，沿轴取节点，由出发的特征与直线交于。我们就可以用下式逼近沿特征方向的导数：， 根据进一步得出： 取充分小，使得位于之间，则时位于左侧，是位于右侧。现设，如右图，在节点和之间作反距离步长线性插值，得， 将代入到，并且用二阶中心差商代替等号右端的二阶微商，就得到逼近的迎风差分格式， 同样，当时，

26、位于右侧，作线性插值得到逼近的迎风差分格式为：， 差分格式的截断误差均为。4.2迎风差分格式的稳定性分析我们首先考虑，当的情况，即对式作稳定性分析。运用方法，令，则式可以改写为， 再以代入上式，得对比等式两端的系数，容易得到增长因子， 然后，将增长因子中的指数函数化作三角函数，再根据三角函数的取值范围，令。于是增长因子又可以写作：， 经计算可得 其中。由于讨论的是一维问题，直接根据条件，该迎风差分格式稳定的充要条件是，成立。而对于实际的问题，应满足令，依然将看作是一个关于的二次抛物线方程，其取值范围是，即该函数的二次项系数大于零（不考虑等于零的情况），则抛物线下凸，所以函数的最大值只能在端点取

27、到。故稳定性要求，从而有， 于是得到了式的稳定性条件。再对此稳定性条件作进一步的讨论，将，代入式，可以得到，由此推出。又因为对流项占优，即系数，当取较大的空间步长时，将分子上的略去，稳定性条件就变为。这使得时间步长的值也容易取得，即很容易满足稳定性条件。同样，当，也对式作相同的稳定性分析，得其稳定性条件为：。 因此，使用迎风差分格式求解对流占优扩散方程，更能够得到稳定的数值解，比起常规的中心差分格式有更大的优越性。具体的数值算例将在最后给出，以说明其良好的性能。 5、迎风加权差分格式我们已经知道迎风格式比起常规的中心差分格式有更大的优越性，陆金甫在其论文中提出，中心差分格式解对流占优扩散方程会

28、出现数值震荡，迎风格式虽然可以避免震荡，但是将会引入过大的人为粘性。那么在中心差分格式和迎风格式的基础上稍做改进，可以得到一种新的差分格式，这就是迎风加权差分格式。这里，我们不单单对其进行中心差分或者迎风差分，而是引入参数对前两种格式进行加权，相当于把前两种格式整合在了一起。5.1加权差分格式的建立我们考虑同样的方程，即式，方程的对流项系数与扩散项系数之间的关系不变，初边值条件也不变。假设对流项系数，建立中心差分格式，即，再建立迎风格式，现在，引入参数，对上述两式中的对流项进行加权，于是便得到了迎风加权差分格式：。 迎风加权差分格式的优点是可以调节参数。事实上，关于的选取很难得到可以遵循的准则

29、，但是，当取某些值得时候，可以得到一些特殊的格式，讨论这些特殊的格式是有意义的。例如，分别是中心差分格式和迎风格式。另外，还可以取作修正中心格式，格式，指数型格式，本文中不讨论这些格式。5.2稳定性分析对式作稳定性分析，其步骤与迎风格式的稳定性分析差不多。先将式改写为下式，其中， ， 于是，我们可以得到以为未知数的增长因子形式：， 经计算可得， 对于实际的问题，要使得差分格式稳定，应满足令，将看作是一个关于的二次抛物线方程，其取值范围是，即要满足稳定性，就要求，从而有， 于是得到了式的稳定性条件。再对此稳定性条件作进一步的讨论，将，代入式，可以得到，由此推出。由于，可将分子上的略去，稳定性条件

30、就变为，。可以看出稳定性条件容易得到满足。6、特征有限差分法求解对流占优扩散方程的另一种改进的差分方法，就是本节介绍的特征有限差分法。特征有限差分法将特征线方法和有限差分方法相结合，给出了一种求解对流占优扩散方程的新的隐式特征差分格式。与迎风格式相比，两种方法主要的差异体现在线性插值上面，这一节主要介绍秦新强提出的一种双线性插值的方法。特征有限差分法的优点是适应性强，且差分格式无条件稳定。同时也适用于变系数方程，数值试验的结果表明在消除数值震荡方面更有效。6.1特征差分格式的建立我们依然考虑同样的方程，即式，方程的对流项系数与扩散项系数之间的关系不变，初边值条件也不变。运用特征线法，可以将改写

31、成如下的形式（其推导过程已经在建立迎风差分格式时详细给出，这里不再赘述）：， 其中，表示与相伴的特征方向，。接下来要对微分算子做离散化处理，而线性插值的区别也体现在这里。离散的方向还是沿着特征线的方向，不同的是，我们利用四个点对作线性插值，以此得到点的值，如图所示，为沿特征线的方向。取空间步长，取节点时间步长，沿轴取节点，然后沿着特征线作差商离散，可得， 而表示点的坐标，要将其表示出来。由图看出，于是，；由，得，并求得，代入，再舍去局部截断误差，则可以得到。 再对扩散项用二阶中心差商来逼近，于是，我们就得到了求解初边值问题的特征差分格式。 6.2双线性插值 作插值所使用的方法较之前有所不同，点

32、的函数值由四个点的函数值通过双线性插值得到。令表示用插值数据求得的双线性插值函数。关于插值函数的构成分两步进行，首先分别以为节点作线性插值，求解在的函数值。再用的值作线性插值，于是就求出了点的函数值，即。具体步骤如下， ， 上面分别得到的是的函数值，接下来再次进行线性插值， 将代入到，就可以得到点的函数值， 于是可以写为 这就是对流占优扩散方程的特征差分格式。而关于稳定性，在本节的开头我们已经提到过，此格式无条件稳定，本文不给出证明。另外，式是一个隐式格式，使用两层上的共四个点，通过双线性插值得到点的函数值。虽然绝对稳定，且计算精度较高，但是计算量非常大。7、数值算例 求解下列对流扩散方程，其

33、中。定解条件为其中为正常数，求解时令。已知方程的精确解为。通过系数的大小可以判断出，这是一个对流占优的问题。要对上述方程进行求解，我们需要用到编程，进而得到其数值解。取空间步长，时间步长。首先可以得到其精确解，如图所示。然后考虑不同的差分法，依次使用中心差分格式、迎风格式、特征有限差分法以及迎风加权格式求解。另外，考虑到迎风加权格式中参数对数值解的影响，我们再分别取进行求解。各数值计算的结果和误差都在下列三维图中给出。图（精确解）图（中心差分）可以从图看出，当对流项占优时，使用中心差分格式出现了明显的数值震荡。图（迎风格式）图（特征有限差分法）图（迎风加权差分格式，）图（迎风加权差分格式，）图

34、（迎风加权差分格式，）图分别对应三种改进的差分法的数值解，可以看出其数值解和精确解比较接近。接下来我们给出各方法数值解与精确解之间的误差图（迎风格式的误差）图（特征有限差分法的误差）图（迎风加权格式的误差，）图（迎风加权格式的误差，）图（迎风加权格式的误差，）结论从给出的数值算例可以看出,当对流项明显占优的时候，采用本文提到的三种改进差分法可以得到较好的近似解，同时在消除数值震荡的方面是有效的。迎风格式采用的是显式方程，计算量较小，能够消除数值震荡，但是数值解的精度不高。特征有限差分格式中含有多个未知的点（即第层上的点），因此属于隐式方程，计算量特别大，从误差图中看出，其得到的数值解的精度较高

35、。而关于迎风加权差分格式，我们从数值算例中看出，其数值解非常接近于精确解；另外，还需要考虑参数的选取对数值解的影响，在本算例中，当时，数值解的精度是最高的。谢辞经过两个月的认真劳作，我的本科毕业论文终于完成了。掩卷沉思，四年的大学读书生涯已倏忽而过，恍惚间四年来的欢乐、痛苦、成长、收获，种种生活情景又一一浮现在眼前。毕业论文能够顺利完成，首先要感谢我的毕业设计指导老师杨艳。论文的字里行间都有杨老师的悉心指导，她为我们答疑解惑，还给我们传授方法。授人以鱼，不如授人以渔。如果没有杨老师的指导，我不能完成的编程。我要感谢大学期间传授我知识的各位老师，是你们为我的大学学习生涯不断添砖加瓦，为我将来研究

36、生阶段的学习打下坚实的基础。我还要感谢关心帮助我的各位同学。没有同学们的帮助，我的毕业论文不会完成的如此顺利。有了你们，我的大学生活才充满乐趣，我的大学生活才完整。最后要感谢各参考文献的作者，你们的研究成果为毕业论文提供了良好的素材，使我幡然顿悟，让我在巨人肩膀上看的更远。区区数言不尽能表达我心中无尽的感谢。祝愿老师们工作顺利，愿同学们前程似锦。再次感谢大学期间我的老师，我的同学。参考文献1，李治平等.偏微分方程数值解M. 北京：人民邮电出版社, 2006.2李荣华. 偏微分方程数值解法M.北京:高等教育出版社,2010.3余德浩，汤华中. 微分方程数值解法M.北京：科学出版社，2004.4秦

37、新强, 马逸尘. 对流占优扩散方程的一种特征差分算法J. 高等学校计算数学学报, 26(1),2004,30-37.5梁栋. 对流扩散方程的一类迎风格式J. 计算数学, 1991(02),133-141.6陆金甫. 关于非定态对流扩散方程的迎风加权差分格式J. 清华大学学报（自然科学版）, 32(3), 1991, 138-148.7秦新强, 李秋芳. 对流占优扩散方程一种新的特征差分算法J. 西安理工大学学报, 19(2),2003.附录（代码）中心差分格式a=1;b=10;L=5;T=5;h=0.5;t=0.02;r=(a*t)/h2;r1=(b*t)/h;u=zeros(L/h+1,T/

38、t+1); time=0.02:t:T;hang=1:L/h+1;u(1,:)=1;for n=1:length(time) for j=2:length(hang)-1 u(j,n+1)=(r-r1)*u(j+1,n)+(1-2*r)*u(j,n)+(r+r1)*u(j-1,n); endendX,Y=meshgrid(0:t:5,0:h:5);Z=u;mesh(X,Y,Z),xlabel(time),ylabel(space),zlabel(u数值解)/迎风格式a=1;b=10;L=5;T=5;h=0.5;t=0.02;r=(a*t)/h2;r1=(b*t)/h;u=zeros(L/h+1

39、,T/t+1);%11行，251列time=0.02:t:T;hang=1:L/h+1;u(1,:)=1;for n=1:length(time) for j=2:length(hang) if j=10 u(j,n+1)=r*u(j+1,n)+(1-2*r-r1)*u(j,n)+(r+r1)*u(j-1,n); else u(j,n+1)=(1-2*r-r1)*u(j,n)+(r+r1)*u(j-1,n); end endendX,Y=meshgrid(0:t:5,0:h:5);Z=u;mesh(X,Y,Z),xlabel(time),ylabel(space),zlabel(u数值解)/迎

40、风加权格式a=1;b=10;L=5;T=5;h=0.5;t=0.02;r=(a*t)/h2;r1=(b*t)/h;q=0.8;u=zeros(L/h+1,T/t+1);%11行，251列time=0.02:t:T;hang=1:L/h+1;u(1,:)=1;for n=1:length(time) for j=2:length(hang) if j=10 u(j,n+1)=(r-r1*(1-q)*u(j+1,n)+(1-2*r-r1*q)*u(j,n)+(r+r1)*u(j-1,n); else u(j,n+1)=(1-2*r-r1*q)*u(j,n)+(r+r1)*u(j-1,n); end

41、 endendX,Y=meshgrid(0:t:5,0:h:5);Z=u;mesh(X,Y,Z),xlabel(time),ylabel(space),zlabel(u数值解)/特征有限差分格式a=1;b=10;L=5;T=5;h=0.5;t=0.02;r=(a*t*(h+b*t)/h;u=zeros(L/h+1,T/t+1);%11行，251列time=0.02:t:T;hang=1:L/h+1;u(1,:)=1;for n=1:length(time) for j=2:length(hang) if j=2 u(j,n+1)=(1/(h+b*t)2-b*h*t)*(r*u(j+1,n)+(

42、h2-2*r)*u(j,n)+(r+b*h*t)*u(j-1,n)+(b*t)2); else if j=10 u(j,n+1)=(1/(h+b*t)2-b*h*t)*(r*u(j+1,n)+(h2-2*r)*u(j,n)+(r+b*h*t)*u(j-1,n)+(b*t)2)*u(j-1,n+1); else u(j,n+1)=(1/(h+b*t)2-b*h*t)*(h2-2*r)*u(j,n)+(r+b*h*t)*u(j-1,n)+(b*t)2)*u(j-1,n+1); end end endendX,Y=meshgrid(0:t:5,0:h:5);Z=u;mesh(X,Y,Z),xlabe

43、l(time),ylabel(space),zlabel(u数值解)/精确解a=1;b=10;L=5;T=5;h=0.5;t=0.02;x=0:0.5:5;n=0:0.02:5;u=zeros(L/h+1,T/t+1);for j=1:length(x) for k=1:length(n) u(j,k)=0.5*(erfc(x(j)-b*n(k)/(2*sqrt(n(k)+exp(b*x(j)*erfc(x(j)+b*n(k)/(2*sqrt(n(k); endendX,Y=meshgrid(0:t:5,0:h:5);Z=u;mesh(X,Y,Z),xlabel(time),ylabel(space),zlabel(u数值解)/31