当前位置:文档之家› 第二章 有限差分基础

第二章 有限差分基础

第二章 有限差分基础
第二章 有限差分基础

第2章 有限差分基础(finite difference method ,FDM )

1.1 偏微分方程的一般形式

()()φφφρρφq x x x u t j j j j +???

? ????Γ

??=??+?? ( 2-1 ) 2.1 网格划分

一般有限差分采用结构化网格划分。即节点对应于当地坐标系统的原点。它的轴同网格线一致。即两个同一族的网格线不相交,且没对网格线对应不同的族。每一个节点可用唯一的一个坐标表示,如(ξ1, ξ2)。网格线能用ξ1=const, ξ2=const 表示。 1D 2D

有限差分法就是要将节点上的偏微分方程用相邻点上的值表示,变成线性代数方程式。

i-1

i

i+1

N

1 N j

j+1 j-1

j 1

1

i-1 i i+1

为流体力学的微分方程的数值求解方法之代表。必要条件:

? 连续领域内的分配有限的网格

?

领域内的函数分布可用网格点上的值代表

1. 计算分子(computational molecule )

5点计算分子 15点计算分子 7点计算分子

这些节点又称为计算分子。方程的个数应与未知数相同,即每个节点有一个方程式。

E

W

N

E

T N

2. T aylor 展开

例如:一维时间变量φ的理论解为φ(t,x),它在离散点上的值为投影(projection )的近似值为:()x i t n ?Λ,φ,

n: 时间的step 数

i:空间的step 数

为了求得此近似解,需对微分方程进行差分近似。利用T aylor 展开可得到几个差分表示形式,仅考虑空间依存问题:

在?x 很小时,位置j ?x 内的物理量φ用φj 来表示,则位置(j+1)?x 上的值φj+1表示为:

()???+???+???+=+i

i i i x x x x 22

2121φφφφ ( 2-2 )

(j-1)?x 上的值φj+1表示为:

()???+???+???-=-i

i i i x x x x 222121φφφφ ( 2-3 ) 2.2 基本差分格式

1. 一阶导数(first derivative )的近似

()x

u or

x

????φρφ

( 2-4 ) i. 向前差分(forward difference ,FDS)

