使用变步长的龙格库塔算法求解常微分方程代码,Java代码库,德仔网

代码语言

代码分类

【Java】使用变步长的龙格库塔算法求解常微分方程

作者:yamorn / 发布于2014/11/13/ 557

本程序使用变步长的龙格库塔算法求解常微分方程f(x,y)=y-2x/y，初值y(0)为1，使用者可以根据自身需要修改要求解的常微分方程


	import java.util.Scanner;
	
	/*四阶龙格-库塔方法，求解常微分方程f(x,y)=y-2x/y*/
	public class variable_stepsize_RungeKutta {
	 public static void variable_stepsize_RungeKutta_method(double Y[],double step,double prec){
	 int n=Y.length-1;
	 double K1=0,K2=0,K3=0,K4=0;
	 double x=0;//用于存放初始横坐标位置
	 double tm_step=0;
	 double temp_Y_1=0;//用于存放前一步长计算得到的值
	 double temp_Y_2=0;//用于存放后一步长计算得到的值
	 int k=0;//用于存放横坐标抵达目标位置的步数
	 double start_v=0;//用于存放每次迭代计算的初始值
	 for(int i=1;i<=n;i++)
	 {
	 x=(i-1)*step;
	 k=1;
	 tm_step=step;
	 temp_Y_1=0;
	 temp_Y_2=0;
	 start_v=Y[i-1];//Y[i-1]作为当前初值
	 K1=start_v-(2*x)/start_v;//K1=f(x(n),y(n))
	 K2=(start_v+(tm_step/2)*K1)-2*(x+tm_step/2)/(start_v+(tm_step/2)*K1);//K2=f(x(n)+h/2,y(n)+h*K1/2)
	 K3=(start_v+(tm_step/2)*K2)-2*(x+tm_step/2)/(start_v+(tm_step/2)*K2);//K3=f(x(n)+h/2,y(n)+h*K2/2)
	 K4=(start_v+tm_step*K3)-2*(x+tm_step)/(start_v+tm_step*K3);//K4=f(x(n)+h,y(n)+step*K3)
	 temp_Y_1=start_v+tm_step*(K1+2*K2+2*K3+K4)/6;
	 tm_step=tm_step/2;//变步长，步长变为原来的1/2
	 k*=2;//步数变为原来的两倍
	 for(int j=1;j<=k;j++){
	 K1=start_v-(2*x)/start_v;//K1=f(x(n),y(n))
	 K2=(start_v+(tm_step/2)*K1)-2*(x+tm_step/2)/(start_v+(tm_step/2)*K1);//K2=f(x(n)+h/2,y(n)+h*K1/2)
	 K3=(start_v+(tm_step/2)*K2)-2*(x+tm_step/2)/(start_v+(tm_step/2)*K2);//K3=f(x(n)+h/2,y(n)+h*K2/2)
	 K4=(start_v+tm_step*K3)-2*(x+tm_step)/(start_v+tm_step*K3);//K4=f(x(n)+h,y(n)+step*K3)
	 start_v=start_v+tm_step*(K1+2*K2+2*K3+K4)/6;
	 x=x+tm_step;
	 }
	 temp_Y_2=start_v;
	
	 while(temp_Y_2-temp_Y_1>prec)//不满足精度要求时步长继续减半
	 {
	 start_v=Y[i-1];
	 x=(i-1)*step;
	 temp_Y_1=temp_Y_2;
	 tm_step=tm_step/2;
	 k*=2;
	 for(int j=1;j<=k;j++){
	 K1=start_v-(2*x)/start_v;//K1=f(x(n),y(n))
	 K2=(start_v+(tm_step/2)*K1)-2*(x+tm_step/2)/(start_v+(tm_step/2)*K1);//K2=f(x(n)+h/2,y(n)+h*K1/2)
	 K3=(start_v+(tm_step/2)*K2)-2*(x+tm_step/2)/(start_v+(tm_step/2)*K2);//K3=f(x(n)+h/2,y(n)+h*K2/2)
	 K4=(start_v+tm_step*K3)-2*(x+tm_step)/(start_v+tm_step*K3);//K4=f(x(n)+h,y(n)+step*K3)
	 start_v=start_v+tm_step*(K1+2*K2+2*K3+K4)/6;
	 x=x+tm_step;
	 }
	 temp_Y_2=start_v;
	
	 }
	 System.out.println("The step size we choose in this phase is:"+tm_step);
	 System.out.println("The value we calculate by this step size is:"+temp_Y_2);
	 Y[i]=temp_Y_2;
	// Y[i]=Y[i-1]+step*(K1+2*K2+2*K3+K4)/6;
	 }
	
	
	
	 }
	 public static void main(String[] args) {
	 Scanner scan=new Scanner(System.in);
	 System.out.println("Input initial y(0):");//Input initial value,when x(0)=0
	 double y_0=scan.nextDouble();
	 System.out.println("Input the number of values:");//Input the number of dots you want
	 int n=scan.nextInt();
	 double Y[]=new double[n+1];
	 Y[0]=y_0;
	 System.out.println("Input the step size:");
	 double step=scan.nextDouble();
	// Runge_Kutta_method(Y, step);
	 System.out.println("Input the precision:");
	 double prec=scan.nextDouble();
	 variable_stepsize_RungeKutta_method(Y, step, prec);//采用变步长的龙格库塔算法求解
	 System.out.println("The value of the dots are:");//输出算法求解结果
	 for(int i=1;i<=n;i++)
	 {
	 System.out.print(Y[i]+" ");
	 if(i==5)
	 {
	 System.out.println();
	 }
	 }
	 System.out.println();
	
	
	
	
	
	 }
	
	
	}

试试其它关键字

同语言下

可能有用的

yamorn贡献的其它代码(2)