function du = satellite3j20 (t, u, GM, RZ)
r =  sqrt(u(1)^2+u(3)^2+u(5)^2);
konst= 1.5*RZ^2*1086e-6/r^2;
gr = -GM/r^3;
du(1) = u(2); 
du(2) = gr*u(1)*(1 - konst*(5*u(5)^2/r^2-1)) ;
du(3) = u(4); 
du(4) = u(3)/u(1)*du(2);
du(5) = u(6); 
du(6) = gr*u(5)*(1+konst*(3-5*u(5)^2/r^2));
endfunction