利用( 1 ) 式,可得到1阶微分的向前差分形式:

)(1x O x

x j j j ?+?-=??+φφφ ( 2-5 ) ii. 向后差分(backward difference ,BDS)

利用( 1 ) 式,可得到1阶微分的向后差分形式

)(1x O x

x j j j ?+?-=??-φφφ ( 2-6 ) iii. 中心差分(central difference ,CDS )

(1)-(2) 得1阶微分的2 次精度中心差分法:

)(2211x O x

x j j j ?+?-=??-+φφφ

( 2-7 ) iv. 上风法、迎风法(upwind difference, UDS ) 与速度有关的微分

()???????<-->--≈??++--;

0if

,;0if ,1111

u x x u u x x u x u i

i i i i i i i φφρφφρφρ ( 2-8 )

2. 二阶导数的近似

i. 中心差分(central difference ,CDS )

利用(j ±1/2)?x 的T ayor 展开,可得过且1阶微分的2次精度的向前向后差分形式:

1

1112

121---+++--=??--=??i i i i i i

i i

i i x x x x x x φ

φφφ

φφ ( 2-9 ) 将上二式相减,得2阶微分的差分方程式中心差分:(?x 相当)

)(22

2

1122x O x

x i i i i ?+?+-=??-+φφφφ ( 2-10 ) 其它还有各种形式。但最常用的是中心差分。对于扩散项常采用

()()

1121

1

1

111121

2

1

2

1

2

1

21-+---+++-+-+---Γ---Γ≈-???

????Γ-??? ????Γ≈?

???????? ?

???Γ??i i i i i i i i i i

i i i i i i i x x x x x x x x x x x x φφφφφφ ( 2-11 )

2.3 数值误差

1. 截断误差(truncation error )

例如:热传递方程

022=??-??t

T

t T α ( 2-12 ) 时间向前差分,空间中心差分得

()

()

02112

1=+-?-

?--++n j n j n

j n

j n j T T T

x t

T T α

( 2-13 )

离散方程同精确方程的误差(i.e. 微分算子L(T)和L ’(T n j ) ) ,称为截断误差。将1

+n j u 进行

时间方向的T ayor 展开和n j u 1±进行

空间方向的T ayor 展开即, () (2122)

21+???+???+=+n

j

n

j n

j n j t u t t u t u u ( 2-14 )

()()() (241)

612144

43

3

322

21n

j

n

j

n

j n

j n

j n j x u x x

u

x x u x x u x

u u ???+

???±???+???±=± ( 2-15 )

因此,此方程的截断误差为

()()()???

?????+???+???-+?+--?-=??-??-++ (12224)

422222

1

1122

n j n j n

j n j n j n

j

n j n

j

n

j x u

x x u x x u u u t

u u x u

x u αα

α ( 2-16 )

截断误差的第一项为O(?t, ?x 2),为时间一次精度,空间二次精度。

2. 舍入误差

计算机的精度引起的误差称为舍入误差( r ound-off error 圆误差)。通常CFD 忽略舍入误差。这是因为,几乎所有的CFD 程序的精度都在3阶以下,舍入精度都比它们下的缘故。舍入精度同计算精度有关。要注意的是,超级计算机的利用。目前由于存储容量的限制,许多使用的是单精度(32bit )的机器。但是使用超级计算机,往往采用特别小的格子,这时需要注意计算机的精度问题。

通常,定常解的收敛以截断误差为3阶为标准。但必须注意单精度同双精度(64bit )的误差相差很大。特别是,在使用广泛应用的一般坐标系上的守恒表达形式的条件下,要达到实质性的收敛解,需采用双精度,可使误差充分减少,但此时如采用单精度,可能误差都达不到3阶。因此,要特别注意收敛的判断,尽可能使用64bit 的计算。

3. 离散误差

截断误差+数值的边界条件?

指在无舍入误差的条件下,差分方程的解(即差分方程的严密解)和微分方程的精确解之差。

i.

耗散误差(dissipative error )

由数值扩散引起的离散误差部分 ii.

分散误差(dispersive error )

由相位误差引起的离散误差。 iii.

延迟相位误差

波的位置在真波后发生 iv.

前进相位误差

波的位置在真波前发生

v.

正确解 耗散误差 分散误差

x

x

x

x

2.4 显式、隐式和半隐式求解格式

对于时间推进问题,有显式、半隐式和隐式三种基本格式。

()00

>=+c cu u x t ( 2-17 )

t u u t u n

i

n i i

?-=

??? ????+1 ( 2-18 ) 空间CDS

1. 显式格式

x u u x u n

i n i ?-=??

? ????-+11 ( 2-19 )

∑=+nb

n

nb

nb n i u a u 1 ( 2-20 ) 2. 隐式格式

x u u x u n i n i ?-=??

? ????+-++1

111 ( 2-21 )

∑=+nb

n i i n nb nb u a u a 1 ( 2-22 ) 3. 半隐式

主要用于非线性问题

x u u x u n i n i ?-=

??

? ????+-+1

1

1 ( 2-23 ) ∑∑=++nb

nb

n

i n n nb n u a u

a

nb nb

1

1 ( 2-24 )

2.5 稳定性条件

1. 相容性

指差分方程接近微分方程的程度。时间方向和空间方向的分割(?t, ?x )变小,截断误差逐渐消失,即差分方程接近原来的微分方程,则称为该差分形式与偏微分方程的相容(consistent )。大部分差分方法满足此条件。如上式中只要满足0→??x t 的条件,则consistent. 2. 稳定性

对于发展性问题,基本可以满足数值的稳定性这一条件。稳定性是指,计算的一步步进行,不管什么原因引起的误差都不会使其成长。 3. 收敛性

指差分解接近真实解的程度。格子变小变细后,其近似解接近于元微分方程解。离散化误差趋近零。即只要使用足够的点,不管要求多小的误差都能达到。注意:不是指收敛于定常流的解问题。

4. Lax 的等同定理(Lax ’s Equivalence Theorem )

收敛的格式一定是稳定而且相容的格式。这对线性的标量的时间发展问题成立。即:

稳定性+相容性为收敛性的必要条件。

即要得到正确的解,必须是用稳定且相容的方法。 5. V on Neumann 稳定性条件

适用于初值问题。基本思想可从小扰动法的角度去考虑。假定某一时刻引入了一个误差矢量,即扰动,如果这一扰动的强度(或振幅)随时间不断增大,则这一格式是不稳定的。反之,扰动的振幅如果是衰减的,或保持不变,格式是稳定的。令振幅比定义为

()()t t t ψψμ/?+= ( 2-25 )

初值稳定性条件: ??1≤μ 其中:中性稳定 1

CFL 条件:Cournat 条件:显式格式的稳定性条件

1≤??x

t

c 2.6 线性标量方程式的物理意义

流体力学的基本现象: 对流和扩散

1. 输运方程式或波动方程式

例1:以对流为代表的典型方程式:输运方程式或波动方程式:为双曲线问题(hyperbolic problem),必须有初期条件。

),0(0∞<<∞->=??+??x t x

u c t u ( 2-26 ) c 为正整数。使用下列初始条件,成为时间发展问题:

)0(),0(21<=≥=x u u x u u ( 2-27 )

它的一种差分格式为:

01

1=?-+?--+x

u u c

t

u u n

j n j n j

n j ( 2-28 )

即:

()

n

j n j n j n j u u x t c u u 11-+-??

?

????-= ( 2-29 ) 当t x c ??=时,n j n j u u 11

-+=。

即在时间步长n 时,1-=j x x 处发生的现象,将在时间步

长为n+1时,在j x x =处发生。实际上,式(11)的解为:

)(),(ct x u t x u -= ( 2-30 )

其解在称之为特性线的ct x =上。只有在t x c ??=时其差分解才与解析解相同。 物理意义

?

波动方程表示上游发生的事以c 的速度向下游发展。 c 为正时:左向右传播;c 为负时:右向左传播;

? 传播速度c 为特性线的斜率; ?

该现象可用上风之概念描述:

式(11)所示形式表示j 位置上的值由j-1上的值即可确定。即为情报传播方向的差分形式。这种差分形式称为上风差分法(up-wind method )。本情况下为向后差分的上风法。这样正确反映情报的传播方向的差分为好的差分方法。如采用差的上风法(应称为下风法),即采用向前差分:

11

=?-+?-++x

u u c

t

u u n

j

n j n j

n j ( 2-31 )

()

n

j n j n j

n j u u x t c u u -??? ????-=++11 ( 2-32 )

当t x c ??=时:

n

j n j n j u u u 1

12++-= ( 2-33 )

实际计算会发现,其解完全不一样,其原因以后再论述。因此,差分化所进行的数值仿真的重要点是必须在物理上包含必要的情报。因此:

? 中心差分因包含着两边的信息,一般认为比较安全; ?

时间微分因注意使用单边的差分,因时间不会倒流。

2. 扩散方程式或热传导方程

例2:扩散方程式(diffusion equation ),又可称为热传导方程(heat equation ), 为抛物性方程(parabolic equation )

022=??-??t u

t u α ( 2-34 )

需要初始条件和边界条件。典型的离散方法为:

()

22

1

11=?+--?--++x u u u t

u u n

j n j n j n

j

n j α

( 2-35 )

()

()

n

j n j n

j n j n j u u u

x t

u u 1

12

1

2-+++-??+=α

( 2-36 )

考虑其简单形式:

()

21

2

=

??x t α

( 2-37 )

则:

()

n j n

j n j

u u u 11121+-++= ( 2-38 )

这方程代表了扩散方程的特性。 物理特性:

? u j 的值由左右两个值的上一时刻的平均值来确定; ? 扩散使现象向平均化方向发展; ?

因此,扩散问题宜用中心差分方法。

差分方法的重要注意点:

遵循物理法则进行离散近似,否则不能保证得到正确解。

完整版有限差分方法概述.doc

有限差分法( Finite Difference Method,简称FDM)是数值方法中最经典的方法,也是计算机数值模拟最早采用的方法,至今仍被广泛运用。该方法将求解域划分为差分网格,用有限个网格节点代替连续的求解域。有限差分法以Taylor级数展开等方法,把控制方程中的导数用网格节点上的函数值的差商代替进行离散,从而建立以网格节点上的值为未知数的代数方程组。该方法是一种直接将微分问题变为代数问题的近似数值解法,数学概念直观,表达简单,是发展较 早且比较成熟的数值方法。对于有限差分格式,从格式的精度来划分,有一阶格式、二阶格式和高阶格式。从差分的空间形式来考虑,可分 为中心格式和逆风格式。考虑时间因子的影响,差分格式还可以分为显格式、隐格式、显隐交替格式等。目前常见的差分格式,主要是上 述几种形式的组合,不同的组合构成不同的差分格式。差分方法主要适用于有结构网格,网格的步长一般根据实际地形的情况和柯朗稳定条件来决定。 构造差分的方法有多种形式,目前主要采用的是泰勒级数展开方法。其基本的差分表达式主要有三种形式:一阶向前差分、一阶向后 差分、一阶中心差分和二阶中心差分等,其中前两种格式为一阶计算精度,后两种格式为二阶计算精度。通过对时间和空间这几种不同差分格式的组合,可以组合成不同的差分计算格式。 下面我们从有限差分方法的基本思想、技术要点、应用步骤三个方面来深入了解一下有限差分方法。 1.基本思想 有限差分算法的基本思想是把连续的定解区域用有限个离散点 构成的网格来代替,这些离散点称作网格的节点;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解。然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解。在采用数值计算方法求解偏微分方程时,再将每一处导数由有限差分近似公式替代,从而把求解偏微分方程的问题转换成求解代数方程的问题,即 所谓的有限差分法。 2.技术要点 如何根据问题的特点将定解区域作网格剖分;如何把原微分

对流扩散方程有限差分方法.

对流扩散方程有限差分方法 求解对流扩散方程的差分格式有很多种,在本节中将介绍以下3种有限差分格式:中心差分格式、Samarskii 格式、Crank-Nicolson 型隐式差分格式。 3.1 中心差分格式 时间导数用向前差商、空间导数用中心差商来逼近,那么就得到了(1)式的中心差分格式]6[ 2 1 11 1122h u u u v h u u a u u n j n j n j n j n j n j n j -+-+++-=-+-τ (3) 若令 h a τ λ=,2h v τ μ=,则(3)式可改写为 )2()(2 111111 n j n j n j n j n j n j n j u u u u u u u -+-+++-+--=μλ (4) 从上式我们看到,在新的时间层1+n 上只包含了一个未知量1 +n j u ,它可以由时间层n 上的值n j u 1-,n j u ,n j u 1+直接计算出来。因此,中心差分格式是求解对 流扩散方程的显示格式。 假定),(t x u 是定解问题的充分光滑的解,将1 +n j u ,n j u 1+,n j u 1-分别在),(n j t x 处 进行Taylor 展开: )(),(),(211ττO t u t x u t x u u n j n j n j n j +??? ?????+==++ )(2),(),(3 22211 h O x u h x u h t x u t x u u n j n j n j n j n j +????????+????????+==++ )(2),(),(3 22211 h O x u h x u h t x u t x u u n j n j n j n j n j +????????+????????-==-- 代入(4)式,有 2 111 1122),(h u u u v h u u a u u t x T n j n j n j n j n j n j n j n j -+-+++---+-= τ )()()(2222 h O v x u v h O a x u a O t u n j n j n j ?-????????-?+????????++????????=τ )()()(222h O v a O x u v x u a t u n j n j n j ?-++????????-??? ?????+????????=τ

有限差分方法计算欧式期权价格

假设当前股票价格为50美元,股票价格波动率sigma=0.3;以该股票为标的资产的欧式看跌期权的执行价格为50美元,期权有效期为5个月;市场上的无风险利率为10%。利用显示差分格式为该期权进行定价。 %%% 显示法求解欧式看跌期权%%% s0=50; %股价 k=50; %执行价 r=0.1; %无风险利率 T=5/12; %存续期 sigma=0.3; %股票波动率 Smax=100; %确定股票价格最大价格 ds=2; %确定股价离散步长 dt=5/1200; %确定时间离散步长 M=round(Smax/ds); %计算股价离散步数,对Smax/ds取整运算 ds=Smax/M; %计算股价离散实际步长 N=round(T/dt); %计算时间离散步数 dt=T/N; %计算时间离散实际步长 matval=zeros(M+1,N+1); vets=linspace(0,Smax,M+1); %将区间[0,Smax]分成M段 veti=0:N; vetj=0:M; %建立偏微分方程边界条件 matval(:,N+1)=max(k-vets,0); matval(1,:)=k*exp(-r*dt*(N-veti)); matval(M+1,:)=0; %确定叠代矩阵系数 a=0.5*dt*(sigma^2*vetj-r).*vetj; b=1-dt*(sigma^2*vetj.^2+r); c=0.5*dt*(sigma^2*vetj+r).*vetj; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%% L=zeros(M-1,M+1); for i=2:M %%建立递推关系 L(i-1,i-1)=a(i); L(i-1,i)=b(i); L(i-1,i+1)=c(i); end for i=N:-1:1 matval(2:M,i)=L*matval(:,i+1); end matval %寻找期权价格进行插值。 Jdown=floor(s0/ds);

有限差分法

班级:通信13-4 姓名: 学号: 指导教师:徐维 成绩: 电子与信息工程学院 信息与通信工程系

求解金属槽的电位分布 1.实验原理 利用有限差分法和matlab软件解决电位在金属槽中的分布。 有限差分法基本思想是把连续的定解区域用有限个离散点构成的网格来代替,这些离散点称作网格的节点;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解.然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解.在采用数值计算方法求解偏微分方程时,若将每一处导数由有限差分近似公式替代,从而把求解偏微分方程的问题转换成求解代数方程的问题。 2.有限差分法 方程的定解问题就是在满足某些定解条件下求微分方程的解。在空间区域的边界上要满足的定解条件称为边值条件。如果问题与时间有关,在初始时刻所要满足的定解条件,称为初值条件。不含时间而只带边值条件的定解问题,称为边值问题。与时间有关而只带初值条件的定解问题,称为初值问题。同时带有两种定解条件的问题,称为初值边值混合问题。 定解问题往往不具有解析解,或者其解析解不易计算。所以要采用可行的数值解法。有限差分方法就是一种数值解法,它的基本思想是先把问题的定义域进行网格剖分,然后在网格点上,按适当的数值微分公式把定解问题中的微商换成差商,从而把原问题离散化为差分格式,进而求出数值解。此外,还要研究差分格式的解的存在性和唯一性、解的求法、解法的数值稳定性、差分格式的解与原定解问题的真解的误差估计、差分格式的解当网格大小趋于零时是否趋于真解(即收敛性),等等。 有限差分方法具有简单、灵活以及通用性强等特点,容易在计算机上实现。 2.1有限差分法原理

有限差分法

利用有限差分法分析电磁场边界问题 在一个电磁系统中,电场和磁场的计算对于完成该系统的有效设计师极端重要的。例如,在系统中,用一种绝缘材料是导体相互隔离是,就要保证电场强度低于绝缘介质的击穿强度。在磁力开关中,所要求的磁场强弱,应能产生足够大的力来驱动开关。在发射系统中进行天线的有效设计时,关于天线周围介质中电磁场分布的知识显然有实质性的意义。 为了分析电磁场,我们可以从问题所涉及的数学公式入手。依据电磁系统的特性,拉普拉斯方程和泊松方程只能适合于描述静态和准静态(低频)运行条件下的情况。但是,在高频应用中,则必须在时域或频域中求解波动方程,以做到准确地预测电场和磁场,在任何情况下,满足边界条件的一个或多个偏微分方程的解,因此,计算电池系统内部和周围的电场和磁场都是必要的。 对电磁场理论而言,计算电磁场可以为其研究提供进行复杂的数值及解析运算的方法,手段和计算结果;而电磁场理论则为计算电磁场问题提供了电磁规律,数学方程,进而验证计算结果。常用的计算电磁场边值问题的方法主要有两大类,其每一类又包含若干种方法,第一类是解析法;第二类是数值法。对于那些具有最简单的边界条件和几何形状规则的(如矩形、圆形等)问题,可用分离变量法和镜像法求电磁场边值问题的解析解(精确解),但是在许多实际问题中往往由于边界条件过于复杂而无法求得解析解。在这种情况下,一般借助于数值法求解电磁场的数值解。 有限差分法,微分方程和积分微分方程数值解的方法。基本思想是把连续的定解区域用有限个离散点构成的网络来代替,这些离散点称作网格的节点;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解。然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解。 差分运算的基本概念: 有限差分法是指用差分来近似取代微分,从而将微分方程离散成为差分方程组。于是求解边值问题即转换成为求解矩阵方程[5]。 对单元函数 ()x f而言,取变量x的一个增量x?=h,则函数()x f的增量可以表示为 ()x f? = ()h x f+-()x f 称为函数()x f 的差分或一阶差分。函数增量还经常表示为 ()x f? = ? ? ? ? ? + 2 h x f - ? ? ? ? ? - 2 h x f

差分方法

一、差分方法 1.1 导数的差分公式 在x 附近对()f x 展开,由泰勒展开公式 ()()()f x h f x f x h '+≈+ 得到前差公式为 ()() ()f x h f x f x h +-'= 同理也可以得到后差公式 ()() ()f x f x h f x h --'= 由后差分公式可以得到二阶导数的差分公式为 2 ()()()2()() ()f x h f x f x h f x f x h f x h h ''+-+-+-''= = 叫中心差分公式。 利用这些公式可以将微分方程写成差分方程。 1.2 热传导方程的差分公式 热传导方程是 2t xx u a u = 可以写成差分形式 2 2 (,)(,)(,)2(,)(,) ()u x t t u x t u x x t u x t u x x t a t x +?-+?-+-?≈?? 即 []2 2 (,)(,)(,)2(,)(,)()t u x t t u x t a u x x t u x t u x x t x ?+?≈+ +?-+-?? 令 ,,0,1,2,...,1x i x t i t i n =?=?=- 上式可以写为(显示格式) []2 2 (,1)(,)(1,)2(,)(1,)()t u i j u i j a u i j u i j u i j x ?+=+ +-+-? 可以证明,上式的稳定条件为 2 2 ()2x t a ??≤,即 221()2t a x ?≤? 稳定且非振荡的条件为

22 1 ()4 t a x ?≤? 截断误差为 2((),)O x t ?? 另一种格式为 2 2 (,)(,)(,)2(,)(,) ()u x t t u x t u x x t t u x t t u x x t t a t x +?-+?+?-+?+-?+?≈?? 即 22 22()()(,1,1)2(,1)(1,1)(,)x x u i j u i j u i j u i j a t a t ????-++--++++=-????? ? 该式称为隐式格式。对任何步长都是恒稳定的。在t ?上取值的唯一限制是,要将截断误差 保持在合理的程度上从而节约计算时间。 截断误差为 2((),)O x t ??。 二、一维热传导方问题 2.1 无限长细杆的热传导 无限长细杆的热传导的定解问题是 2(,0)()t xx u a u u x x ??=? =? 利用Fourier 变换求得问题的解是 2 2()4(,)()x a t u x t d ξ?ξξ--+∞ -∞?? =???? 其中取初始温度分布如下: 1,01()0,0,1x x x x ?≤≤?=? <>? 这是在区间0—1之间高度为1的一个矩形脉冲,于是得到 2 (,)u x t ξ=? 可以用图1所示的瀑布图来表示稳定随时间与空间的变化。 从图中可以看到,在开始时,温度分布是原点附近的一个脉冲状得分布,随着时间的增加,热量向两边传播,形成一个平缓的波包,不难想象如果时间足够长,最终杆上的温度会全

有限差分法

有限差分法 有限差分法(Finite Differential Method, FDM ) 什么是有限差分法 有限差分法是指用泰勒技术展开式将变量的导数写成变量,在不同时间或空间点值的差分形式的方法。 按时间步长和空间步长将时间和空间区域剖分成若干网格,用未知函数在网格结(节)点上的值所构成的差分近似代替所用偏微分方程中出现的各阶导数,从而把表示变量连续变化关系的偏微分方程离散为有限个代数方程,然后解此线性代数方程组,以求出溶质在各网格结(节)点上不同时刻的浓度。 有限差分法的基本步骤 (1)剖分渗流区,确定离散点。将所研究的水动力弥散区域按某种几何形状(如矩形、任意多边形等)剖分成网络系统。 (2)建立水动力弥散问题的差分方程组。 (3)求解差分方程组。采用各种迭代法,如点逐次超松驰方法(SOR)、线逐次超松驰方法(LSOR)、迭代的交替方向隐式方法(IADI)及强隐式方法(SID)等。 (1) 现在分别对时间(从0时刻到到期日)和股票价格(S max )为可达到的足够高的股票价格) 进行分割,即\triangle S=S_{max}/M,\triangle T/N,这样就分别有N+1个时间段和M+1个股票价格,建立如图(所示的坐标方格,将定解区域网格化,坐标方格上的点(i,j )对应时刻 和股票价格 ,用变量f i ,j 表示(i,j )点的期权价格。 2.建立差分格式 (1)内含的有限差分方法 其步骤可分为以下几步:

(1)求前向差分近似:(2) 后向差分格式:(3) 将(2),(3)式平均可更加对称地求出的近似,即 (4) (2)求用前向差分近似: (5) (3)求 (6) (4)将(4),(5),(6)式代入(1)式可得到内含有限差分公式: + b j f i,j?c j f i,j + 1 = f i + 1,j(7) a j f i,j? 1 其中: i=0,1,…,N-1。j=0,1…,M-1 针对看跌期权和看涨期权可分别求出方程的边界条件: 看跌期权: 看涨期权: (5)利用边界条件和(7)式可以给出M-1个联立方程组: + b j f N? 1,j + c j f N? 1,j + 1j=1,2…,M-1 a j f N? 1,j? 1

有限差分法及其应用

有限差分法及其应用 1有限差分法简介 有限差分法(FDM)是计算机数值模拟最早采用的方法,至今仍被广泛运用。该方程将解域划分为差分网格,用有限个网络节点代替连续的求解域。有限差分法通过泰勒级数展开等方法,把控制方程中的导数用网格节点上的函数值得差商代替进行离散,从而建立以网格节点上的值为未知数的代数方程组。该方法是一种直接将微分问题变为代数问题的近似值解法,数学概念直观,表达简单,是发展较早且比较成熟的数值方法。 2有限差分法的数学基础 有限差分法的数学基础是用差分代替微分,用差商代替微商而用差商代替微商的意义是用函数在某区域内的平均变化率来代替函数的真是变化率。而根据泰勒级数展开可以看出,用差商代替微商必然会带来阶段误差,相应的用差分方程代替微分方程也会带来误差,因此,在应用有限差分法进行计算的时候,必须注意差分方程的形式,建立方法及由此产生的误差。 3有限差分解题基本步骤 有限差分法的主要解题步骤如下: 1)建立微分方程 根据问题的性质选择计算区域,建立微分方程式,写出初始条件和边界条件。 2)构建差分格式 首先对求解域进行离散化,确定计算节点,选择网格布局,差分形式和步长;然后以有限差分代替无线微分,以差商代替微商,以差分方程代替微分方程及边界条件。 3)求解差分方程 差分方程通常是一组数量较多的线性代数方程,其求解方法主要包括两种:精确法和近似法。其中精确法又称直接发,主要包括矩阵法,高斯消元法及主元素消元法等;近似法又称间接法,以迭代法为主,主要包括直接迭代法,间接迭代法以及超松弛迭代法。4)精度分析和检验 对所得到的数值进行精度与收敛性分析和检验。 4商用有限差分软件简介 商用有限差分软件主要包括FLAC、UDEC/3DEC和PFC程序,其中,FLAC是一个基于显式有限差分法的连续介质程序,主要用来进行土质、岩石和其他材料的三维结构受力特性模拟和塑性流动分析;UDEC/3DEC是针对岩体不连续问题开发,用于模拟非连续介质在静,动态载荷作用下的反应;PFC是利用显式差分算法和离散元理论开发的微、细观力学程序,它是从介质的基本粒子结构的角度考虑介质的基本力学特性,并认为给定介质在不同应力条件下的基本特征主要取决于粒子之间接粗状态的变化,适用于研究粒状集合体的破裂和破裂发展问题,以及颗粒的流动(大位移)问题。

