常微分方程差分解法.入门.多解法

毕业论文

题 目 抛物型方程的差分解法

学 院

专 业 信息与计算科学

班 级 学 生 王丹丹

学 号

指导教师 王宣欣

二〇一二 年 五 月二十五日

摘 要

偏微分方程的数值解法在数值分析中占有重要的地位,很多科学技术问题的数值计算包括了偏微分方程的数值解问题【1】。近三十多年来,数值解法的理论和方法都有了很大的发展,而且在各个科学技术的领域中应用也愈来愈广泛。本文的研究主要集中在依赖于时间的问题,借助于简单的常系数扩散方程,介绍抛物型方程的差分解法。本文以基本概念和基本方法为主,同时结合算例实现算法。

第一部分介绍偏微分方程及差分解法的基本概念,引入本文的研究对象——常系u2u数扩散方程:a2,xR,t0 tx

第二部分介绍上述方程的几种差分格式及每种格式的相容性、收敛性与稳定性。 第三部分通过算例检验每种差分格式的可行性。

关键词:偏微分方程;抛物型;差分格式;收敛性;稳定性;算例

ABSTRACT The numerical solution of partial differential equation holds an important role in numerical analysis .Many problems of compution in the field of science and techology include the numerical solution of partial differential equation. For more than 30 years, the theory and method of the numerical computation made a great development and its applications in various fields of science and technology are more and more widely. This paper focuses on the problems based on time. I will use object-constant diffusion equation to introduces the finite difference method of parabolic equation. This paper mainly focus on the basic concept ,basic method and simple numerical example.

The first part of this paper introduces partial differential equations and basic concepts of finite difference method.I will introduce the object-constant diffusion equation for the

u2ufirst time. a2,xR,t0 tx

The second part of this paper introduces several difference schemes of the above equation and their compatibility ,convergence and stability.

The third part tests the accuracy of each scheme.

Key words:partial differential equation;parabolic;difference scheme;convergence;stability;application

目 录

摘要..............................................................………………..................I ABSTRACT………………………………………………………………………...…II 目录…………………………………………………………………………………..……III 1前言……….…………………………………………….….……………………….....….1 2基本概念和定理………………………………………………………………………….2

2.1抛物型方程的基本概念.....……………….………….…………………................2

2.1.1偏微分方程的定义……………………………………………………….…2

2.1.2抛物型方程的定义…………………………………………………….……2

2.1.3初边值条件的定义…………………………………………………….……3

2.2 差分方法的基本思想……………………………………………………..………3

2.3网格剖分...................................................................................................................4

2.4截断误差的基本概念...............................................................................................5

2.5相容性的基本概念...................................................................................................7

2.6收敛性的基本概念...................................................................................................7

2.7稳定性的基本概念...................................................................................................8

2.7.1判断稳定性的直接法…………………………………………………...…..8

2.7.2判断稳定性的Fourier方法……………………………………………....….9

3常系数扩散方程的差分格式及其相容性、收敛性和稳定性分析…….........................12

3.1向前差分格式…………………………………………………............................12

3.2向后差分格式........................................................................................................13

3.3 Crank-Nicolson格式..............................................................................................14

3.4 Richardson格式.....................................................................................................16

4差分解法的应用...............................................................................................................18

结论......................................................................................................................................25

参考文献..................................................... .................. .....................................................26 致谢......................………………….……………………..…….…………...…………….27 附录……………………………………………………………………………………...……………..28

1前言

微积分方程这门学科产生于十八世纪,欧拉在他的著作中最早提出了弦振动的二阶方程,随后不久,法国数学家达朗贝尔也在他的著作《论动力学》中提出了特殊的偏微分方程[2]。偏微分方程得到迅速发展是在十九世纪,那时候,数学物理问题的研究繁荣起来了,许多数学家都对数学物理问题的解决做出了贡献。对于偏微分方程的求解,虽然具有明确表达式的解析解是很好的结果,但是能求出解析解的情况却十分有限,即使是最简单的双变量二阶线性常系数偏微分方程,也往往难以得到解析解。这是因为方程的解除了取决于方程本身的复杂度外,还要考虑到边界条件的复杂性[3]。很复杂的二阶偏微分方程,也许因为边界条件的简单性存在很简单的解析解,但是如果边界条件稍微复杂,就算是二阶常微分方程也没有解析解。

从目前的研究现状来看,偏微分方程数值解的理论和方法都日趋成熟,很多学术论文都在力求寻找更为精确且性质良好的求解方法。而且在实际问题中常会遇到多个自变量,非线性的方程或方程组;它们还可能是混合型的偏微分方程(如机翼的跨声速绕流),其解包含着各种间断(如激波间断、按触间断等)[4]。非线性问题的差分法求解是十分困难的。抛物型方程的数值解法目前有傅里叶算法(SSPE)、有限差分法(FDM)、有限元法(FEM)等,每种方法都有自己的适用范围,虽然SSPE 的效率高,但本文将选择使用更容易处理阻抗边界条件的FDM。由于FDM对网格间距要求足够小,计算效率很依赖计算机硬件速度,21 世纪前,大多是FDM的理论推导和误差分析,直到2007 年国际上才出现公开发表的使用FDM求解抛物型方程的实验并得出简单模型的计算结果。随着电子计算机的发展,在解决各种非线性问题中,FDM法得到了很快的发展,并且出现了许多新的思想和方法,如守恒差分格式,时间相关法、分步法等。

本文将从基本概念和基本方法入手,通过简单的常系数扩散方程,介绍抛物型方程的差分解法及其简单的实际应用,起到初步介绍偏微分方程数值解法的目的,希望有助于初学者了解相关基本知识,培养进一步学习的兴趣。

2基本概念和定理

2.1抛物型方程的基本概念

2.1.1偏微分方程的定义

偏微分方程在科学研究和工程技术中常常会出现,比如核反应和核爆炸过程的数学模型、飞行器设计过程中的空气动力学问题等等。

定义2.1 含有未知函数u(x1,x2,…,xn,t)的偏导数的方程称为偏微分方程。如果方程中对于未知函数和它的所有偏导数都是线性的,这样的方程称为线性偏微分方程,否则称它为非线性偏微分方程。方程中出现的偏导数的最高阶是方程的阶。

下面举出几个典型的偏微分方程:

(1)Laplace方程[5]

2u2u2u…+20 22x1x2xn

其中uu(x),该方程称为调和方程。在力学、电学常常遇到的势函数满足这个方程。

(2)对流扩散方程

uuu2u2u2uabFk(222) xytxyt

其中uu(x,y,t)表示流场中某种物质的浓度,(a,b)是流速。

(3)波动方程

22u2u2u2ua(2…+2)F(x,t) t2x1x22xn

其中uu(x),而F(x,t)给定。在一些声学、光学和力学的波动问题中常常出现这类方程。

2.1.2抛物型方程的定义

定义2.2 设uu(x1,x2,…,xn),其中xn可以是时间变量t,二阶拟线性方程指

n2uuabcuf,其中aij,b,c和f可以与x1,x2,…,xn有关,也可以ijixxxi,j1i1ijin

与u和u有关。 xi

n2uubicuf,设矩阵A[aij]是定义2.3 对于二阶拟线性方程aijxxxi,j1i1ijin

一个nn的矩阵,如果A的特征值至少有一个为零,则该方程称为抛物型方程。

考虑常系数扩散方程 u2ua2,xR,t0tx (2.1)

其中u是扩散过程中某种物质的浓度,a是扩散系数。显然它是二阶线性方程,其中

a01Aa11a,a12a21a220,它的矩阵为A,的特征值为,20,所以100a

方程(2.1)是一个抛物型方程。

在下文中我们将以该方程为基本模型讨论适用于抛物型方程定解问题的几种差分格式。

2.1.3初边值条件的定义

定义2.4 对于偏微分方程我们都是在一些特定条件下求方程的解,这样的条件称为定解条件。如果在R的某个区域内求解方程,在的边界上给出u的条件称n

为边界条件。在超平面tt0给出的条件称为初始条件。给出了方程和定解条件,就

构成了一个定解问题。

按照定解条件的不同给法,可将方程(2.1)的定解问题分为两类。

第一类,初值问题:空间变量x的变化范围是x。求具有所需次数偏微分的函数u(x,y),满足(2.1)和初始条件

u(x,0)(x),x

第二类,初边值问题:空间变量x的变化范围是axb。求具有所需次数偏微分的函数u(x,y),满足(2.1)和初始条件

u(x,0)g(x),x

及边界条件

u(a,t)1(t),u(b,t)2(t)

2.2 差分方法的基本思想

