MATLAB改进欧拉法与四阶龙格-库塔求解一阶常微分方程
姓名:樊元君 学号:02 日期:.11.06
一、实验目的
掌握改进欧拉法与四阶龙格-库塔求解一阶常微分方程的初值问题。二、实验内容
分别写出改进欧拉法与四阶龙格-库塔求解的算法,编写程序上机调试出结果,要求所编程序适用于任何一阶常微分方程的数值解问题,即能解决这一类问题,而不是某一个问题。
实验中以下列数据验证程序的正确性。
求步长h=0.25。实验注意事项
的精确解为,通过调整步长,观察结果的精度的变化function [] = GJOL(h,x0,y0,X,Y)
format long
h=input('h=');
x0=input('x0=');
y0=input('y0=');
disp('输入的范围是:');
X=input('X=');Y=input('Y=');
n=round((Y-X)/h);
i=1;x1=0;yp=0;yc=0;
for i=1:1:n
x1=x0+h;
yp=y0+h*(-x0*(y0)^2);%yp=y0+h*(y0-2*x0/y0);%
yc=y0+h*(-x1*(yp)^2);%yc=y0+h*(yp-2*x1/yp);%
y1=(yp+yc)/2;
x0=x1;y0=y1;
y=2/(1+x0^2);%y=sqrt(1+2*x0);%
fprintf('结果=%.3f,%.8f,%.8f\n',x1,y1,y);
end
end
四阶龙格库塔程序源代码:
function [] = LGKT(h,x0,y0,X,Y)
format long
h=input('h=');
x0=input('x0=');
y0=input('y0=');
disp('输入的范围是:');
X=input('X=');Y=input('Y=');
n=round((Y-X)/h);
i=1;x1=0;k1=0;k2=0;k3=0;k4=0;
for i=1:1:n
x1=x0+h;
k1=-x0*y0^2;%k1=y0-2*x0/y0;%
k2=(-(x0+h/2)*(y0+h/2*k1)^2);%k2=(y0+h/2*k1)-2*(x0+h/2)/(y0+h/2*k1);%
k3=(-(x0+h/2)*(y0+h/2*k2)^2);%k3=(y0+h/2*k2)-2*(x0+h/2)/(y0+h/2*k2);%
k4=(-(x1)*(y0+h*k3)^2);%k4=(y0+h*k3)-2*(x1)/(y0+h*k3);%
y1=y0+h/6*(k1+2*k2+2*k3+k4);%y1=y0+h/6*(k1+2*k2+2*k3+k4);%
x0=x1;y0=y1;
y=2/(1+x0^2);%y=sqrt(1+2*x0);%
fprintf('结果=%.3f,%.7f,%.7f\n',x1,y1,y);
end
end
五、运行结果:
改进欧拉格式结果:
四阶龙格库塔结果:
步长分别为:0.25和0.1时,不同结果显示验证了步长减少,对于精度的提高起到很大作用,有效数字位数明显增加。
六、实验小结:
通过这次实验学习,首先第一点对改进欧拉格式和四阶龙格库塔的原理推导有了深入的理解,改进欧拉格式采用(预报+校正)模式得到较精确的原函数数值解;而四阶龙格库塔则采用多预报几个点的斜率值,采用加权平均作为平均斜率的近似值的思想达到更高精度的数值解,二阶龙格库塔的特例就是改进后的欧拉格式。
七、思考题如何对四阶龙格-库塔法进行改进,以保证结果的精度。
,如果某次计算结果精度,,反复将步长减半,直到,这时的可作为结果。
(2)步长加倍。对于给定精度,如果某次计算结果精度,,反复将步长加倍,直到,这时的可作为结果。
6