有限差分法

有限差分法有限差分法 finite difference method 微分方程和积分微分方程数值解的方法。基本思想是把连续的定解区域用有限个离散点构成的网格来代替,这些离散点称作网格的节点;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解。然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解。 有限差分法的主要内容包括:如何根据问题的特点将定解区域作网格剖分;如何把原微分方程离散化为差分方程组以及如何解此代数方程组。此外为了保证计算过程的可行和计算结果的正确,还需从理论上分析差分方程组的性态,包括解的唯一性、存在性和差分格式的相容性、收敛性和稳定性。对于一个微分方程建立的各种差分格式,为了有实用意义,一个基本要求是它们能够任意逼近微分方程,这就是相容性要求。另外,一个差分格式是否有用,最终要看差分方程的精确解能否任意逼近微分方程的解,这就是收敛性的概念。此外,还有一个重要的概念必须考虑,即差分格式的稳定性。因为差分格式的计算过程是逐层推进的,在计算第n+1层的近似值时要用到第n层的近似值,直到与初始值有关。前面各层若有舍入误差,必然影响到后面各层的值,如果误差的影响越来越大,以致差分格式的精确解的面貌完全被掩盖,这种格式是不稳定的,相反如果误差的传播是可以控制的,就认为格式是稳定的。只有在这种情形,差分格式在实际计算中的近似解才可能任意逼近差分方程的精确解。关于差分格式的构造一般有以下3种方法。最常用的方法是数值微分法,比如用差商代替微商等。另一方法叫积分插值法,因为在实际问题中得出的微分方程常常反映物理上的某种守恒原理,一般可以通过积分形式来表示。此外还可以用待定系数法构造一些精度较高的差分格式。 有限差分方法(FDM)是计算机数值模拟最早采用的方法,至今仍被广泛

