有限差分法求解偏微分方程

南京理工大学

课程考核论文

课程名称：高等数值分析

论文题目：有限差分法求解偏微分方程姓名：罗晨

学号：

成绩：

有限差分法求解偏微分方程

一、主要内容

1.有限差分法求解偏微分方程，偏微分方程如一般形式的一维抛物线型方程：

22(,)()u u

f x t t x

αα∂∂-=∂∂其中为常数

具体求解的偏微分方程如下：

22001

(,0)sin()(0,)(1,)00

u u x t x u x x u t u t t π⎧∂∂-=≤≤⎪∂∂⎪⎪⎪

=⎨⎪⎪

==≥⎪⎪⎩

2.推导五种差分格式、截断误差并分析其稳定性；

3.编写MATLAB 程序实现五种差分格式对偏微分方程的求解及误差分析；

4.结论及完成本次实验报告的感想。

二、推导几种差分格式的过程：

有限差分法（finite-difference methods ）是一种数值方法通过有限个微分方程近似求导从而寻求微分方程的近似解。有限差分法的基本思想是把连续的定解区域用有限个离散点构成的网格来代替；把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；把原方程和定解条件中的微商用差商来近似，积分用积分和来近似，于是原微分方程和定解条件就近似地代之以代数方程组，即有限差分方程组，解此方程组就可以得到原问题在离散点上的近似解。

推导差分方程的过程中需要用到的泰勒展开公式如下：

()2100000000()()()()()()()......()(())

1!2!!

n n n f x f x f x f x f x x x x x x x o x x n +'''=+-+-++-+-

(2-1)

求解区域的网格划分步长参数如下：

11k k k k

t t x x h τ

++-=⎧⎨

-=⎩ (2-2)

古典显格式

2.1.1 古典显格式的推导

由泰勒展开公式将(,)u x t 对时间展开得

2,(,)(,)(

)()(())i i k i k k k u

u x t u x t t t o t t t

∂=+-+-∂ (2-3)

当1k t t +=时有

21,112,(,)(,)(

)()(())(,)()()

i k i k i k k k k k i k i k u

u x t u x t t t o t t t

u

u x t o t

ττ+++∂=+-+-∂∂=+⋅+∂ (2-4)

得到对时间的一阶偏导数

1,(,)(,)()=()i k i k i k u x t u x t u

o t ττ+-∂+∂ (2-5) 由泰勒展开公式将(,)u x t 对位置展开得

223,,21(,)(,)()()()()(())2!k i k i k i i k i i u u

u x t u x t x x x x o x x x x

∂∂=+-+-+-∂∂

(2-6)

当11i i x x x x +-==和时，代入式(2-6)得

2231,1,1122

231,1,1121(,)(,)()()()()(())2!1(,)(,)()()()()(())

2!i k i k i k i i i k i i i i i k i k i k i i i k i i i i

u u

u x t u x t x x x x o x x x x u u u x t u x t x x x x o x x x x ++++----⎧∂∂=+-+-+-⎪⎪∂∂⎨∂∂⎪=+-+-+-⎪∂∂⎩

(2-7)

因为1k k x x h +-=，代入上式得

2231,,22

231,,21(,)(,)()()()2!1(,)(,)()()()

2!i k i k i k i k i k i k i k i k

u u

u x t u x t h h o h x x

u u u x t u x t h h o h x x +-⎧∂∂=+⋅+⋅+⎪⎪∂∂⎨∂∂⎪=-⋅+⋅+⎪∂∂⎩

(2-8)

得到对位置的二阶偏导数

22

11,22

(,)2(,)(,)()()i k i k i k i k u x t u x t u x t u o h x h

+--+∂=+∂ (2-9)

将式(2-5)、(2-9)代入一般形式的抛物线型偏微分方程得

21112(,)(,)

(,)2(,)(,)(,)()i k i k i k i k i k i k u x t u x t u x t u x t u x t f x t o h h αττ

++---+⎡⎤-=++⎢⎥⎣⎦

(2-10

)

为了方便我们可以将式(2-10)写成

1112

2k k

k k k k i i i i i i u u u u u f h ατ

++-⎡⎤--+-=⎢⎥⎣⎦

(2-11) ()11

12

2k k k k k k i i i i i i u u u

u u f h

τα

τ++---

-+=

(2-12)

最后得到古典显格式的差分格式为

()111(12)k k k k k i i i i i u ra u r u u f ατ++-=-+++

(2-13)

2

r h τ

=

其中，古典显格式的差分格式的截断误差是2()o h τ+。

2.1.2 古典显格式稳定性分析

古典显格式(2-13)写成矩阵形式为

()112k k k h h

h u ra I raC u f τ+=-++⎡⎤⎣⎦ (2-14)

12212,(,,......,,)k k k k k

h N N r u u u u u h τ

--=

=其中。

(1)(1)

01

010*********N N C -⨯-⎡⎤⎢⎥⎢⎥

⎢⎥=⎢⎥⎢⎥⎢⎥⎣⎦L L M

M L 上面的C 矩阵的特征值是：2cos()

1,2,......,1C j h j N λπ==-

()12H ra I raC =-+

()()()2

12=122cos()

121cos()14sin 1,2,......,1

2H j C ra ra ra ra j h ra j h j h

ra j N λλπππ=-+-+=--=-=-

(2-15)