实验五-常微分方程数值解.doc

实验五常微分方程数值解————————————————————————————————作者：————————————————————————————————日期：个人收集整理勿做商业用途个人收集整理勿做商业用途个人收集整理勿做商业用途实验五常微分方程数值解欧拉法算法说明对于xi,i=0，1,2,…，n,取步长h为定值时,有h=xi+1-xi，EURLER法的计算公式为：yi+1=yi+h*f(xi,yi）；i=0,1，2,…,n。程序中主要符号说明a为x的下界,b为x的上界，h为步长，n为循环次数（即为x的数值点数减一），x0、y0为循环初始值，x1、y1为输出值。计算流程框图开始读入数据x0，y0，h，nFORI=1TONX0+H=>X1,Y0+H*F(X0,Y0)=Y1输出X1，Y1X1=>X0Y1=>Y0NEXT结束程序清单symsabhnx0y0x1y1；a=0;b=0.1；n=5;x0=0;y0=1;h=（b—a)/n；forj=1:nx1=x0+h;y1=y0-h＊0。9*y0/(1+2＊x0);x1y1x0=x1;y0=y1；end算例及输入数据说明算例：求解初值问题当x=0，0.02,0.04，…，0.10时的数值解。输入数据说明:a=0；a为x的下界b=0。1；b为x的上界n=5;n为循环次数,即为x的数值点数减一x0=0；x（0）的值y0=1；y（0)的值程序运行结果及结果分析运行结果：x1=0.0200y1=0.9820x1=0.0400y1=0.9650x1=0。0600y1=0。9489x1=0.0800y1=0。9337x1=0。1000y1=0.9192结果分析:欧拉法计算简单，但计算效率并不高，计算精度很低，局部截断误差较大.改进欧拉法算法说明对于xi,i=0,1,2,…，n，取步长h为定值时，有h=xi+1-xi，EURLER法的计算公式为:yp=yi+h＊f(xi，yi）yc=yi+h*f(xi+1，yp）yi+1=（yp+yc)/2；i=0,1,2，…，n。程序中主要符号说明a为x的下界，b为x的上界，h为步长，n为循环次数（即为x的数值点数减一)，x0、y0为循环初始值，yp、yc为运算中间值，x1、y1为输出值.计算流程框图开始读入数据x0，y0，h，nFORI=1TONX0+H=>X1Y0+H*F(X0,Y0)=>YPY0+H*F(X0,Y0)=>YC(YP+YC)/2=>Y1输出X1，Y1X1=>X0Y1=>Y0NEXT结束程序清单symsabhnx0y0ypycx1y1；a=0;b=0。1;n=5；x0=0;y0=1;yp=1;yc=1；h=(b—a)/n;forj=1：nx1=x0+h;yp=y0—h＊0。9*y0/(1+2*x0)；yc=y0—h＊0.9＊yp/（1+2*x1);y1=（yp+yc)/2;x1y1x0=x1;y0=y1;end算例及输入数据说明算例：求解初值问题当x=0，0。02,0.04，…，0。10时的数值解。输入数据说明:a=0；a为x的下界b=0.1；b为x的上界n=5；n为循环次数，即为x的数值点数减一x0=0;x(0)的值y0=1;y（0）的值yp=1;由y(0)的值决定yc=1;由y（0）的值决定程序运行结果及结果分析x1=0.0200y1=0.9825x1=0。0400y1=0。9660x1=0.0600y1=0.9503x1=0。0800y1=0。9354x1=0。1000y1=0。9212结果分析:计算过程比欧拉法较复杂，但改进欧拉法先用欧拉法求出预报值,再利用公式求出校正值，局部截断误差比欧拉法低了一阶，较大程度地提高了计算精度.龙格库塔法算法说明对于xi，i=0，1，2，…，n，取步长h为定值时，有h=xi+1-xi,EURLER法的计算公式为：yi+1＝yi＋h＊(K1+2＊K2+2*K3+K4）/6K1=f(xi,yi)K2=f(xi+h/2,yi+h*K1/2）K3=f（xi+h/2,yi+h*K2/2）K4=f(xi+h,yi+h*K3);i=0,1，2，…，n。程序中主要符号说明h为步长,n为循环次数（即为x的数值点数减一)，x0、y0为循环初始值,k1，k2，k3,k4为运算中间值，x1、y1为输出值。计算流程框图开始读入数据x0，y0，h，nFORI=1TONX0+H=>X1F(X0,Y0)=>K1F(X0+H/2,Y0+K1*H/2)=>K2F(X0+H/2,Y0+K2*H/2)=>K3F(X0+H,Y0+K3*H)=>K4Y0+(K1+2K2+2K3+K4)=>Y1输出X1，Y1X1=>X0Y1=>Y0NEXT结束程序清单symshnx0y0x1y1k1k2k3k4;x0=0;y0=1；