% Takes one time step of size dt using the rungekutta method
% starting from @arg v and where the function F(v,t) gives the
% time derivatives. coeffs is a struct holding the RK coefficients
% for the specific method.
% Also returns the stage approximations (V) and stage rates (K).
function [v, V, K] = rungekutta(v, t , dt, F, coeffs)
    % Compute the intermediate stages k
    K = zeros(length(v), coeffs.s);
    V = zeros(length(v), coeffs.s);
    for i = 1:coeffs.s
        u = v;
        for j = 1:i-1
            u = u + dt*coeffs.a(i,j)*K(:,j);
        end
        V(:,i) = u;
        K(:,i) = F(u,t+coeffs.c(i)*dt);
    end
    % Compute the updated solution as a linear combination
    % of the intermediate stages.
    for i = 1:coeffs.s
        v = v + dt*coeffs.b(i)*k(:,i);
    end
end