2nd order coupled differential equations Working on a mechanical problem I obtain a system of two coupled differential equations (Eq.1 and Eq.2) in r and theta using polar coordinates system.
I have summarized everything on the figure attached with this post and I hope that it is quite explicit. Would you have informations on the numerical method to be used to solve such a system (using Scilab or Matlab or whatever).
Many thanks to all.

 A: Here is a solution.
First I post the equations system (I corrected 2 mistakes).

Now, I propose a Scilab code (equivalent to Matlab) to solve it by creating a first order system of differential equations.

t0 =  0; // initial time (s)
tmax = 30; // solving max time (s)
tinc = 0.005; // time increment (s)

t = t0:tinc:tmax;

// Differential equations system

function du = equadiff(t,u)

  // d for derivative with respect to time t
  // u(1) = theta
  // u(2) = d(theta)/dt
  // u(3) = r
  // u(4) = d(r)/dt
  // du(1) = d(theta)/dt
  // du(2) = d²(theta)dt²
  // du(3) = d(r)/dt
  // du(4) = d²(r)/dt²

  F1=1300;
  F2=1000;
  l=13;
  k=68000;
  menv = 112.14;
  m0 = 56.07;
  g=9.81;

  m = (menv - m0) * (1-abs(cos(u(1)))) + m0;
  du(1) = u(2);
  du(2) = (F2*sin(u(1)) + (F1-m*g).*cos(u(1)) - 2*m.*u(4).*u(2)) ./ (m.*u(3));
  du(3) = u(4);
  du(4) = ((F1-m*g).*sin(u(1)) - F2*cos(u(1)) - k*(u(3)-l) + m.*u(3).*u(2).*u(2)) ./ m;

endfunction

// Integration

u0 = [0,0,13,0]';
[u] = ode(u0,t0,t,equadiff);

// Plot theta for example
plot(t,u(1,:));