差分方法又称为有限差分方法或网格法,是求偏微分方程定解问题的数值解中应

用最广泛的方法之一。它的基本思想是:先对求解区域作网格剖分,将自变量的连续变化区域用有限离散点(网格点)集代替;将问题中出现的连续变量的函数用定义在网格点上离散变量的函数代替;通过用网格点上函数的差商代替导数,将含连续变量的偏微分方程定解问题化成只含有限个未知数的代数方程组(称为差分格式)。如果差分格式有解,且当网格无限变小时其解收敛于原微分方程定解问题的解,则差分格式的解就作为原问题的近似解(数值解)。

衡量一个差分格式是否经济实用,由多方面的因素决定,主要有:

1、计算简单。显格式无需解方程组,计算较隐格式简单。但有些隐格式左端的系数矩阵式与n无关的三对角矩阵,相当于求解三对角矩阵系数的方程组,用消元法也是可取的。

2、收敛性和收敛速度。当网格比固定,步长趋于零时,差分解un

j应收敛到真解,

并希望有尽可能快的收敛速度。

3、稳定性。由于初始数据有误差,并不可避免的有舍入误差,因此必然要关心误差传递时是无限增长还是可以被控制。

因此,用差分方法求偏微分方程定解问题一般需要解决以下问题:

1、选取网格;

2、对微分方程及定解条件选择差分近似,列出差分格式;

3、求解差分格式;

4、讨论差分格式解对于微分方程解的收敛性及误差估计。

下面我们就按照这样的步骤来进行讨论。

2.3网格剖分的概念

用有限差分方法求解偏微分方程问题必须把连续问题进行离散化,因此首先要对求解区域进行网格剖分[6]。

在求解区域的平面画出两族平行于坐标轴的直线,把平面分成矩形网格。这样的直线称为网格线,其交点称为网格点或节点[7]。平行于t轴的直线可以是等距的,可设距离为x0,记为h,称为空间步长。平行于x轴的直线则大多是不等距的,设为t0,记为,称为时间步长。

下面用例子来说明不同区域的剖分。

例2.1 双曲型方程和抛物型方程的初值问题,求解区域为:

1(x,t)|x,t0

两族网格线可以写为:

xxjjxjh,j0,1,2,…,

ttnntn,n0,1,2,…

1的网格剖分见图2.1.

图2.1

例2.2双曲型方程和抛物型方程的初边值问题,设其求解区域是

2(x,t)|0xl,t0

这个区域的网格由平行于t轴的直线族

xxj,j0,1,…,J

与平行于x轴的直线族

ttn,n0,1,2,… l所构成,其中xjxjh,xh;tnntn。 j

2的网格剖分见图

2.2

图 2.2

2.4有限差分格式的截断误差

为叙述方便,先引入下面的差分记号:

向前差分

tv(x,t)v(x,tt)v(x,t),

xv(x,t)v(xx,t)v(x,t)

向后差分

tv(x,t)v(x,t)v(x,tt),

xv(x,t)v(x,t)v(xx,t).

中心差分

1

2 11xv(x,t)v(xx,t)v(xx,t).22tv(x,t)v(x,tt)v(x,tt),12

二阶中心差分

x2u(x,t)v(xx,t)2v(x,t)v(xx,t)

那么对于扩散方程(2.1)的解,关于t的向前差分的Taylor级数展开为

u12u13u2(x,t)(x,t)3…. (2.2) tu(x,t)u(x,t)u(x,t)(x,t) 23t2t6t

对变量x进行Taylor级数展开为

2u14u2(x,t)h4…. u(x,t)2(x,t)h4x12x2

x (2.3)

定义2.5 用微分方程的解u(xj,tn)来替代差分格式中的全部近似解un

j,这样得到

的方程两边的差就是截断误差[8]。如果一个差分格式的截断误差TO(p)O(hq),则称差分格式对是p阶精度,对h是q阶精度。

例2.3 考虑下列差分格式的截断误差

1nunujj

a1nnun2uujjj1

h20, (2.4)

假定u(x,t)是充分光滑的,将u(xj,tn)带入上式的左侧,利用(2.2)式,(2.3)式有

u2u12ua4u2T(x,t)(a2)h…24tx2t12x

12ua4u2h…242t12x

这样就得到差分格式(2.4)的截断误差为O()O(h2)。事实上,对于不在边界上的任何一点(x,t),可以定义(2.4)式的截断误差T(x,t)为

11

T(x,t)tu(x,t)a2x2u(x,t),

h

由截断误差的定义及以上例子可知,只要网格剖分得很细,即和h很小,那么

偏微分方程的解近似地满足相应的差分方程。其实,一个有限差分格式的截断误差表示了用u(xj,tn)(偏微分方程的解)代替unj(差分方程的解)的差分方程与偏微分方程在点(xj,tn)上的差。

2.5有限差分格式的相容性

从偏微分方程建立差分方程时,总是要求当0,h0时差分方程能与微分方程充分“接近”,这就导致了差分方程的一个基本特征,差分格式的相容性。

考虑一般性的问题,设L为微分算子,初值问题可以叙述为

Lu0,

 (2.5)

u(x,0)g(x).其差分格式可以统一写成

1

Lhun un (2.6) jj,

其中Lh是一个依赖于和h的线性算子,称为差分算子,它把定义在第n层上的函数

n1层上的函数un1。 un变换到定义在第jj

设(2.6)式为(2.5)式的差分格式,则相应的截断误差应是 1

T(xj,tn)(Lhu(xj,tn)u(xj,tn)),

定义2.6 设u(x,t)是定解问题(2.5)的充分光滑解,(2.6)式为求解问题(2.5)所得的差分格式,如果,当,h0时有

T(xj,tn)0,

则称差分格式(2.6)与定解问题(2.5)是相容的。由此可知,相容性表达了微分方程与差分方程间的关系。

2.6有限差分格式的收敛性

构造一个差分格式,它是否能在实际中应用,应该考虑当时间步长和空间步长h

无限缩小时,差分格式的解能否逼近微分方程问题的解,这就是差分格式的收敛性问题。

定义2.7 设u(x,t)是偏微分方程的解,unj是逼近这个微分方程的差分格式的“真解”,所谓真解是指在求解差分格式的过程中忽略了各种类型的误差。如果当时间步

n

长和空间步长h趋向于0时,enju(xj,tn)uj0,我们称差分格式是收敛的。

应该注意的是,收敛性和相容性是两个完全不同的概念。对于一个相容的差分格式,从定义判断其是否收敛是不容易的,通常需要寻求一些判别差分格式的收敛准则。不加证明地,我们给出如下定理。

定理2.1(Lax等价定理) 给定一个适定的线性初值问题以及与其相容的差分格式,则差分格式的稳定性是差分格式收敛性的充分必要条件

在应用中,差分格式的相容性是容易验证的,只要使截断误差趋于0就可以了。有了Lax等价定理,我们可以着重于差分格式的稳定性的讨论,一般不再讨论收敛性问题。差分格式一旦具有稳定性,就可以用差分格式计算出偏微分方程的近似解来。

2.7有限差分格式的稳定性

2.7.1判断稳定性的直接法

1

利用有限差分格式进行计算时是按时间逐层推进的,计算第n1层上的值unj要1

用到第n层上的计算结果,而前一步的舍入误差必然会影响到unj的值。所以要分析

这种误差传播的情况,我们希望误差的影响不会越来越大以至于掩盖差分格式的解的面貌,这就是所谓稳定性问题。

差分格式可以统一写成(2.6)式,重复应用该式,有

n0

unjLhuj

'

(2.6 )

为了度量误差及其他应用,引入范数

2||un||h(un)hj

j

现在给出差分格式(2.6)的稳定性描述。

12

0nn

定义2.8 设u0 有一个误差,则就有一个误差ujjjj。如果存在一个正常数K,

使得当0,nT时,一致地有

n

||K|0| | | (2.7) ||

则称差分格式(2.6)是稳定的。

2.7.2判断稳定性的Fourier方法

按差分格式稳定性的定义来直接验证其稳定性往往比较复杂,对于现行常系数偏微分方程初值问题可以用Fourier变换来进行求解和研究。Fourier积分和Fourier变换是很多数学分支用到的有利工具,本文用于差分方法的分析。

1、Fourier变换的定义:

从Fourier级数出发进行推导可以得到Fourier积分,所谓Fourier积分公式,是指定义在(,)上的函数v(x)的一个关系式,设

11i(x)

。令v(x)v()eddv()

22Fourier变换[9]。





|v(x)|dx,有关系式





v(x)e

ix

dx,则v()称为v(x)的