隐式有限差分编程

#include #include #define N 11 void main() { int i,k,j; long double x=1000,T=1000,u=0.001; long double h[N][N]={0},f[N]={0},Z[N]={0}; long double g[N]={0},m[N]={0},q[N]={0}; long double a,b,c,t; for(i=0,k=0;k0) /*解三对角线方程的算法第一步*/ { printf("请输入时间步长t的值:"); scanf("%lf",&t); if(t<=0) { printf(" 输入有误! \n"); exit(0); } printf("\n"); a=T; c=a; b=-(2*T+(u*x*x)/t); g[1]=c/b; m[1]=a; q[1]=b; for(i=2;i0;i--) {

for(j=1;j

二维问题的有限差分方法

西北农林科技大学实习报告 学院:理学院 专业年级:信计061 姓名:袁金龙 学号:15206012 课程:微分方程数值解 报告日期:2008-12-3 实习二、二维问题的有限差分方法 一) 实习问题: 二维经典初边值问题: 2 22 2,01(,0),01(0,)(1,)0,01x u u te t t t u x x x u t u t t ???=+<≤?????=<

研究有限差分格式稳定性的其他方法 - 报告

2015 年秋季学期研究生课程考核 (读书报告、研究报告) 考核科目:偏微分方程数值解法 学生所在院(系):理学院数学系 学生所在学科:数学 学生姓名:H i t e r 学号:1X S012000 学生类别: 考核结果阅卷人

