O problema do Pêndulo Invertido é um clássico da teoria de controle. Existem muitas versões possíveis mas vamos analisar o seguinte modêlo:
que está detalhado na página da Wikipédia. As equações do movimento são: $$ \begin{gather}\left ( M + m \right ) \ddot x - m l \ddot \theta \cos \theta + m l \dot \theta^2 \sin \theta = F \\ (-g \sin \theta - \ddot x \cos \theta + l \ddot \theta) = 0\end{gather}$$ Pensaremos na força $F$ (aceleração do carrinho) como o parâmetro de controle. O sistema é uma equação diferencial de segunda ordem que pode ser colocado na forma explicita como:
$$ \begin{equation*} \begin{array}{l} { \begin{bmatrix} (M + m) && - m l \cos \theta \\ -\cos \theta && l \end{bmatrix} } \begin{bmatrix} \ddot x \\ \ddot \theta \end{bmatrix} = { \begin{bmatrix} F - m l \dot \theta^2 \sin \theta \\ g \sin \theta \end{bmatrix}} \\ \\ A \ddot u = B \text{, onde: } u = \begin{bmatrix} x \\ \theta \end{bmatrix} \end{array} \end{equation*}$$ Resolvendo este sistema temos
$$ \begin{gather} A^{-1} = \frac{1}{\det(A)}adj(A) \\ det A = (M + m)l - ml \cos^2 \theta = l[(M + m) -m\cos^2 \theta] = l(M + (1-\cos^2 \theta)m) = l(M + m \sin^2 \theta) \\ adj(A) = \begin{bmatrix} l & ml \cos \theta \\ \cos \theta & (M + m) \end{bmatrix} \end{gather}$$daí: $$ \begin{gather} \begin{bmatrix} \ddot x \\ \ddot \theta \end{bmatrix} = \frac{1}{l(M + m \sin^2 \theta)} \begin{bmatrix} l & ml \cos \theta \\ \cos \theta & (M + m) \end{bmatrix}\begin{bmatrix} F - m l \dot{\theta}^2 \sin \theta \\ g \sin \theta \end{bmatrix} \end{gather}$$
O sistema acima pode ser colocado na forma de um sistema de primeira ordem $$\begin{equation*} \begin{array}{l} \dot x = v \\ \dot v = \frac{Fl - ml^2 w^2 \sin \theta + mlg \cos \theta \sin \theta}{l(M + m \sin^2 \theta) } \\ \dot \theta = w \\ \dot w = \frac{ F \cos \theta - ml w^2 \cos \theta \sin \theta + (M + m)g \sin \theta }{l(M + m \sin^2 \theta) } \end{array} \end{equation*}$$ Tomando inicialmente $F=0$, o ponto $(x_0,0,0,0)$ é um ponto de equilíbrio do sistema.
O lado direito do sistema acima pode ser pensando como uma função $G: \mathbb{R}^5 \to \mathbb{R}^4$, ou seja $$G(x,v,\theta,w,F) = (G_1(x,v,\theta,w,F), G_2(x,v,\theta,w,F), G_3(x,v,\theta,w,F), G_4(x,v,\theta,w,F))$$ com $$ \begin{gather} G_1(x,v,\theta,w,F) = v \\ G_2(x,v,\theta,w,F) = \frac{Fl - ml^2 w^2 \sin \theta + mlg \cos \theta \sin \theta}{l(M + m \sin^2 \theta) } \\[0.5mm] G_3(x,v,\theta,w,F) = w \\[0.5mm] G_4(x,v,\theta,w,F) = \frac{ F \cos \theta - ml w^2 \cos \theta \sin \theta + (M + m)g \sin \theta }{l(M + m \sin^2 \theta) } \end{gather}$$
Vamos denotar $\mathbf{x} = (x,v,\theta,w)$. O sistema linearizado em torno do ponto de equilíbrio $(x_0,0,0,0,0)$ é o sistema: $$ \dot{\mathbf{y}} = \frac{\partial{G}}{\partial{\mathbf{x}}}\mathbf{y} + \frac{\partial{G}}{\partial{F}}U$$ onde Em nosso caso: $$ \frac{\partial{G}}{\partial{\mathbf{x}}} = \begin{pmatrix}0 & 1 & 0 & 0 \\ 0 & 0 & \frac{mg}{M} & 0 \\ 0 & 0 & 0 & 1 \\ 0 & 0 & \frac{(M+m)g}{lM} & 0\end{pmatrix} $$ e $$ \frac{\partial{G}}{\partial{F}} = \begin{pmatrix} 0 \\ \frac{1}{M} \\ 0 \\ \frac{1}{lM} \end{pmatrix} $$
%matplotlib inline
from matplotlib import pyplot as plt
import numpy as np
from scipy.integrate import odeint
def rhs(y,t):
return [y[1], (F*l - m*l**2*y[3]**2*np.sin(y[2])+m*l*g*np.cos(y[2])*np.sin(y[2]))/(l*(M+m*np.sin(y[2])**2)),y[3], (F*np.cos(y[2]) -m*l*y[3]**2*np.sin(y[2])*np.cos(y[2])+(M+m)*g*np.sin(y[2]))/(l*(M+m*np.sin(y[2])**2))]
#parametros
M=5
m=1
l=1
g=9.8
F=0
# malha de tempo
t = np.arange(0,200,0.1)
# condição inicial
p0=[1,0.1,0,-0.1]
sol = odeint(rhs,p0,t)
theta = sol[:,2]
omega = sol[:,3]
plt.plot(t,theta)
plt.grid()
plt.plot(t,omega)
plt.grid()
plt.plot(theta,omega)
plt.grid()