function val=funval(z)
x=z(1); y=z(2); v=z(3); w=z(4);
val=[v;w;-x/(x^2+y^2)^(3/2);-y/(x^2+y^2)^(3/2)];