研究有限差分格式稳定性的其他方法 摘要 偏微分方程的求解一直是大家比较关心的一个问题,而有限差分格式则是求解偏微分方程时常用并且有效的一个方法。因此,研究有限差分格式的性质就显得尤为重要。在课上我们已经跟着老师学习了运用Fourier方法研究有限差分格式的稳定性,但是在很多研究有限差分格式稳定性的问题中仅仅会用Fourier方法是不够的,所以在本篇论文中,将会介绍其他三种常用的研究有限差分格式稳定性的方法,分别是:Hirt启示型方法、直接方法(或称矩阵方法)和能量不等式方法。 关键字:偏微分方程;有限差分格式;稳定性 Abstract The solution of partial differential equations has been more concerned with a problem, and the finite difference scheme is a common and effective method for solving partial differential equations. Therefore, it is very important to study the character of the finite difference scheme. We have followed the teacher to learn the use of Fourier method of finite difference scheme stability, but in a lot of research on the stability of finite difference scheme is only used Fourier method is not enough, so in this paper, will introduce the other three kinds of commonly used in the study of finite difference scheme stability method, respectively is: Hirt enlightenment method, direct method (or matrix method) and energy inequality method. Key words: partial differential equation; finite difference scheme; stability 1 前言 微分方程的定解问题就是在满足某些定解条件下求微分方程的解。在空间区域的边界上要满足的定解条件称为边值条件。如果问题与时间有关,在初始时刻所要满足的定解条件,称为初值条件。不含时间而只带边值条件的定解问题,称为边值问题。与时间有关而只带初值条件的定解问题,称为初值问题。同时带有两种定解条件的问题,称为初值边值混合问题。定解问题往往不具有解析解,或者其解析解不易计算。所以要采用可行的数值解法。有限差分方法就是一种数值解法,它的基本思想是先把问题的定义域进行网格剖分,然后在网格点上,按适当的数值微分公式把定解问题中的微商换成差商,从而把原问题离散化为差分格式,进而求出数值解。此外,还要研究差分格式的解的存在性和唯一性、解的求法、解法的数值稳定性、差分格式的解与原定解问题的真解的误差估计、差分格式的解当网格大小趋于零时是否趋于真解(即收敛性),等等。有限差分方法具有简单、灵活以及通用性强等特点,容易在计算机上实现。在课上我们已经跟着老师学习了运用Fourier方法研究有限差分格式的稳定性,但是在很多研究有限差分格式稳定性的问题中仅仅会用Fourier方法是不够的,所以在本篇论文中,将会介绍其他三种常用的研究有限差分格式稳定性的方法,分别是:Hirt 启示型方法、直接方法和能量不等式方法。 2 Hirt启示性方法 2.1 方法概述 Hirt启示性方法是一种近似分析方法。主要是把差分格式在某确定点上作泰勒级数近似

