function dy = lorentzf(t,y)
global qm E B betam
dy=zeros(6,1);
dy(1)=y(4);
dy(2)=y(5);
dy(3)=y(6);
v=sqrt(power(y(4),2)+power(y(5),2)+power(y(6),2));
dy(4)=qm*(y(5)*B(3)-y(6)*B(2)+E(1))-betam*v*y(4);
dy(5)=qm*(y(6)*B(1)-y(4)*B(3)+E(2))-betam*v*y(5);
dy(6)=qm*(y(4)*B(2)-y(5)*B(1)+E(3))-betam*v*y(6);
end

