package odend; import Jama.*; /** * solver using classical 4th order Runge-Kutta method */ public class ODESolverRK4 extends ODESolver { /** sets ode and initial values */ public ODESolverRK4(ODE ode) { super(ode); } /** * computes solution at next time step<br> * using the classical 4th order Runke-Kutta scheme */ public void nextStep(double h) { Matrix k1 = ode.f(x, t); double t1 = t + h / 2.0; Matrix x1 = x.plus(k1.times(h / 2.0)); // x1 = x + h*k1/2.0 Matrix k2 = ode.f(x1, t1); double t2 = t1; Matrix x2 = x.plus(k2.times(h / 2.0)); // x2 = x + h*k2/2.0 Matrix k3 = ode.f(x2, t2); double t3 = t + h; Matrix x3 = x.plus(k3.times(h)); // x3 = x + k3*h Matrix k4 = ode.f(x3, t3); // x += (k1 + 2.0*k2 + 2.0*k3 + k4)*h/6.0 x.plusEquals(k1.times(h / 6.0)); x.plusEquals(k2.times(h / 3.0)); x.plusEquals(k3.times(h / 3.0)); x.plusEquals(k4.times(h / 6.0)); t += h; } }