Solution to the inverted pendulum problem (system of non-linear ODE)

Question

A classic control systems problem is that of the translational inverted pendulum, as shown below. By applying a force $f(t)$ to the car (e.g. through a DC motor attached to the car), the goal is to adequately adjust $x(t)$ and $\theta(t)$ such that the rod is mostly vertical (i.e. $\theta = 0$). Let $M_p$ be the mass of the pendulum (actually, the ball), $M_c$ the mass of the car, $L$ the length of the rod, and $g$ the magnitude of the acceleration of gravity.

After some assumptions (i.e. massless rod, zero moment of inertia for the ball, no friction between wheels of the car and the ground), but without any approximations, the system of two second-order non-linear ODEs is

$\cos{[\theta(t)]} \ddot x(t) - L \ddot \theta(t) + g \sin{[\theta(t)]} = 0 \tag*{}$

$-M_p L \cos{[\theta(t)]} \ddot \theta(t) + M_p L \sin{[\theta(t)]} \dot \theta^2(t) + (M_c + M_p) \ddot x(t) = f(t) \tag*{}$

Using Mathematica, I'd like to get the particular solution for this problem, given the initial conditions $\theta(0) = \theta_0$, $\dot \theta(0) = \theta_1$, $x(0) = x_0$ and $\dot x(0) = x_1$. Now, I don't need a symbolic solution; we can assume all initials conditions to be zero except $\theta(0) = \pi/10$, and letting $M_p = 0.1$ (in kg), $M_c = 10$, $L = 0.03$ (in meters), $g = 9.81$ (in meters per second squared), and $f(t) = 0$. In such case, after defining the variables, this is my code for the ODEs:

DSolve[{Cos[θ[t]] x''[t] - L θ''[t] + g Sin[θ[t]] == 0, -Mp L Cos[θ[t]] θ''[t] + Mp L Sin[θ[t]] (θ'[t])^2 + (Mc + Mp) x''[t] == 0 , θ[0] == θ0, θ'[0] == θ1, x[0] == x0, x'[0] == x1}, {θ[t], x[t]}, t]

However, Mathematica outputs just what I typed, without solving:

I only want the expression of the unknown functions $\theta(t)$ and $x(t)$. Do I have to use NDSolve? If so, how? I know this system is unstable but is controllable and observable.

Answer 1

Do I have to use NDSolve?

Yes

If so, how?

Something like

ClearAll[θ,t,g,x];
L   =  1;
Mp  =  2;
Mc  =  3;
θ0  = 5 Degree;
x1  = 0;
x0  = 4;
θ1  = 0;
g   = 9.81;

ode1 = Cos[θ[t]] x''[t] - L θ''[t] + g Sin[θ[t]] == 0;
ode2 = - Mp L Cos[θ[t]] θ''[t] + Mp L Sin[θ[t]] (θ'[t])^2 + (Mc + Mp) x''[t] == 0;
ic  = {θ[0] == θ0, θ'[0] == θ1, x[0] == x0, x'[0] == x1};

sol = NDSolve[{ode1,ode2,ic},{θ, x}, {t,0,5}]

See help on NDSolve on how to do processing on result, i.e. plot, etc..