多期双重差分法,政策实施时间不同的处理方法

多期双重差分法,政策实施时间不同的处理方法 今天,计量经济圈主要给圈友引荐一些平时在咱们社群问得比较多的问题——多期双重差分法和一些要点。我们想检验修建地铁对城市环境污染的影响,那么我们想到的是使用DID方法来得到因果关系。但是,我们有疑惑的地方是,各个城市修地铁的时间有先有后,而标准的双重差分方法一般要求t为同一时间点,比如20xx年。 对于这个问题,我们可以采用多期DID方法,将所有还没有修建地铁的城市作为控制组,把已经修建地铁的城市作为处理组,即使最终所有城市都修建了地铁,我们也可以把还没有修建地铁之时的城市作为控制组。 简单点讲,就是每个修建地铁的城市的DID交互项在数据中显示的不一样,因为DID交互项是两个虚拟变量的乘积:treated(是不是修建了地铁)和time(修建地铁的时间)。 这个DID的交互项等于1的情况是,这个城市在具体某年修建了地铁,而对于在修建地铁之前的年份,这个城市的DID 交互项等于0。这就表明,我们在多期DID使用中不再有统一的政策实施年份,而是允许每个城市都有自己的政策实施年份。 这样是不是有助于解决我们遇到的大部分问题。对于那些压根到目前为止都没有地铁的城市,那他的DID(自然不用说)

