程序求助
用C语言编写的程序,结果和matlab差很多画相图的程序#include <stdio.h>
#include <math.h> /* definitions for math functions */
/*--- function prototypes: ------------------------------------------*/
double max(double,double);
double min(double,double);
/*-------------------------------------------------------------------*/
#define NDIM 7
#define NPARAM 1
void
function(double *v, double time, double *param, double *deriv){
double x1, x2, x3, x4, x5, x6, x7, R;
x1 = v; x2 = v; x3 = v; x4 =v; x5 =v; x6 =v; x7=v;
R= param;
deriv = -16*x1+R*(-x2*x3+3*x3*x4+8*x6*x7)-4;
deriv = -10*x2+R*(2.8*x1*x3-1.2*x3*x5+7.2*x4*x6);
deriv = -2*x3-R*(6*x1*x2+10*x1*x4+2*x2*x5+14*x4*x7);
deriv = -26*x4+R*(-(14/13)*x1*x3-(4/13)*x2*x6+(38/13)*x3*x7);
deriv = -8*x5+R*(4*x1*x6+2*x2*x3);
deriv = -8*x6-R*(4*x1*x5+12*x1*x7+8*x2*x4);
deriv = -40*x7-R*(0.8*x1*x6+1.2*x3*x4);
}
void
runge_kutta (double *x,double *t, double *tau, double *param, double *xout){
double F1, F2, F3, F4, xtemp;
double t_half, half_tau, t_full;
int i;
half_tau = 0.5 * *tau;
function (x, *t, param, F1);
t_half = *t + half_tau;
for (i = 0;i < NDIM; i++) xtemp = x + half_tau*F1;
function (xtemp, t_half, param, F2);
for (i = 0;i < NDIM; i++) xtemp = x + half_tau*F2;
function (xtemp, t_half, param, F3);
t_full = *t + *tau;
for (i = 0;i < NDIM; i++) xtemp = x + *tau * F3;
function (xtemp, t_full, param, F4);
for (i = 0;i < NDIM; i++)
xout = x + *tau/6 *(F1 + F4 + 2*(F2 + F3));
}
int
main(int argc, char **argv) {
double x={-0.25,0.0,-0.1,0.0,0.0,0.0,0.0}, param={44.5},xout;
double time=0.0,tau=0.0001;
int i,j;
FILE *fout;
fout=fopen("phase.dat","w");
for (i=0;i<500000;i++){
runge_kutta(x,&time,&tau,param,xout);
for (j=0;j<NDIM;j++){
x=xout;
}
time=time+tau;
runge_kutta(x,&time,&tau,param,xout);
i=i+1;
if (i==10000){
time=0.0;
tau=0.001;
}
if (i>150000){
fprintf(fout,"%f %f %f %f\n",time,x,x,x);
//printf("%f %f %f %f\n",time,x,x,x);
}
}
fclose(fout);
}
对了，变步长的也用了，结果也不对,不知道到底是怎么回事.请高人指点
E-mail:F_inzaghi@126.com
页:
[1]