微分方程数值解

1、数学与计算科学学院实 验 报 告实验项目名称 常微分方程数值解 所属课程名称 微分方程数值解法 实 验 类 型 验证 实 验 日 期 班 级 学 号 姓 名 成 绩 一、实验概述：【实验目的】掌握求解常微分方程的欧拉法，Runge-Kutta方法（二阶的预估校正法，经典的四阶R-K法）【实验原理】1，欧拉法的格式如下： yn+1=yn+hfxn,yn y0=y(x0)其中n从0开始取值，fxn,yn为微分方程中的dy/dx的值，而y0是我们根据题目设定的初值。2，预估校正法格式如下：yn+10=yn+hfxn,yn yn+1=yn+h2fxn,yn+f(xn+1,yn+10)y0=f(x0)

2、其中第一个式子为预报格式，第二个式子为校正格式，第三个式子也是给定的初值。3，经典四阶Runge-Kutta格式，亦即R-K法格式如下：yn+1=yn+h6（k1+2k2+2k3+k4）k1=f(xn,yn)k2=f(xn+h2,yn+hk12)k3=f(xn+h2,yn+hk22)k4=f(xn+h,yn+hk3)【实验环境】Microsoft Visual C+6.0二、实验内容：【实验方案】用欧拉法，预估校正法, 经典的四阶龙格库塔方法求解初值问题 dy/dx=，初值y（0）=1； 将计算结果与精确解为比较在区间0,1上分别取步长h=0.1; 0.05时进行计算。对三个算法的收敛性进行分

3、析比较，【实验过程】（实验步骤、记录、数据、分析）1，部分算法代码：(1)欧拉算法代码void ruler(double xn+1,double yn+1,double h)/欧拉法 int j;for(j=0;j<n;j+) yj+1=yj+h*f(xj,yj); (2)预步校正法代码void yugu_jz(double xn+1,double yn+1,double h)/预估校正法double yyn+12;int j;yy01=y0; for(j=0;j<n;j+) yyj+10=yyj1+h*f(xj,yyj1); yyj+11=yyj1+h/2*(f(xj,yyj1)

4、+f(xj,yyj+10);/预估校正法的解 yj+1=yyj+11; (3)Runge_Kutta法代码void Runge_Kutta(double xn+1,double yn+1,double h)/Runge_Kutta法double k1,k2,k3,k4;int i;for(i=0;i<n;i+)k1=f(xi,yi);k2=f(xi+h/2,yi+h/2*k1);k3=f(xi+h/2,yi+h/2*k2);k4=f(xi+h,yi+h*k3);yi+1=yi+h/6*(k1+2*k2+2*k3+k4);(4)目标函数及解析解double f(double a,doubl

5、e b)/目标函数double fun;fun=a*exp(-a)-b;return fun;void precise(double xn+1,double yn+1)/精确解int j;for(j=0;j<n;j+) yj+1=(pow(xj+1,2)+2)*exp(-xj+1)/2; 2，结果分析：根据下面的实验结果，可以得出：当步长取0.1时，三种算法收敛性的关系为：欧拉法<预步校正法<Runge_Kutta法，相应的预步校正法精确度为普遍在10-2，略高于欧拉法，Runge_Kutta法的精确度达到10-7，远远优于欧拉法和预步校正法。当步长取0.05时，三种算法收敛

6、性的关系不变，但三种方法的精确度均有提高，相对于步长0.1时，Runge_Kutta法收敛速度高于另外两种，可得其误差项阶数高于欧拉法和预步校正法，即结果均与理论相符。【实验结论】（结果）当步长为0.1时，结果及误差如下：当步长为0.05时，结果及误差如下：【实验小结】（收获体会） 通过本次实验，我熟悉了欧拉法，预步校正法以及Runge_Kutta法的算法思想，还有它们之间的收敛性差异，精确度差异，同时对也提高了自己的上机编程能力，感受到了理论与实践相结合的乐趣。三、指导教师评语及成绩：评 语评语等级优良中及格不及格1.实验报告按时完成,字迹清楚,文字叙述流畅,逻辑性强2.实验方案设计合理3.

7、实验过程（实验步骤详细,记录完整,数据合理,分析透彻）4实验结论正确. 成 绩： 指导教师签名： 批阅日期：附录1：源 程 序#include<stdio.h>#include<math.h>#define n 10void precise(double xn+1,double yn+1);/精确解void ruler(double xn+1,double yn+1,double h);/欧拉法void yugu_jz(double xn+1,double yn+1,double h);/预估校正法void Runge_Kutta(double xn+1,double

8、yn+1,double h);/Runge_Kutta法double f(double a,double b);/定义函数void main() printf("*该程序解决微分方程y=x*e(-x)-y的数值解问题*n"); printf("*区间为0,1，步长h=0.1nn"); double xn+1,y4n+1; double h; int i,k; h=0.1;/步长赋值 /初值赋值 x0=0.0; for(i=0;i<4;i+) yi0=1; for(i=1;i<n+1;i+) xi=x0+h*i; precise(x,y0);

9、ruler(x,y1,h); yugu_jz(x,y2,h); Runge_Kutta(x,y3,h); printf(" x 精确解 欧拉法的解 预估校正法 Runge_Kuttan"); for(k=0;k<n+1;k+)/打印结果 printf(" %-3.2lf %-12.8lf %-12.8lf %-12.8lf %-12.8lfn",xk,y0k,y1k,y2k,y3k); printf(" x 精确解 欧拉法误差 预校法误差 R_K法误差n"); for(k=0;k<n+1;k+)/打印误差 printf(&

10、quot; %-3.2lf %-12.8lf %e %e %en",xk,y0k,fabs(y1k-y0k),fabs(y2k-y0k),fabs(y3k-y0k); void precise(double xn+1,double yn+1)/精确解int j;for(j=0;j<n;j+) yj+1=(pow(xj+1,2)+2)*exp(-xj+1)/2; void ruler(double xn+1,double yn+1,double h)/欧拉法 int j;for(j=0;j<n;j+) yj+1=yj+h*f(xj,yj); void yugu_jz(dou

11、ble xn+1,double yn+1,double h)/预估校正法double yyn+12;int j;yy01=y0; for(j=0;j<n;j+) yyj+10=yyj1+h*f(xj,yyj1); yyj+11=yyj1+h/2*(f(xj,yyj1)+f(xj,yyj+10);/预估校正法的解 yj+1=yyj+11; void Runge_Kutta(double xn+1,double yn+1,double h)/Runge_Kutta法double k1,k2,k3,k4;int i;for(i=0;i<n;i+)k1=f(xi,yi);k2=f(xi+h

12、/2,yi+h/2*k1);k3=f(xi+h/2,yi+h/2*k2);k4=f(xi+h,yi+h*k3);yi+1=yi+h/6*(k1+2*k2+2*k3+k4);double f(double a,double b)/目标函数double fun;fun=a*exp(-a)-b;return fun;附录2：实验报告填写说明 1实验项目名称：要求与实验教学大纲一致。2实验目的：目的要明确，要抓住重点，符合实验教学大纲要求。3实验原理：简要说明本实验项目所涉及的理论知识。4实验环境：实验用的软、硬件环境。5实验方案（思路、步骤和方法等）：这是实验报告极其重要的内容。概括整个实验过程。对于验证性实验，要写明依据何种原理、操作方法进行实验，要写明需要经过哪几个步骤来实现其操作。对于设计性和综合性实验，在上述内容基础上还应该画