2、线性常系数方程的Fourier方法

由于解un为了应用Fourier方法讨论方程(2.1)j及初值g(xi)只在网格点上有意义,的第一类初值问题,需要扩充这些函数的定义域,使它们在整个实轴R上都有定义。令

U(x,tn)unj,xj

hhxxj22hh

(x)g(xj),xjxxj

22

这样,U(x,tn),(x)对任意xR都有定义。

例2.4 以差分格式 为例讨论Fourier方法。 (2.9)式可以写为



1ununjj

a

n1

unjuj

h2

0 (2.9)

U(x,nt1)U(xn,t)a[U(nx,t)

U(xnh,t)]

(2.10)

ik(xh)

U(k,t)edkn



对(2.10)式两边用Fourier积分,可以得到

U(k,t)edkn1

ikx

U(k,t)edkan



ikx

U(k,t)edkn



ikx

消去eikx可得

U(k,tn1)[1a(1e

ikh

)]U(k,tn)

以上方法推广到一般形式的差分格式可以得到

U(k,tn1)G(,k)U(k,tn)

式中G(,k)为增长因子,显然差分格式(2.9)的增长因子为G(,k)[1a(1eikh)]。

不加证明地,我们给出以下结论:

差分格式稳定的充分必要条件是存在常数00,K0使得当0,nT,kR时,有

||G(,k)n||K

上式也称为von Neumann条件[10]。

3、线性常系数方程组的Fourier方法

稳定性的概念及相关Fourier方法的推导可以推广到线性常系数方程组。这在下面我们讨论Richardson差分格式时需要用到。

一般的差分方程组可以写为

u

nj

n1j

al

A(x,)T

a

j

l

a

u

nj

其中uRp,Aa(x,)Rpp。由于hg(),即h和满足一定关系,故在Aa(x,)中仅标出。令

C(xj,)

al

A(x,)T

a

j

l

a

则有

u

n1j

(2.11) C(xj,)u

nj

其中C(xj,)称为差分算子,上式为一个差分格式。 由(2.11)式可以得出

u[C(xj,)]u

如果C(xj,)不依赖于xj,即为常系数差分方程组,则可利用Fourier积分得到





nj

n

0j

U(k,tn1)G(,k)U(k,tn),





U(k,tn1)G(,k)U(k,tn),



其中U(k,tn)Rp,G(,k)Rpp,称G(,k)为增长矩阵。

类似于差分方程的情况,我们有如下结论:

差分格式(2.11)稳定的充要条件是存在常数00,K0使得当0,nT,

kR时,有

||G(,k)n||K

定理2.2 差分格式(2.11)稳定的必要条件是当0,nT,kR时,有

|j(G(,k))|1M,j1,2,…,p,

其中j(G(,k)n)表示G(,k)的特征值,M为常数。

定理2.3[11] 当G(,k)只有一个元素时,von Neumann条件是差分格式稳定的充要条件。

对于Fourier方法,还有很多定理和推论可以帮助我们更加容易的判断不同性质的差分格式的稳定性,本文仅给出了需要用到的基本定理。

3常系数扩散方程的差分格式及其相容性、收敛性和稳定性分析

用Taylor级数展开方法求解偏微分方程实际上等价于用差商来近似微商得到相应的差分格式,除此之外也可以用积分的方法。下面用Taylor级数展开的方法来讨论扩散方程初值问题

u2u

