Leapfrog/Midpoint ODE Method - Incorrect

I simply don't understand why Leapfrog isn't working correctly, Euler does. The algorithm for Leapfrog should be:

y[n+1] = y[n-1] + 2*h*f(x[n],y[n])

x[n+1] = x[n] + h;

Following page 33 of: Here

Please help!

Leapfrog.m

function [ matrix ] = LeapfrogMethod( fun, initX, initY,...
    steplength, maximum, opt_print )
%LEAPFROGMETHOD Solves the differential equation supplied

% Define variables
f = fun;
x0 = initX;
y0 = initY;
h = steplength;
n = maximum;

% Set up matrix
range = n/h;
matrix = zeros(range,2);

% Calculate first value
x1 = x0 + h;
euler = EulerMethod(f,x0,y0,h,h,0);
y1 = euler(1,2);
matrix(1,1) = x1;
matrix(1,2) = y1;
if (~exist('opt_print', 'var'))
    fprintf('x = %3g, y = %3g\n',matrix(1,1),matrix(1,2))
end

% Method
for J = 2:range
    matrix(J,2) = y0 + 2*h*f(x1,y1);
    matrix(J,1) = x1 + h;
    if (~exist('opt_print', 'var'))
        fprintf('x = %3g, y = %3g\n',matrix(J,1),matrix(J,2))
    end
    x1 = matrix(J,1);
    y0 = y1;
    y1 = matrix(J,2);
end
%
end

EulerMethod.m

function [ matrix ] = EulerMethod( fun, initX, initY,...
    steplength, maximum, opt_print )
%EULERMETHOD Solves the differential equation supplied

% Define variables
f = fun;
x0 = initX;
y0 = initY;
h = steplength;
n = maximum;

% Set up matrix
range = n/h;
matrix = zeros(range,2);

% Method
for J = 1:range
    matrix(J,2) = y0 + h*f(x0,y0);
    matrix(J,1) = x0 + h;
    if (~exist('opt_print', 'var'))
        fprintf('x = %3g, y = %3g\n',matrix(J,1),matrix(J,2))
    end
    x0 = matrix(J,1);
    y0 = matrix(J,2);
end
%
end

[Patddy]

Apologies, it was unstable and I was later asked why it was unstable. Nothing wrong with the code.

Answer 1

[Varun Bhaskar]

Hi,

Can you elaborate on the question?

Thanks

[Patddy]

Apologies, it was unstable and I was later asked why it was unstable. Nothing wrong with the code.