数值分析常微分数值解的求法C语言.docx
《数值分析常微分数值解的求法C语言.docx》由会员分享,可在线阅读,更多相关《数值分析常微分数值解的求法C语言.docx(13页珍藏版)》请在冰点文库上搜索。
数值分析常微分数值解的求法C语言
本科生课程设计报告
实习课程
数值分析
学院名称
管理科学学院
专业名称
信息与计算科学
学生姓名
学生学号
指导教师
实验地点
实验成绩
二〇一六年六月二〇一六年六
摘要
常微分方程数值解法是计算数学的一个分支.是解常微分方程各类定解问题的数值方法.现有的解析方法只能用于求解一些特殊类型的定解问题,实用上许多很有价值的常微分方程的解不能用初等函数来表示,常常需要求其数值解.所谓数值解,是指在求解区间内一系列离散点处给出真解的近似值.这就促成了数值方法的产生与发展.?
关键词:
数值解法;常微分
1.实验内容和要求
常微分方程初值问题
有精确解
。
要求:
分别取步长h=0.1,0.01,0.001,采用改进的Euler方法、4阶经典龙格-库塔R-K方法和4阶Adams预测-校正方法计算初值问题。
以表格形式列出10个等距节点上的计算值和精确值,并比较他们的计算精确度。
其中多步法需要的初值由经典R-K法提供。
运算时,取足以表示计算精度的有效位。
2.算法说明
2.1函数及变量说明
表1函数及变量定义
函数/变量名
类型
说明
f
double
要求的微分方程
p_f
double
微分方程原函数
Euler,Euler_Pro
K_R,Adams
double
欧拉,改进欧拉,经典K_R,4阶Adams预测-校正方法
Main
Void
主函数
1、欧拉方法:
(i=0,1,2,3,......n-1)
(其中a为初值)
2、改进欧拉方法:
(i=0,1,2......n-1)(其中a为初值)
3、经典K-R方法:
4、4阶adams预测-校正方法
预测:
校正:
Adsms内插外插公式联合使用称为Adams预测-校正系统,利用外插公式计算预测,用内插公式进行校正。
计算时需要注意以下两点:
1、外插公式为显式,内插公式为隐式。
故用内插外插公式时需要进行迭代。
2、这种预测-校正法是四步法,计算Yn+1时,不但用到前一步信息,而且要用到更前三步信息
,
,因此它不是自动开始的,实际计算时必须借助某种单步法,用Runge-Kutta格式为它提供初始值
,依据局部截断误差公式得:
进一步将Adams预测-校正系统加工成下列方案:
运用上述方案计算
时,要先一步的信息
和更前一步的信息
因此在计算机之前必须给出初值
和
,
可用其他单步法来计算,
则一般令他为零。
3.源程序
#include
#include
#include
#include
//微分方程
doublef(doublex,doubley)
{
return(-y+cos(2*x)-2*sin(2*x)+2*x*exp(-x));
//returnx*pow(y,-2)*2.0/3.0;
}
//原函数
doublep_f(doublex)
{
returnx*x*exp(-x)+cos(2*x);
//returnpow(1.0+pow(x,2),1.0/3.0);
}
//求精确解
double*Exact(doublea,doubleb,doubleh)
{
inti;
doublec=(b-a)/h;
double*y=newdouble[c+1];
for(i=0;i<=c;i++)
y[i]=p_f(a+i*h);
returny;
}
//欧拉方法
double*Euler(doublea,doubleb,doubleh,doubley0)
{
inti;
doublec=(b-a)/h;
double*y=newdouble[c+1];
y[0]=y0;
for(i=1;i<=c;i++)
{
y[i]=y[i-1]+h*f(a+(i-1)*h,y[i-1]);
//printf("%f\n",f(a+(i-1)*h,y[i-1]));
}
returny;
}
//改进欧拉方法
double*Euler_Pro(doublea,doubleb,doubleh,doubley0)
{
inti;
doublec=(b-a)/h,yb;
double*y=newdouble[c+1];
y[0]=y0;
for(i=1;i<=c;i++)
{
yb=y[i-1]+h*f(a+(i-1)*h,y[i-1]);
y[i]=y[i-1]+0.5*h*(f(a+(i-1)*h,y[i-1])+f(a+i*h,yb));
}
returny;
}
//经典K-R方法
double*K_R(doublea,doubleb,doubleh,doubley0)
{
doubleK1,K2,K3,K4,x;
inti;
doublec=(b-a)/h;
double*y=newdouble[c+1];
y[0]=y0;
for(i=1;i<=c;i++)
{
x=a+(i-1)*h;
K1=f(x,y[i-1]);
K2=f(x+0.5*h,y[i-1]+0.5*h*K1);
K3=f(x+0.5*h,y[i-1]+0.5*h*K2);
K4=f(x+h,y[i-1]+h*K3);
y[i]=y[i-1]+h*(K1+2*K2+2*K3+K4)/6;
}
returny;
}
//4阶Adams预测-校正方法
double**Adams(doublea,doubleb,doubleh,doubley0)
{
inti;
double*y;
doublecount=(b-a)/h;
doubledy[100]={0},x[4]={0};
doublep_0=0,p_1=0,c=0;
double**r=newdouble*[count+1];//动态初始化,储存预测值和校值
for(i=0;i<=count;i++)
r[i]=newdouble[2];
if(r==NULL)
{
printf("内存分配失败\n");
exit(0);
}
for(i=0;i{
r[i][0]=0;
r[i][1]=0;
}
y=K_R(a,b,h,y0);
for(i=0;i<4;i++)
{
x[i]=a+h*i;
dy[i]=f(x[i],y[i]);
r[i][0]=y[i];
}
i=3;
while(x[3]
{
x[3]=x[3]+h;
p_1=y[3]+h*(55*dy[3]-59*dy[2]+37*dy[1]-9*dy[0])/24;//预估
c=p_1+251*(c-p_0)/270;//修正
r[i][0]=c;
c=f(x[3],c);//求f
c=y[3]+h*(9*c+19*dy[3]-5*dy[2]+dy[1])/24;//校正
y[3]=c-19*(c-p_1)/270;//修正
r[i][1]=y[3];
dy[0]=dy[1];
dy[1]=dy[2];
dy[2]=dy[3];
dy[3]=f(x[3],y[3]);
p_0=p_1;
i++;
}
returnr;
}
voidmain()
{
inti,selet;
doublea=0,b=2,h=0.001,y0=1;
while
(1)
{
printf("请输入下限,上限,步长,初值:
(用空格隔开)\n");
scanf("%lf%lf%lf%lf",&a,&b,&h,&y0);
doublec=(b-a)/h;
double*yx,*y;
yx=Exact(a,b,h);
while
(1)
{
printf("1、欧拉2、改进欧拉3、经典K-R4、4阶Adams预测-校正5、修改数据6、退出\n");
scanf("%d",&selet);
system("cls");
if(selet==1)
{
y=Euler(a,b,h,y0);
printf("欧拉方法:
\n");
}
if(selet==2)
{
y=Euler_Pro(a,b,h,y0);
printf("改进欧拉方法:
\n");
}
if(selet==3)
{
y=K_R(a,b,h,y0);
printf("经典K-R方法(4阶龙格-库塔方法):
\n");
}
if(selet==4)
{
double**r;
r=Adams(a,b,h,y0);
printf("4阶Adams预测-校正方法:
\n");
printf("x值预估值校正值误差\n");
printf("%0.3f%lf--\n",0,r[0][0]);
printf("%0.3f%lf--\n",0.1,r[1][0]);
printf("%0.3f%lf--\n",0.2,r[2][0]);
for(i=3;i<=10;i++)
printf("%0.3f%f%f%0.20f\n",a+i*h,r[i][0],r[i][1],fabs(r[i][0]-r[i][1]));
}
if(selet==1||selet==2||selet==3)
{
printf("x值计算值准确值误差\n");
for(i=0;i<=10;i++)
printf("%0.3f%lf%lf%0.20f\n",(a+i)*h,y[i],yx[i],fabs(yx[i]-y[i]));
}
if(selet==5)
break;
if(selet==6)
exit
(1);
continue;
}
}
}
4.实验结果
1、改进欧拉方法
步长
x值
0.1
0.01
0.001
精确解
2、经典K-R方法
步长
x值
0.1
0.01
0.001
精确解
3、4阶adams预测-校正法
步长
x值
0.1
0.01
0.001
精确解
5.对比分析
欧拉方法
改进欧拉方法
经典K-R方法
4阶adams预测-校正方法
0.1
0.01
0.001
0.000000000000005088
通过表格可以看出:
1,在步长相同的时候,结果准确度经典K-R>改进欧拉>欧拉方法>4阶adams预测校正方法
2,同一个方法的情况下,步长越小结果准确度越高;
参考文献:
[1]朱建新,李有法,数值计算方法(第三版),高等教育出版社;
[2]Adam预测-校正系统
学生学习心得
学生(签名):
2016年6月17日
诚信承诺
本人郑重声明所呈交的课程报告是本人在指导教师指导下进行的研究工作及取得的研究成果。
据我所知,除了文中特别加以标注的地方外,论文中不包含其他人已经发表或撰写过的研究成果。
与我一同工作的同学对本文研究所做的贡献均已在报告中作了明确的说明并表示谢意。
学生(签名):
任课
教师
评语
成绩评定:
任课教师(签名):
年月日