Adams预测-校正系统

2012-2013（1）专业课程实践论文
Adams 预测-校正系统
高强强，0818180212，R数学08-2班
王岳，0818180108，R数学08-1班
一、算法理论
显示格式和隐式格式都具有四阶精度，两种格式可匹配成下列Adams 预测-校正
系统：
Adams 预测 - 校正
预测 ()3-2-1-n 19-3759f -5524
n n n n n f f f h y y ++=+, ()111,+++=n n n y x f f ;
校正 ()
2-1-115-19924n n n n n n f f f f h y y +++=++ ; ()111,+++=n n n y x f f .
Adams 内插外插公式联合使用称作 Adams 预测 - 校正系统，用外插公式计算
预测，用内插公式进行校正。

计算时需要注意以下两点：
1 、外插公式为显式，内插公式为隐式。

故用内插公式计算需要进行迭代。

2 、这种预测 - 校正法是四步法，计算1+n y 时，不但用到前一步的信息n x ，n y ，
而且要用到更前三步的信息 1-n f ,2-n f ,3-n f ，因此它不是自开始的，实际计算时必
须借助于某种单步法，用Runge-Kutta 格式为它提供开始值 1y ,2y ,3y 。

依据局部截断误差公式得
()()()n n n x y h y x y 5511720
251-≈++, ()()()n n n x y h y x y 5511720
19-≈
++. 进一步将Adams 预测-校正系统加工成下列方案：
预测 ()3-2-1-19-3759-5524
n n n n n n f f f f h y p ++=+, 改进 ()n n n n p c p m -72025111+=++, 计算 ()111,+++=n n n m x f m ,
校正 ()2-1-115-19924n n n n n n f f f m h y c ++'+
=++, 改进 ()1111
-72019-++++=n n n n p c c y , 计算 ()111,+++='n n n
y x f y ; 运用上述方案计算1+n y 时，要用到先一步的信息n y ，n
y '，n n c p -和更前一步的
信息1-n y ，因此在计算机之前必须给出初值1y 和11-c p ，1y 可用其他单步法来计
算，11-c p 则一般令它等于零。

开始
输入区间（a,b),步长h,初始值
x[0],y[0]
调用RK函数求出y[1],y[2,y[3]]
计算dy[1],dy[2],dy[3]
开始迭代
p_0=0
x[3]<b
y
预估值p_1,修正值c,并打印预估值，调用函数f求出dy[i]校正值c，求出修正值y[i],输出y[i],p_0=p_1
x[3]=x[3]+h
j结束
#include<stdio.h>
#include<math.h>
double f(double x,double y)
{
double z;
z=-y+cos(2*x)-2*sin(2*x)+2*x*exp(-x);
return(z);
}
void RK(double x[4],double y[4],double h) //用经典4阶R-K法计算前三个点，以得到较好的初值
{
int i;
double k[4];
for (i=1;i<4;i++)
{
k[0]=h*f(x[i-1],y[i-1]);
k[1]=h*f(x[i-1]+h/2,y[i-1]+k[0]/2);
k[2]=h*f(x[i-1]+h/2,y[i-1]+k[1]/2);
k[3]=h*f(x[i-1]+h,y[i-1]+k[2]);
y[i]=y[i-1]+(k[0]+2*k[1]+2*k[2]+k[3])/6;
x[i]=x[i-1]+h;
printf("y(%lf)=%.15lf\n",x[i],y[i]);
}
}
main()
{
int i;
double a,b;
double h;
double c=0;
double p_0=0,p_1;
double x[4],y[4],dy[100];
double t;
printf("请输入区间（a,b）:\n");
scanf("%lf,%lf",&a,&b);
printf("请输入步长h:\n");
scanf("%lf",&h);
printf("请输入初值（x_0,y_0）:\n");
scanf("%lf,%lf",&x,&y); //计时开始
RK(x,y,h); //使用R-K法开始，计算前三个点
for (i=0;i<4;i++)
dy[i]=f(x[i],y[i]);
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; //修正printf("预测值y(%lf)=%.15lf\n",x[3],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; //修正printf("校正值y(%lf)=%.15lf\n",x[3],y[3]); //输出y //为利用前一次计算结果并节省存储空间，赋值
dy[0]=dy[1];
dy[1]=dy[2];
dy[2]=dy[3];
dy[3]=f(x[3],y[3]);
p_0=p_1;
}
}
四、算法实现
利用Adams 预测 - 校正系统，求解初值问题。

例1.求解初值问题
⎩
⎨⎧=++=1 y(0))exp(-**2)*sin(2*2-)*cos(2y - y'x x x x 解：运行程序
(1) 显示出 请输入求解区间()b a ,：输入b a ,的值分别为0,1,回车。

(2) 显示出 请输入步长h ：输入h 的值为0.1，回车。

(3) 显示出 请输入初值()0_,0_y x ：输入0_,0_y x 的值分别为0，1，回车。

(4) 显示结果：y(0.100000)=0.989114945734689
y(0.200000)=0.953810207163844
y(0.300000)=0.892009234874488
y(0.400000)=0.803953952833856
y(0.500000)=0.691930591306842
y(0.600000)=0.559925160841388
y(0.700000)=0.413289029050956
y(0.800000)=0.258366269650948
y(0.900000)=0.102115037648165
y(1.000000)=-0.048270950301088 y(1.100000)=-0.185729663537347
例2. 求解初值问题 ⎩⎨⎧==0
y(0) y - 1 y' 解：运行程序
(1) 显示出 请输入求解区间()b a ,：输入b a ,的值分别为0,1,回车。

(2) 显示出 请输入步长h ：输入h 的值为0.1，回车。

(3) 显示出 请输入初值()0_,0_y x ：输入0_,0_y x 的值分别为0，0，回车。

(4) 显示结果：y(0.100000)=0.095162500000000
y(0.200000)=0.181269098593750
y(0.300000)=0.259181577998822
y(0.400000)=0.329679857927070
y(0.500000)=0.393469225089630
y(0.600000)=0.451188240636286
y(0.700000)=0.503414567069739
y(0.800000)=0.550670902638676
y(0.900000)=0.593430205058103
y(1.000000)=0.632120423225642
y(1.100000)=0.667128781599515。