classdef Rungekutta4SecondOrder < time.Timestepper
    properties
        F
        dt
        t, n
        v, v_t
    end


    methods
        % Create a time stepper for
        %   v_tt = F(t,v,v_t),  v(t0) = v0, v_t(t0) = v0t
        % with step size dt, by rewriting on first order form
        function obj = Rungekutta4SecondOrder(F, dt, t0, v0, v0t)
            obj.F = F;
            obj.dt = dt;
            obj.t = t0;
            obj.n = 0;

            obj.v = v0;
            obj.v_t = v0t;
        end

        function [v,t] = getV(obj)
            v = obj.v
            t = obj.t;
        end

        function [vt,t] = getVt(obj)
            vt = obj.v_t;
            t = obj.t;
        end

        function obj = step(obj)
            t = obj.t;
            v = obj.v;
            v_t = obj.v_t;
            dt = obj.dt;

            k1 = obj.F(t,          v,                            v_t);
            k2 = obj.F(t + 1/2*dt, v + 1/2*dt*v_t,               v_t + 1/2*dt*k1);
            k3 = obj.F(t + 1/2*dt, v + 1/2*dt*v_t + 1/4*dt^2*k1, v_t + 1/2*dt*k2);
            k4 = obj.F(t + dt,     v + dt*v_t + 1/2*dt^2*k2,     v_t + dt*k3);

            obj.v = v + dt*v_t + dt^2*(1/6)*(k1 + k2 + k3);
            obj.v_t = v_t + dt*(1/6)*(k1 + 2*k2 + 2*k3 + k4);
            obj.t = obj.t + obj.dt;
            obj.n = obj.n + 1;
        end
    end

    methods (Static)
        function k = getTimeStep(lambda)
            k = rk4.get_rk4_time_step(lambda);
        end
    end
end