clear
clf
global m q
q=10;
m=2;
ii=1;
%intial position and velocity
x0=2;
y0=0;
z0=1;
vx0=5;
vy0=0;
vz0=0;
y0=[x0;y0;z0;vx0;vy0;vz0];

t=0;
tmax=50;
dt=.1;
yn=y0;

%arrays
tt=zeros(tmax/dt,1);
xt=zeros(tmax/dt,1);
yt=zeros(tmax/dt,1);
zt=zeros(tmax/dt,1);

while (t<tmax)
    x=yn(1);
    y=yn(2);
    z=yn(3);
    tt(ii)=t;
    xt(ii)=x;
    yt(ii)=y;
    zt(ii)=z;
   % E = Efield(x,y,z);
    %E=[1,0,0];
    %B = Bfield(x,y,z);
    %B=[0,1,0];
    f1 = emode(t     ,yn);
    f2 = emode(t+dt/2,yn+(dt/2).*f1);
    f3 = emode(t+dt/2,yn+(dt/2).*f2);
    f4 = emode(t+dt  ,yn+ dt   .*f3);
    yn = yn+dt/6.*(f1+2*f2+2*f3+f4);
    t = t + dt;
   
    plot3(xt(1:ii),yt(1:ii),zt(1:ii))
    saveas(gcf, ['electromag',sprintf('%03d',ii),'.png'])
    ii=ii+1;
    drawnow
    
end