dydxdz?f(x,y,z),y(x0)?y0
(1)
?g(x,y,z),z(x0)?z0dx
1、欧拉法的计算公式? yn?1?yn?hf(xn,yn,zn) zn?1?zn?hg(xn,yn,zn)
y0?y(x0),z0?z0(x0), n?0,1,2,? (2) 2、改进的欧拉法的计算公式?
对n=0,1,2,?,计算??
y(0)n?1?yn?hf(xn,yn,zn) (0)
zn?1?zn?hg(xn,yn,zn)
y(k?1)(k)(k)n?1?yn?h2[f(xn,yn,zn)?f(xn?1,yn?1,zn?1)] (3)
z(k?1)?hn?12[f(x(k)(k)
?znn,yn,zn)?f(xn?1,yn?1,zn?1)]
k?0,1,2,?
当连续两次迭代结果之差的绝对值小于给定的精度ε,即 y(k?1)n?1?y(k)n?1?ε z(k?1)n?1?z(k)n?1?ε
时,取y(k?1)n?1?yn?1,z(k?1)n?1?zn?1,然后转入下一步运算。?
3、四阶标准龙格—库塔公式??
yn?1?yn?16(k1?2k2?2k3?k4) z1n?1?zn?6(m1?2m2?2m3?m4)
k1?hf(xn,yn,zn) m1?hg(xn,yn,zn)
k12?hf(xhn?2,yn?k12,zn?m2) (4)
mk12?hg(xhn?2,yn?2,zn?m12)
21
k3?hf(xhn?2,yk2n?2,zn?m22)
mk2m23?hg(xn?h2,yn?2,zn?2)
k4?hf(xn?h,yn?k3,zn?m3) m4?hg(xn?h,yn?k3,zn?m3)
n?0,1,2,?
4、四阶亚当斯显式公式?? h
yn?1?yn?24(55fn?59fn?1?37fn?2?9fn?3)
zhn?1?zn?24(55gn?59gn?1?37gn?2?9gn?3) (5)
n?4,5,6,? ? 其中 ?
? fn?i?f(xn?i,yn?,zn?i)
gn?i?g(xn?i,yn?,zn?i) i?0,1,2,3 ??需要注意,前几步的值由其它单步法求得。? 7.5.2二阶方程?
设有二阶微分方程初值问题??
y???g(x,y,y?) y(x0)?y0,y?(x0)?y?0 (6)
令z=y′,则(6)化为一阶微分方程初值问题?
dydx?z
dzdx?g(x,y,z) (7)
y(x0)?y0,z(x0)?y?0
(7)式中
dydx?z?f(x,y,z),于是应用四阶龙格—库塔公式(4),有??k1?hzn?y?n
22
?令
k2?h(zn?m12m22)?hy?n?hm12hm23
2k3?h(zn?)?hy?n?
k4?h(zn?m3)?hy?n?hm??则??
yn?1?yn?16(k1?2k2?2k3?k4)?yn?hy?n?16(m1?m2?m3)
又因为??
zn?y?n?1,zn?1?y?n?1
zn?1?zn?16(m1?2m2?2m3?m4)
?所以?
y?n?1?zy?n?16(m1?2m2?2m3?m4)
于是得初值问题(6)的计算公式为?
yn?1?yn?hy?n?1616(m1?m2?m3)
y?n?1?y?n?(m1?2m2?2m3?m4)m1?hg(xn,yn,y?n)
m2?hg(xn?h2h2,yn?121212hy?n,y?n?141212m1)
12m2)m3?hg(xn?,yn?hy?n?hm1,y?n?
m4?hg(xn?h,yn?hy?n?hm2,y?n?m3)
第六节 边值问题的数值解法
?在具体求解微分方程时,需要附加某种定解条件。微分方程与定解条件一起组成定解问题。对于高阶微分方程,定解条件通常有两种给法,第一种是给出积分 曲线在初始时刻的性态,称这类条件为初始条件,相对应的定解问题称为初值问题;第二种是给出积分曲线在首末两端的性态,称这类条件为边值条件,相
23
对应的定解问题为边值问题 。 ?本节以二阶微分方程??
y″=f(x,y,y′),a≤x≤b??
为例讨论边值问题,其边值条件 可分为三类:?
第一边值条件??
y(a)=a, y(b)=β
??第二边值条件?
? y′(a)=a, y′(b)= β??
第三边值条件??
y′(a)-a0y(a)=a , y′(b)- β0y(b)=β??
其中 ??a0≥0,?β0≥0,?a0 +β0>0.??
本节介绍解线性方程边值问题的差分方法以及 适用于非线性方程边值问题的试射法
。?
7.6.1解线性方程边值问题的差分方法?
二阶线性方程的一般形式为?
? y″+p(x)y′+q(x)y=f(x), a≤x≤b? (1)
首先将区间[a,b]进行等距划分,分点
? xi?a?ih, i=0,1,2,?,n? 其中?
h?b?an
一般地称?x0,?a与xn=b为边界点,称x1,x2,?,xn?1为内部节点。
其次,在各个xi节点上,将y′,y″用差商近似表示。这里要求有相同阶数的截断误差,以保证精度协调。对于内部节点,二阶导数用二阶中心差商表示,得
yi?1?2yi?yi?1h2?y??(xi)?O(h), i=1,2,?,n-1
2一 阶导数用一阶中心差商表示,得??
24
yi?1?yi?12h?y?(xi)?O(h), i=1,2,?,n-1
2? 假设?yi?y(xi),,则
y??(xi)?yi?1?2yi?yi?1h2,
y?(xi)?yi?1?yi?12h,于是得方程(1)的差分方程
yi?1?2yi?yi?1h2+
piyi?1?yi?1qiy?fi+
2h
i=1,2,?,n-1 (2)
其中pi?p(xi),qi?q(xi),fi?f(xi)?.?
将(2)整理,可以写成下列形式
aiyi?1?biyi?ciyi?1?di? ? i=1,2,?,n-1 (3)
其中 ??
ai?1?1212hpi2bi??2?hqici?1?2hpi di?hfi
(3)是含有n+1个未知数yi (i=0,1,2,?,n)的线性方程组,方程的个数为n-1。要使方程组(3)有唯一解,还需要有两上边值条件补充两个方程。对于第一边值条件,直接就得到两个方程?
y0?a yn??
于是得到第一边值问题的差分方程组?
?y1??d1?a1a2?????y2d2????b1c1???????an?2?a2b2? (4) y??d?n?2??n?2??y????n?1??dn?1?cn?1??c2 25