function [T,Y]= leapfrog2(f,t,y)
% LEAPFROG2 Second-order Leapfrog solver
% Input values are column vector of times t 
%                  and column vector initial conditions y(6)= (x,y,z,vx,vy,vz)
% Output values are column vectors time T
% and Y(size(t,1),6)=(dx/dt,dy/dt,dz/dt,dvx/dt,dvy/dt,dvz/dt)
% First step is initialized using Euler's method

nn=size(t,1);
%Set initial conditions
T(1)=t(1);
Y(:,1)=y;

%Perform a single Euler step to initialize
i=2;
T(i)=t(i);
dt=T(i)-T(i-1);
dydt=lorentz(T(i),Y(:,i-1));
Y(:,i)=Y(:,i-1) + dt*dydt(:);

%Perfrom Leapfrog steps
for i= 3:nn
  T(i)=t(i);
  dt=T(i)-T(i-2);
  dydt=lorentz(T(i),Y(:,i-1));
  Y(:,i)=Y(:,i-2) + dt*dydt(:);
end

%Transpose matrices for output
Y=Y';
T=T';