大神们看下这个程序怎么改
#include<iostream>#include<iomanip>
using namespace std;
void RK4(double(*f)(double t,double x,double y,double z),double(*g)(double t,double x,double y,double z),double(*k)(double t,double x,double y,double z),double initial,double resu,double h)
{
double f1,f2,f3,f4,g1,g2,g3,g4,k1,k2,k3,k4,t0,x0,y0,z0,x1,y1,z1;
t0=initial;x0=initial;y0=initial;z0=initial;
f1=f(t0,x0,y0,z0); g1=g(t0,x0,y0,z0); k1=k(t0,x0,y0,z0);
f2=f((t0+h/2),(x0+h*f1/2),(y0+h*g1/2),(z0+h*k1/2));g2=g((t0+h/2),(x0+h*f1/2),(y0+h*g1/2),(z0+h*k1/2)); k2=k((t0+h/2),(x0+h*f1/2),(z0+h*k1/2));
f3=f((t0+h/2),(x0+h*f2/2),(y0+h*g2/2),(z0+h*k2/2));g3=g((t0+h/2),(x0+h*f2/2),(y0+h*g2/2),(z0+h*k2/2)); k3=k((t0+h/2),(x0+h*f2/2),(z0+h*k2/2));
f4=f((t0+h),(x0+h*f3),(y0+h*g3),(z0+h*k3)); g4=((t0+h),(x0+h*f3),(y0+h*g3),(z0+h*k3)); k4=k((t0+h),(x0+h*f3),(z0+h*k3));
x1=x0+h*(f1+2*f2+2*f3+f4)/6;
y1=y0+h*(g1+2*g2+2*g3+g4)/6;
z1=z0+h*(k1+2*k2+2*k3+k4)/6;
resu=t0+h;resu=x1;resu=y1;resu=z1;
}
int main()
{
double f(double t,double x,double y,double z);
double g(double t,double x,double y,double z);
double k(double t,double x,double y,double z);
double initial,resu;
double t,step;
double a,b,H;
int i;
cout<<"输入所求微分方程组的初值t0,x0,y0,z0:";
cin>>initial>>initial>>initial>>initial;
cout<<"输入所求微分方程组的微分区间:";
cin>>a>>b;
cout<<"输入所求微分方程组所分解子区间的个数step:";
cin>>step;
cout<<setiosflags(ios::right)<<setiosflags<<(ios::fixed)<<setiosflags(10);
H=(b-a)/step;
cout<<initial<<setw(18)<<initial<<setw(18)<<initial<<setw(18)<<initial<<endl;
for(i=0;i<step;i++)
{
RK4(f,g,k,initial,resu,H);
cout<<resu<<setw(20)<<resu<<setw(20)<<resu<<setw(20)<<resu<<endl;
initial=resu;initial=resu;initial=resu;initial=resu;
}
return(0);
}
double f(double t,double x,double y,double z)
{
double dx;
dx=-10*x+10*y;
return(dx);
}
double g(double t,double x,double y,double z)
{
double dy;
dy=30*x-y-x*z;
return(dy);
}
double k(double t,double x,double y,double z)
{
double dz;
dz=-(8/3)*z+x*y;
return(dz);
}
页:
[1]