a2,xR,t0tx (3.1)

)g(x)x, R u(x,0 (3.2)

的差分格式。

3.1向前差分格式

假定偏微分方程初值问题的解u(x,t)是充分光滑的。由Taylor级数展开有

u(xj,tn1)u(xj,tn)un

[]jO(),

t

u(xj,tn1)u(xj,tn1)un2

[]O(),j2t

u(xj1,tn)u(xj,tn)un

[]jO(h),hx

u(xj,tn)u(xj1,tn)[u]nO(h),

j

hx

u(xj1,tn)u(xj1,tn)[u]nO(h2),

j

2hx

u(x,t)2u(x,t)u(x,t)2

uj1njnj1n[2]njO(h2).

2hx  (3.3)

1、差分格式的推导过程

利用(3.3)中第一式和第六式有:

u(xj,tn1)u(xj,tn)

a

u(xj1,tn)2u(xj,tn)u(xj1,tn)

h2

u2un

[a2]jO(h2)tx

如果u是(3.1)式的光滑解,即u是满足

u2u

a2tx

的光滑函数,那么,方程(3.1)可以用如下差分方程来近似:

1ununjj

a

nn

unj12ujuj1

h2

0

, (3.4)

这就是微分方程(3.1)的向前差分格式。它所联系的网点分布为:

2、截断误差

在前面的分析中,我们曾经利用这种差分格式来讨论截断误差的概念,得知其截断误差为O()O(h2)。 3、稳定性分析

用Fourier方法来讨论其稳定性,容易求出(3.1)式的增长因子是

G(,k)14asin2

kh

2,

1

,那么有|G(,k)|1,即满足von Neumann条件,所以向前

h22

1

差分格式的稳定性条件是a。这种情况的差分格式是条件稳定的。

2

其中

,如果a

4、收敛性分析

由Lax定理可知,当差分格式稳定时一定是收敛的,所以当a格式收敛。

1

时,向前差分2

3.2向后差分格式

1、差分格式的推导过程

1

上面构造的差分格式是显式的,即在时间层tn1上的每一个unj可以独立地根据在

时间层tn上的值unj得出。但是如果采用

u(xj,tn)u(xj,tn1)

[

un

]jO()t,

和(3.3)中的六式,就可以得到扩散方程(3.1)的另一个差分格式

n1

unjuj

nn

unj12ujuj1

a

h

2

0

(3.5)

这是微分方程(3.1)的向后差分格式。它所联系的网点分布为

:

向后差分格式与向前差分格式不同,在新时间层n上包含了3个未知量ununj1,j,

nn

unj1,因此不能由uj1直接计算出uj,这种差分格式称为隐式格式。

2、截断误差

考虑其截断误差

122uah24u

T(x,t)tu(x,t)a2xu(x,t)(x,)(,t)24

h2t12x,

1

其中(t,t),(xh,xh)。因此其截断误差为O()O(h2)。 3、稳定性分析

先把差分格式变形为

1n11n

aunaunj1(12a)ujj1uj

h

式的增长因子为

其中

2

,令ujvje

nhnjki

,并把它带入上面方程并消去公因子e

ikjh

,容易求出(3.5)

G(,k)

1kh

14asin2

2

由于a0,所以对任何网格比都有|G(,k)|1。由定理2.3知,差分格式(3.5)是稳定的。 4、收敛性分析

由Lax定理可知,当差分格式稳定时一定是收敛的,所以差分格式(3.5)收敛。

3.3 Crank-Nicloson格式

1、差分格式的推导过程

将向前差分格式与向后差分格式作算术平均,即得Crank-Nicolson格式:

1nunujj

a1n11nnn

[(ununj12ujj1)(uj12ujuj1)]0 (3.6)2

2h

[12]

Crank-Nicolson格式又叫六点对称格式,因为这种格式所联系的网点分布为

:

Crank-Nicolson格式也是隐式格式。 2、截断误差 令

Lhu

nj

1ununjj

n1n1n!nnn

auj12ujuj1uj12ujuj1[] 222hh

将截断误差

T(x,t)Lhu(xj,tn)[Lu]nj

在(xj,t

n

1)(t2

n

1

2

1

(n))处展开,则得

2

T(x,t)O(2h2)

因此其截断误差为O(2h2)。 3、稳定性分析

仍然用Fourier方法。先把差分格式变形为

nnn1n11

aunaunj12(1a)ujauj1auj12(1a)ujj1

nikjhikjh

令un,带入上面方程并消去公因子e,容易求出(3.6)式的增长因子为 jvje

kh

G(,k)

12asin2

2

12asin2

由于a0,显然对任何网格比都有|G(,k)|1,由定理2.3知,差分格式(3.6)是稳定的。 4、收敛性分析

由Lax定理可知,当差分格式稳定时一定是收敛的,所以差分格式(3.6)收敛。

3.4 Richardson格式

1、差分格式的推导过程

利用(3.3)式中的第二式和第六式可以得到Richardson格式:

11

ununjj

2

a

nn

unj12ujuj1

h2

0

(3.7)

从(3.7)式可以看出,这种格式联系到三个时间层,称为三层格式。一般地,一个多于二层的差分格式成为多层差分格式。它所联系的网点分布为:

2、截断误差

这个格式是为了提高差分方程的截断误差所设计的,利用上述方法我们可以很容易的求出其截断误差为O(2h2)。单从截断误差这一角度来考虑,这个格式是比较好的,但是以后的分析可以看到,它是没有实用价值的。 3、稳定性分析

首先将(3.7)式变形为

11nn

unun2a(unjjj12ujuj1)

讨论多层差分格式,一般先化成与其等价的二层差分方程组。Richardson差分方程的等价二层差分方程组为

n1nnnnujvj2a(uj12ujuj1),

n1n

vjuj.

nT

如果令u[un,vjj],则上面的方程组可以写成

n1j

nj

u

nj

n

2a004an

uj1000n2a

uj000n

uj1 0

设uveikjh,将其带入上式并消去公因子eikjh,可以得

济南大学毕业论文

v

n1

2kh8asin2

1

1

vn 0

因此增长矩阵为

2kh8asin

G(,k)2

1

1 0

其特征值为

1,2

1

kh224kh2

4asin(116asin)

22

2

1kh224kh2

14asin(116asin)

22

2

则有

kh

2

由定理2.2可知Richardson格式是不稳定的。这种情况的差分格式称为绝对不稳定的。

|1|14asin2

4、Richardson格式的改进

1953年Du Fort和Frankle对Richardson格式进行了修改,提出了如下的DuFort-Frankle格式[13]

1n1

unujj

2

a

n1n1n

un(uu)uj1jjj1

h2

0

这是一个无条件稳定的差分格式,但是其相容性是有条件的。在这里对DuFort-Frankle格式我们不做过多讨论。

事实证明,我们无法构造出无条件相容和无条件稳定的显示格式。

济南大学毕业论文

当x0.2时,数值解与解析解比较如表4.4

表4.4

当t0.07时,数值解与解析解比较如图4.2:

图4.2

从结果来看,向后差分格式是稳定的,它的最大优点在于稳定性、收敛性对于一切网格比都成立,网格比增大意味着工作量减少,所以隐式格式当1时的工作量比显示格式要少,然而从截断误差的角度来看,h也不能太大,否则太大,误差也就太大,所以更高精度的格式也是十分需要的。 (3)Crank-Nicloson格式为

1ununjj

an1n1n1nnn

[(u2uu)(u2uuj1jj1j1jj1)]0 22h

写成便于计算的格式为

1n11nnnununj16ujj1uj12ujuj1

(4.4)

初边值条件的格式为

结 论

首先本文给出了偏微分方程及抛物型方程的基本概念,其次给出了差分格式的基本概念,重点介绍了差分格式的收敛性及稳定性,随后针对常系数扩散方程给出了各种差分格式的推导过程及其各种性质,最后举出数值算例对给出的各种差分格式进行应用,以观察其可行性。结果证明,向前差分格式、向后差分格式、Crank-Nicloson格式都是稳定的,或者说在一定条件下是稳定的,Richardson格式是绝对不稳定的。

数值计算主要研究如何利用计算机更好的解决各种数学问题,包括连续系统离散化和离散形方程的求解,在对具体问题进行分析时,数值计算比解析解的方法更方便。几乎所有的偏微分方程都能求出数值解,数值解虽然是近似值,但对于许多方程,近似解可以足够精确,满足应用需求,编程易实现,可以分析误差的界[15]。因此,对于这个问题的研究是有很强的实用性的。

参 考 文 献

[1]黎阳.抛物型方程的有限差分解法及其在复杂磁场中的应用[D].武汉:武汉理工大学,2010:11-12.

[2]李德元,抛物型方程差分方法引论[M].北京:科学出版社,1982:80-92. [3]关治,陆金甫.数值分析基础[M].北京:高等教育出版社,1979:200–260.

[4]葛坤.一类抛物型方程的数值解法研究与应用[D].南京:南京林业大学,2011:11-26.

[5]Per Brinch Hansen.Numerical Solution of Laplace's Equation[J].Electrical Engineering and Computer Science Technical Reports,1992(4):34-46.

[6]孙雪莉.抛物型方程任意偶数阶精度的两层显示格式[D].西安:西安建筑科技大学,2009:10-17. [7] 张宝林,谷同祥,莫则尧.数值并行计算原理与方法[M].北京:国防工业出版社,1999:42–43. [8]杨丹丹.求解抛物型方程的一种有限差分并行格式[D].吉林:吉林大学数学研究所,2010:10-15. [9]陆金甫,关治.偏微分方程数值解[M].北京:清华大学出版社,2004:28-34.

[10]林鹏程.解双抛物型方程的高精度稳定差分格式[J].福州大学学报,2004(4):120–121. [11]张耀科,陆如枫.初边值问题稳定性分析的一个标注[J].工科数学,1994(1):64-78. [12]许秀娟.两类抛物型方程的有限差分解法[D].哈尔滨:哈尔滨工业大学,2009:8-21.

[13]C. Fuhse,T. Salditt.Finite-difference field calculations for two-dimensionally confined x-ray waveguides[J].Applied Optics ,July 2006,45(19): 507–511.

[14]王爱锋,曲小钢.解抛物型方程的一种隐式差分格式[J].纺织高校基础科学学报,2010(4):32-37.

[15]Adel A.S.Almarashi.Approximation Solution of Fractional Partial Differential Equation By Neural Network[J].Development Of Mathematics,2011(10):98-126.

致 谢

历时将近两个月的时间终于将这篇论文写完,在论文的写作过程中遇到了无数的困难和障碍,都在同学和老师的帮助下度过了。尤其要强烈感谢王宣欣老师,她对我进行了无私的指导和帮助,不厌其烦的帮助进行论文的修改和改进。另外,在校图书馆查找资料的时候,图书馆的老师也给我提供了很多方面的支持与帮助。在此向帮助和指导过我的各位老师表示最中心的感谢!

感谢这篇论文所涉及到的各位学者。本文引用了数位学者的研究文献,如果没有各位学者的研究成果的帮助和启发,我将很难完成本篇论文的写作。

感谢我的同学和朋友,在我写论文的过程中给予我了很多你问素材,还在论文的撰写和排版灯过程中提供热情的帮助。

由于我的学术水平有限,所写论文难免有不足之处,恳请各位老师和学友批评和指正!

27 - -

附录

追赶法程序如下:

#include

#define N 9

main()

{

double A[N][N]={{6,-1,0,0,0,0,0,0,0},{-1,6,-1,0,0,0,0,0,0},{0,-1,6,-1,0,0,0,0,0},{0,0,-1,6,-1,0,0,0,0},{0 ,0,0,-1,6,-1,0,0,0},{0,0,0,0,-1,6,-1,0,0},{0,0,0,0,0,-1,6,-1,0},{0,0,0,0,0,0,-1,6,-1},{0,0,0,0,0,0,0,-1,6}};b

[N]={0.898922,1.709851,2.353408,2.766597,2.908972,2.766597,2.353408,1.709851,0.898922},m[9];

}

28 - -int i; A[0][1]=A[0][1]/A[0][0]; for(i=1;i=0;i--) b[i]=b[i]-A[i][i+1]*b[i+1]; for(i=0;i

毕业论文

题 目 抛物型方程的差分解法

学 院

专 业 信息与计算科学

班 级 学 生 王丹丹

学 号

指导教师 王宣欣

二〇一二 年 五 月二十五日

摘 要

偏微分方程的数值解法在数值分析中占有重要的地位,很多科学技术问题的数值计算包括了偏微分方程的数值解问题【1】。近三十多年来,数值解法的理论和方法都有了很大的发展,而且在各个科学技术的领域中应用也愈来愈广泛。本文的研究主要集中在依赖于时间的问题,借助于简单的常系数扩散方程,介绍抛物型方程的差分解法。本文以基本概念和基本方法为主,同时结合算例实现算法。

第一部分介绍偏微分方程及差分解法的基本概念,引入本文的研究对象——常系u2u数扩散方程:a2,xR,t0 tx

第二部分介绍上述方程的几种差分格式及每种格式的相容性、收敛性与稳定性。 第三部分通过算例检验每种差分格式的可行性。

关键词:偏微分方程;抛物型;差分格式;收敛性;稳定性;算例

ABSTRACT The numerical solution of partial differential equation holds an important role in numerical analysis .Many problems of compution in the field of science and techology include the numerical solution of partial differential equation. For more than 30 years, the theory and method of the numerical computation made a great development and its applications in various fields of science and technology are more and more widely. This paper focuses on the problems based on time. I will use object-constant diffusion equation to introduces the finite difference method of parabolic equation. This paper mainly focus on the basic concept ,basic method and simple numerical example.

The first part of this paper introduces partial differential equations and basic concepts of finite difference method.I will introduce the object-constant diffusion equation for the

u2ufirst time. a2,xR,t0 tx

The second part of this paper introduces several difference schemes of the above equation and their compatibility ,convergence and stability.

The third part tests the accuracy of each scheme.

Key words:partial differential equation;parabolic;difference scheme;convergence;stability;application

目 录

摘要..............................................................………………..................I ABSTRACT………………………………………………………………………...…II 目录…………………………………………………………………………………..……III 1前言……….…………………………………………….….……………………….....….1 2基本概念和定理………………………………………………………………………….2

2.1抛物型方程的基本概念.....……………….………….…………………................2

2.1.1偏微分方程的定义……………………………………………………….…2

2.1.2抛物型方程的定义…………………………………………………….……2

2.1.3初边值条件的定义…………………………………………………….……3

2.2 差分方法的基本思想……………………………………………………..………3

2.3网格剖分...................................................................................................................4

2.4截断误差的基本概念...............................................................................................5

2.5相容性的基本概念...................................................................................................7

2.6收敛性的基本概念...................................................................................................7

2.7稳定性的基本概念...................................................................................................8

2.7.1判断稳定性的直接法…………………………………………………...…..8

2.7.2判断稳定性的Fourier方法……………………………………………....….9

3常系数扩散方程的差分格式及其相容性、收敛性和稳定性分析…….........................12

3.1向前差分格式…………………………………………………............................12

3.2向后差分格式........................................................................................................13

3.3 Crank-Nicolson格式..............................................................................................14

3.4 Richardson格式.....................................................................................................16

4差分解法的应用...............................................................................................................18

结论......................................................................................................................................25

参考文献..................................................... .................. .....................................................26 致谢......................………………….……………………..…….…………...…………….27 附录……………………………………………………………………………………...……………..28

1前言

微积分方程这门学科产生于十八世纪,欧拉在他的著作中最早提出了弦振动的二阶方程,随后不久,法国数学家达朗贝尔也在他的著作《论动力学》中提出了特殊的偏微分方程[2]。偏微分方程得到迅速发展是在十九世纪,那时候,数学物理问题的研究繁荣起来了,许多数学家都对数学物理问题的解决做出了贡献。对于偏微分方程的求解,虽然具有明确表达式的解析解是很好的结果,但是能求出解析解的情况却十分有限,即使是最简单的双变量二阶线性常系数偏微分方程,也往往难以得到解析解。这是因为方程的解除了取决于方程本身的复杂度外,还要考虑到边界条件的复杂性[3]。很复杂的二阶偏微分方程,也许因为边界条件的简单性存在很简单的解析解,但是如果边界条件稍微复杂,就算是二阶常微分方程也没有解析解。

从目前的研究现状来看,偏微分方程数值解的理论和方法都日趋成熟,很多学术论文都在力求寻找更为精确且性质良好的求解方法。而且在实际问题中常会遇到多个自变量,非线性的方程或方程组;它们还可能是混合型的偏微分方程(如机翼的跨声速绕流),其解包含着各种间断(如激波间断、按触间断等)[4]。非线性问题的差分法求解是十分困难的。抛物型方程的数值解法目前有傅里叶算法(SSPE)、有限差分法(FDM)、有限元法(FEM)等,每种方法都有自己的适用范围,虽然SSPE 的效率高,但本文将选择使用更容易处理阻抗边界条件的FDM。由于FDM对网格间距要求足够小,计算效率很依赖计算机硬件速度,21 世纪前,大多是FDM的理论推导和误差分析,直到2007 年国际上才出现公开发表的使用FDM求解抛物型方程的实验并得出简单模型的计算结果。随着电子计算机的发展,在解决各种非线性问题中,FDM法得到了很快的发展,并且出现了许多新的思想和方法,如守恒差分格式,时间相关法、分步法等。

本文将从基本概念和基本方法入手,通过简单的常系数扩散方程,介绍抛物型方程的差分解法及其简单的实际应用,起到初步介绍偏微分方程数值解法的目的,希望有助于初学者了解相关基本知识,培养进一步学习的兴趣。

2基本概念和定理

2.1抛物型方程的基本概念

2.1.1偏微分方程的定义

偏微分方程在科学研究和工程技术中常常会出现,比如核反应和核爆炸过程的数学模型、飞行器设计过程中的空气动力学问题等等。

定义2.1 含有未知函数u(x1,x2,…,xn,t)的偏导数的方程称为偏微分方程。如果方程中对于未知函数和它的所有偏导数都是线性的,这样的方程称为线性偏微分方程,否则称它为非线性偏微分方程。方程中出现的偏导数的最高阶是方程的阶。

下面举出几个典型的偏微分方程:

(1)Laplace方程[5]

2u2u2u…+20 22x1x2xn

其中uu(x),该方程称为调和方程。在力学、电学常常遇到的势函数满足这个方程。

(2)对流扩散方程

uuu2u2u2uabFk(222) xytxyt

其中uu(x,y,t)表示流场中某种物质的浓度,(a,b)是流速。

(3)波动方程

22u2u2u2ua(2…+2)F(x,t) t2x1x22xn

其中uu(x),而F(x,t)给定。在一些声学、光学和力学的波动问题中常常出现这类方程。

2.1.2抛物型方程的定义

定义2.2 设uu(x1,x2,…,xn),其中xn可以是时间变量t,二阶拟线性方程指

n2uuabcuf,其中aij,b,c和f可以与x1,x2,…,xn有关,也可以ijixxxi,j1i1ijin

与u和u有关。 xi

n2uubicuf,设矩阵A[aij]是定义2.3 对于二阶拟线性方程aijxxxi,j1i1ijin

一个nn的矩阵,如果A的特征值至少有一个为零,则该方程称为抛物型方程。

考虑常系数扩散方程 u2ua2,xR,t0tx (2.1)

其中u是扩散过程中某种物质的浓度,a是扩散系数。显然它是二阶线性方程,其中

a01Aa11a,a12a21a220,它的矩阵为A,的特征值为,20,所以100a

方程(2.1)是一个抛物型方程。

在下文中我们将以该方程为基本模型讨论适用于抛物型方程定解问题的几种差分格式。

2.1.3初边值条件的定义

定义2.4 对于偏微分方程我们都是在一些特定条件下求方程的解,这样的条件称为定解条件。如果在R的某个区域内求解方程,在的边界上给出u的条件称n

为边界条件。在超平面tt0给出的条件称为初始条件。给出了方程和定解条件,就

构成了一个定解问题。

按照定解条件的不同给法,可将方程(2.1)的定解问题分为两类。

第一类,初值问题:空间变量x的变化范围是x。求具有所需次数偏微分的函数u(x,y),满足(2.1)和初始条件

u(x,0)(x),x

第二类,初边值问题:空间变量x的变化范围是axb。求具有所需次数偏微分的函数u(x,y),满足(2.1)和初始条件

u(x,0)g(x),x

及边界条件

u(a,t)1(t),u(b,t)2(t)

2.2 差分方法的基本思想

差分方法又称为有限差分方法或网格法,是求偏微分方程定解问题的数值解中应

用最广泛的方法之一。它的基本思想是:先对求解区域作网格剖分,将自变量的连续变化区域用有限离散点(网格点)集代替;将问题中出现的连续变量的函数用定义在网格点上离散变量的函数代替;通过用网格点上函数的差商代替导数,将含连续变量的偏微分方程定解问题化成只含有限个未知数的代数方程组(称为差分格式)。如果差分格式有解,且当网格无限变小时其解收敛于原微分方程定解问题的解,则差分格式的解就作为原问题的近似解(数值解)。

衡量一个差分格式是否经济实用,由多方面的因素决定,主要有:

1、计算简单。显格式无需解方程组,计算较隐格式简单。但有些隐格式左端的系数矩阵式与n无关的三对角矩阵,相当于求解三对角矩阵系数的方程组,用消元法也是可取的。

2、收敛性和收敛速度。当网格比固定,步长趋于零时,差分解un

j应收敛到真解,

并希望有尽可能快的收敛速度。

3、稳定性。由于初始数据有误差,并不可避免的有舍入误差,因此必然要关心误差传递时是无限增长还是可以被控制。

因此,用差分方法求偏微分方程定解问题一般需要解决以下问题:

1、选取网格;

2、对微分方程及定解条件选择差分近似,列出差分格式;

3、求解差分格式;

4、讨论差分格式解对于微分方程解的收敛性及误差估计。

下面我们就按照这样的步骤来进行讨论。

2.3网格剖分的概念

用有限差分方法求解偏微分方程问题必须把连续问题进行离散化,因此首先要对求解区域进行网格剖分[6]。

在求解区域的平面画出两族平行于坐标轴的直线,把平面分成矩形网格。这样的直线称为网格线,其交点称为网格点或节点[7]。平行于t轴的直线可以是等距的,可设距离为x0,记为h,称为空间步长。平行于x轴的直线则大多是不等距的,设为t0,记为,称为时间步长。

下面用例子来说明不同区域的剖分。

例2.1 双曲型方程和抛物型方程的初值问题,求解区域为:

1(x,t)|x,t0

两族网格线可以写为:

xxjjxjh,j0,1,2,…,

ttnntn,n0,1,2,…

1的网格剖分见图2.1.

图2.1

例2.2双曲型方程和抛物型方程的初边值问题,设其求解区域是

2(x,t)|0xl,t0

这个区域的网格由平行于t轴的直线族

xxj,j0,1,…,J

与平行于x轴的直线族

ttn,n0,1,2,… l所构成,其中xjxjh,xh;tnntn。 j

2的网格剖分见图

2.2

图 2.2

2.4有限差分格式的截断误差

为叙述方便,先引入下面的差分记号:

向前差分

tv(x,t)v(x,tt)v(x,t),

xv(x,t)v(xx,t)v(x,t)

向后差分

tv(x,t)v(x,t)v(x,tt),

xv(x,t)v(x,t)v(xx,t).

中心差分

1

2 11xv(x,t)v(xx,t)v(xx,t).22tv(x,t)v(x,tt)v(x,tt),12

二阶中心差分

x2u(x,t)v(xx,t)2v(x,t)v(xx,t)

那么对于扩散方程(2.1)的解,关于t的向前差分的Taylor级数展开为

u12u13u2(x,t)(x,t)3…. (2.2) tu(x,t)u(x,t)u(x,t)(x,t) 23t2t6t

对变量x进行Taylor级数展开为

2u14u2(x,t)h4…. u(x,t)2(x,t)h4x12x2

x (2.3)

定义2.5 用微分方程的解u(xj,tn)来替代差分格式中的全部近似解un

j,这样得到

的方程两边的差就是截断误差[8]。如果一个差分格式的截断误差TO(p)O(hq),则称差分格式对是p阶精度,对h是q阶精度。

例2.3 考虑下列差分格式的截断误差

1nunujj

a1nnun2uujjj1

h20, (2.4)

假定u(x,t)是充分光滑的,将u(xj,tn)带入上式的左侧,利用(2.2)式,(2.3)式有

u2u12ua4u2T(x,t)(a2)h…24tx2t12x

12ua4u2h…242t12x

这样就得到差分格式(2.4)的截断误差为O()O(h2)。事实上,对于不在边界上的任何一点(x,t),可以定义(2.4)式的截断误差T(x,t)为

11

T(x,t)tu(x,t)a2x2u(x,t),

h

由截断误差的定义及以上例子可知,只要网格剖分得很细,即和h很小,那么

偏微分方程的解近似地满足相应的差分方程。其实,一个有限差分格式的截断误差表示了用u(xj,tn)(偏微分方程的解)代替unj(差分方程的解)的差分方程与偏微分方程在点(xj,tn)上的差。

2.5有限差分格式的相容性

从偏微分方程建立差分方程时,总是要求当0,h0时差分方程能与微分方程充分“接近”,这就导致了差分方程的一个基本特征,差分格式的相容性。

考虑一般性的问题,设L为微分算子,初值问题可以叙述为

Lu0,

 (2.5)

u(x,0)g(x).其差分格式可以统一写成

1

Lhun un (2.6) jj,

其中Lh是一个依赖于和h的线性算子,称为差分算子,它把定义在第n层上的函数

n1层上的函数un1。 un变换到定义在第jj

设(2.6)式为(2.5)式的差分格式,则相应的截断误差应是 1

T(xj,tn)(Lhu(xj,tn)u(xj,tn)),

定义2.6 设u(x,t)是定解问题(2.5)的充分光滑解,(2.6)式为求解问题(2.5)所得的差分格式,如果,当,h0时有

T(xj,tn)0,

则称差分格式(2.6)与定解问题(2.5)是相容的。由此可知,相容性表达了微分方程与差分方程间的关系。

2.6有限差分格式的收敛性

构造一个差分格式,它是否能在实际中应用,应该考虑当时间步长和空间步长h

无限缩小时,差分格式的解能否逼近微分方程问题的解,这就是差分格式的收敛性问题。

定义2.7 设u(x,t)是偏微分方程的解,unj是逼近这个微分方程的差分格式的“真解”,所谓真解是指在求解差分格式的过程中忽略了各种类型的误差。如果当时间步

n

长和空间步长h趋向于0时,enju(xj,tn)uj0,我们称差分格式是收敛的。

应该注意的是,收敛性和相容性是两个完全不同的概念。对于一个相容的差分格式,从定义判断其是否收敛是不容易的,通常需要寻求一些判别差分格式的收敛准则。不加证明地,我们给出如下定理。

定理2.1(Lax等价定理) 给定一个适定的线性初值问题以及与其相容的差分格式,则差分格式的稳定性是差分格式收敛性的充分必要条件

在应用中,差分格式的相容性是容易验证的,只要使截断误差趋于0就可以了。有了Lax等价定理,我们可以着重于差分格式的稳定性的讨论,一般不再讨论收敛性问题。差分格式一旦具有稳定性,就可以用差分格式计算出偏微分方程的近似解来。

2.7有限差分格式的稳定性

2.7.1判断稳定性的直接法

1

利用有限差分格式进行计算时是按时间逐层推进的,计算第n1层上的值unj要1

用到第n层上的计算结果,而前一步的舍入误差必然会影响到unj的值。所以要分析

这种误差传播的情况,我们希望误差的影响不会越来越大以至于掩盖差分格式的解的面貌,这就是所谓稳定性问题。

差分格式可以统一写成(2.6)式,重复应用该式,有

n0

unjLhuj

'

(2.6 )

为了度量误差及其他应用,引入范数

2||un||h(un)hj

j

现在给出差分格式(2.6)的稳定性描述。

12

0nn

定义2.8 设u0 有一个误差,则就有一个误差ujjjj。如果存在一个正常数K,

使得当0,nT时,一致地有

n

||K|0| | | (2.7) ||

则称差分格式(2.6)是稳定的。

2.7.2判断稳定性的Fourier方法

按差分格式稳定性的定义来直接验证其稳定性往往比较复杂,对于现行常系数偏微分方程初值问题可以用Fourier变换来进行求解和研究。Fourier积分和Fourier变换是很多数学分支用到的有利工具,本文用于差分方法的分析。

1、Fourier变换的定义:

从Fourier级数出发进行推导可以得到Fourier积分,所谓Fourier积分公式,是指定义在(,)上的函数v(x)的一个关系式,设

11i(x)

。令v(x)v()eddv()

22Fourier变换[9]。





|v(x)|dx,有关系式





v(x)e

ix

dx,则v()称为v(x)的

2、线性常系数方程的Fourier方法

由于解un为了应用Fourier方法讨论方程(2.1)j及初值g(xi)只在网格点上有意义,的第一类初值问题,需要扩充这些函数的定义域,使它们在整个实轴R上都有定义。令

U(x,tn)unj,xj

hhxxj22hh

(x)g(xj),xjxxj

22

这样,U(x,tn),(x)对任意xR都有定义。

例2.4 以差分格式 为例讨论Fourier方法。 (2.9)式可以写为



1ununjj

a

n1

unjuj

h2

0 (2.9)

U(x,nt1)U(xn,t)a[U(nx,t)

U(xnh,t)]

(2.10)

ik(xh)

U(k,t)edkn



对(2.10)式两边用Fourier积分,可以得到

U(k,t)edkn1

ikx

U(k,t)edkan



ikx

U(k,t)edkn



ikx

消去eikx可得

U(k,tn1)[1a(1e

ikh

)]U(k,tn)

以上方法推广到一般形式的差分格式可以得到

U(k,tn1)G(,k)U(k,tn)

式中G(,k)为增长因子,显然差分格式(2.9)的增长因子为G(,k)[1a(1eikh)]。

不加证明地,我们给出以下结论:

差分格式稳定的充分必要条件是存在常数00,K0使得当0,nT,kR时,有

||G(,k)n||K

上式也称为von Neumann条件[10]。

3、线性常系数方程组的Fourier方法

稳定性的概念及相关Fourier方法的推导可以推广到线性常系数方程组。这在下面我们讨论Richardson差分格式时需要用到。

一般的差分方程组可以写为

u

nj

n1j

al

A(x,)T

a

j

l

a

u

nj

其中uRp,Aa(x,)Rpp。由于hg(),即h和满足一定关系,故在Aa(x,)中仅标出。令

C(xj,)

al

A(x,)T

a

j

l

a

则有

u

n1j

(2.11) C(xj,)u

nj

其中C(xj,)称为差分算子,上式为一个差分格式。 由(2.11)式可以得出

u[C(xj,)]u

如果C(xj,)不依赖于xj,即为常系数差分方程组,则可利用Fourier积分得到





nj

n

0j

U(k,tn1)G(,k)U(k,tn),





U(k,tn1)G(,k)U(k,tn),



其中U(k,tn)Rp,G(,k)Rpp,称G(,k)为增长矩阵。

类似于差分方程的情况,我们有如下结论:

差分格式(2.11)稳定的充要条件是存在常数00,K0使得当0,nT,

kR时,有

||G(,k)n||K

定理2.2 差分格式(2.11)稳定的必要条件是当0,nT,kR时,有

|j(G(,k))|1M,j1,2,…,p,

其中j(G(,k)n)表示G(,k)的特征值,M为常数。

定理2.3[11] 当G(,k)只有一个元素时,von Neumann条件是差分格式稳定的充要条件。

对于Fourier方法,还有很多定理和推论可以帮助我们更加容易的判断不同性质的差分格式的稳定性,本文仅给出了需要用到的基本定理。

3常系数扩散方程的差分格式及其相容性、收敛性和稳定性分析

用Taylor级数展开方法求解偏微分方程实际上等价于用差商来近似微商得到相应的差分格式,除此之外也可以用积分的方法。下面用Taylor级数展开的方法来讨论扩散方程初值问题

u2u

a2,xR,t0tx (3.1)

)g(x)x, R u(x,0 (3.2)

的差分格式。

3.1向前差分格式

假定偏微分方程初值问题的解u(x,t)是充分光滑的。由Taylor级数展开有

u(xj,tn1)u(xj,tn)un

[]jO(),

t

u(xj,tn1)u(xj,tn1)un2

[]O(),j2t

u(xj1,tn)u(xj,tn)un

[]jO(h),hx

u(xj,tn)u(xj1,tn)[u]nO(h),

j

hx

u(xj1,tn)u(xj1,tn)[u]nO(h2),

j

2hx

u(x,t)2u(x,t)u(x,t)2

uj1njnj1n[2]njO(h2).

2hx  (3.3)

1、差分格式的推导过程

利用(3.3)中第一式和第六式有:

u(xj,tn1)u(xj,tn)

a

u(xj1,tn)2u(xj,tn)u(xj1,tn)

h2

u2un

[a2]jO(h2)tx

如果u是(3.1)式的光滑解,即u是满足

u2u

a2tx

的光滑函数,那么,方程(3.1)可以用如下差分方程来近似:

1ununjj

a

nn

unj12ujuj1

h2

0

, (3.4)

这就是微分方程(3.1)的向前差分格式。它所联系的网点分布为:

2、截断误差

在前面的分析中,我们曾经利用这种差分格式来讨论截断误差的概念,得知其截断误差为O()O(h2)。 3、稳定性分析

用Fourier方法来讨论其稳定性,容易求出(3.1)式的增长因子是

G(,k)14asin2

kh

2,

1

,那么有|G(,k)|1,即满足von Neumann条件,所以向前

h22

1

差分格式的稳定性条件是a。这种情况的差分格式是条件稳定的。

2

其中

,如果a

4、收敛性分析

由Lax定理可知,当差分格式稳定时一定是收敛的,所以当a格式收敛。

1

时,向前差分2

3.2向后差分格式

1、差分格式的推导过程

1

上面构造的差分格式是显式的,即在时间层tn1上的每一个unj可以独立地根据在

时间层tn上的值unj得出。但是如果采用

u(xj,tn)u(xj,tn1)

[

un

]jO()t,

和(3.3)中的六式,就可以得到扩散方程(3.1)的另一个差分格式

n1

unjuj

nn

unj12ujuj1

a

h

2

0

(3.5)

这是微分方程(3.1)的向后差分格式。它所联系的网点分布为

:

向后差分格式与向前差分格式不同,在新时间层n上包含了3个未知量ununj1,j,

nn

unj1,因此不能由uj1直接计算出uj,这种差分格式称为隐式格式。

2、截断误差

考虑其截断误差

122uah24u

T(x,t)tu(x,t)a2xu(x,t)(x,)(,t)24

h2t12x,

1

其中(t,t),(xh,xh)。因此其截断误差为O()O(h2)。 3、稳定性分析

先把差分格式变形为

1n11n

aunaunj1(12a)ujj1uj

h

式的增长因子为

其中

2

,令ujvje

nhnjki

,并把它带入上面方程并消去公因子e

ikjh

,容易求出(3.5)

G(,k)

1kh

14asin2

2

由于a0,所以对任何网格比都有|G(,k)|1。由定理2.3知,差分格式(3.5)是稳定的。 4、收敛性分析

由Lax定理可知,当差分格式稳定时一定是收敛的,所以差分格式(3.5)收敛。

3.3 Crank-Nicloson格式

1、差分格式的推导过程

将向前差分格式与向后差分格式作算术平均,即得Crank-Nicolson格式:

1nunujj

a1n11nnn

[(ununj12ujj1)(uj12ujuj1)]0 (3.6)2

2h

[12]

Crank-Nicolson格式又叫六点对称格式,因为这种格式所联系的网点分布为

:

Crank-Nicolson格式也是隐式格式。 2、截断误差 令

Lhu

nj

1ununjj

n1n1n!nnn

auj12ujuj1uj12ujuj1[] 222hh

将截断误差

T(x,t)Lhu(xj,tn)[Lu]nj

在(xj,t

n

1)(t2

n

1

2

1

(n))处展开,则得

2

T(x,t)O(2h2)

因此其截断误差为O(2h2)。 3、稳定性分析

仍然用Fourier方法。先把差分格式变形为

nnn1n11

aunaunj12(1a)ujauj1auj12(1a)ujj1

nikjhikjh

令un,带入上面方程并消去公因子e,容易求出(3.6)式的增长因子为 jvje

kh

G(,k)

12asin2

2

12asin2

由于a0,显然对任何网格比都有|G(,k)|1,由定理2.3知,差分格式(3.6)是稳定的。 4、收敛性分析

由Lax定理可知,当差分格式稳定时一定是收敛的,所以差分格式(3.6)收敛。

3.4 Richardson格式

1、差分格式的推导过程

利用(3.3)式中的第二式和第六式可以得到Richardson格式:

11

ununjj

2

a

nn

unj12ujuj1

h2

0

(3.7)

从(3.7)式可以看出,这种格式联系到三个时间层,称为三层格式。一般地,一个多于二层的差分格式成为多层差分格式。它所联系的网点分布为:

2、截断误差

这个格式是为了提高差分方程的截断误差所设计的,利用上述方法我们可以很容易的求出其截断误差为O(2h2)。单从截断误差这一角度来考虑,这个格式是比较好的,但是以后的分析可以看到,它是没有实用价值的。 3、稳定性分析

首先将(3.7)式变形为

11nn

unun2a(unjjj12ujuj1)

讨论多层差分格式,一般先化成与其等价的二层差分方程组。Richardson差分方程的等价二层差分方程组为

n1nnnnujvj2a(uj12ujuj1),

n1n

vjuj.

nT

如果令u[un,vjj],则上面的方程组可以写成

n1j

nj

u

nj

n

2a004an

uj1000n2a

uj000n

uj1 0

设uveikjh,将其带入上式并消去公因子eikjh,可以得

济南大学毕业论文

v

n1

2kh8asin2

1

1

vn 0

因此增长矩阵为

2kh8asin

G(,k)2

1

1 0

其特征值为

1,2

1

kh224kh2

4asin(116asin)

22

2

1kh224kh2

14asin(116asin)

22

2

则有

kh

2

由定理2.2可知Richardson格式是不稳定的。这种情况的差分格式称为绝对不稳定的。

|1|14asin2

4、Richardson格式的改进

1953年Du Fort和Frankle对Richardson格式进行了修改,提出了如下的DuFort-Frankle格式[13]

1n1

unujj

2

a

n1n1n

un(uu)uj1jjj1

h2

0

这是一个无条件稳定的差分格式,但是其相容性是有条件的。在这里对DuFort-Frankle格式我们不做过多讨论。

事实证明,我们无法构造出无条件相容和无条件稳定的显示格式。

济南大学毕业论文

当x0.2时,数值解与解析解比较如表4.4

表4.4

当t0.07时,数值解与解析解比较如图4.2:

图4.2

从结果来看,向后差分格式是稳定的,它的最大优点在于稳定性、收敛性对于一切网格比都成立,网格比增大意味着工作量减少,所以隐式格式当1时的工作量比显示格式要少,然而从截断误差的角度来看,h也不能太大,否则太大,误差也就太大,所以更高精度的格式也是十分需要的。 (3)Crank-Nicloson格式为

1ununjj

an1n1n1nnn

[(u2uu)(u2uuj1jj1j1jj1)]0 22h

写成便于计算的格式为

1n11nnnununj16ujj1uj12ujuj1

(4.4)

初边值条件的格式为

结 论

首先本文给出了偏微分方程及抛物型方程的基本概念,其次给出了差分格式的基本概念,重点介绍了差分格式的收敛性及稳定性,随后针对常系数扩散方程给出了各种差分格式的推导过程及其各种性质,最后举出数值算例对给出的各种差分格式进行应用,以观察其可行性。结果证明,向前差分格式、向后差分格式、Crank-Nicloson格式都是稳定的,或者说在一定条件下是稳定的,Richardson格式是绝对不稳定的。

数值计算主要研究如何利用计算机更好的解决各种数学问题,包括连续系统离散化和离散形方程的求解,在对具体问题进行分析时,数值计算比解析解的方法更方便。几乎所有的偏微分方程都能求出数值解,数值解虽然是近似值,但对于许多方程,近似解可以足够精确,满足应用需求,编程易实现,可以分析误差的界[15]。因此,对于这个问题的研究是有很强的实用性的。

参 考 文 献

[1]黎阳.抛物型方程的有限差分解法及其在复杂磁场中的应用[D].武汉:武汉理工大学,2010:11-12.

[2]李德元,抛物型方程差分方法引论[M].北京:科学出版社,1982:80-92. [3]关治,陆金甫.数值分析基础[M].北京:高等教育出版社,1979:200–260.

[4]葛坤.一类抛物型方程的数值解法研究与应用[D].南京:南京林业大学,2011:11-26.

[5]Per Brinch Hansen.Numerical Solution of Laplace's Equation[J].Electrical Engineering and Computer Science Technical Reports,1992(4):34-46.

[6]孙雪莉.抛物型方程任意偶数阶精度的两层显示格式[D].西安:西安建筑科技大学,2009:10-17. [7] 张宝林,谷同祥,莫则尧.数值并行计算原理与方法[M].北京:国防工业出版社,1999:42–43. [8]杨丹丹.求解抛物型方程的一种有限差分并行格式[D].吉林:吉林大学数学研究所,2010:10-15. [9]陆金甫,关治.偏微分方程数值解[M].北京:清华大学出版社,2004:28-34.

[10]林鹏程.解双抛物型方程的高精度稳定差分格式[J].福州大学学报,2004(4):120–121. [11]张耀科,陆如枫.初边值问题稳定性分析的一个标注[J].工科数学,1994(1):64-78. [12]许秀娟.两类抛物型方程的有限差分解法[D].哈尔滨:哈尔滨工业大学,2009:8-21.

[13]C. Fuhse,T. Salditt.Finite-difference field calculations for two-dimensionally confined x-ray waveguides[J].Applied Optics ,July 2006,45(19): 507–511.

[14]王爱锋,曲小钢.解抛物型方程的一种隐式差分格式[J].纺织高校基础科学学报,2010(4):32-37.

[15]Adel A.S.Almarashi.Approximation Solution of Fractional Partial Differential Equation By Neural Network[J].Development Of Mathematics,2011(10):98-126.

致 谢

历时将近两个月的时间终于将这篇论文写完,在论文的写作过程中遇到了无数的困难和障碍,都在同学和老师的帮助下度过了。尤其要强烈感谢王宣欣老师,她对我进行了无私的指导和帮助,不厌其烦的帮助进行论文的修改和改进。另外,在校图书馆查找资料的时候,图书馆的老师也给我提供了很多方面的支持与帮助。在此向帮助和指导过我的各位老师表示最中心的感谢!

感谢这篇论文所涉及到的各位学者。本文引用了数位学者的研究文献,如果没有各位学者的研究成果的帮助和启发,我将很难完成本篇论文的写作。

感谢我的同学和朋友,在我写论文的过程中给予我了很多你问素材,还在论文的撰写和排版灯过程中提供热情的帮助。

由于我的学术水平有限,所写论文难免有不足之处,恳请各位老师和学友批评和指正!

27 - -

附录

追赶法程序如下:

#include

#define N 9

main()

{

double A[N][N]={{6,-1,0,0,0,0,0,0,0},{-1,6,-1,0,0,0,0,0,0},{0,-1,6,-1,0,0,0,0,0},{0,0,-1,6,-1,0,0,0,0},{0 ,0,0,-1,6,-1,0,0,0},{0,0,0,0,-1,6,-1,0,0},{0,0,0,0,0,-1,6,-1,0},{0,0,0,0,0,0,-1,6,-1},{0,0,0,0,0,0,0,-1,6}};b

[N]={0.898922,1.709851,2.353408,2.766597,2.908972,2.766597,2.353408,1.709851,0.898922},m[9];

}

28 - -int i; A[0][1]=A[0][1]/A[0][0]; for(i=1;i=0;i--) b[i]=b[i]-A[i][i+1]*b[i+1]; for(i=0;i


相关文章

  • 常系数线性差分方程的微分解法
  • [***********][***********][***********]77生的整体素质也起到了重要的作用:符合国家教育部提出的; 重在参与:重在普及:扩大受益面 >以竞赛促教改 撰写论文的能力? 查阅资料的能力:分析问题? 解 ...查看


  • 二阶变系数齐次微分方程通解的求法
  • ((高等数学研究1:?,($$" 一类二阶变系数齐次微分方程通解的求法 李永利! 摘要!!桑改莲!(河南质量工程职业学院基础部!河南平顶山!!"#$$$)作为文[%]两种情形的统一推广,给出一类二阶变系数线性微分方程的通 ...查看


  • 10-345二阶线性微分方程
  • 10.3 可降阶的高阶微分方程 一. y ( n ) = f ( x ) 型的微分方程 二. y ′′ = f ( x , y ′) 型的微分方程 三. y ′′ = f ( y , y ′) 型的微分方程 第一种类型 y( n) = f ...查看


  • 差分方程及其应用
  • 差分方程及其应用 在经济与管理及其它实际问题中,许多数据都是以等间隔时间周期统计的.例如,银行中的定期存款是按所设定的时间等间隔计息,外贸出口额按月统计,国民收入按年统计,产品的产量按月统计等等.这些量是变量,通常称这类变量为离散型变量.描 ...查看


  • 考研讲义数三经济部分
  • 第十三章 微积分在经济学中的经济应用 (数三) <考试要求> 1. 掌握导数的经济意义(含边际与弹性的概念). 2. 了解差分与差分方程及其通解与特解等概念. 3. 掌握一阶常系数线性差分方程的求解方法. 4. 会应用一阶差分方 ...查看


  • 差分方程模型
  • 第7章 差分方程模型 在第三章中我们看到,动态连续模型用微分方程方法建立,与此相应,当时间变化离散后,可以用差分方程方法建立.有些实际问题既可建立连续模型,又可建立离散模型.本章将主要讨论差分方程模型. 7.1市场经济中的蛛网模型 在自由贸 ...查看


  • 水力瞬变特征线法和隐式差分法的对比分析
  • '*(' 油气储运 () ) *年 ! ! " !!!!!! " 设计计算 !! " ! ! ! ! ! ! " 蒋仕章 摘 要严重得多- 水力瞬变特征线法和隐式差分法 的对比分析 蒋仕章 # 蒲家宁 ...查看


  • 基于有限差分法的微分方程离散化求解
  • 基于有限差分法的微分方程离散化求解 [摘要] 目前偏微分方程数值求解的方法主要有两种,即有限差分法和有限元方法.本文论述了基于有限差分法的微分方程求解,离散化过程,并对结果进行了分析. [关键词] 有限差分法 离散化 数值模拟 1. 前言 ...查看


  • 常微分方程常用数值解法
  • 第一章 绪论 1.1 引言 常微分方程是现代数学的一个重要分支,是人们解决各种实际问题的有效工具.微分方程的理论和方法从17世纪末开始发展起来,很快成了研究自然现象的强有力工具,在17到18世纪,在力学.天文.科学技术.物理中,就已借助微分 ...查看


热门内容