// This program uses the fourth-order Runge-Kutta method to solve a first order differential equation
// of the form dy/dx = f(x,y) over an interval of x. Here this is symbolized by (d/dx)y = der. 

// Import all classes needed for program
import java.util.Vector;
import java.lang.Double;
import java.lang.Math;
   
// Define class RungeKutta  
public class RungeKutta{
    
    // The method RK1 will step through the independent variable x to solve a first 
    // order differential equation over an interval xStart <= x <= End. The initial 
    // condition for Y is y = yStart. The method requires initial x and y values, 
    // stepsize h, and final x-value End. The variables xOutv and yOutv
    // are the output "arrays". These output "arrays" are called Vectors (java.util.Vector class). 
    // Vectors are act as arrays that can dynamically add elements through method *.addElement().  
    public Vector RK1 (double xStart, double yStart, double h, double End, Vector xOutv, Vector yOutv) {

        // Assign values. Recall that all variables must have an initial assignment
        double x = xStart;
        double y = yStart;
        double end = End;
        double yout = 0.0;
        
        // Place initial values in Vectors.
        // Add xStart and yStart values to Vectors xOutv and yOutv. Since the Vector
        // method addElement(Object) 
        // (i.e., the argument of addElement must be an object) is used, 
        // the double values x+h and yout must be made an instance 
        // (i.e., c and d) of the class Double. 
        // The Double class wraps a value of the primitive type double in an object. 
        // An object of type Double contains a single field whose type is double. 
        Double a = new Double(xStart);
        Double b = new Double(yStart);
        xOutv.addElement(a);
        yOutv.addElement(b);
       
        // Check for endpoint. Begin loop over values of indepedent variable x
        do {
     
            // Call Runge-Kutta integration method
            yout = fork(x, y, h, yout); 
            
            // Add x and y values to Vectors. Since the Vector method addElement(Object) (i.e.,
            // the argument of addElement must be an object) is used, the double values 
            // x+h and yout must be made an instance (i.e., c and d) of the class Double. 
            // The Double class wraps a value of the primitive type double in an object. 
            // An object of type Double contains a single field whose type is double. 
            Double c = new Double(x+h);
            Double d = new Double(yout);
            xOutv.addElement(c);
            yOutv.addElement(d);
                   
            // Prepare for next step
            x = x+h;
            y = yout;

        }
        while(x <= end);
            
        // Print to screen Vector xOutv
        System.out.println(xOutv);
        System.out.println(yOutv);
        
        // At or passed endpoint so exit loop and return
        System.out.println ("Integration is complete");
        
        // This statement is required even though we are not returning the value of xOutv to 
        // anywhere. Not sure why.
        return xOutv;
    }       
    
    // Takes one fourth-order Runge-Kutta step of stepsize h
    private double fork (double x, double y, double h, double yout){
 
        // Initialize and define local variables. Note: These values for k will be 
        // reassigned before they are used.  
        double k1 = 1.0000;
        double k2 = 1.0000;
        double k3 = 1.0000;
        double k4 = 1.0000;
        double der = y;
        
        // Perform calculations. The method deriv(x, y, der) is called and 
        // returns the value for der in its place in the equation.
        k1 = h*deriv(x,y,der);
        k2 = h*deriv(x+0.5*h,y+0.5*k1,der);
        k3 = h*deriv(x+0.5*h,y+0.5*k2,der);
        k4 = h*deriv(x+h,y+k3,der);
        
        // Set output
        double youtTerm =(k1+2*k2+2*k3+k4);
   
        // Calculate yout. Important: yout, y, and k1,...,k4 are of type double. If (1/6)*() 
        // was used we would get 0.0 in return, since 1/6 is int/int which equals int.
        yout = y+((k1+2*k2+2*k3+k4)/6);
        
        // Return value of yout to the "RK1" method
        return yout;
    } 
    // Method deriv: contains derivative of function to be integrated. The word "double"
    // comes before "deriv" since this is the data type being returned to the fork method.
    private double deriv (double x, double y, double der){
        
        // Derivative of differential equation (right-hand side) goes here
        // For this example, dy/dx = x + y + 2
        der = x + y + 2; 
        
        // Return the value of the derivative to the "fork" method
        return der;
    }
    
    // Main method
    public static void main (String args[]){
        
        // Vector() is a class so we must create a new instance of this class
        Vector xOutv = new Vector();
        Vector yOutv = new Vector();
        
        // Create a new instance of the class RungeKutta()
        RungeKutta integrate = new RungeKutta();
        
        System.out.println("Calling RungeKutta" );
        
        // The instance "integrate" calls the RK1 method
        // Initial x is 0.0, initial y is 0.0, intial stepsize is 1.0, last x is 5.0
        integrate.RK1(0.0, 0.0, 1.0, 5.0, xOutv, yOutv);        
    }  
}