4次ルンゲクッタ法によるロジスティック方程式のシミュレーションの実装例です

RK4_KogisticEquation.c

/* 4次ルンゲクッタ法によるロジスティック方程式の解の計算 */
#include <stdio.h>

#define K 1000.0		// 個体数のしきい値
#define r 0.9			// 1個体の増加率
#define N_LOOP 10000	// ループ回数

// 積分する関数g(t,N)
double g(double t,double N_t){
	return r*(K-N_t)*N_t/K;
}

int main(void){
	int k;				// ループカウンタ
	double dt = 1e-3;	// 刻み幅Δt
	double N = 1.0;		// 個体数t（初期値を代入しておく)
	double k1,k2,k3,k4;	// 計算に使う諸変数
	double a,m,b;		// t_k,t_k+0.5,t_k+1を格納する
	
	for( k = 0 ; k < N_LOOP ; k++ ){
		// t_k,t_k+0.5,t_k+1を計算
		a = k*dt;
		m = (k+0.5)*dt;
		b = (k+1.0)*dt;
		
		// k1,...,k4の値を計算
		k1 = g(a,N);
		k2 = g(m,N+k1*dt*0.5);
		k3 = g(m,N+k2*dt*0.5);
		k4 = g(b,N+k3*dt);
		
		// 積分して代入
		N = N + (dt/6.0)*(k1+2.0*k2+2.0*k3+k4);
		
		// 値を表示
		printf("%lf,%e\n",b,N);
	}
	
	return 0;
}