#include #include #include "nrutil.h" FILE *fpdata1; void rk4(double y[], double dydx[], int n, double x, double h, double yout[], void (*derivs)(double x, double yt[], double dyt[])) // nはベクトル(配列)次数 { int i; double xh,hh,h6,*dym,*dyt,*yt; dym=vector(n); // dymの先頭ポインタを代入 dyt=vector(n); // dytの先頭ポインタを代入 yt=vector(n); // ytの先頭ポインタを代入 hh=h*0.5; h6=h/6.0; xh=x+hh; // 4次のルンゲクッタ (*derivs)(x, y, dydx); // 2004/08/05 追加, 第1回目計算 for (i=0;i