就是等于0,因为他的treated始终是为0,属于我们的控制组样本。注意,现在就是一个普通的xtreg回归,但是这里有些地方需要注意。第一,我们平时经常看到的 treated+time+treated*time+协变量的标准DID组合已经不见了,现在只剩下了treated*time这个DID交互项和协变量了。第二,我们尽量控制一下城市的个体效应和时间效应,来消除那些会影响DID交互项估计的不可观测因素和时间效应。下面这个多期DID模型就是如此的,αt是时间效应,βi 是城市效应,Xit是随着时间变动的协变量,BC*After就属于咱们感兴趣的DID估计量。 第三,这里面的treated(就是BC)虚拟变量当然可以灵活地替换为其他连续变量,比如,我们不仅对是否修建地铁对环境影响感兴趣,更是对修建地铁的里程对环境影响感兴趣。我们可以把BC替换成地铁的里程(length),然后我们的准DID 交互项就是length*After。这种DID设置的灵活性让这种方法有很大的适用性。 如果有时候我们不知道处理组具体怎么选择,那该如何设计方法呢?比如我们想要研究一下,美国政府对那些破产的按揭房(金融危机之后的事情)兴起了一个维护修理的政策举动,那这些房子就不至于破败不堪而影响了周围房子的价格。此时,我们就想看看这个政策举动对周围房子的价格的影响,但我们并不知道到底多远的距离才叫“周围”。

(完整版)有限差分方法概述

有限差分法(Finite Difference Method,简称FDM)是数值方法中最经典的方法,也是计算机数值模拟最早采用的方法,至今仍被广泛运用。该方法将求解域划分为差分网格,用有限个网格节点代替连续的求解域。有限差分法以Taylor级数展开等方法,把控制方程中的导数用网格节点上的函数值的差商代替进行离散,从而建立以网格节点上的值为未知数的代数方程组。该方法是一种直接将微分问题变为代数问题的近似数值解法,数学概念直观,表达简单,是发展较早且比较成熟的数值方法。对于有限差分格式,从格式的精度来划分,有一阶格式、二阶格式和高阶格式。从差分的空间形式来考虑,可分为中心格式和逆风格式。考虑时间因子的影响,差分格式还可以分为显格式、隐格式、显隐交替格式等。目前常见的差分格式,主要是上述几种形式的组合,不同的组合构成不同的差分格式。差分方法主要适用于有结构网格,网格的步长一般根据实际地形的情况和柯朗稳定条件来决定。 构造差分的方法有多种形式,目前主要采用的是泰勒级数展开方法。其基本的差分表达式主要有三种形式:一阶向前差分、一阶向后差分、一阶中心差分和二阶中心差分等,其中前两种格式为一阶计算精度,后两种格式为二阶计算精度。通过对时间和空间这几种不同差分格式的组合,可以组合成不同的差分计算格式。 下面我们从有限差分方法的基本思想、技术要点、应用步骤三个方面来深入了解一下有限差分方法。 1.基本思想 有限差分算法的基本思想是把连续的定解区域用有限个离散点构成的网格来代替,这些离散点称作网格的节点;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解。然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解。在采用数值计算方法求解偏微分方程时,再将每一处导数由有限差分近似公式替代,从而把求解偏微分方程的问题转换成求解代数方程的问题,即所谓的有限差分法。 2.技术要点 如何根据问题的特点将定解区域作网格剖分;如何把原微分

小振幅波有限差分方法数值模拟

小振幅波有限差分方法数值模拟 摘要:文章研究了二维水槽下线性势流方程小振幅波的数值解,采用crank-nicolson隐式有限差分方法对线性势流方程组进行离散.通过坐标变换将随不规则物理区域变换成一个固定的正方形计算区域,运用交错网格对计算区域进行分割。数值模拟自由面小振幅波的波高,计算出自由面波高的误差值.数值结果表明数值解很好地吻合解析解以及小振幅波具有驻波的线性特征和波节不动性。abstract: the model of small-amplitude water wave has been developed based on potential flow theory. there are numerical experiments of water wave motion undertaken in a two dimensional tank. crank-nicolson’s implicit scheme is adopted to discretize the equations about the linear potential flow. firstly, coordinate transformation is introduced to map the time-dependent irregular physical domain to a time-independent fixed square calculation domain. secondly, a staggered mesh system is applied to divide calculation domain. the numerical results agree well with analytic result and small-amplitude wave have standing wave linear features and node immobility. 关键词:小振幅波;crank-nicolson隐式有限差分方法;交错网格 key words: small-amplitude wave;crank-nicolson implicit

时域有限差分法(FDTD算法)的基本原理及仿真

时域有限差分法(FDTD算法)的基本原理及仿真

