我想计算洛伦茨系统使用龙格库塔方法，但我找不到我的代码有错误的地方。当我运行它时，系统进入静态点，我应该获得一个蝴蝶(洛伦兹吸引子)。我认为这是Runge Kutta方法部分的'for'循环中的一些东西。这里是我的代码龙格库塔在C为洛伦兹方程

#include

#include

#include

double f(int,double,double []);

double s,r,b;

FILE *output;

int main()

{

output=fopen("lorenzdata.dat","w");

int j,N=3;

double x[2],dt,y[2],K1[2],K2[2],K3[2],K4[2],ti,t,i;

printf("\n Introduce parameters s, r and b sepparated by a space:");

scanf("%lf %lf %lf",&s,&r,&b);

printf("\n Introduce initial conditions t0,x0,y0 and z0:");

scanf("%lf %lf %lf %lf",&ti,&x[0],&x[1],&x[2]);

printf("\n Introduce time step:");

scanf("%lf",&dt);

printf("\n Specify time to compute in the equations:");

scanf("%lf",&t);

/* Loop para Runge Kutta 4 */

do

{

/* printf("%9.4f %9.4f %9.4f %9.4f \n",ti,x[0],x[1],x[2]); */

i++;

for(j=0;j

{

K1[j] = f(j,ti,x);

}

for(j=0;j

{

y[j] = x[j]+(K1[j]/2)*dt;

}

for(j=0;j

{

K2[j] = f(j,ti+dt/2,y);

}

for(j=0;j

{

y[j] = x[j]+(K2[j]/2)*dt;

}

for(j=0;j

{

K3[j] = f(j,ti+dt/2,y);

}

for(j=0;j

{

y[j] = x[j]+(K3[j])*dt;

}

for(j=0;j

{

K4[j] = f(j,ti+dt,y);

}

for(j=0;j

{

x[j] += dt*((K1[j]/6)+(K2[j]/3)+(K3[j]/3)+(K4[j]/6));

}

ti +=dt;

fprintf(output,"%9.4f %9.4f %9.4f \n",x[0],x[1],x[2]);

}while(i*dt <= t);

fclose(output);

return 0;

}

/* Definimos la funcion */

double f(int m,double h,double x[])

{

if(m==0)

{

return s*(x[1]-x[0]);

}

else if(m==1)

{

return x[0]*(r-x[2])-x[1];

}

else if(m==2)

{

return x[0]*x[1]-b*x[2];

}

}

在此先感谢

编辑

我在Linux上重新写的代码(以更简化的方式)，并且它运行正常，看来我在Windows上编辑了一些东西(也许是编码)，当我编译代码时，它抛出了很多infinites值。不知道为什么，我花了很多时间才注意到这一点。 感谢您的帮助

+0

双h作为函数f的输入点是什么？它似乎没有被使用。有必要吗？ –

2012-02-20 00:28:40

+0

mm nop，我只是为了做一个更一般化的代码，以防万一diff方程式明确地取决于t –

2012-02-20 03:05:04