郑州大学研究生课程（2012-2013学年第一学期）数值分析
Numerical Analysis
习题课
第八章常微分方程数值解法
待求解的问题：一阶常微分方程的初值问题/* Initial-Value Problem */:
?????=∈=0
)(]
,[),(y a y b a x y x f dx dy
解的存在唯一性（“常微分方程”理论）：只要f (x , y )
在[a , b ] ×R 1 上连续，且关于y 满足Lipschitz 条件，即存在与x , y 无关的常数L 使
对任意定义在[a , b ] 上的y 1(x ) 和y 2(x ) 都成立，则上述IVP 存在唯一解。
1212|(,)(,)|||
f x y f x y L y y ?≤?一、要点回顾
§8.2 欧拉(Euler)法
通常取(常数),则Euler 法的计算格式
h h x x i i i ==?+1??
?=+=+)
(),(001x y y y x hf y y i i i i i =0,1,…,n ( 8.2 )
§8.2 欧拉(Euler)法(1) 用差商近似导数
))
(,()()()()(1n n n n n n x y x hf x y x y h x y x y +=′+≈+??
?=+=+)
(),(01a y y y x hf y y n n n n 差分方程初值问题向前Euler 方法h
x y x y x y n n n )
()()(1?≈
′+))
(,()
()(1n n n n x y x f h
x y x y ≈?+))(,()(n n n x y x f x y =′
§8.2 欧拉(Euler)法若用向后差商近似导数，即
))
(,()()(111++++≈n n n n x y x hf x y x y ??
?=+=+++)
(),(0111a y y y x hf y y n n n n h
x y x y x y n n n )()()(11?≈
′++向后Euler 方法
))(,()
()(111+++≈?n n n n x y x f h x y x y ))
(,()(111+++=′n n n x y x f x y
§8.2 欧拉(Euler)法（2）用数值积分方法
∫
+=?+1
))(,()()(1n n
x x n n dx
x y x f x y x y ∫
∫
++=′1
1
))(,()(n n
n n
x x x x dx
x y x f dx x y 1111(,)(,()), (,)(,()),
()()n n n n n n n n f x y f x y x f x y f x y x y x y y x y ++++≈≈≈≈分别用左矩形和右矩形公式，即
代替上式右端的积分,并注意 ,分别得到
1111(,) (,)
n n n n n n n n y y h f x y y y h f x y ++++=+=+，。向前欧拉公式和向后欧拉公式：
§8.2 欧拉(Euler)法若对积分用梯形公式，则得
[]
))(,())(,(2
)()(111++++≈?n n n n n n x y x f x y x f h
x y x y []?????
=++≈+++
)(),(),(20111a y y y x f y x f h y y n n n n n n 梯形欧拉公式
§8.2 欧拉(Euler)法
欧拉方法的收敛性
112
12 () Taylor , (,),()()()()()
2 ()(,())() (1)
2
n n n n n n n n n n n n n n y x x x x h y x y x h y x hy x y h
y x hf x y x y ξξξ+++?∈′′′=+=++′′=++将在点展开%1
(),(8.2) () (,()) (2)
n n
n
n
n y x y y x h f x y x ++ 假定已知准确值 利用欧拉公式，定义
§8.2 欧拉(Euler)法
局部截断误差
%()111
12
() ()()(,()) (). 2n n n n n n n n T y x y y x y x hf x y x h
y ξ++++=?=?+′′=2
2
21 max |()|,|||()|().
22
a x
b n n M y x h h T y M O h ξ≤≤+′′=′′=≤=令则
称为局部截断误差
§8.2 欧拉(Euler)法欧拉方法的收敛性
定义
若给定方法的局部截断误差满足
则称该方法是P 阶的，或称为具有P 阶精度。
1
1||(),
p n T O h
++=
§8.2 欧拉(Euler)法
整体截断误差
111().
n n n e y x y +++=?记112121112 ,,(), (),, () ,(),,,.
n n n n n n y y y y y x y x y x y x e y y y +++L L L 因为计算 时用到的,,是 的近似值每步产生的误差会累积到计算的误差中因此 与 ,, 都有关称整体截断误差为%%%11111111
1
1
|||()||()||| ||||. (3)
n n n n n n n n n n e y x y y x y y y T y y ++++++++++=?≤?+?≤+?
§8.2 欧拉(Euler)法
欧拉方法的收敛性
%%()1111,(8.2)|||()
(,())
(,)| |() ||(,())(,)| |() ||() | () (1)n n n n n n n n n n n n n n n n n n n n y y y y y x h f x y x y h f x y y x y h f x y x f x y y x y hL y x y hL ++++?=+?+≤?+?≤?+?≤+对应用欧拉公式得
李普希兹条件||. (4)n e
由此知，当
max ||,(4)(3)k k
T T =记 将代入得
[]112
12
1
01
1
2
||(1)||
(1)(1)|| (1)(1)|| (1)(1)(1) (1)
||
(1)1(1)1
()
(1)1 (n n n n n
n n n e T hL e T hL T hL e T hL T hL e T hL T hL T hL T hL e hL hL T O h hL hL
O +??+++≤++≤++++≤++++≤≤++++++++++?+?≤=
+?≤L L ).
h 10,||0, n h e +→→欧拉公式是一阶有收敛的。
§8.2 欧拉(Euler)法
向后欧拉公式
#隐式欧拉法或向后欧拉法
#/* implicit Euler method or backward Euler method*/
11()()
()n n n y x y x y x h
++?′≈
x n +1点向后差商近似导数
111111()()()
()()(,)
n n n n n n n n n n y x y x hy x y x y y x y y h f x y ++++++′≈+↑≈≈=+代入隐式或后退欧拉公式
§8.2 欧拉(Euler)法
向后欧拉公式
由于未知数y n +1同时出现在等式的两边，故称为隐式/* implicit */欧拉公式，而前者称为显式/* explicit */ 欧拉公式。隐式公式不能直接求解，一般需要用Euler 显式公式得到初值，然后用Euler 隐式公式迭代求解。因此隐式公式较显式公式计算复杂，但稳定性好（后面分析）。隐式欧拉公式中的未知数y n +1可通过以下迭代法求解：
0)1(1)()
1
11(,)(,)n n n n k k n n n n y y h f x y y
y h f x y +++++?=+?=+?（
§8.2 欧拉(Euler)法
向后欧拉公式
1(,) 0,1,...n n n n y y h f x y n +=+=比较欧拉显式公式和隐式公式及其局部截断误差2
31112()()()h n n n n T y x y y x O h +++′′=?=+%显
式公式
111
(,)
n n n n y y h f x y +++=+隐式
公式
2311
12()()()h n n n n T y x y y x O h +++′′=?=?+%
若将这两种方法进行算术平均，即可消除误差
的主要部分/*leading term*/而获得更高的精度,称为梯形法
#梯形公式/* trapezoid formula */
—显、隐式两种算法的平均
111[(,)(,)]
2
n n n n n n h
y y f x y f x y +++=++R y x
y
O h
++
+
=?=
例8.2.3对初值问题??
?==+′1
)0(0
y y y 证明用梯形公式求得的近似解为n
n h h y ???
?
??+2?2=并证明当步长h →0时,y n 收敛于精确解证明: 解初值问题的梯形公式为
x
e
?)]
,(),([1+1+1++2
+=n n n n n n y x f y x f h
y y y y x f ?=),(∵][211++??+=n n n n y y h
y y ∴整理成显式n n y h h y ??
?
???+2?2=1+反复迭代,得到
1
+2?3
1?2
1
+??
?
???+2?2==??????+2?2=??????+2?2=??????+2?2=y h h y h h y h h y h h y n n n n n ...10=y 22n
x
n
h y e h ????
=→??+??
∵
∴
梯形公式
单步
好
隐
1阶
欧拉隐式公式
单步
差
显
1阶
欧拉显式公式
步数
稳
定
性
显
隐
精
度
局部截断误差
公式
()
3
(3)
3n
h
y x
()
2
(2)
2n
h
y x
()
2
(2)
2n
h
y x
?
欧拉法小结
单步
好
隐
2阶
§8.3 改进欧拉(Euler)方法
先用欧拉公式(8.2)求出一个初步的近似值,称为预测值, 它的精度不高, 再用梯形公式对它校正一次,即迭代一次,求得y n+1,称为校正值, 这种预测-校正方法称为改进的欧拉公式：
1
n y +[]??
???++=+=++++校正预测 ),(),(2
),(1111n n n n n n n n n n y x f y x f h
y y y x hf y y 称为Euler 公式与梯形公式的预测—校正系统。
《计算机数学基础》数值部分第五单元辅导 14 常微分方程的数值解法 一.重点内容 1. 欧拉公式： ）心知1）a 儿+1 =儿 + hfg ，儿） m 1、 伙=0丄2,…川一 1） I 无=x Q +kh 局部截断误差是0（*）。 2. 改进欧拉公式： 预报一校正公式： 预报值 _v*+1 =儿+ hf （x k ,儿） - h - 校正值 y M = y k +-［f （x kt y k ） + /（x A+1, y M ）］ 即 儿+1 =儿+ ￡ "（忑'儿）+心+「儿+ hfg ,儿））］ 或表成平均的形式： 儿=儿+ hfg ,儿） '儿=儿+"（无+】，儿） +K ） 改进欧拉法的局部截断误差是0（2） 3. 龙格一库塔法 二阶龙格一库塔法的局部截断误差是0（爪） 三阶龙格一库塔法的局部截断误差是0（护） 四阶龙格F 塔法公式：儿计=儿+ 2（匕+ 2心+ 2? + ?） 四阶龙格一库塔法的局部截断误差是0（爪）。 二实例 y' = — y — xv f2（0 < x < 0.6） 例1用欧拉法解初值问题｛ ' ? -取步长/匸02计算过程保留 b （o ）= 1 4位小数。 解/i=0.2. f （x ）= —y —xy 2<,首先建立欧拉迭代格式 y*+i =儿+ hf g,y k ） = y k -hy k -hx k y ； =0.2 儿（4 - x k y k ）（k = 0,1,2） K 2=f(x n +^h, yk+-hK\)t gg+舟人，>'n +y/?A3）；
当k=0, xi=0.2 时，已知x（）=0,y（）=l,有 y(0?2)今i=0?2X l(4-0X 1)=0.8000 当k=\. M=0?4时，已知“=0?2」尸0?8,有 y(0?4)今2=0.2 X 0.8X(4-0.2X0.8)=0.614 4 当k=2, xs=0.6 时，已知x2=0.4,y2=0.6144,有 y(0?6)今3=0.2 X0.6144X (4-0.4 X 0.4613)=0.8000 「J, ,2 ?_ ZX 例2用欧拉预报一校正公式求解初值问题\y + v +V sinx=,取步长/?=0.2,计算 .y ⑴=1 y(0.2),y(0.4)的近似值，计算过程保留5位小数。 解步长力=0.2,此时/(x,y)=—y—fsiiu 欧拉预报一校正公式为： 预报值兀I = y k + hfg y k) - I J_ 校正值)3=儿+尹(忑，儿)+ fg,儿+1)] 有迭代格式] 预报值儿+] = y k 4-h(-y k -y； sin x k) =y k (0?8-0?2儿sin x k) < h 、—— 2 校止值y如]=儿 +尸[(一片一力sinxJ + LN+i-yl sin.v I+1)] ——?> =儿(°?9一0?1儿sin心)一0?1(儿+| +y；j sin心利) 、"M=0.別=1」)=1 时，Xj=1.2> 有 儿=yo(°?8-O?2yo sinx0) = 1 x (0.8-02x lsin 1) = 0.63171 y(1.2) ?= lx(0.9-0.1xlxsinl)-0.1(0.63171+0.631712sinl.2) = 0.71549 当 T xi=1.2, yi=0.71549 时,x2=1.4,有 y2 =儿(0.8-0?2儿sinXj) = 0.71549x(0.8-02x0.71549sinl.2) =0.47697 y(14) z y2 = 0.71549x(0.9-0.1x0.71549xsin 1.2)-0.1(0.47697+ 0.476972 sin 1.4) =0.52608 V = 8 — 3y 例3写出用四阶龙格一库塔法求解初值问题^ ‘的计算公式，取步长/匸0.2计 b(0) = 2 算y(0.4)的近似值。讣算过程保留4位小数。 解此处.心,刃=8 —3”四阶龙格一库塔法公式为 艰=儿 + % + 2? + 2勺 + ?) 1 h, y n+ y/?A3)： 本例计算公式为: 0 2 呱严儿+三(32?+2?+心
2014级硕士研究生数值分析上机实习（第一次） 姓名：学号：学院： 实习题目：分别用二分法和Newton迭代法求方程x3■ 2x210x-20=0的根.实习目的：掌握两种解法，体会两种解法的收敛速度. 实习要求：用C程序语言编程上机进行计算，精确到8位有效数字. 报告内容： 1.确定实根的个数以及所在区间 2.将最后两次计算结果填入下表（保留8位数字）： 3.实习过程中遇到哪些问题？如何解决？有何心得体会?
4.两种解法的计算程序（此页写不下时可以加页）：
2014级硕士研究生数值分析上机实习（第二次）姓名：学号：学院： 实习题目：计算8阶三对角矩阵A=tridiag（0.235, 1.274, 0.235）的行列式.实习目的：掌握计算行列式的方法. 实习要求：首先选择一种算法，然后用C程序语言编程上机进行计算.报告内容： 1.简单描述所采用的算法: 2?计算结果： A 3.实习过程中遇到哪些问题？如何解决？有何心得体会?
4.写出C语言计算程序（此页写不下时可以加页）：
2014级硕士研究生数值分析上机实习（第三次） 姓名：学号：学院： 分别用Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组实习题目： 2lx + 9.8y+ 3.4z= 6.7 <2.7x + 1.8y+ 7.2z= 2.4 8.6x + 1.5y + 3.4z = 1.9 实习目的：感受两种迭代法的收敛速度. 首先构造收敛的Jacobi迭代法和Gauss-Seidel迭代法，然后用实习要求： C程序语言编程上机进行求解，初始值均取为0,精确到4位小 数. 报告内容： 1.写出收敛的Jacobi迭代法和Gauss-Seidel迭代法:
2009年春季工学硕士研究生学位课程（数值分析）真题试卷 (总分：28.00，做题时间：90分钟) 一、填空题(总题数：6，分数：12.00) 1.填空题请完成下列各题，在各题的空处填入恰当的答案。（分数： 2.00） __________________________________________________________________________________________ 解析： 2.已知x=0．045，y=2．013_____ （分数：2.00） __________________________________________________________________________________________ 正确答案：(正确答案：0．902×10 -4) 解析： 3.已知矩阵1 =______，‖A‖ 2 =______． （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 4.设函数f(x)=2x 3 -x+1，则f(x)以x 0 =-1，x 1 =0，x 2 =1为插值节点的二次插值多项式为______．（分数：2.00） __________________________________________________________________________________________ 正确答案：(正确答案：x+1) 解析： 5.设函数f(x)∈C 2 [x 0 -h，x 0 +h]，h＞0，则 （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 6.______，该公式的代数精度为_____． （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 二、计算题(总题数：2，分数：4.00) 7.(0，+∞)内实根的分布情况，并用迭代法求出该方程在(0，+∞)内的全部实根，精确至3位有效数字． （分数：2.00） __________________________________________________________________________________________ 正确答案：(正确答案：设，显然f(x)=0在(2，+∞)内无根．在(0，2]内，f"(x)=cosx-
i.常微分方程初值问题数值解法 常微分方程初值问题的真解可以看成是从给定初始点出发的一条连续曲线。差分法是常微分方程初值问题的主要数值解法，其目的是得到若干个离散点来逼近这条解曲线。有两个基本途径。一个是用离散点上的差商近似替代微商。另一个是先对微分方程积分得到积分方程，再利用离散点作数值积分。 i.1 常微分方程差分法 考虑常微分方程初值问题：求函数()u t 满足 (,), 0du f t u t T dt =<≤ （i.1a ） 0(0)u u = (i.1b) 其中(,)f t u 是定义在区域G : 0t T ≤≤, u <∞上的连续函数，0u 和T 是给定的常数。我们假设(,)f t u 对u 满足Lipschitz 条件，即存在常数L 使得 121212(,)(,), [0,]; ,(,)f t u f t u L u u t T u u -≤-?∈∈-∞∞ （i.2） 这一条件保证了（i.1）的解是适定的，即存在，唯一，而且连续依赖于初值0u 。 通常情况下，（i.1）的精确解不可能用简单的解析表达式给出，只能求近似解。本章讨论常微分方程最常用的近似数值解法-差分方法。先来讨论最简单的Euler 法。为此，首先将求解区域[0,]T 离散化为若干个离散点： 0110N N t t t t T -=<< <<= （i.3） 其中n t hn =，0h >称为步长。 在微积分课程中我们熟知，微商（即导数）是差商的极限。反过来，差商就是微商的近似。在0t t =处，在（i.1a ）中用向前差商 10()()u t u t h -代替微商du dt ，便得 10000()()(,())u t u t hf t u t ε=++ 如果忽略误差项0ε，再换个记号，用i u 代替()i u t 便得到 1000(,)u u hf t u -= 一般地，我们有 1Euler (,), 0,1, ,1n n n n u u hf t u n N +=+=-方法： （i.4） 从(i.1b) 给出的初始值0u 出发，由上式可以依次算出1,,N t t 上的差分解1,,N u u 。
2005~2006学年第一学期硕士研究生期末考试试题（A 卷） 科目名称：数值分析 学生所在院： 学号： 姓名： 注意：所有的答题内容必须答在答题纸上，凡答在试题或草稿纸上的一律无效。 一、（15分）设求方程 0cos 2312=+-x x 根的迭代法 k k x x cos 3 2 41+=+ (1) 证明对R x ∈?0,均有*lim x x k k =∞ →,其中*x 为方程的根. (2) 此迭代法收敛阶是多少? 证明你的结论. 二、（12分）讨论分别用Jacobi 迭代法和Gauss-Seidel 迭代法求解下列方程组的收敛性。 ??? ??=++-=++=-+. 022,1, 122321 321321x x x x x x x x x 三、（8分）若矩阵??? ? ? ??=a a a a A 000002，说明对任意实数0≠a ，方程组b AX =都是非病态的。（范数用∞?） 四、（ 求)(x f 的Hermite 插值多项式)(3x H ，并给出截断误差)()()(3x H x f x R -=。 五、（10分）在某个低温过程中，函数 y 依赖于温度x （℃）的试验数据
为 已知经验公式的形式为 2bx ax y += ，试用最小二乘法求出 a ，b 。 六、（12分）确定常数 a ，b 的值，使积分 [ ] dx x b ax b a I 2 1 1 2 ),(?--+= 取得最小值。 七、（14分）已知Legendre(勒让德)正交多项式)(x L n 有递推关系式： ?? ? ? ???=+-++===-+),2,1()(1)(112)()(, 1)(1110 n x L n n x xL n n x L x x L x L n n n 试确定两点的高斯—勒让德（G —L ）求积公式 ? -+≈1 1 2211)()()(x f A x f A dx x f 的求积系数和节点，并用此公式近似计算积分 ?=2 11 dx e I x 八、（14分）对于下面求解常微分方程初值问题 ?????==0 0)() ,(y x y y x f dx dy 的单步法： ??? ? ??? ++==++=+) ,() ,()2 121(1 21211 hk y h x f k y x f k k k h y y n n n n n n
研究生《数值分析》教学大纲 课程名称：数值分析 课程编号：S061005 课程学时：64 学时 课程学分： 4 适用专业：工科硕士生 课程性质：学位课 先修课程：高等数学，线性代数，计算方法，Matlab语言及程序设计 一、课程目的与要求 “数值分析”课是理工科各专业硕士研究生的学位课程。主要介绍用计算机解决数学问题的数值计算方法及其理论。内容新颖，起点较高，并加强了数值试验和程序设计环节。通过本课程的学习，使学生熟练掌握各种常用的数值算法的构造原理和过程分析，提高算法设计和理论分析能力，并且能够根据数学模型，提出相应的数值计算方法编制程序在计算机上算出结果。力求使学生掌握应用数值计算方法解决实际问题的常用技巧。 二、教学内容、重点和难点及学时安排： 第一章? 数值计算与误差分析( 4学时) 介绍数值分析的研究对象与特点，算法分析与误差分析的主要内容。 第一节数值问题与数值方法 第二节数值计算的误差分析 第三节数学软件工具----MATLAB 语言简介 重点：误差分析 第二章? 矩阵分析基础( 10学时) 建立线性空间、赋范线性空间、内积空间的概念，为学习以后各章打好基础。矩阵分解是解决数值代数问题的常用方法，掌握矩阵的三角分解、正交分解、奇异值分解，并能够编写算法程序。 第一节? 矩阵代数基础
第二节? 线性空间 第三节? 赋范线性空间 第四节? 内积空间和内积空间中的正交系 第五节矩阵的三角分解 第六节矩阵的正交分解 第七节矩阵的奇异值分解 难点：内积空间中的正交系。矩阵的正交分解。 重点：范数，施密特(Schmidt) 正交化过程，正交多项式，矩阵的三角分解, 矩阵的正交分解。 第三章? 线性代数方程组的数值方法( 12学时) 了解研究求解线性代数方程组的数值方法分类及直接法的应用范围。高斯消元法是解线性代数方程组的最常用的直接法，也是其它类型直接法的基础。在此方法基础上加以改进，可得选主元的高斯消元法、按比例增减的高斯消元法，其数值稳定性更高。掌握用列主元高斯消元法解线性方程组及计算矩阵的行列式及逆，并且能编写算法程序。掌握矩阵的直接三角分解法：列主元LU 分解，Cholesky分解。了解三对角方程组的追赶法的分解形式及数值稳定性的充分条件。掌握矩阵条件数的定义，并能利用条件数判别方程组是否病态以及对方程组的直接方法的误差进行估计。 迭代解法是求解大型稀疏方程组的常用解法。熟练掌握雅可比迭代法、高斯- 塞德尔迭代法及SOR 方法的计算分量形式、矩阵形式，并能在计算机上编出三种方法的程序用于解决实际问题。了解极小化方法：最速下降法、共轭斜量法。迭代法的收敛性分析是研究解线性代数方程组的迭代法时必须考虑的问题。对于上述常用的迭代法，须掌握其收敛的条件。而对一般的迭代法，掌握其收敛性分析的基本方法和主要结果有助于进一步探究新的迭代法。 第一节求解线性代数方程组的基本定理 第二节高斯消元法及其计算机实现 第三节矩阵分解法求解线性代数方程组 第三节? 误差分析和解的精度改进 第四节? 大型稀疏方程组的迭代法 第五节? 极小化方法 难点：列主元高斯消元法，直接矩阵三角分解。迭代法的收敛性，雅可比迭代法，高斯-塞德尔迭代法，SOR 迭代法。
昆明理工大学2010级硕士研究生考试试卷 （注：考试时间150分钟；所有答案，包括填空题答案一律答在答题纸上，否则不予记分。） 一、 填空（每空2分，共24分） 1．近似数490.00的有效数字有 位，其相对误差限为 。 2．设7 4 ()431f x x x x =+++，则017[2,2,......2]f = ，018 [2,2,......2]f = 。 3．设4()2,[1,1]f x x x =∈-，()f x 的三次最佳一致逼近多项式为 。 4．1234A ??=??-??，1A = ，A ∞= ，2A = 。 5．210121012A -????=-????-?? ，其条件数2()Cond A = 。 6．2101202A a a ????=?????? ，为使分解T A L L =?成立（L 是对角线元素为正的下三角阵），a 的取 值范围应是 。 7．给定方程组121 122 ,x ax b a ax x b -=?? -+=?为实数。当a 满足 且02ω 时，SOR 迭代法收敛。 8．对于初值问题/ 2 100()2,(0)1y y x x y =--+=，要使用欧拉法求解的数值计算稳定，应限定步长h 的范围是 。 二、 推导计算 （15分）
（小数点后至少保留5位）。（15分） 3．确定高斯型求积公式 01 1010 ()()(),(0,1)f x d x A f x A f x x x ≈+ ∈? 的节点01,x x 及积分系数01,A A 。（15分） 三、 证明 1． 在线性方程组AX b =中，111a a A a a a a ?? ??=?????? 。证明当112a - 时高斯-塞德尔法 收敛，而雅可比法只在11 22 a - 时才收敛。 （10分） 2． 给定初值02 0, x a ≠以及迭代公式 1(2) ,(0,1,2...., 0) k k k x x a x k a +=-=≠ 证明该迭代公式是二阶收敛的。（7分） 3． 试证明线性二步法 212(1)[(3)(31)]4 n n n n n h y b y by b f b f ++++--=+++ 当1b ≠-时，方法是二阶，当1b =-时，方法是三阶的。（14分）
实验4 常微分方程的数值解 【实验目的】 1．掌握用MATLAB软件求微分方程初值问题数值解的方法； 2．通过实例用微分方程模型解决简化的实际问题； 3．了解欧拉方法和龙格-库塔方法的基本思想和计算公式，及稳定性等概念。 【实验内容】 题3 小型火箭初始重量为1400kg，其中包括1080kg燃料。火箭竖直向上发射时燃料燃烧率为18kg/s，由此产生32000N的推力，火箭引擎在燃料用尽时关闭。设火箭上升时空气阻力正比于速度的平方，比例系数为0.4kg/m，求引擎关闭瞬间火箭的高度、速度、加速度，及火箭到达最高点的时的高度和加速度，并画出高度、速度、加速度随时间变化的图形。 模型及其求解 火箭在上升的过程可分为两个阶段，在全过程中假设重力加速度始终保持不变，g=9.8m/s2。 在第一个过程中，火箭通过燃烧燃料产生向上的推力，同时它还受到自身重力（包括自重和该时刻剩余燃料的重量）以及与速度平方成正比的空气阻力的作用，根据牛顿第二定律，三个力的合力产生加速度，方向竖直向上。因此有如下二式： a=dv/dt=(F-mg-0.4v2)/m=(32000-0.4v2)/(1400-18t)-9.8 dh/dt=v 又知初始时刻t=0，v=0，h=0。记x(1)=h，x(2)=v，根据MATLAB 可以求出0到60秒内火箭的速度、高度、加速度随时间的变化情况。程序如下： function [ dx ] = rocket( t,x ) a=[(32000-0.4*x(2)^2)/(1400-18*t)]-9.8; dx=[x(2);a]; end ts=0:1:60;
x0=[0,0]; [t,x]=ode45(@rocket,ts,x0); h=x(:,1); v=x(:,2); a=[(32000-0.4*(v.^2))./(1400-18*t)]-9.8; [t,h,v,a]; 数据如下： t h v a 0 0 0 13.06 1.00 6.57 13.19 13.30 2.00 26.44 26.58 1 3.45 3.00 59.76 40.06 13.50 4.00 106.57 53.54 13.43 5.00 16 6.79 66.89 13.26 6.00 240.27 80.02 12.99 7.00 326.72 92.83 12.61 8.00 425.79 105.22 12.15 9.00 536.99 117.11 11.62 10.00 659.80 128.43 11.02 11.00 793.63 139.14 10.38 12.00 937.85 149.18 9.71 13.00 1091.79 158.55 9.02 14.00 1254.71 167.23 8.33 15.00 1425.93 175.22 7.65 16.00 1604.83 182.55 6.99 17.00 1790.78 189.22 6.36 18.00 1983.13 195.27 5.76 19.00 2181.24 200.75 5.21 20.00 2384.47 205.70 4.69 21.00 2592.36 210.18 4.22 22.00 2804.52 214.19 3.79 23.00 3020.56 217.79 3.41 24.00 3240.08 221.01 3.07 25.00 3462.65 223.92 2.77 26.00 3687.88 226.56 2.50 27.00 3915.58 228.97 2.27
郑州大学研究生课程（2012-2013学年第一学期）数值分析 Numerical Analysis 习题课 第八章常微分方程数值解法
待求解的问题：一阶常微分方程的初值问题/* Initial-Value Problem */: ?????=∈=0 )(] ,[),(y a y b a x y x f dx dy 解的存在唯一性（“常微分方程”理论）：只要f (x , y ) 在[a , b ] ×R 1 上连续，且关于y 满足Lipschitz 条件，即存在与x , y 无关的常数L 使 对任意定义在[a , b ] 上的y 1(x ) 和y 2(x ) 都成立，则上述IVP 存在唯一解。 1212|(,)(,)||| f x y f x y L y y ?≤?一、要点回顾
§8.2 欧拉(Euler)法 通常取(常数),则Euler 法的计算格式 h h x x i i i ==?+1?? ?=+=+) (),(001x y y y x hf y y i i i i i =0,1,…,n ( 8.2 )
§8.2 欧拉(Euler)法(1) 用差商近似导数 )) (,()()()()(1n n n n n n x y x hf x y x y h x y x y +=′+≈+?? ?=+=+) (),(01a y y y x hf y y n n n n 差分方程初值问题向前Euler 方法h x y x y x y n n n ) ()()(1?≈ ′+)) (,() ()(1n n n n x y x f h x y x y ≈?+))(,()(n n n x y x f x y =′
数值分析（研究生，2008-12-15） 1.（10分）求函数???≤≤++<≤-+=1 0,101,1sin )(2x x x x x x f 在区间[-1，1]上的最佳平方逼近式 x e a x a a x 210)(++=φ。 2．（15分）利用乘幂法计算下列矩阵的主特征值和相应的特征向量 ???? ??????----110141012，初始向量为T x ]0,0,1[0=（要求结果有三位有效数字）。同时计算该矩阵的1-条件数和谱条件数。
3.（15分）已知函数x x f sin )(=在36.0,3 4.0,32.0210===x x x 处的值分别为352274.0,333487.0,314567.0210===y y y 。用Lagrange 插值多项式对3167.0=x 的函数值进行近似计算，并估计近似计算的误差界。
4.（15分）用Newton 迭代法求方程0ln 2=+x x 在区间（0，2 π）内的解，选择你认为合适的初始点，计算方程的根，使得近似解具有四位有效数字。请从理论上估计达到所需精度所需的迭代次数。
5.（15分）用Gauss-Seidel 迭代法解方程组 ?????? ????-=????????????????????---542834*********x x x 取初始近似向量0[0,0,0]T x =，估计达到4位有效数字需要的迭代次数，并实际计算之。就该具体问题分析计算过程中总的乘除法计算量。
6. （10分）应用拟牛顿法解非线性方程组 ?????=-+=-+. 12,2322112221x x x x x x 取T x ]1,0[)0(= ，终止容限210-=ε。 7.（10分） 求解矛盾方程组 ???????=++=++=++=++2 32328.12221321321 321321x x x x x x x x x x x x
2011年秋季工学硕士研究生学位课程（数值分析）真题试卷B (总分：28.00，做题时间：90分钟) 一、填空题(总题数：6，分数：12.00) 1.填空题请完成下列各题，在各题的空处填入恰当的答案。（分数： 2.00） __________________________________________________________________________________________ 解析： 2.设｜x｜>>1______ （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 3.求积分∫ a b f（x）dx的两点Gauss公式为______ （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 4.设∞ =______，‖A‖ 2 =______． （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 5.给定f(x)=x 4，以0为三重节点，2为二重节点的f(x)的Hermite插值多项式为______． （分数：2.00） __________________________________________________________________________________________ 正确答案：(正确答案：x 4) 解析： 6.己知差分格式r≤______时，该差分格式在L ∞范数下是稳定的． （分数：2.00） __________________________________________________________________________________________ 正确答案：() 解析： 二、计算题(总题数：2，分数：4.00) 7.给定方程lnx-x 2+4=0，分析该方程存在几个根，并用迭代法求此方程的最大根，精确至3位有效数字．（分数：2.00） __________________________________________________________________________________________ 正确答案：(正确答案：令f(x)=lnx-x 2 +4，则f"(x)= -2x,当x= 时，f"(x)=0. 注意到 f(0．01)=-0．6053＜0,f(1)=3＞0,f(3)=-3．9014＜0,而当时，f"(x)＞0，当时，f"(x)＜
数值分析（研究生，2008-12-15） （ 分）求函数???≤≤++<≤-+=1 0,101,1sin )(2x x x x x x f 在区间?? ， 上的最佳平方逼近式 x e a x a a x 210)(++=φ。 ．（ 分）利用乘幂法计算下列矩阵的主特征值和相应的特征向量 ???? ??????----110141012，初始向量为T x ]0,0,1[0=（要求结果有三位有效数字）。同时计算该矩阵的 条件数和谱条件数。
（ 分）已知函数x x f sin )(=在36.0,34.0,32.0210===x x x 处的值分别为352274.0,333487.0,314567.0210===y y y 。用????????插值多项式对3167.0=x 的函数值进行近似计算，并估计近似计算的误差界。
（ 分）用??????迭代法求方程0ln 2=+x x 在区间（ ，2 π）内的解，选择你认为合适的初始点，计算方程的根，使得近似解具有四位有效数字。请从理论上估计达到所需精度所需的迭代次数。
?（ 分）用??◆????????●迭代法解方程组 ?????? ????-=????????????????????---542834*********x x x 取初始近似向量0[0,0,0]T x =，估计达到 位有效数字需要的迭代次数，并实际计算之。就该具体问题分析计算过程中总的乘除法计算量。
? （ 分）应用拟牛顿法解非线性方程组 ?????=-+=-+. 12,2322112221x x x x x x 取T x ]1,0[)0(= ，终止容限210 -=ε。 （ 分） 求解矛盾方程组 ???????=++=++=++=++2 32328 .12221 321321321321x x x x x x x x x x x x
硕士研究生 《数值分析》练习题 一、判断题 1、用Newton 切线法求解非线性线性方程可以任选初值。 （ ） 2、求解非线性线性方程，Newton 切线法比弦截法迭代次数多。 （ ） 3、若n n A R ?∈非奇异，用Jacobi 迭代法求解线性方程组Ax b =必收敛。（ ） 4、Lagrange 插值法与Newton 插值法得到同一个插值多项式。 （ ） 二、填空题 1、近似数 3.14108937a =关 于π具 位有效数字。 2、双点弦截法具有 阶收敛速度。 3、求方程x x e =根的单点弦截法迭代公式是 。 4、设2112A ?? = ? ?? ? ，则()A ρ= 。 5、若(),0,1,2,3i l x i =是以01231,3,,x x x x ==为插值节点的Lagrange 插值基函数，则()()3 3012i i i x l =-=∑ 。 6、由下数据表确定的代数插值多项式的不超过 次。 7、若()8754321f x x x x =+-+，则差商[]0,1,2,,8f = 。 8、拟合三点()()()0,1,1,3,2,2A B C 的 直线是y = 。 三、分析与计算题 1、设()14,2,3515T A x -??==-?? -?? ，求∞=,2,1,,p x A p p 和()1A cond 。
2、1001012,20253A x -???? ? ? == ? ? ? ?-???? ，试计算p p x A ,，p=1,2,∞，和1)(A c o n d 。 3、线性方程组,0Ax b b =≠，用Jacobi 迭代法是否收敛，为什么？其中 122111221A -?? ?=-- ? ?--?? 。 4、线性方程组,0Ax b b =≠，用Jacobi 迭代法是否收敛，为什么？其中 2-11=11111-2A ?? ???? ???? 。 5、已知函数表如下: ⑴ ()111.75ln11.75L ≈、估计截断误差并说明结果有几位有效数字； ⑵ ()211.75ln11.75N ≈、估计截断误差并说明结果有几位有效数字。 6、已知函数表 如下: ⑴用Lagrange 插值法求ln 0.55的近似值()10.55N 、估计截断误差并说明结果的有效数字； ⑵用Newton 插值法求ln 0.55的近似值()20.55N 、估计截断误差并说明结果的有效数字。 7、已知数据如下，求满足条件的Hermite 插值多项式。
第一章绪论 1.1 引言 常微分方程是现代数学的一个重要分支，是人们解决各种实际问题的有效工具。微分方程的理论和方法从17世纪末开始发展起来，很快成了研究自然现象的强有力工具，在17到18世纪，在力学、天文、科学技术、物理中，就已借助微分方程取得了巨大的成就。1864年Leverrer根据这个方程预见了海王星的存在，并确定出海王星在天空中的位置。现在，常微分方程在许多方面获得了日新月异的应用。这些应用也为常微分方程的进一步发展提供了新的问题，促使人们对微分方程进行更深入的研究，以便适应科学技术飞速发展的需要。 研究常微分方程常用数值解是数学工作者的一项基本的且重要的工作。在国内外众多数学家的不懈努力，使此学科基本上形成了一套完美的体系。微分方程的首要问题是如何求一个给定方程的通解或特解。到目前为止，人们已经对许多微分方程得出了求解的一般方法。由于在生产实际和科学研究中所遇到的微分方程问题比较复杂，使这些问题的解即使能求出解析表达式，也往往因计算量太大而难于求出，而对于一些典型的微分方程则可以运用基本方法求出其解析解，并可以根据初值问题的条件把其中的任意常数确定下来。 由于求通解存在许多困难，人们就开始研究带某种定解条件的特解。首先是Cauchy对微分方程初始解的存在惟一性进行了研究。目前解的存在惟一性、延拓性、大范围的存在性以及解对初始解和参数的延续性和可微性等理论问题都已发展成熟。与此同时，人们开始采取各种近似方法来求微分方程的特解，例如求微分方程数值解的Euler折线法、Runge-Kutta法等，可以求得若干个点上微分方程的近似解。最后，由于当代高科技的发展为数学的广泛应用和深入研究提供了更好的手段。用计算机结合Matlab软件求方程的精确解、近似解，对解的性态进行图示和定性、稳定性研究都十分方便有效。 本章先介绍常微分的一般概念、导出微分方程的一些典型例子及求解微分方程的思路分析。从而得到常微分方程的常用数值解法。
研究生《数值分析》课程作业（二） 姓名： 学号： 专业： 1、据如下函数值表，建立二次的Lagrange 插值多项式及Newton 插值多项式。 20012222()()()()()()() (1)(2)(0)(2)(-0)(1)59 3143 (01)(02)(10)(12(20)(21)22 L x f x l x f x l x f x l x x x x x x x x x =++-----=? +?+?=-+------解： 二次 l agr ange插值 ） Newton 插值多项式： 200100120122()()[,](-)[,,](-)(-) 5559 32(0)(0)(1)32()3 2222 N x f x f x x x x f x x x x x x x x x x x x x x x =++=-?-+--=-+-=-+ ()y f x =2、已知单调连续函数在如下采样点处的函数值 *()0[2,4],f x x =求方程在内根的近似值使误差尽可能小。 解：1 ()()y f x x f y -==解： 对的反函数进行二次插值
1110201122012010210122021(0)(0)(0)(0)(0)(0) (0)() ()() ()()()()()() (0 2.25)(05)(03)(05)(03)(0 2.25) 2 3.54( 3 2.25)(35)(2.253)(2.255)(53)(5 2.25) y y y y y y L f y f y f y y y y y y y y y y y y y ---------=++--------+-+-=? +?+? ----+-+- 2.945 ≈()(1)01(1)1()[,]()(,),()[,],() ()()()() (1)! ,n n n n n n n n f x a b f x a b a x x x b L x x a b f R x f x L x x n a b x ξωξ+++≤<<<≤∈=-=+∈ 3、证明：设在上连续，在内存在，节点是满足拉格朗日插值条件的多项式，则对任何插值余项 这里（）且依赖于。 0110101(0,1,,)()()0()()()()()()()()[,]()()()()()()() (),,,(k n n k n n n n n n x k n R x R x R x K x x x x x x x K x x K x x x a b t f t L t K x t x t x t x t x x x x t ωφφφ+===---==----- 证由条件知节点是的零点，即。于是其中是与有关的待定函数。 现把看成上的固定点，作函数 根据插值条件和余项定义，知在点及处均为零。故明：1111)[,]2()[,]1()()[,]()(,)(,),()()(1)!()0 ()()(,),(1)! n n n n a b n t a b n t t a b n t a b a b f n K x f K x a b x n φφφφξφξξξξ++++'+'''+∈=-+==∈+（） （） （）（）在上有个零点，根据罗尔定理，在内至少有个零点。对再应用罗尔定理，可知在内至少 有个零点。依次类推，在上至少有一个零点，记为 使 于是 ， 且依赖于于是得到插值余项。 证毕。 44、试用数据表建立不超过次的埃尔米特插值多项式。 解：（用重节点的均差表建立埃尔米特多项式）
2014级硕士研究生试卷 科目： 数值分析 考试时间： 出题教师： 集体 考生姓名： 专业： 学号： 不予计分；可带计算器。 一、 填空题（每空2分，共30分） 1．设14.30=x 是准确值21.30=* x 的近似值，则近似值x 有 位有效数字，近 似值x 的相对误差为 。 2．函数)(x f 过点(0,1), (1,3)和(2,9)，对应的基函数分别为)(),(),(210x l x l x l ，过这三个节点的二次拉格朗日插值多项式为 ，余项为 。 3. 已知0)1(,3)1(,0)2(=-==f f f ，二阶均差]1,1,2[-f = 。 4．方程012 3 =--x x 在5.10 =x 附近有个根，构造不动点迭代收敛的格式 为 ,若用牛顿法迭代求根，其收敛阶是 。 5．设???? ? ??=2021012a a A ，为了使A 可分解成T LL A =，其中L 是对角元素为正的下三角矩阵， 则a 的取值范围 。 6． 设????? ??-----=232221413A ，??? ? ? ??-=111x ，则∞||||Ax ，1||||A = ， 2||||A = 。 7．设U L D A --=,b Ax =的Gauss-Seidel 迭代的矩阵形式b Ux Lx Dx k k k ++=++)()1() 1(， 其迭代矩阵为 ，该迭代格式收敛的充要条件__________________。 8．求解一阶常微分方程初值问题?? ???=<<-=1)0(1 0,2' y x y x y y ，取步长1.0=h 的Euler 法公式为 ，其截断误差的首项为 。
“传热流动的数值分析”2015年大作业 1. 2维条件下的无粘、不可压缩流体通过出口和入口流过箱体，具体情况如图所示，求该箱体内的流线情况，所有单位为厘米。 （1）流线方程为：22220x y ψψ ??+=?? 使用Gauss-Seidel 线迭代，0.1x y ?=?=，误差0.005ξ=，结果输出中，包括在y=0,1,2,3,4,5 处的所有X 处对应的流函数值。 （2）设出口处纵向速度V =0，试采用PSOR 方法，0.1x y ?=?=，计算在x=0,1,2,3,4,5 处的所有Y 处对应的流函数值，以及不同的松弛系数和迭代次数的关系曲线（至少三个系数）。 答：（1）该问题为稳态问题，流线方程为椭圆型方程，在求解方程时，首先对方程在计算域内进行离散。计算域为：{}(,)05,05x y x y Ω=≤≤≤≤，在离散时，0.1x y h ?=?=?=，因此可以得到流线方程的差分方程为： 1,,1,,1,,1 22 220i j i j i j i j i j i j h h ψψψψψψ+-+--+-++=?? （1） 整理后可得： 1,1,,1,1 ,4 i j i j i j i j i j ψψψψψ+-+-+++= （2） 在本题中，采用Gauss-Seidel 线迭代方法进行求解，扫描方向选为自左向右，此时有 111,1,11,1,1,4 n n n n i j i j i j i j n i j ψψψψψ ++++--+++++= （3） 由于是线推进，因此在当前线方向求解时，之前的扫描线上的参数已经得到更新，所以 方程可改写为： 11,1,11 1,1,,111444n n i j i j n n n i j i j i j ψψψψψ+-++++-++-+-= （4） 其中1 1,n i j ψ+-看着当前迭代层中的已知变量。