时域有限差分法(FDTD 算法) 时域有限差分法是1966年K.S.Yee 发表在AP 上的一篇论文建立起来的,后被称为Yee 网格空间离散方式。这种方法通过将Maxwell 旋度方程转化为有限差分式而直接在时域求解, 通过建立时间离散的递进序列, 在相互交织的网格空间中交替计算电场和磁场。 FDTD 算法的基本思想是把带时间变量的Maxwell 旋度方程转化为差分形式,模拟出电子脉冲和理想导体作用的时域响应。需要考虑的三点是差分格式、解的稳定性、吸收边界条件。有限差分通常采用的步骤是:采用一定的网格划分方式离散化场域;对场内的偏微分方程及各种边界条件进行差分离散化处理,建立差分格式,得到差分方程组;结合选定的代数方程组的解法,编制程序,求边值问题的数值解。 1.FDTD 的基本原理 FDTD 方法由Maxwell 旋度方程的微分形式出发,利用二阶精度的中心差分近似,直接将微分运算转换为差分运算,这样达到了在一定体积内和一段时间上对连续电磁场数据的抽样压缩。 Maxwell 方程的旋度方程组为: E E H σε +??=??t H H E m t σμ-??-=?? (1) 在直角坐标系中,(1)式可化为如下六个标量方程: ???????????+??=??-??+??=??-??+??=??-??z z x y y y z x x x y z E t E y H x H E t E x H z H E t E z H y H σεσεσε,????? ??? ??? -??-=??-??-??-=??-??-??-=??-??z m z x y y m y z x x m x y z H t H y E x E H t H x E z E H t H z E y E σμσμσμ (2) 上面的六个偏微分方程是FDTD 算法的基础。 Yee 首先在空间上建立矩形差分网格,在时刻t n ?时刻,F(x,y,z)可以写成 ),,(),,,(),,,(k j i F t n z k y j x i F t z y x F n =????= (3) 用中心差分取二阶精度: 对空间离散: ()[] 2 ),,21(),,21() ,,,(x O x k j i F k j i F x t z y x F n n x i x ?+?--+≈???= ()[] 2 ),21,(),21,() ,,,(y O y k j i F k j i F y t z y x F n n y j y ?+?--+≈???= ()[] 2 )21,,()21,,() ,,,(z O z k j i F k j i F z t z y x F n n z k z ?+?--+≈???=

有限差分,有限元,有限体积等离散方法的区别介绍

有限差分,有限元,有限体积等等离散方法的区别介绍 一、区域离散化 所谓区域离散化,实质上就是用一组有限个离散的点来代替原来连续的空间。实施过程是;把所计算的区域划分成许多互不重迭的子区域,确定每个子区域的节点位置及该节点所代表的控制容积。节点:需要求解的未知物理量的几何位置;控制容积:应用控制方程或守恒定律的最小几何单位。一般把节点看成是控制容积的代表。控制容积和子区域并不总是重合的。在区域离散化过程开始时,由一系列与坐标轴相应的直线或曲线簇所划分出来的小区域称为子区域。网格是离散的基础,网格节点是离散化物理量的存储位置。 大家都知道,常用的离散化方法有:有限差分法,有限元法,有限体积法。 1. 有限差分法是数值解法中最经典的方法。它是将求解区域划分为差分网格,用有限个网格节点代替连续的求解域,然后将偏微分方程(控制方程)的导数用差商代替,推导出含有离散点上有限个未知数的差分方程组。这种方法发展比较早,比较成熟,较多用于求解双曲线和抛物线型问题。用它求解边界条件复杂、尤其是椭圆型问题不如有限元法或有限体积法方便。 2. 有限元法是将一个连续的求解域任意分成适当形状的许多微小单元,并于各小单元分片构造插值函数,然后根据极值原理(变分或加权余量法),将问题的控制方程转化为所有单元上的有限元方程,把总体的极值作为各单元极值之和,即将局部单元总体合成,形成嵌入了指定边界条件的代数方程组,求解该方程组就得到各节点上待求的函数值。对椭圆型问题有更好的适应性。有限元法求解的速度较有限差分法和有限体积法慢,在商用CFD软件中应用并不广泛。目前的商用CFD软件中,FIDAP采用的是有限元法。 3. 有限体积法又称为控制体积法,是将计算区域划分为网格,并使每个网格点周围有一个互不重复的控制体积,将待解的微分方程对每个控制体积积分,从而得到一组离散方程。其中的未知数十网格节点上的因变量。子域法加离散,就是有限体积法的基本方法。就离散方法而言,有限体积法可视作有限元法和有限差分法的中间产物。 4. 有限分析法:同有限差分法一样,用一系列网格线将区域离散,所不同的是每个节点与相邻8个邻点组成。在计算单元中把控制方程中的非线形项局部线形化,并对该单元上未知函数的变化型线作出假设,把所选定型线表达式中的系数和常数项用单元边界节点上未知的变量值来表示,这样该单元内的被求问题就转化为第一类边界条件下的一个定解问题,可以找出分析解;然后利用这一分析解,得出该单元中点及边界上8个邻点上未知值间的代数方程,此即为单元中点的离散方程。两种离散方法外节点法:节点在子域的四角,先定节点位置而计算相应的界面内节点法:节点在子域中心,子域与控制容积重合。计算时先定界面后算出节点位置。 5. 边界元法(Boundary Element Method,BEM)上面四种方法都必须对整个区域作离散化处理,用分布在整个区域上的有限个节点上函数的近似值来代替连续问题的解。在边界元方法中应用格林函数公式,并通过选择适当的权函数把空间求解域的偏微分方程转换成为其边界上的积分方程,它把求解区中任一点的求解变量(如温度)与边界条件联系了起来。通过离散化处理,由积分方程导出边界节点上未知值的代数方程。解出边界上的未知值后就可以利用边界积分方程来获得内部任一点的被求函数之值。边界元法的最大优点是,可以使求解问题的空间维数降低一阶,从而使计算工作量及所需计算机容量大大减小。边界元法推广应用的一个最大限制是,需要已知所求解偏微分方程的格林函数基本解。虽然对不少偏微分方程这种基本解业已找出,但对Navier-Stoles方程这样的非线性偏微分方程,至今尚未找到其基本解。目前的一种处理方式是,把Navier-Stokes方程中的非线性项看作是扩散方程的源项并通过迭代的方式来求解,但一般只能获得Re较低情形的解。最近文献中采用高阶涡量—流函数方程的边界元方法,已使获得顶盖驱动流稳定解的Re高达10000。 格子—Boltzmann方法(Lattice-Boltzmann method,LBM)格子—Boltzmann方法是基于分子运动论的一种模拟流体流动的数值方法。在上述各种数值方法中,把本质上是离散的介质先假定是连续的,在此基础上建立起了N-S方程,然后又再把它离散化。在LBM中不再基于连续介质的假设,而是把流体看成是许多只有质量没有体积的微粒所组成,这些微粒可以向空间若干个方向任意运动。通过其质量、动量守恒的

相关主题
文本预览
相关文档 最新文档