#include<iostream>
#include<math.h>
using namespace std;
double fax(double x,double y)
{
return y-2*x/y;
}
double Euler(double h,double x,double y)
{
return y+h*fax(x,y);
}
double RungKutta(double h,double x,double y)
{
double k1,k2,k3,k4;
k1=fax(x,y);
k2=fax(x+h/2,y+k1*h/2);
k3=fax(x+h/2,y+k2*h/2);
k4=fax(x+h,y+h*k3);
return y+(k1+2*k2+2*k3+k4)*h/6;
}
double progressEuler(double h,double x,double y)
{
double y1;
y1=y+h*fax(x,y);
return y+(fax(x,y)+fax(x+h,y1))*h/2;
}
double accurate(double x)
{
double y=pow(1+2*x,1.0/2);
return y;
}
int main()
{
double h,x,r,y;
double n1=0,n2=0,n3=0;
cout<<"please input x0,y,h,r:"<<endl;
cin>>x>>y>>h>>r;
n1=n2=n3=y;
cout<<"x\t y(4阶龙格)\t y(改进)\t y(Euler)\t y(精确)"<<endl;
for(x;x<r-h;x=x+h)
{
n1=RungKutta(h,x,n1);
n2=progressEuler(h,x,n2);
n3=Euler(h,x,n3);
cout<<x+h<<"\t"<<n1<<" \t"<<n2<<" \t"<<n3<<" \t"<<accurate(x)<<endl;
}
return 0;
}
