数值计算方法实验4

1实验报告学院（系）名称：姓名学号专业班级实验项目常微分方程的数值解法课程名称数值计算方法课程代码实验时间实验地点批改意见：成绩教师签字：实验目的：编写四阶龙格—库塔方法算法程序，上机调试通过。实验环境：硬件环境：IBM－PC或兼容机软件环境：Windows操作系统编程语言：C语言实验内容：一．用四阶龙格－库塔公式求解常微分方程二．用经典龙格－库塔公式求解初值问题取步长h=0.1，进行计算实验步骤：一．四阶龙格-库塔法步骤1.在给定区间[xn,xn+1]内，寻找四个考查点，以便求得四个斜率值，即2.根据预报－校正系统以及加权平均的思想，可得如下的原始四阶龙格－库塔法计算公式3.常用的四阶公式4.截断误差为O(h5)龙格－库塔四阶法的特点1.龙格－库塔四阶公式精度较高，可满足一般工程计算的要求2.每次计算yn+1时，只用到前一步的计算结果yn，因此在已知y0的条件下，可以自动的进行计算3.可以在计算过程中随时改变步长h4.缺点是每前进一步需要多次调用函数f(x,y)，工作量较大，并且误差不容易估计2.0,6.20.21)0.2(hxyyxdxdy1)0(0)0(10,4'23'zyxzyzzyy43211110krqxkqpxkpxkxrnqnpnn对应斜率为对应斜率为对应斜率为对应斜率为))(,())(,(),(),()(3214213121443322111ckbkakrhyrhxfkvkukqhyqhxfkphkyphxfkyxfkkkkkhyynnnnnnnnnn342312143211,2,212,21),()22(6hkyhxfkkhyhxfkkhyhxfkyxfkkkkkhyynnnnnnnnnn2二．经典龙格-库塔法步骤（两个未知函数）思想与四阶龙格-库塔法相同，公式如下：经典龙格-库塔法（两个未知函数）的流程图与龙格库塔法相同，只是输入部分增加Z0，将循环部分的公式进行替换，添加Z的函数和L1到L4即可，输出部分增加Zi即可。实验结果：一．四阶龙格-库塔法输入想x，y的初值和步长h以及循环次数n后，输出x1=2.2y1=1.356505x2=2.4y2=1.661361x3=2.6y3=1.939104二．经典龙格-库塔法（两个未知函数）输入想x，y，z的初值和步长h以及循环次数n后，输出x1=0.1y1=0.247867z1=1.152704x6=0.6y6=6.505331z6=7.054143x2=0.2y2=0.632872z2=1.451603x7=0.7y7=10.859663z7=11.356249x3=0.3y3=1.246186z3=1.987004x8=0.8y8=18.024563z8=18.473892x4=0.4y4=2.237883z4=2.908204x9=0.9y9=29.823739z9=30.230309x5=0.5y5=3.855161z5=4.461692x10=1y10=49.263345z10=49.631225四阶龙格-库塔法流程图：334334223223112112114321143211,,,,21,2,2121,2,212,2,212,2,21),,(),,()22(6)22(6hLzhkyhxgLhLzhkyhxfkhLzkhyhxgLhLzkhyhxfkLhzkhyhxgLLhzkhyhxfkzyxgLzyxfkLLLLhzzkkkkhyynnnnnnnnnnnnnnnnnnnnnnnnnnnn3附录（源程序及运行结果）：一．四阶龙格-库塔法#includestdio.h#includemath.hdoublef(doublex,doubley){return(x/y);}voidmain(){doublex0,y0,h,k1,k2,k3,k4,x1,y1;intn=0,N;printf(请输入x和y的初值：);scanf(%lf,%lf,&x0,&y0);printf(请输入步长h：);scanf(%lf,&h);printf(请输入n的值：);scanf(%d,&N);while(n!=N){x1=x0+h;k1=f(x0,y0);k2=f(x0+h/2,y0+h*k1/2);k3=f(x0+h/2,y0+h*k2/2);k4=f(x1,y0+h*k3);y1=y0+h*(k1+2*k2+2*k3+k4)/6;printf(x%d=%lf,y%d=%lf\n,n+1,x1,n+1,y1);n=n+1;x0=x1;y0=y1;}}运行结果：二．经典龙格-库塔法（两个未知函数）#includestdio.h#includemath.hdoublef(doublex,doubley,doublez){return(3*y+2*z);}doubleg(doublex,doubley,doublez){return(4*y+z);}voidmain(){doublex0,y0,z0,h,k1,k2,k3,k4,l1,l2,l3,l4,x1,y1,z1;intn=0,N;printf(请输入x,y,z的初值：);scanf(%lf,%lf,%lf,&x0,&y0,&z0);printf(请输入步长h：);scanf(%lf,&h);printf(请输入n的值：);scanf(%d,&N);while(n!=N){x1=x0+h;k1=f(x0,y0,z0);l1=g(x0,y0,z0);k2=f(x0+h/2,y0+h*k1/2,z0+h*l1/2);l2=g(x0+h/2,y0+h*k1/2,z0+h*l1/2);k3=f(x0+h/2,y0+h*k2/2,z0+h*l2/2);l3=g(x0+h/2,y0+h*k2/2,z0+h*l2/2);k4=f(x1,y0+h*k3,z0+h*l3);l4=g(x1,y0+h*k3,z0+h*l3);y1=y0+h*(k1+2*k2+2*k3+k4)/6;z1=z0+h*(l1+2*l2+2*l3+l4)/6;printf(x%d=%lf,y%d=%lf,z%d=%lf\n,n+1,x1,n+1,y1,n+1,z1);n=n+1;x0=x1;y0=y1;z0=z1;}}运行结果：