配列を用いたC言語プログラミングについて
以下のルンゲクッタ法を用いたプログラムに配列などを使いさらに短くしたいのですが
どのような方法が有りますか？
#include <stdio.h>
#include <math.h>
double f1(double t1,double w,double x,double y,double z);
double f2(double t1,double w,double x,double y,double z);
double f3(double t1,double w,double x,double y,double z);
double f4(double t1,double w,double x,double y,double z); //箱Aの関数
double g1(double t1,double a,double b,double c,double d);
double g2(double t1,double a,double b,double c,double d);
double g3(double t1,double a,double b,double c,double d);
double g4(double t1,double a,double b,double c,double d); //箱Bの関数
int main(void)
{
double t1,w,x,y,z,a,b,c,d,dt,t1max,t2max,lam,gam,lat,dw,dx,dy,dz,da,db,dc,dd ;
double k1[4],k2[4],k3[4],k4[4],l1[4],l2[4],l3[4],l4[4] ; ///宣言
t1 = 0.0;
dt = 0.3;
t1max = 40.0; //時間初期値
w = 200.0;
x = 40.0;
y = 30.0;
z = 30.0; ///箱A初期値（w：感受性人口、x:潜伏人口、y:感染人口、z:隔離人口）
a = 20.0;
b = 8.0;
c = 12.0;
d = 10.0; ///箱B初期値（a:感受性人口、b:潜伏人口、c:感染人口,d:隔離人口）
for(t1=0.0;t1<=t1max;t1+=dt) {
k1[0]=dt*f1(t1,w,x,y,z);
k1[1]=dt*f2(t1,w,x,y,z);
k1[2]=dt*f3(t1,w,x,y,z);
k1[3]=dt*f4(t1,w,x,y,z);
k2[0]=dt*f1(t1+dt/2.0,w+k1[0]/2.0,x+k1[1]/2.0,y+k1[2]/2.0,z+k1[3]/2.0);
k2[1]=dt*f2(t1+dt/2.0,w+k1[0]/2.0,x+k1[1]/2.0,y+k1[2]/2.0,z+k1[3]/2.0);
k2[2]=dt*f3(t1+dt/2.0,w+k1[0]/2.0,x+k1[1]/2.0,y+k1[2]/2.0,z+k1[3]/2.0);
k2[3]=dt*f4(t1+dt/2.0,w+k1[0]/2.0,x+k1[1]/2.0,y+k1[2]/2.0,z+k1[3]/2.0);
k3[0]=dt*f1(t1+dt/2.0,w+k2[0]/2.0,x+k2[1]/2.0,y+k2[2]/2.0,z+k2[3]/2.0);
k3[1]=dt*f2(t1+dt/2.0,w+k2[0]/2.0,x+k2[1]/2.0,y+k2[2]/2.0,z+k2[3]/2.0);
k3[2]=dt*f3(t1+dt/2.0,w+k2[0]/2.0,x+k2[1]/2.0,y+k2[2]/2.0,z+k2[3]/2.0);
k3[3]=dt*f4(t1+dt/2.0,w+k2[0]/2.0,x+k2[1]/2.0,y+k2[2]/2.0,z+k2[3]/2.0);
k4[0]=dt*f1(t1+dt,w+k3[0],x+k3[1],y+k3[2],z+k3[3]);
k4[1]=dt*f2(t1+dt,w+k3[0],x+k3[1],y+k3[2],z+k3[3]);
k4[2]=dt*f3(t1+dt,w+k3[0],x+k3[1],y+k3[2],z+k3[3]);
k4[3]=dt*f4(t1+dt,w+k3[0],x+k3[1],y+k3[2],z+k3[3]); ///箱Aルンゲクッタ
l1[0]=dt*g1(t1,a,b,c,d);
l1[1]=dt*g2(t1,a,b,c,d);
l1[2]=dt*g3(t1,a,b,c,d);
l1[3]=dt*g4(t1,a,b,c,d);
l2[0]=dt*g1(t1+dt/2.0,a+l1[0]/2.0,b+l1[1]/2.0,c+l1[2]/2.0,d+l1[3]/2.0);
l2[1]=dt*g2(t1+dt/2.0,a+l1[0]/2.0,b+l1[1]/2.0,c+l1[2]/2.0,d+l1[3]/2.0);
l2[2]=dt*g3(t1+dt/2.0,a+l1[0]/2.0,b+l1[1]/2.0,c+l1[2]/2.0,d+l1[3]/2.0);
l2[3]=dt*g4(t1+dt/2.0,a+l1[0]/2.0,b+l1[1]/2.0,c+l1[2]/2.0,d+l1[3]/2.0);
l3[0]=dt*g1(t1+dt/2.0,a+l2[0]/2.0,b+l2[1]/2.0,c+l2[2]/2.0,d+l2[3]/2.0);
l3[1]=dt*g2(t1+dt/2.0,a+l2[0]/2.0,b+l2[1]/2.0,c+l2[2]/2.0,d+l2[3]/2.0);
l3[2]=dt*g3(t1+dt/2.0,a+l2[0]/2.0,b+l2[1]/2.0,c+l2[2]/2.0,d+l2[3]/2.0);
l3[3]=dt*g4(t1+dt/2.0,a+l2[0]/2.0,b+l2[1]/2.0,c+l2[2]/2.0,d+l2[3]/2.0);
l4[0]=dt*g1(t1+dt,a+l3[0],b+l3[1],c+l3[2],d+l3[3]);
l4[1]=dt*g2(t1+dt,a+l3[0],b+l3[1],c+l3[2],d+l3[3]);
l4[2]=dt*g3(t1+dt,a+l3[0],b+l3[1],c+l3[2],d+l3[3]);
l4[3]=dt*g4(t1+dt,a+l3[0],b+l3[1],c+l3[2],d+l3[3]); ///箱Bルンゲクッタ
w=w+((k1[0]+2.0*k2[0]+2.0*k3[0]+k4[0])/6.0);
x=x+((k1[1]+2.0*k2[1]+2.0*k3[1]+k4[1])/6.0);
y=y+((k1[2]+2.0*k2[2]+2.0*k3[2]+k4[2])/6.0);
z=z+((k1[3]+2.0*k2[3]+2.0*k3[3]+k4[3])/6.0);
a=a+((l1[0]+2.0*l2[0]+2.0*l3[0]+l4[0])/6.0);
b=b+((l1[1]+2.0*l2[1]+2.0*l3[1]+l4[1])/6.0);
c=c+((l1[2]+2.0*l2[2]+2.0*l3[2]+l4[2])/6.0);
d=d+((l1[3]+2.0*l2[3]+2.0*l3[3]+l4[3])/6.0);
}
return 0;
}
