function dx = RR(t,x)
dx = zeros(4,1);
dx(1) = x(2);
dx(2) = -3/4*(12*sin(x(1))^2*sin(x(3))*cos(x(3))*x(2)^2+12*cos(x(1))*cos(x(3))^2*sin(x(1))*x(2)^2-120*sin(x(1))*sin(x(3))*cos(x(3))+sin(x(1))*cos(x(3))*x(4)^2+3*sin(x(1))*cos(x(3))*x(4)^2*sin(x(3))^2+3*sin(x(1))*cos(x(3))^3*x(4)^2-12*sin(x(1))*sin(x(3))^2*cos(x(1))*x(2)^2-cos(x(1))*sin(x(3))*x(4)^2-3*cos(x(1))*sin(x(3))^3*x(4)^2+60*cos(x(1))+60*cos(x(1))*cos(x(3))^2+180*sin(x(3))^2*cos(x(1))-3*cos(x(1))*sin(x(3))*x(4)^2*cos(x(3))^2-12*cos(x(1))^2*cos(x(3))*sin(x(3))*x(2)^2)/(3*cos(x(1))^2+1+3*sin(x(1))^2+3*cos(x(3))^2+9*cos(x(3))^2*sin(x(1))^2+9*sin(x(3))^2*cos(x(1))^2+3*sin(x(3))^2-18*cos(x(1))*cos(x(3))*sin(x(1))*sin(x(3)));
dx(3) = x(4);
dx(4) = 3/(3*cos(x(1))^2+1+3*sin(x(1))^2+3*cos(x(3))^2+9*cos(x(3))^2*sin(x(1))^2+9*sin(x(3))^2*cos(x(1))^2+3*sin(x(3))^2-18*cos(x(1))*cos(x(3))*sin(x(1))*sin(x(3)))*(60*cos(x(3))*cos(x(1))^2-40*cos(x(3))-120*cos(x(3))*sin(x(1))^2+12*cos(x(3))*sin(x(1))*x(2)^2*cos(x(1))^2+4*cos(x(3))*sin(x(1))*x(2)^2+12*cos(x(3))*sin(x(1))^3*x(2)^2-12*sin(x(3))*cos(x(1))^3*x(2)^2-4*sin(x(3))*cos(x(1))*x(2)^2-12*sin(x(3))*cos(x(1))*x(2)^2*sin(x(1))^2-3*cos(x(1))^2*cos(x(3))*sin(x(3))*x(4)^2+3*cos(x(1))*cos(x(3))^2*sin(x(1))*x(4)^2+180*sin(x(1))*sin(x(3))*cos(x(1))-3*sin(x(1))*sin(x(3))^2*cos(x(1))*x(4)^2+3*sin(x(1))^2*sin(x(3))*cos(x(3))*x(4)^2);