四阶龙格库塔法的JAVA编程

Java语言的二三四阶龙格库塔程序如下:

public class D {

static double[] x=new double[6]

static double[] y=new double[6]

public static double function(double a,double b)

{

return b-2*a/b

}

/*n表示几等分，n+1表示他输出的个数*/

public static void RungeKutta(double y0,double a,double b,int n,int style)

{

double h=(b-a)/n,k1,k2,k3,k4

int i

x[0]=a

y[0]=y0

switch(style)

{

case 2:

for(i=0i<ni++)

{

x[i+1]=x[i]+h

k1=function(x[i],y[i])

k2=function(x[i]+h/2,y[i]+h*k1/2)

y[i+1]=y[i]+h*k2

}

break

case 3:

for(i=0i<ni++)

{

x[i+1]=x[i]+h

k1=function(x[i],y[i])

k2=function(x[i]+h/2,y[i]+h*k1/2)

k3=function(x[i]+h,y[i]-h*k1+2*h*k2)

y[i+1]=y[i]+h*(k1+4*k2+k3)/6

}

break

case 4:

for(i=0i<ni++)

{

x[i+1]=x[i]+h

k1=function(x[i],y[i])

k2=function(x[i]+h/2,y[i]+h*k1/2)

k3=function(x[i]+h/2,y[i]+h*k2/2)

k4=function(x[i]+h,y[i]+h*k3)

y[i+1]=y[i]+h*(k1+2*k2+2*k3+k4)/6

}

break

}

}

public static void main(String[] args) {

//例子求y'=y-2*x/y(0<x<1)y0=1

System.out.println("用二阶龙格-库塔方法")

RungeKutta(1,0,1,5,2)

for(int i=0i<6i++)

System.out.printf("x[%d]=%f,y[%d]=%f\n",i,x[i],i,y[i])

System.out.println("用三阶龙格-库塔方法")

RungeKutta(1,0,1,5,3)

for(int i=0i<6i++)

System.out.printf("x[%d]=%f,y[%d]=%f\n",i,x[i],i,y[i])

System.out.println("用四阶龙格-库塔方法")

RungeKutta(1,0,1,5,4)

for(int i=0i<6i++){

System.out.printf("x[%d]=%f,y[%d]=%f\n",i,x[i],i,y[i])

}

}

}

运行结果:

用二阶龙格-库塔方法

x[0]=0.000000,y[0]=1.000000

x[1]=0.200000,y[1]=1.183636

x[2]=0.400000,y[2]=1.342656

x[3]=0.600000,y[3]=1.485014

x[4]=0.800000,y[4]=1.615225

x[5]=1.000000,y[5]=1.736182

用三阶龙格-库塔方法

x[0]=0.000000,y[0]=1.000000

x[1]=0.200000,y[1]=1.183244

x[2]=0.400000,y[2]=1.341729

x[3]=0.600000,y[3]=1.483408

x[4]=0.800000,y[4]=1.612727

x[5]=1.000000,y[5]=1.732472

用四阶龙格-库塔方法

x[0]=0.000000,y[0]=1.000000

x[1]=0.200000,y[1]=1.183229

x[2]=0.400000,y[2]=1.341667

x[3]=0.600000,y[3]=1.483281

x[4]=0.800000,y[4]=1.612514

x[5]=1.000000,y[5]=1.732142

四阶龙格库塔法精度为4，属于单步递推法，单步递推法的基本思想是从（x(i),y(i)）点出发，以某一斜率沿直线达到（x(i+1),y(i+1)）点，从上述定义可以看出，龙格库塔实质上是求一阶微分方程，对其进行排列后就可以进入Matlab进行编程。

但是如果将一阶导看作变量，则二阶导也不过是这个变量的一阶导而已，对于下述二阶方程：基本思想如下令位移为q的一阶导，即位移的一阶导（速度）为q的二阶导。

是用一个简单的迭代式来画分叉混沌图的，从这个例子里，不难学到如何如何画分叉图，其原理不难，两个循环来在2维平面作图，一个if筛选迭代的结果。

还是类似上面那个例子，这里给出不一样的代码实现方法。没看懂的可以再看一遍，里面主要的还是两个for循环，下面这个代码是用第三个for循环来实现对迭代结果的筛选的，这里不必纠结细节，功能实现用if或者for都可以。

依次类推，如果在区间[xi,xi+1]内多预估几个点上的斜率值K1、K2、……Km，并用他们的加权平均数作为平均斜率K*的近似值，显然能构造出具有很高精度的高阶计算公式。经数学推导、求解，可以得出四阶龙格－库塔公式，也就是在工程中应用广泛的经典龙格－库塔算法这样就完成了matlab的编程。

---