第七章常微分方程数值解 本章介绍求解微分方程数值解的基本思想和方法。 许多科学和技术问题的数学模型,常常归结为一个 含有自变量、未知函数和它的一阶导数和高阶导数 的方程,称为常微分方程。它是描述运动、变化规 律的重要数学方法之 分为两类: 1.初值问题,即给出未知函数及导数在初始点的值 的问题; 2.边值问题,即给出未知函数及(或)它的某些导 数在区间两个端点的值的问题。 ●常徼分方程初值问题 这里仅讨论一阶微分方程初值问题的数值解 法,即,所要讨论的问题为 y'=f(x,y),a<x≤b,|y|k∞ y(a)=yo, y0∈R
第七章 常微分方程数值解 本章介绍求解微分方程数值解的基本思想和方法。 许多科学和技术问题的数学模型，常常归结为一个 含有自变量、未知函数和它的一阶导数和高阶导数 的方程，称为常微分方程。它是描述运动、变化规 律的重要数学方法之一。 分为两类： 1.初值问题，即给出未知函数及导数在初始点的值 的问题； 2.边值问题，即给出未知函数及（或）它的某些导 数在区间两个端点的值的问题。 ⚫常微分方程初值问题 这里仅讨论一阶微分方程初值问题的数值解 法，即，所要讨论的问题为： = = ( ) , , ( , ), , | | y a y0 y0 R y f x y a x b y
其中∫为x,y的已知函数,y为给定的初值。 ●初值问题的适定性(其解是否唯一存在) 将区域:[小]×R记为G,即G=[ab]×R 设f:G→>R为连续映射,若存在常数L>0使得不 等式|f(x,y)-f(xy2)≤Ly-y2 对一切(x,y1)(x,y2)∈G都成立,则称f(x,y)在G 上关于y满足 Lipschitz条件,而式中的常数L称 为 Lipschitz常数。 补充定理:初值问题 y’=f(x,y),a<x≤b,|yk<∞ y(a)=yo ∈R 当f(x,y)在G上连续,且关于y满足 Lipschitz 条件,则其解存在且唯一. 本章,总是假设所讨论的问题满足解的存在唯一性 定理条件。 ●数值解法
其中 f 为 x y, 的已知函数， 0 y 为给定的初值。 ⚫ 初值问题的适定性(其解是否唯一存在) 将区域： a b R , 记为 G ，即 G a b R = , 设 f G R : → 为连续映射，若存在常数 L 0 使得不 等式 1 2 1 2 f x y f x y L y y ( , ) ( , ) − − 对一切 1 2 ( , ),( , ) x y x y G 都成立，则称 f x y ( , ) 在 G 上关于 y 满足 Lipschitz 条件，而式中的常数 L 称 为 Lipschitz 常数。 补充定理：初值问题 = = ( ) , , ( , ), , | | y a y0 y0 R y f x y a x b y ， 当 f (x, y) 在 G 上连续，且关于 y 满足 Lipschitz 条件，则其解存在且唯一. 本章，总是假设所讨论的问题满足解的存在唯一性 定理条件。 ⚫数值解法
寻求微分方程初值间题的解y(x)在一系列离散点 q=X<x1<x2<…<xx=b 上的近似值y,y,y2,…,yN的方法 n}:问题的数值解。 数值解所满足的离散方程统称为差分格式. 步长:h=x1-x1,一般取定步长h=h §1 Euler方法 ● Euler公式 将初值问题 y’=f(x,y),a<x≤b,|yk<∞ ∈R 的求解区间[a,bN等分, 分点xn=a+mh,(n=0,2…,N), b 其中 h N
寻求微分方程初值问题的解 y(x) 在一系列离散点 a = x0 x1 x2 x N = b 上的近似值 , , , , N y y y y 0 1 2 的方法. yn : 问题的数值解。 数值解所满足的离散方程统称为差分格式. 步长： i = i − i−1 h x x ，一般取定步长 hi = h §1 Euler 方法 ⚫ Euler 公式 将初值问题 = = ( ) , , ( , ), , | | y a y0 y0 R y f x y a x b y 的求解区间 [a,b] N 等分， 分点 x a nh, (n 0,1,2, ,N) n = + = , 其中 N b a h − =
将y(xn+1)在x,进行二阶 Tay lor展开 (n+d=y(n)+y(n(xn+1-xn)+2on(xm+1-xn) y(n)+hy(rn)+y(sn) 5n E(rn, xn+p y(n)+hf(n,y(n))+y(sn) 2 p y(+i=y(rn)+hf(n,y())+y(sn) h 记:7=y(n),截去T得近似计算公式: y(n+day(n)+hf(xn,y(n)) 可得y(xn)的近似值yn的递推公式 yn+=yn+hf(en, yn), n=0, 1, .,N-1 (*1) 上式称为 Euler公式 初值问题给定了y0,代入(*1),依次解出数值解y
将 1 ( ) n y x + 在 n x 进行二阶 Taylor 展开： 2 1 1 1 ''( ) ( ) ( ) '( )( ) ( ) 2! n n n n n n n n y y x y x y x x x x x + + + = + − + − 2 ( ) '( ) ''( ) 2 n n n h = + + y x hy x y 1 ( , ) n n n x x + 2 ( ) ( , ( )) ''( ) 2 n n n n h = + + y x hf x y x y 即 1 ( ) n y x + 2 ( ) ( , ( )) ''( ) 2 n n n n h = + + y x hf x y x y 记： 2 1 ''( ) 2 n h T y = ，截去 T1 得近似计算公式： ( ) ( ) ( , ( )) n 1 n n n y x y x + hf x y x + . 可得 ( ) n y x 的近似值 n y 的递推公式 yn+1 = yn + hf (xn , yn ), n = 0,1, , N −1 (*1) 上式称为 Euler 公式。 初值问题给定了 0 y ，代入(*1)，依次解出数值解 1 y ， ， N y
7称为 Euler方法的局部截断误差。 例1用 Euler公式解初值问题 2x y=y 0<x≤1 y(0)=1 解:取h=0.1, Euler公式的具体形式为 n+1 yn+hf(n,yn)=yn+h(n-n) 其中xn=mh=0.1n,(n=0,1,…10), 已知y=1,则有 2. )=1+0.1 2x 0.2 y2=yi+h(y )=1.1+0.1(1.1-x,)=1.191818 依次计算可得 y3,y4,y5,y6,y,y8,y9,少
T1 称为 Euler 方法的局部截断误差。 例 1 用 Euler 公式解初值问题 = = − (0) 1 , 0 1 2 y x y x y y 解：取 h = 0.1 ，Euler 公式的具体形式为 1 2 ( , ) ( ) n n n n n n n n x y y hf x y y h y y + = + = + − 其中 x = nh = 0.1n, n = 0,1, ,10) n ( ， 已知 0 y =1 ，则有 0 1 0 0 0 2 ( ) 1 0.1 1.1 x y y h y y = + − = + = 1 2 1 1 1 2 0.2 ( ) 1.1 0.1(1.1 ) 1.191818 1.1 x y y h y y = + − = + − = … … … 依次计算可得 3 4 5 6 7 8 9 10 y , y , y , y , y , y , y , y
其部分结果见下表 数值准确解误差 解yny(x) 021.1918181.1832160.008602 041.3582131.3416410016572 0.61.5089661.4832400025726 0.8 6497831.6124520.037331 1.01.7847701.7320510.052719 其中,准确解为:y=√1+2x。 可见 Euler方法的计算结果精度不太高。 几何意义:用折线近似代替方程的真解曲线。故 该方法常常也称为欧拉折线法。见P161图7-1 隐式(后退的) Euler公式 将y(xn)在xn+1进行 Taylor展开: (xn)=y(xn+1)+y(xn+1)(xn-xn+1)+ n+1 n+1 n+1 n+1) y"(5n+1)
其部分结果见下表 n x 数值 解 n y 准确解 ( ) n y x 误差 1.0 0.8 0.6 0.4 0.2 1.784770 1.649783 1.508966 1.358213 1.191818 1.732051 1.612452 1.483240 1.341641 1.183216 0.052719 0.037331 0.025726 0.016572 0.008602 其中，准确解为： y = 1+ 2x 。 可见 Euler 方法的计算结果精度不太高。 几何意义：用折线近似代替方程的真解曲线。故 该方法常常也称为欧拉折线法。见 P161 图 7－1。 ⚫ 隐式(后退的)Euler 公式 将 ( ) n y x 在 n 1 x + 进行 Taylor 展开： 1 2 1 1 1 1 ''( ) ( ) ( ) '( )( ) ( ) 2! n n n n n n n n y y x y x y x x x x x + = + − + − + + + + 2 1 1 1 ( ) '( ) ''( ) 2 n n n h y x hy x y = − + + + +
n+1∈(xn,xn+1),令mn=5n+1,则,mn∈(xn,xn+1) y(x b2b2 ,)=y(xn+1)-h(xn+1,y(xn+1)+y(mn) 即y(xn+1)=y(xn)+hf(xn+1y(xn+1) y"(n) h2 令T y"(7n),截去72,则 y(n+1=y(xn)+hf(xn+1,y (n+D)) 得到数值求解公式: yn+1=yn+ hf(rn+1,yn+1 (*2) (*2)式称为隐式 Euler公式。即要求出yn+t必需要 解方程,这样的公式叫做隐式公式,而 Euler公式 是显式公式。 T2称为隐式 Euler方法的局部截断误差 ●方法的阶 1.定义 若局部截断误差(将准确解y(x)代入公式的
1 1 ( , ) n n n x x + + ，令 n n = +1 ，则， 1 ( , ) n n n x x + ( ) n y x 2 1 1 1 ( ) ( , ( )) ''( ) 2 n n n n h = − + y x hf x y x y + + + 即 1 ( ) n y x + 2 1 1 ( ) ( , ( )) ''( ) 2 n n n n h = + − y x hf x y x y + + 令 2 2 ''( ) 2 n h T y = − ，截去 T2 ，则 1 1 1 ( ) ( ) ( , ( )) n n n n y x y x hf x y x + + + = + 得到数值求解公式： 1 1 1 ( , ) n n n n y y hf x y + + + = + (*2) (*2)式称为隐式 Euler 公式。即要求出 n 1 y + 必需要 解方程，这样的公式叫做隐式公式，而 Euler 公式 是显式公式。 T2 称为隐式 Euler 方法的局部截断误差 ⚫ 方法的阶 1．定义: 若局部截断误差（将准确解 y x( ) 代入公式的
左、右两端,其左端与右端之差) T=y(xn1)-yn+1=O(h+), 则称该数值方法具有p阶精度。p越大,精度越高, 数值方法越好。 ① Euler方法的精度 y(n+D=y(n)+hf(n,y(n))+y(sn) T=y(n+D-ym+1=y(n+1)-Ly(rn)+hf(rn,y(n)) h 2 y(Sn n5n+1 即:71=O(h2) ∴ Euler方法具有一阶精度。 ②后退的 Euler方法的精度 . y(n+1=y(n)+ hf(n+1, y(n +D-y(nn) 2=y(n++=y(n+i-Ly(n)+ hf(xm+,y(xn+D))
左、右两端，其左端与右端之差） 1 1 1 ( ) ( ) p T y x y O h n n + = − = + + ， 则称该数值方法具有 p 阶精度。 p 越大，精度越高， 数值方法越好。 ①Euler 方法的精度 1 ( ) n y x + 2 ( ) ( , ( )) ''( ) 2 n n n n h = + + y x hf x y x y 1 1 1 1 ( ) ( ) [ ( ) ( , ( ))] T y x y y x y x hf x y x = − = − + n n n n n n + + + 2 ''( ) 2 n h = y 1 ( , ) n n n x x + 即： 2 1T O h = ( ) Euler 方法具有一阶精度 。 ②后退的 Euler 方法的精度 1 ( ) n y x + 2 1 1 ( ) ( , ( )) ''( ) 2 n n n n h = + − y x hf x y x y + + 2 1 1 1 1 1 ( ) ( ) [ ( ) ( , ( ))] T y x y y x y x hf x y x = − = − + n n n n n n + + + + +
h y(n)7n∈(xn2xn+1) 即: O(h2) 后退的 Euler方法也是具有一阶精度 ●改进的 Euler方法 h 考虑到 Euler方法的局部截断误差为y(n), 而后退的 Euler方法的局部截断误差为 h2 3ym),其中n,mn的不一定相同。 我们将 Euler公式:yn+1=yn+bf(xn,yn) 和后退的 Euler公式:yn+1=yn+h(xmn+1,ymn+1) 两端求平均值: m4l=n+I(r,, y,)+f( 1, n-1(*3) (*3)式称为平均(梯形)公式
2 ''( ) 2 n h = − y 1 ( , ) n n n x x + 即： 2 2 T O h = ( ) 后退的 Euler 方法也是具有一阶精度 ⚫ 改进的 Euler 方法 考虑到 Euler 方法的局部截断误差为 2 ''( ) 2 n h y ， 而 后 退 的 Euler 方 法 的 局 部 截 断 误 差 为 2 ''( ) 2 n h − y ，其中 n ， n 的不一定相同。 我们将 Euler 公式： 1 ( , ) n n n n y y hf x y + = + 和后退的 Euler 公式： 1 1 1 ( , ) n n n n y y hf x y + + + = + 两端求平均值： 1 1 1 ( , ) ( , ) 2 n n n n n n h y y f x y f x y + + + = + + (*3) (*3)式称为平均（梯形）公式
其局部截断误差为: T=y(n+1-v(n)+elf(xn, y(rn))+f(n+1,y(n+DI h y(n+1-y(n-y(rn)-y(+) 2 将y(xn+1)及y(xn+1)均在点xn处三阶 Taylor展开 1;(xn)6 y(n+=y(n)+hy'(xn)+oy h )+-y"(xn) y(n+=y(rn)+ hy(n)+y(n)+y →T=0(h 这种方法具有二阶精度 平均(梯形)公式 Yn+=yn+alf(n,yn)+f(n+i,]n+) 为隐式公式,一般用迭代法求解,迭代初值由 Euler 公式提供,只迭代一次即得如下预测一校正型公式 预测yn+1=yn2+bf(xn2yn) 校正yn1=yn+=[f(xnyn)+f(xn+13n+1)
其局部截断误差为： 1 1 1 ( ) { ( ) [ ( , ( )) ( , ( ))]} 2 n n n n n n h T y x y x f x y x f x y x = − + + + + + 1 1 ( ) ( ) '( ) '( ) 2 2 n n n n h h y x y x y x y x = − − − + + 将 1 ( ) n y x + 及 1 ( ) n y x + 均在点 n x 处三阶 Taylor 展开 2 3 1 ( ) ( ) ( ) ( ) '''( ) 2 6 n n n n n h h y x y x hy x y x y x + = + + + 2 3 (4) 1 ( ) ( ) ( ) ( ) ( ) 2 6 n n n n n h h y x y x hy x y x y x + = + + + T = 3 O h( ) 这种方法具有二阶精度 平均（梯形）公式 1 1 1 ( , ) ( , ) 2 n n n n n n h y y f x y f x y + + + = + + 为隐式公式，一般用迭代法求解，迭代初值由 Euler 公式提供，只迭代一次即得如下预测－校正型公式 1 1 1 1 ( , ) [ ( , ) ( , )] 2 n n n n n n n n n n y y hf x y h y y f x y f x y + + + + = + = + + 预测 校正