数值分析常微分数值解的求法C语言Word下载.docx

上传人:b****3 文档编号:17049278 上传时间:2022-11-28 格式:DOCX 页数:12 大小:45.34KB
下载 相关 举报
数值分析常微分数值解的求法C语言Word下载.docx_第1页
第1页 / 共12页
数值分析常微分数值解的求法C语言Word下载.docx_第2页
第2页 / 共12页
数值分析常微分数值解的求法C语言Word下载.docx_第3页
第3页 / 共12页
数值分析常微分数值解的求法C语言Word下载.docx_第4页
第4页 / 共12页
数值分析常微分数值解的求法C语言Word下载.docx_第5页
第5页 / 共12页
点击查看更多>>
下载资源
资源描述

数值分析常微分数值解的求法C语言Word下载.docx

《数值分析常微分数值解的求法C语言Word下载.docx》由会员分享,可在线阅读,更多相关《数值分析常微分数值解的求法C语言Word下载.docx(12页珍藏版)》请在冰豆网上搜索。

数值分析常微分数值解的求法C语言Word下载.docx

以表格形式列出10个等距节点上的计算值和精确值,并比较他们的计算精确度。

其中多步法需要的初值由经典R-K法提供。

运算时,取足以表示计算精度的有效位。

2.算法说明

函数及变量说明

表1函数及变量定义

函数/变量名

类型

说明

f

double

要求的微分方程

p_f

微分方程原函数

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<

>

//微分方程

doublef(doublex,doubley)

{

return(-y+cos(2*x)-2*sin(2*x)+2*x*exp(-x));

//returnx*pow(y,-2)*;

}

//原函数

doublep_f(doublex)

returnx*x*exp(-x)+cos(2*x);

//returnpow+pow(x,2),;

//求精确解

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)

y[0]=y0;

for(i=1;

{

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]));

}

//改进欧拉方法

double*Euler_Pro(doublea,doubleb,doubleh,doubley0)

doublec=(b-a)/h,yb;

yb=y[i-1]+h*f(a+(i-1)*h,y[i-1]);

y[i]=y[i-1]+*h*(f(a+(i-1)*h,y[i-1])+f(a+i*h,yb));

//经典K-R方法

double*K_R(doublea,doubleb,doubleh,doubley0)

doubleK1,K2,K3,K4,x;

doublec=(b-a)/h;

x=a+(i-1)*h;

K1=f(x,y[i-1]);

K2=f(x+*h,y[i-1]+*h*K1);

K3=f(x+*h,y[i-1]+*h*K2);

K4=f(x+h,y[i-1]+h*K3);

y[i]=y[i-1]+h*(K1+2*K2+2*K3+K4)/6;

//4阶Adams预测-校正方法

double**Adams(doublea,doubleb,doubleh,doubley0)

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];

//动态初始化,储存预测值和校值

=count;

r[i]=newdouble[2];

if(r==NULL)

printf("

内存分配失败\n"

);

exit(0);

count;

r[i][0]=0;

r[i][1]=0;

y=K_R(a,b,h,y0);

4;

x[i]=a+h*i;

dy[i]=f(x[i],y[i]);

r[i][0]=y[i];

i=3;

while(x[3]<

b)

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=,y0=1;

while

(1)

请输入下限,上限,步长,初值:

(用空格隔开)\n"

scanf("

%lf%lf%lf%lf"

&

a,&

b,&

h,&

y0);

double*yx,*y;

yx=Exact(a,b,h);

1、欧拉2、改进欧拉3、经典K-R4、4阶Adams预测-校正5、修改数据6、退出\n"

%d"

selet);

system("

cls"

if(selet==1)

y=Euler(a,b,h,y0);

欧拉方法:

\n"

if(selet==2)

y=Euler_Pro(a,b,h,y0);

改进欧拉方法:

if(selet==3)

经典K-R方法(4阶龙格-库塔方法):

if(selet==4)

double**r;

r=Adams(a,b,h,y0);

4阶Adams预测-校正方法:

x值预估值校正值误差\n"

%%lf--\n"

0,r[0][0]);

,r[1][0]);

,r[2][0]);

for(i=3;

=10;

%%f%f%\n"

a+i*h,r[i][0],r[i][1],fabs(r[i][0]-r[i][1]));

if(selet==1||selet==2||selet==3)

x值计算值准确值误差\n"

%%lf%lf%\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值

精确解

2、经典K-R方法

3、4阶adams预测-校正法

5.对比分析

欧拉方法

改进欧拉方法

经典K-R方法

4阶adams预测-校正方法

通过表格可以看出:

1,在步长相同的时候,结果准确度经典K-R>

改进欧拉>

欧拉方法>

4阶adams预测校正方法

2,同一个方法的情况下,步长越小结果准确度越高;

参考文献:

[1]朱建新,李有法,数值计算方法(第三版),高等教育出版社;

[2]Adam预测-校正系统

学生学习心得

学生(签名):

2016年6月17日

诚信承诺

本人郑重声明所呈交的课程报告是本人在指导教师指导下进行的研究工作及取得的研究成果。

据我所知,除了文中特别加以标注的地方外,论文中不包含其他人已经发表或撰写过的研究成果。

与我一同工作的同学对本文研究所做的贡献均已在报告中作了明确的说明并表示谢意。

任课

教师

评语

成绩评定:

任课教师(签名):

年月日

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 医药卫生 > 基础医学

copyright@ 2008-2022 冰豆网网站版权所有

经营许可证编号:鄂ICP备2022015515号-1