## Wednesday, 23 July 2014

### ODEs with Java using the Apache Commons Math library

It drives me nuts when I struggle to get my head around something just to, when I eventually do, discover that it is not really as complicated as the documentation makes it out to be. Simple examples are what I need. Here is an example of how to implement the Lorenz System. Not that I really know much about the Lorenz System but it is always a simple system to use when you want to try out a new way of solving differential equations. More information about the Lorenz System can be found in Wikipedia.

The first thing to do is to create a class containing the equations. The Lorenz system consists of the following three equations:

Doing this using the Apache Commons Math library means that we have to create a class that implements the FirstOrderDifferentialEquations interface. This interface requires two methods to be implemented, which are computeDerivatives and getDimension:

 ``` 01 import org.apache.commons.math3.exception.DimensionMismatchException; 02 import org.apache.commons.math3.exception.MaxCountExceededException; 03 import org.apache.commons.math3.ode.FirstOrderDifferentialEquations; 04  05 public class Lorenz implements FirstOrderDifferentialEquations { 06  07         double sigma = 10.0; 08         double beta = 8 / 3; 09         double rho = 28.0; 10         double X, Y, Z; 11  12         @Override 13         public void computeDerivatives(double t, double[] y, double[] yDot) 14                         throws MaxCountExceededException, ``` ```                           DimensionMismatchException { 15                 X = y; 16                 Y = y; 17                 Z = y; 18                 yDot = sigma * (Y - X); 19                 yDot = X * (rho - Z) - Y; 20                 yDot = (X * Y) - (beta * Z); 21         } 22  23         @Override 24         public int getDimension() { 25                 // TODO Auto-generated method stub 26                 return 3; 27         } 28 }``` Java2html

The next step is to write the code that will call the integrator. There are various integrators. For this application we'll use the ClassicalRungeKuttaIntegrator. You will also need a step handler to capture each step of the integration:

 ``` 01 import org.apache.commons.math3.ode.FirstOrderDifferentialEquations; 02 import org.apache.commons.math3.ode.nonstiff.ClassicalRungeKuttaIntegrator; 03 import org.apache.commons.math3.ode.sampling.StepHandler; 04  05 public class LorenzMain { 06  07         public static void main(String[] args) { 08  09                 ClassicalRungeKuttaIntegrator rk4 = ``` ```                           new ClassicalRungeKuttaIntegrator(0.01); 10                 FirstOrderDifferentialEquations ode = new Lorenz(); 11                 double[] y = new double[] { 10.0, -2.0, 50 }; // initial state 12                 StepHandler stepHandler = new ``` ```                           WriteToFileStepHandler("lorenz.csv"); 13                 rk4.addStepHandler(stepHandler); 14                 rk4.integrate(ode, // equations 15                                 0.0, // start time 16                                 y, // initial conditions 17                                 100.0, // end time 18                                 y); // result 19                 for (int i = 0; i < y.length; i++) { 20                         System.out.println(y[i]); 21                 } 22         } 23 }``` Java2html

The step handler has to implement the StepHandler class. For the step handler below I created a constructor that takes a filename so that the results can be written to a csv file. Obviously you can do whatever you want with this data:

 ``` 01 import java.io.File; 02 import java.io.PrintWriter; 03 import java.util.ArrayList; 04  05 import org.apache.commons.math3.exception.MaxCountExceededException; 06 import org.apache.commons.math3.ode.sampling.StepHandler; 07 import org.apache.commons.math3.ode.sampling.StepInterpolator; 08  09  10 public class WriteToFileStepHandler implements StepHandler { 11         ArrayList steps = new ArrayList(); 12         String filename; 13          14         public WriteToFileStepHandler(String filename) { 15                 this.filename = filename; 16         } 17         @Override 18         public void handleStep(StepInterpolator interpolator, boolean isLast) 19                         throws MaxCountExceededException { 20                 double t = interpolator.getCurrentTime(); 21                 double[] y = interpolator.getInterpolatedState(); 22                         String onestep = "" + t; 23                         for (int i = 0; i < y.length; i++) { 24                                 onestep += ", " + y[i]; 25                         } 26                         steps.add(onestep); 27                 if (isLast) { 28                         try { 29                                 PrintWriter writer = ``` `30                                 new PrintWriter(new File(filename), ` ```                                   "UTF-8");                                     31                                 for (String step : steps) { 32                                         writer.println(step); 33                                 } 34                                 writer.close(); 35                         } catch (Exception e) { 36                         } 37                         ; 38                 } 39         } 40  41         @Override 42         public void init(double t0, double[] y0, double t) { 43                 // TODO Auto-generated method stub 44  45         } 46  47 }``` Java2html
Last but not the least I read the file with R and created the graph: