计算方法4_常微分方程数值解法.ppt

对象 思路 计算过程： 4.1 Euler 公式 ?思想: 用向前差商近似代替微商. 几何意义 欧拉法（续） ?思想: 用向后差商近似代替微商. 欧拉法（续） 欧拉法（续） 截断误差 Def4.1 设y(xn) 是(4.1)式的精确解，yn是(4.2)式欧拉法得到的近似解，称y(xn)-yn为欧拉法的整体截断误差. 二步欧拉法的局部截断误差 Def4.5 假设yn=y(xn) , yn－1=y(xn－1),称Rn+1=y(xn+1)-yn+1为二步欧拉法的局部截断误差. 数值积分法 ? 梯形公式:将显示欧拉公式,隐式欧拉公式平均可得 从而y(xn+1)=yn+1+h fy(xn+1,η)[y(xn+1)-yn+1] /2 +O(h3) ∴y(xn+1)-yn+1 = h fy(xn+1,η)[y(xn+1)-yn+1] /2 +O(h3) ∴y(xn+1)-yn+1=O(h3)/[1-hfy(xn+1,η)/2]=O(h3) ∴梯形公式的截断误差为O(h3)，其精度是2阶。 欧拉公式的比较 4.2 改进的Euler法 4.3龙格—库塔方法 为了构造函数φ使得(4.5)式成为高阶方法，Taylor 令 Runge-Kutta法的思想：在［xn,xn+1］内多预报几个点的ki值并用其加权平均值作为k*近似值。而构造出具有更高精度的计算公式。 二阶龙格—库塔方法 上述方程组有四个未知量，只有三个方程，有无穷多组解。 取任意一组解便得一种二阶龙－库公式。 要使其具有三阶精度，必须使局部截断误差为O(h4) 类似可推出四阶龙格-库塔公式，常用的有 例：经典Runge-Kutta法 Step 1: 输入a,b,y0 及N Step 2: (b-a)/N=h,a=x,y0=y Step 3: 输出 (x,y) Step 4: For I=1 T0 N hf(x,y)=k1 hf(x+h/2,y+ k1/2)= k2 hf(x+h/2,y+k2/2)=k3 hf(x+h,y＋k3)=k4 y+(k1+2k2+2k3+k4)/6=y x+h=x 输出(x,y) Step 5 : 停止 自适应龙格－库塔法 用户提出问题I ： 定理：对于问题I若用P阶龙格－库塔法计算y(xn+1)在步长折半前后的近似值分别为yn+1(h), yn+1(h/2)则有误差公式 4.5 收敛与稳定性 显式单步法： yn+1=yn+hφ(xn,yn,h) φ(x,y,h)为增量函数，它依赖于f,仅是xn,yn, h的函数 改进Euler法的收敛性 4.5.2 稳定性 Euler法的稳定性 Euler法：yn+1=yn+hf(xn,yn) 考察模型方程 y’=λy，(λ0) 即yn+1=(1+hλ)yn 隐式Euler法稳定性 隐式Euler法：yn+1=yn+hf(xn+1,yn+1) 对于模型方程 y’=λy(λ0) 则yn+1=yn+hλyn+1 = yn+1=yn/(1-hλ) 梯形公式稳定性 4.6 多步法 多步法的截断误差 4.7 常微分方程组和高阶微分方程的数值解法 4.8 常微分方程边值问题的数值解法 打靶法（Shooting Method） 有限差分法(Finite-Difference Methods) 当c1=c2=1/2, a2=b21=1时二阶Runge-Kutta公式为 yn+1=yn+k1/2+k2/2 k1=hf(xn,yn) k2=hf(xn+h,yn+k1) 此即改进Euler法 取c2=0 ，c2=1，a2=1/2，b21=1/2 yn+1=yn+k2 k1=hf(xn,yn) k2=hf(xn+h/2,yn+k1/2) 此为中点法或变形的 Euler公式 三阶龙格－库塔法是用三个值k1,k2,k3的加权平均来近似k*， 即有 yn+1=yn+c1k1+c2k2+c3k3 k1=hf(xn,yn) k2=hf(xn+a2h,yn+b21k1) k3=hf(xn+a3h,yn+b31k1+b32k2) 类似二阶龙格－库塔法的推导，c1,c2,c3,a2,a3,b21,b31,b32应满足 c1+c2+c3=1 a2=b21 a3=b31+b32 c2a2+c3a3=1/2 c2a22+c3a32=1/3 c3a32=1/6 由该方程组任意解可得三阶龙格-库塔公式 例：Kutta公式 kn+1=yn+(k1+4k2+k3)/6 k1=hf(xn,yn) k2=hf(xn+h/2,yn+k1/2) k3=hf(xn+h,yn-k1+2k2) yn+1=y