完美版第10章习题解答数值分析.docx
《完美版第10章习题解答数值分析.docx》由会员分享，可在线阅读，更多相关《完美版第10章习题解答数值分析.docx（12页珍藏版）》请在冰豆网上搜索。
完美版第10章习题解答数值分析
第十章习题解答
1、用Euler方法及改进的Euler方法求解初值问题
y二x-yx[0,1]
Iy（0）=2
取h=0.1，并将计算结果与精确值相比较。
解：
f（x,y）=X-y,由Euler公式及改进的Euler方法，代入h=0.1，有
果如下
n=0
1
2
3
4
5
6
7
8
9
10
xn
=0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
Yn
=2
1.8000
1.6300
1.4870
1.3683
1.2715
1.1944
1.1350
1.0915
1.0623
1.0461
Yn
=2
1.8150
1.6571
1.5237
1.4124
1.3212
1.2482
1.1916
1.1499
1.1217
1.1056
y
2
1.8145
1.6562
1.5225
1.4110
1.3佃6
1.2464
1.1898
1.1480
1.1197
1.1036
yn为Euler方法的结果，yn为改进的Euler方法的结果，y为精确解。
2、用梯形公式求解初值问题
t2
3、试用Euler公式计算积分fetdt在点x=0.5,1,1.5,2的近似值。
x2x2
解:
f（x,y）二2xe由Euler公式得y“yn2*0.5x“en，计算可得
n=0
1
2
3
4
x；=0
0.5
1
1.5
2
yn=0
0.6420
2.0011
6.7450
34.0441
4、定初值问题
r•
y
工siny
X-
X。
y（xo）=y。
试用Taylor展开法导出一个三阶的显式公式。
I
解：
由Taylor公式，并代入y=siny可得
y（Xnh^yXn）yx（nh）y2；n）h2y:
;h30h（4）
yn1ynsinynh皿h2口项助治
5、已知初值问题
x[0,1]
R-K方法求解此问题，取步长为h二0.5.
y二x2x-yy（0H0试分别用改进的Euler方法和四阶
2
解：
f（x,y）=Xx-y,由改进的Euler方法和四阶R-K方法，代入h=0・5，
广—
yn卄yn+hf（Xn,yn）
改进的Euler方法h-
[yn卄yn+y（f（Xn,yn）+f（Xn*,y诂））
四阶R-K算法yn1二ynh（k「2k?
2kskq）
6
ki=f（Xn,yn）
11
k2=f（Xn中；h,yn+：
hkj
〈22
11
k3=f（Xn+：
h,yn+：
hk2）
22
*4=f（Xn+h,yn+hk3）
计算可得
n
=0
1
2
Xn
=0
0・5
1
yn
=0
0.1875
0.7109
yn
=0
0.1439
0.6329
其中yn为改进的Euler方法的结果,
yn
为四阶R-K
方法的结果。
6、试证明对任意的参数a，以下Runge-Kutta公式是一个二阶公式，并导出其数值稳定条件。
h
yn1二ynh（k2k3
k=f（xn,yn）
k2=f（Xn+ah,yn+hki）
K=f（Xn+（1-a）h,yn+（1-a）hki）
证明：
将k2,k3做二元Taylor展开
k2二f（Xn,yn）ahfx'（Xn,yn）hk1fy'（Xn,yn）O（h、
''2
k^f（Xn,yn）（1-a）hfx（Xn,yn）（1-a）hk1fy（Xn,yn）O（h）
代入得
h''2
yn1二yn?
（2f（Xn,yn）hfx（Xn,%）人匕彳丫仕“』“）O（h））
再将y（xnj在点xn展开
y（XnG=y（Xn）+y（Xn）h+^2X^h2+O（h3，式中y'（Xn）=fXnyn）
代入后有
y（Xn）二fxfyf
II
fx+ffy23
y（Xn1）=y（Xn）fh-hO（h）
故En^=y（焉出）一y.勺=O3h）即对任意的参数a，公式是二阶公式。
下面讨论公式的数值稳定条件：
取模型方程「二爼y，将f（x,y）二纭y代入ki得到
ki=f（Xn』n）=^yn
“k2=f（Xn+ah,%+hki）=A（y^hk）=x（1十hh）y“
k3=f（Xn+（1-a）h,yn+（1-a）hki）=丸（1+丸（1一a）h）y
h
再代入yn.1二yn•-（k2k3）得到
a2
yn1Tn[1h（1-2）伸）]
7、试证明以下Runge-Kutta公式是一个三阶公式，并导出其数值稳定条件。
1yn1二yn（2k13k24k3）
9
k1二hf（Xn,yn）
11
k2二hf（Xn-h,yn
22
33
k3二hf（Xn—h,yn—k?
）
44
证明：
将k2*3做三元Taylor展开
22
k^h（fhfx'如fy'fXX”呼fXy''牛fyy”）0（h'））
222424
22
3'3'19h"9hk?
”9k?
”3
k^h（f^hfX[k2fy[（需fXX*fXy晴fyy）0（h'））
代入得
丄h
9h
1
9f
1
32""
2"3
yn1=yn（9f
9
+
2
fx
2
fy
-（hfxx2hffxyk
【2fyy）0（h））
h2
1
h2
1
h2
nn11
2"4
二ynhf
2
fx
+
2
ffy
+
6
（f+2ff+ff
xxxyxy
fffyy）0（h）
再将y（XnJ在点x
n展开
y（Xn1）=y（Xn）y（Xn）h^^X^h2^^X^h30（h4）
式中y（Xn）二f（Xn,yn）,y'（Xn）二f；f；f
'"''''''2''
y（Xn）=fxx+2fxyf+fxfyf+ffyy代入后有
fX+ffy2h''''''2''4
y（Xn1）=y（xn）fh—2!
—h莎（fxx2fxyffxfyfffyy）O（h）
4故Enni=y（焉申）—y,申=Oh）即公式是三阶公式。
下面讨论公式的数值稳定条件：
取模型方程■=y,将f（x,y）八y代入ki得到
再代入yn厂ynT（2k13k24ka）得到
yn1=yn[1h如）2讣）3]
26
p.1213
是亠1叩h-（h）2-（h）3|n26
1213
绝对稳定区为|1h-（h）-（h）3卜：
1
26
用Euler方法求解下列问题，从数值稳定性条件考虑，对步长应做什么限制?
解：
由Euler公式yn1二y“hf（Xn"n）
式变为yn1=（4-6h）yn，两式相减得到yn1-yn1=（4-6h）（yn-yn）
1
故当|1-61<|即0£h<—时，
（1）用Euler方法求解是数值稳定的。
3
10xnyn
对
（2）ynd=ynh（1哼），若第n步和第n+1步分别有误差，则上
1十Xn
式变为yn冷=yn…10%；yn，两式相减得到
1+Xn
10x2
故当|1与h|：
：
：
1,即0：
：
：
h时，
（2）用Euler方法求解是数值稳定的。
1十Xn5
9、用二阶的Adams预估校正公式求解初值问题
yt_y
ly（0）=
0
取步长h—0・2求出x—1.
0寸的近似值,
表头用改进的Euler方法（保留小数点
后三位）。
解：
计算结果如下：
n=0
1
2
345
Xn=0
0.2
0.4
0.60.81.0
yn=0
0.180
0.329
0.4510.5510.632
10、对初值问题
y二-10y
y（x°）=y°
用以下二阶R-K方法求解，并导出其绝对稳定域。
1
yn1二yn尹1k2
』k1=hf（Xn,yn）
匕=hf（Xn+h,yn+kJ
解：
f（x,y）——10y，代入二阶r-k方法可得ynr=yn-50hyn，若第n步和
——P
yn1-yn1=（4-50h）（yn-yn），于是*1=|1~50h|，其绝对稳定区为
|1-50h|1。
11、试推导三阶的Adams显式和隐式，并写出其带修正的预估校正公式。
解：
三阶Adams显式公式的推导：
Xn_1
由公式y（Xn1）=y（Xn）亠If（X,y）dX，取Xn/,Xn」，Xn作差值节点，得
（X—XnJL）（X—Xn）
f（x,y）二P2（x）R2（x）—-f（X2,y（Xn4）
（Xn/—Xn」）（Xn/—Xn）
（x-Xn/）（X-Xn）
（X-Xn」）（X-Xn/）
（xnJ_xn_2）（xnJ_xn）（x^_XnJ）（xn
f（3）（）
（X-X-t）（X-Xn」）（X-Xn）
3!
f（xn-1，y（Xn」））
将P?
（X）代入公式，并作代换X=Xn，th可得
h
y（Xnl）i（Xn）12（23fn「6fn」5f®
Xn+9
（Xn?
Xn）
局部截断误差为En1二R2（x）dx9h4y⑷（）
冷24
三阶Adams隐式公式的推导：
Xn1
由公式y（Xn1）=y（Xn）•.f（x,y）dx，取Xn4,Xn,Xn'1作差值节点，得
Xn
f（X,y）二P2（X）R2（X）二（XXn1）（XXn）f（Xn」,y（Xn」））
（X-Xn4）（X-Xn1）
（Xn_i-Xn卅）（Xn_t-Xn）
（X-Xn/）（X-Xn）
f（Xn,y（Xn））+f（Xn^y（XnG）
.f（3）（）
（Xn-Xn4）（Xn~Xni）（Xn1-X._|）（Xn1-X.）
3!
（X-XnJ（X-Xn）（X-Xn1）3!
将P2（X）代入公式，并作代换X二Xn•th可得
h
y（Xn1）=y（Xn）—（5fn「8—fnJ
12
Xn+1
局部截断误差为Enj=R2（x）dx—h4y⑷（）-（xn/,xn1）
24
Xn
带修正的预估校正公式:
表头：
三阶RK公式
h
P:
卩“勺二丫区）（23fn-16fn15fn2）
12
M:
9
mn1=Pn1（Cn-Pn）（第一次修正值取0）
10
E:
fn1=f（Xn1,mn1）
C:
h
Cn1二yn（5f（xnmn1）8fn-fn」）
12
M:
1
yn1=Cn1_10（Cn1_Pn1）
E:
fn1二f（Xn1,Yn1），为下一步计算『n2用
升级会员 