如何在下面的程序中加一个庞加莱截面呢？？谢谢
#include <stdio.h>#include <math.h>
#define N 10000
#define h 0.01
double f(double x,double y)
{
return(10*(y-x));
}
double p(double x,double y,double z ,double c)
{
return(c*x-x*z-y);
}
double q(double x,double y,double z)
{
return(x*y-8/3*z);
}
doublerunge_kutta(double x0,double y0,double z0)
{
double x,y,z,c;
double d1,d2,d3,d4,k1,k2,k3,k4,f1,f2,f3,f4;
int i,j;
x=x0;
y=y0;
z=z0;
FILE *fp;
fp=fopen("lorenz分岔4.txt","w");
for(j=1;j<=500;j++)
{
c=1+j;
for(i=1;i<=N;i++)
{
d1=f(x,y);
k1=p(x,y,z,c);
f1=q(x,y,z);
d2=f(x+h*d1/2,y+h*k1/2);
k2=p(x+h*d1/2,y+h*k1/2,z+h*f1/2,c);
f2=q(x+h*d1/2,y+h*k1/2,z+h*f1/2);
d3=f(x+h*d2/2,y+h*k2/2);
k3=p(x+h*d2/2,y+h*k2/2,z+h*f2/2,c);
f3=q(x+h*d2/2,y+h*k2/2,z+h*f2/2);
d4=f(x+h*d3,y+h*k3);
k4=p(x+h*d3,y+h*k3,z+h*f3,c);
f4=q(x+h*d3,y+h*k3,z+h*f3);
x=x+h*(d1+2*d2+2*d3+d4)/6;
y=y+h*(k1+2*k2+2*k3+k4)/6;
z=z+h*(f1+2*f2+2*f3+f4)/6;
if(i>9500)
{
fprintf(fp," %8f, %8f \n",c,x);
printf("%8f ,%8f \n",c,x);
}
}
}
}
main()
{
double x0=-10,y0=0,z0=20;
runge_kutta(x0,y0,z0);
}
这么大的一个论坛居然没有会的吗？
页:
[1]