La peonza

La peonza (II)

Momento angular

El momento angular de un sólido rígido que gira con velocidad angular ω alrededor de un eje instantáneo que pasa por el punto fijo O es

$\vec{L} = \sum_{i} m_{i} ({\vec{r}}_{i} \times \vec{v_{i}}) = \sum_{i} m_{i} (\vec{r_{i}} \times (\vec{ω} \times \vec{r_{i}}))$

Teniendo en cuenta el resultado del triple producto vectorial

$\vec{a} \times (\vec{b} \times \vec{c}) = (\vec{a} \cdot \vec{c}) \vec{b} - (\vec{a} \cdot \vec{b}) \vec{c}$

Obtenemos

$\begin{array}{l} \sum_{i} m_{i} ((\vec{r_{i}} \cdot \vec{r_{i}}) \vec{ω} - (\vec{r_{i}} \cdot \vec{ω}) \vec{r_{i}}) = \sum_{i} m_{i} (r_{i}^{2} \vec{ω} - (\vec{r_{i}} \cdot \vec{ω}) \vec{r_{i}}) = \\ {\begin{cases} L_{x} = \sum_{i} m_{i} ((x_{i}^{2} + y_{i}^{2} + z_{i}^{2}) ω_{x} - (x_{i} ω_{x} + y_{i} ω_{y} + z_{i} ω_{z}) x_{i}) \\ L_{y} = \sum_{i} m_{i} ((x_{i}^{2} + y_{i}^{2} + z_{i}^{2}) ω_{y} - (x_{i} ω_{x} + y_{i} ω_{y} + z_{i} ω_{z}) y_{i}) \\ L_{z} = \sum_{i} m_{i} ((x_{i}^{2} + y_{i}^{2} + z_{i}^{2}) ω_{z} - (x_{i} ω_{x} + y_{i} ω_{y} + z_{i} ω_{z}) z_{i}) \end{cases} \\ {\begin{cases} L_{x} = ω_{x} (\sum_{i} m_{i} (y_{i}^{2} + z_{i}^{2})) - ω_{y} (\sum_{i} m_{i} x_{i} y_{i}) - ω_{z} (\sum_{i} m_{i} x_{i} z_{i}) \\ L_{y} = - ω_{x} (\sum_{i} m_{i} x_{i} y_{i}) + ω_{y} (\sum_{i} m_{i} (x_{i}^{2} + z_{i}^{2})) - ω_{z} (\sum_{i} m_{i} y_{i} z_{i}) \\ L_{z} = - ω_{x} (\sum_{i} m_{i} x_{i} z_{i}) - ω_{y} (\sum_{i} m_{i} y_{i} z_{i}) + ω_{z} (\sum_{i} m_{i} (x_{i}^{2} + z_{i}^{2})) \end{cases} \end{array}$

En forma matricial

$\begin{array}{l} (\begin{matrix} L_{x} \\ L_{y} \\ L_{z} \end{matrix}) = (\begin{matrix} I_{x x} & I_{x y} & I_{x z} \\ I_{x y} & I_{y y} & I_{y z} \\ I_{x z} & I_{y z} & I_{z z} \end{matrix}) (\begin{matrix} ω_{x} \\ ω_{y} \\ ω_{z} \end{matrix}) \\ \vec{L} = I \vec{ω} \end{array}$

I se denomina tensor de inercia y es una matriz simétrica

En general, el vector momento angular $\vec{L}$ , no es paralelo al vector velocidad angular, $\vec{ω}$

Para una distribución continua de masa

$\begin{array}{l} I_{x y} = - \int_{V} ρ (x, y, z) x y \cdot d V \\ I_{x x} = \int_{V} ρ (x, y, z) (y^{2} + z^{2}) d V \end{array}$

ρ(x,y,z) es la densidad que podría depnder de las coordendas del punto considerado del sólido. Estudiaremos sólidos homogéneos en los que la densidad ρ es constante.

Energía cinética

La energía cinética del sistema de partículas es

$E_{k} = \frac{1}{2} \sum_{i} m_{i} v_{i}^{2} = \frac{1}{2} \sum_{i} m_{i} (\vec{ω} \times \vec{r_{i}}) \cdot \vec{v_{i}}$

Teniendo en cuenta la propiedad del producto mixto

$\vec{a} \cdot (\vec{b} \times \vec{c}) = (\vec{a} \times \vec{b}) \vec{c}$

Obtenemos

$\begin{array}{l} E_{k} = \frac{1}{2} \sum_{i} m_{i} \vec{v_{i}} \cdot (\vec{ω} \times \vec{r_{i}}) = \frac{1}{2} \sum_{i} m_{i} \cdot \vec{r_{i}} (\vec{v_{i}} \times \vec{ω}) = \frac{1}{2} \vec{ω} \sum_{i} m_{i} (\vec{r_{i}} \times \vec{v_{i}}) \\ = \frac{1}{2} \vec{ω} \cdot \vec{L} = \frac{1}{2} \vec{ω} \cdot (I \vec{ω}) = \frac{1}{2} (\begin{matrix} ω_{x} & ω_{y} & ω_{z} \end{matrix}) (\begin{matrix} I_{x x} & I_{x y} & I_{x z} \\ I_{x y} & I_{y y} & I_{y z} \\ I_{x z} & I_{y z} & I_{z z} \end{matrix}) (\begin{matrix} ω_{x} \\ ω_{y} \\ ω_{z} \end{matrix}) = \\ \frac{1}{2} (\begin{matrix} ω_{x} & ω_{y} & ω_{z} \end{matrix}) (\begin{matrix} I_{x x} ω_{x} + I_{x y} ω_{y} + I_{x z} ω_{z} \\ I_{x y} ω_{x} + I_{y y} ω_{y} + I_{y z} ω_{z} \\ I_{x z} ω_{x} + I_{y z} ω_{y} + I_{z z} ω_{z} \end{matrix}) = \\ \frac{1}{2} (I_{x x} ω_{x}^{2} + + I_{y y} ω_{y}^{2} + I_{z z} ω_{z}^{2} + 2 I_{x y} ω_{x} ω_{y} + 2 I_{x z} ω_{x} ω_{z} + 2 I_{y z} ω_{y} ω_{z}) \end{array}$

Los elementos del tensor de inercia dependen del origen O y de la orientación de los ejes coordenados. Es posible enocontar unos ejes, denominados principales de inercia para los cuales solamente los elementos diagonales del tensor son no nulos

$I = (\begin{matrix} I_{1} & 0 & 0 \\ 0 & I_{2} & 0 \\ 0 & 0 & I_{3} \end{matrix})$

La energía cinética y el momento angular serían

$\begin{array}{l} E_{k} = \frac{1}{2} (I_{1} ω_{1}^{2} + I_{2} ω_{2}^{2} + I_{3} ω_{3}^{2}) \\ {\begin{cases} L_{1} = I_{1} ω_{1} \\ L_{2} = I_{2} ω_{2} \\ L_{3} = I_{2} ω_{3} \end{cases} \end{array}$

El momento angular y la velocidad angular tendrán la misma dirección a lo largo de un eje principal de inercia

Tensor de inercia

Calculamos el tensor de inercia de una placa y de un cono macizo (versión simple de una peonza)

Placa homogénea cuadrada

Sea una placa cuadrada homogénea de densidad σ kg/m², de lado a situada en el plano XY

Los elementos del tensor de inercia son

$\begin{array}{l} I_{x x} = \int_{0}^{a} \int_{0}^{a} σ y^{2} d x \cdot d y = σ \int_{0}^{a} (\int_{0}^{a} y^{2} d y) d x = \frac{m}{a^{2}} \frac{a^{3}}{3} a = \frac{1}{3} m a^{2} \\ I_{y y} = \int_{0}^{a} \int_{0}^{a} σ x^{2} d x \cdot d y = σ \int_{0}^{a} (\int_{0}^{a} x^{2} d x) d y = \frac{m}{a^{2}} \frac{a^{3}}{3} a = \frac{1}{3} m a^{2} \\ I_{z z} = \int_{0}^{a} \int_{0}^{a} σ (x^{2} + y^{2}) d x \cdot d y = σ {\int_{0}^{a} (\int_{0}^{a} x^{2} d x) d y + \int_{0}^{a} (\int_{0}^{a} y^{2} d y) d x =} \frac{m}{a^{2}} (\frac{a^{3}}{3} a + \frac{a^{3}}{3} a) = \frac{1}{3} m a^{2} \\ I_{x y} = - \int_{0}^{a} \int_{0}^{a} σ x y \cdot d x \cdot d y = - σ (\int_{0}^{a} x d x) (\int_{0}^{a} y d y) = - \frac{m}{a^{2}} \frac{a^{2}}{2} \frac{a^{2}}{2} = - \frac{1}{4} m a^{2} \\ I_{x z} = 0 \\ I_{y z} = 0 \end{array}$

El tensor de inercia es

$I = (\begin{matrix} \frac{1}{3} m a^{2} & - \frac{1}{4} m a^{2} & 0 \\ - \frac{1}{4} m a^{2} & \frac{1}{3} m a^{2} & 0 \\ 0 & 0 & \frac{2}{3} m a^{2} \end{matrix})$

Cono macizo homogéneo

Sea un cono macizo, homogéneo de densidad ρ kg/m³. de altura h y radio de l abase R. El ángulo del vértice es 2α, con tanα=R/h

Para dibujar el cono, se ha empleado el código

theta=pi/3; %angulo del cono
r=linspace(0,1,30);
phi=linspace(0,2*pi,30);
[r,phi]=meshgrid(r,phi);
x=r.*cos(phi)*sin(theta);
y=r.*sin(phi)*sin(theta);
z=r*cos(theta);
hold on
surfl(x,y,z) %cono
surfl(x,y,cos(theta)*ones(length(r))) %tapa
shading interp
colormap(gray);
hold off
axis off
view(20,30)

El elemento de volumen en coordenadas cilíndricas es dV=(r·dφ)·dr·dz, tal como se muestra en a figura. Tenemos que calcular las integrales

$\begin{array}{l} I_{x x} = \int_{V} ρ (y^{2} + z^{2}) d V \\ I_{y y} = \int_{V} ρ (x^{2} + z^{2}) d V = I_{x x} \\ I_{z z} = \int_{V} ρ (x^{2} + y^{2}) d V = \int_{V} ρ r^{2} \cdot d V \end{array}$

El cono es un sólido de revolución, por simetría, I_xx=I_yy

$\begin{array}{l} I_{x x} = \int_{V} ρ (y^{2} + z^{2}) d V = ρ \int_{0}^{2 π} d φ \int_{0}^{R} 2 π r \cdot d r \int_{h r / R}^{h} (r^{2} \sin^{2} φ + z^{2}) d z \\ ρ \int_{0}^{2 π} d φ \int_{0}^{R} r \cdot (r^{2} \sin^{2} φ (h - \frac{h r}{R}) + \frac{h^{3}}{3} - \frac{1}{3} \frac{h^{3} r^{3}}{R^{3}}) d r = \\ ρ \int_{0}^{2 π} d φ \int_{0}^{R} (\sin^{2} φ (h r^{3} - \frac{h r^{4}}{R}) + \frac{h^{3} r}{3} - \frac{1}{3} \frac{h^{3} r^{4}}{R^{3}}) d r = \\ ρ \int_{0}^{2 π} {\sin^{2} φ (h \frac{R^{4}}{4} - \frac{h R^{4}}{5}) + \frac{h^{3} R^{2}}{6} - \frac{1}{3} \frac{h^{3} R^{2}}{5}} d φ = \\ ρ \int_{0}^{2 π} {\sin^{2} φ (\frac{h R^{4}}{20}) + \frac{h^{3} R^{2}}{10}} d φ = ρ {(\frac{h R^{4}}{20}) \int_{0}^{2 π} \sin^{2} φ \cdot d φ + \frac{h^{3} R^{2}}{10} \int_{0}^{2 π} d φ} = \end{array}$

Teniendo en cuenta el resultado de la integral, llegamos a

$\begin{array}{l} \int \sin^{2} φ \cdot d φ = \frac{1}{2} (φ - \frac{1}{2} \sin (2 φ)) \\ I_{x x} = ρ \frac{h R^{2}}{5} π (h^{2} + \frac{R^{2}}{4}) \end{array}$

Calculamos I_zz

$\begin{array}{l} I_{z z} = \int_{V} ρ (x^{2} + y^{2}) d V = \int_{V} ρ r^{2} \cdot d V \\ = ρ \int_{0}^{2 π} d φ \int_{0}^{R} r^{2} \cdot r d r \int_{h r / R}^{h} d z = ρ \int_{0}^{2 π} d φ \int_{0}^{R} r^{3} \cdot (h - \frac{h r}{R}) d r = \\ ρ \int_{0}^{2 π} d φ \int_{0}^{R} (h r^{3} - \frac{h r^{4}}{R}) d r = ρ \frac{h R^{4}}{20} \int_{0}^{2 π} d φ = ρ \frac{h R^{4}}{20} 2 π \\ I_{z z} = ρ \frac{h R^{4}}{10} π \end{array}$

Calculamos los elementos fuera de la diagonal principal, I_xy

$\begin{array}{l} I_{x y} = - \int_{V} ρ x y d V = - ρ \int_{0}^{2 π} \cos φ \cdot \sin φ \cdot d φ \int_{0}^{R} r \cdot r^{2} d r \int_{h r / R}^{h} d z = \\ - ρ \int_{0}^{2 π} \cos φ \cdot \sin φ \cdot d φ \int_{0}^{R} r^{3} (h \frac{r}{R} - h) d r = \\ - ρ \frac{h R^{4}}{20} \int_{0}^{2 π} \cos φ \cdot \sin φ \cdot d φ = - ρ \frac{h R^{4}}{20} \frac{1}{2} \int_{0}^{2 π} \sin (2 φ) d φ = 0 \end{array}$

Del mismo modo I_xz=I_yz=0

$\begin{array}{l} I_{x z} = - \int_{V} ρ x z d V = - ρ \int_{0}^{2 π} \cos φ \cdot d φ \int_{0}^{R} r \cdot r d r \int_{h r / R}^{h} z d z = \\ - ρ \int_{0}^{2 π} \cos φ \cdot d φ \int_{0}^{R} r^{2} \frac{1}{2} ({(h \frac{r}{R})}^{2} - h^{2}) d r = \\ ρ \frac{h^{2} R^{3}}{15} \int_{0}^{2 π} \cos φ \cdot d φ = ρ \frac{h^{2} R^{3}}{15} (\sin (2 π) - \sin (0)) = 0 \end{array}$

El tensor de inercia es

$I = (\begin{matrix} ρ \frac{h R^{2}}{5} π (h^{2} + \frac{R^{2}}{4}) & 0 & 0 \\ 0 & ρ \frac{h R^{2}}{5} π (h^{2} + \frac{R^{2}}{4}) & 0 \\ 0 & 0 & ρ \frac{h R^{4}}{10} π \end{matrix})$

Angulos de Euler

Utilizamos el concepto de componente de un vector $\vec{a}$

$\begin{array}{l} \vec{a} = a_{x} \hat{i} + a_{y} \hat{j} + a_{z} \hat{k} \\ \vec{a} = (\vec{a} \cdot \hat{i}) \hat{i} + (\vec{a} \cdot \hat{j}) \hat{j} + (\vec{a} \cdot \hat{k}) \hat{k} \end{array}$

El sistema de referencia inercial fijo en el espacio es x, y z

Primera rotación alrededor del eje z un ángulo φ

El sistema de referencia rotado un ángulo φ es X, Y, Z=z

$\begin{array}{l} \hat{i} = (\hat{i} \cdot \hat{I}) \hat{I} + (\hat{i} \cdot \hat{J}) \hat{J} + (\hat{i} \cdot \hat{K}) \hat{K} = \cos φ \hat{I} - \sin φ \hat{J} \\ \hat{j} = (\hat{j} \cdot \hat{I}) \hat{I} + (\hat{j} \cdot \hat{J}) \hat{J} + (\hat{j} \cdot \hat{K}) \hat{K} = \sin φ \hat{I} + \cos φ \hat{J} \\ \hat{k} = (\hat{k} \cdot \hat{I}) \hat{I} + (\hat{k} \cdot \hat{J}) \hat{J} + (\hat{k} \cdot \hat{K}) \hat{K} = \hat{K} \end{array}$

Segunda rotación alrededor del nuevo eje X un ángulo θ

El sistema de referencia rotado un ángulo θ es X'=X, Y', Z'

$\begin{array}{l} \hat{I} = (\hat{I} \cdot \hat{I'}) \hat{I'} + (\hat{I} \cdot \hat{J'}) \hat{J'} + (\hat{I} \cdot \hat{K'}) \hat{K'} = \hat{I'} \\ \hat{J} = (\hat{J} \cdot \hat{I'}) \hat{I'} + (\hat{J} \cdot \hat{J'}) \hat{J'} + (\hat{J} \cdot \hat{K'}) \hat{K'} = \cos θ \hat{J'} - \sin θ \hat{K'} \\ \hat{K} = (\hat{K} \cdot \hat{I'}) \hat{I'} + (\hat{K} \cdot \hat{J'}) \hat{J'} + (\hat{K} \cdot \hat{K'}) \hat{K'} = \sin θ \hat{J'} + \cos θ \hat{K'} \end{array}$

Tercera rotación alrededor del nuevo eje Z' un ángulo ψ

El sistema de referencia rotado un ángulo ψ es x', y', z'=Z'

$\begin{array}{l} \hat{I'} = (\hat{I'} \cdot \hat{i'}) \hat{i'} + (\hat{I'} \cdot \hat{j'}) \hat{j'} + (\hat{I'} \cdot \hat{k'}) \hat{k'} = \cos ψ \hat{i'} - \sin ψ \hat{j'} \\ \hat{J'} = (\hat{J'} \cdot \hat{i'}) \hat{i'} + (\hat{J'} \cdot \hat{j'}) \hat{j'} + (\hat{J'} \cdot \hat{k'}) \hat{k'} = \sin ψ \hat{i'} + \cos ψ \hat{j'} \\ \hat{K'} = (\hat{K'} \cdot \hat{i'}) \hat{i'} + (\hat{K'} \cdot \hat{j'}) \hat{j'} + (\hat{K'} \cdot \hat{k'}) \hat{k'} = \hat{k'} \end{array}$

Relacionamos el sistema de referencia inercial (x,y,z) y el sistema de referencia (x',y',z') ligado al sólido rígido

$\begin{array}{l} {\begin{cases} \hat{I'} = \cos ψ \hat{i'} - \sin ψ \hat{j'} \\ \hat{J'} = \sin γ ψ + \cos ψ \hat{j'} \\ \hat{K'} = \hat{k'} \end{cases} {\begin{cases} \hat{I} = \hat{I'} \\ \hat{J} = \cos θ \hat{J'} - \sin θ \hat{K'} \\ \hat{K} = \sin θ \hat{J'} + \cos θ \hat{K'} \end{cases} {\begin{cases} \hat{i} = \cos φ \hat{I} - \sin φ \hat{J} \\ \hat{j} = \sin φ \hat{I} + \cos φ \hat{J} \\ \hat{k} = \hat{K} \end{cases} \\ \hat{i} = \cos φ \hat{I} - \sin φ \hat{J} = \cos φ \hat{I'} - \sin φ (\cos θ \hat{J'} - \sin θ \hat{K'}) = \cos φ \hat{I'} - \sin φ \cos θ \hat{J'} + \sin φ \sin θ \hat{K'} = \\ \cos φ (\cos ψ \hat{i'} - \sin ψ \hat{j'}) - \sin φ \cos θ (\sin ψ \hat{i'} + \cos ψ \hat{j'}) + \sin φ \sin θ \hat{k'} \\ = (\cos φ \cos ψ - \sin φ \cos θ \sin ψ) \hat{i'} - (\cos φ \sin ψ + \sin φ \cos θ \cos ψ) \hat{j'} + \sin φ \sin θ \hat{k'} \\ \hat{j} = \sin φ \hat{I} + \cos φ \hat{J} = \sin φ \hat{I'} + \cos φ (\cos θ \hat{J'} - \sin θ \hat{K'}) = \sin φ \hat{I'} + \cos φ \cos θ \hat{J'} - \cos φ \sin θ \hat{K'} = \\ \sin φ (\cos ψ \hat{i'} - \sin ψ \hat{j'}) + \cos φ \cos θ (\sin ψ \hat{i'} + \cos ψ \hat{j'}) - \cos φ \sin θ \hat{k'} = \\ = (\sin φ \cos ψ + \cos φ \cos β \sin ψ) \hat{i'} + (- \sin φ \sin ψ + \cos φ \cos θ \cos ψ) \hat{j'} - \cos φ \sin θ \hat{k'} \\ \hat{k} = \sin θ \hat{J'} + \cos θ \hat{K'} = \sin θ (\sin ψ \hat{i'} + \cos ψ \hat{j'}) + \cos θ \hat{k'} = (\sin θ \sin ψ) \hat{i'} + (\sin θ \cos ψ) \hat{j'} + \cos θ \hat{k'} \end{array}$

El resultado es

${\begin{cases} \hat{i} = (\cos φ \cos ψ - \sin φ \cos θ \sin ψ) \hat{i'} - (\cos φ \sin ψ + \sin φ \cos θ \cos ψ) \hat{j'} + \sin φ \sin θ \hat{k'} \\ \hat{j} = (\sin φ \cos ψ + \cos φ \cos θ \sin ψ) \hat{i'} + (- \sin φ \sin ψ + \cos φ \cos θ \cos ψ) \hat{j'} - \cos φ \sin θ \hat{k'} \\ \hat{k} = (\sin θ \sin ψ) \hat{i'} + (\sin θ \cos ψ) \hat{j'} + \cos θ \hat{k'} \end{cases}$

En forma matricial

$(\begin{matrix} \hat{i} \\ \hat{j} \\ \hat{k} \end{matrix}) = (\begin{matrix} \cos φ & - \sin φ & 0 \\ \sin φ & \cos φ & 0 \\ 0 & 0 & 1 \end{matrix}) (\begin{matrix} 1 & 0 & 0 \\ 0 & \cos θ & - \sin θ \\ 0 & \sin θ & \cos θ \end{matrix}) (\begin{matrix} \cos ψ & - \sin ψ & 0 \\ \sin ψ & \cos ψ & 0 \\ 0 & 0 & 1 \end{matrix}) (\begin{matrix} \hat{i'} \\ \hat{j'} \\ \hat{k'} \end{matrix})$

Velocidad angular

La dirección de la velocidad angular es la del eje de rotación, sentido, la regla del sacacorchos

El eje Z=z, es el primer eje de rotación con velocidad angular dφ/dt
El eje X, es el segundo eje de rotación con velocidad angular dθ/dt
El eje z'=Z', es el tercer eje de rotación con velocidad angular dψ/dt

$\begin{array}{l} \vec{ω} = \frac{d φ}{d t} \hat{k} + \frac{d θ}{d t} \hat{I} + \frac{d ψ}{d t} \hat{k'} = \frac{d φ}{d t} (\sin θ \hat{J'} + \cos θ \hat{K'}) + \frac{d θ}{d t} \hat{I'} + \frac{d ψ}{d t} \hat{k'} = \\ \frac{d φ}{d t} (\sin θ (\sin ψ \hat{i'} + \cos ψ \hat{j'}) + \cos θ \hat{k'}) + \frac{d θ}{d t} (\cos ψ \hat{i'} - \sin ψ \hat{j'}) + \frac{d ψ}{d t} \hat{k'} = \\ \vec{ω} = (\frac{d φ}{d t} \sin θ \sin ψ + \frac{d θ}{d t} \cos ψ) \hat{i'} + (\frac{d φ}{d t} \sin θ \cos ψ - \frac{d θ}{d t} \sin ψ) \hat{j'} + (\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t}) \hat{k'} \end{array}$

El resultado es

${\begin{cases} ω_{1} = \frac{d φ}{d t} \sin θ \sin ψ + \frac{d θ}{d t} \cos ψ \\ ω_{2} = \frac{d φ}{d t} \sin θ \cos ψ - \frac{d θ}{d t} \sin ψ \\ ω_{3} = \frac{d φ}{d t} \cos θ + \frac{d ψ}{d t} \end{cases}$

La energía cinética se rotación de la peonza simétrica de momentos de inercia, I₁=I₂ e I₃

$\begin{array}{l} E_{k} = \frac{1}{2} (I_{1} ω_{1}^{2} + I_{2} ω_{2}^{2} + I_{3} ω_{3}^{2}) = \\ \frac{1}{2} I_{1} {(\frac{d φ}{d t} \sin θ \sin ψ + \frac{d θ}{d t} \cos ψ)}^{2} + \frac{1}{2} I_{2} {(\frac{d φ}{d t} \sin θ \cos ψ - \frac{d θ}{d t} \sin ψ)}^{2} + \frac{1}{2} I_{3} {(\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t})}^{2} \\ E_{k} = \frac{1}{2} I_{1} ({(\frac{d φ}{d t})}^{2} \sin^{2} θ + {(\frac{d θ}{d t})}^{2}) + \frac{1}{2} I_{3} {(\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t})}^{2} \end{array}$

Ecuaciones del movimiento

La energía potencial de la peonza respecto del punto de apoyo es mghcosθ. Siendo h la distancia entre el centro de masas y el punto de apoyo

$L = \frac{1}{2} I_{1} ({(\frac{d φ}{d t})}^{2} \sin^{2} θ + {(\frac{d θ}{d t})}^{2}) + \frac{1}{2} I_{3} {(\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t})}^{2} - m g h \cos θ$

La lagrangiana depende de tres coordenadas φ, θ y ψ, tenemos tres ecuaciones del movimiento

${\begin{cases} \frac{d}{d t} \frac{\partial L}{\partial \dot{φ}} - \frac{\partial L}{\partial φ} = 0 \\ \frac{d}{d t} \frac{\partial L}{\partial \dot{ψ}} - \frac{\partial L}{\partial ψ} = 0 \\ \frac{d}{d t} \frac{\partial L}{\partial \dot{θ}} - \frac{\partial L}{\partial θ} = 0 \end{cases}$

La primera ecuación φ, describe el movimiento de precesión y la tercera ecuación en θ, el movimiento de nutación

Constantes del movimiento

La lagrangiana no depende del ángulo φ, se conserva la cantidad P_φ

$\begin{array}{l} p_{φ} = \frac{\partial L}{\partial \dot{φ}} = cte \\ p_{φ} = I_{1} (\frac{d φ}{d t}) \sin^{2} θ + I_{3} (\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t}) \cos θ \\ p_{φ} = (I_{1} \sin^{2} θ + I_{3} \cos^{2} θ) \frac{d φ}{d t} + I_{3} \frac{d ψ}{d t} \cos θ \end{array}$

La lagrangiana no depende del ángulo ψ, se conserva la cantidad P_ψ

$\begin{array}{l} p_{ψ} = \frac{\partial L}{\partial \dot{ψ}} = cte \\ p_{ψ} = I_{3} (\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t}) \end{array}$

Tenemos un sistema de dos ecuaciones, despejamos dφ/dt y dψ/dt

$\begin{array}{l} \frac{d φ}{d t} = \frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ} \\ \frac{d ψ}{d t} = \frac{p_{ψ}}{I_{3}} - \frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ} \cos θ \end{array}$

Si p_φ<p_ψ, la velocidad angular dφ/dt se anula para el ángulo cosθ_φ=p_φ/p_ψ. Si θ<θ_φ, entonces dφ/dt<0, en caso contrario, dφ/dt>0

La energía es constante

$E = \frac{1}{2} I_{1} ({(\frac{d φ}{d t})}^{2} \sin^{2} θ + {(\frac{d θ}{d t})}^{2}) + \frac{1}{2} I_{3} {(\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t})}^{2} + m g h \cos θ$

El segundo término es una constante del movimiento

$\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t} = \frac{p_{φ} - p_{γ} \cos θ}{I_{1} \sin^{2} θ} \cos θ + \frac{p_{ψ}}{I_{3}} - \frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ} \cos θ = \frac{p_{ψ}}{I_{3}}$

El resultado es

$E = \frac{1}{2} I_{1} ({(\frac{d φ}{d t})}^{2} \sin^{2} θ + {(\frac{d θ}{d t})}^{2}) + \frac{1}{2} I_{3} {(\frac{p_{ψ}}{I_{3}})}^{2} + m g h \cos θ$

La energía E se puede escribir en términos de las otras dos constantes del movimiento, p_φ y p_ψ

$\begin{array}{l} E = \frac{1}{2} I_{1} {(\frac{d θ}{d t})}^{2} + \frac{1}{2} I_{1} {(\frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ})}^{2} \sin^{2} θ + \frac{p_{ψ}^{2}}{2 I_{3}} + m g h \cos θ \\ E = \frac{1}{2} I_{1} {(\frac{d θ}{d t})}^{2} + \frac{1}{2} \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1} \sin^{2} θ} + \frac{p_{ψ}^{2}}{2 I_{3}} + m g h \cos θ \end{array}$

La energía E' es también constante

$\begin{array}{l} E' = E - \frac{p_{ψ}^{2}}{2 I_{3}} \\ E' = \frac{1}{2} I_{1} {(\frac{d θ}{d t})}^{2} + E_{p} (θ) \end{array}$

Definimos una energía potencial efectiva

$E_{p} (θ) = \frac{1}{2} \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1} \sin^{2} θ} + m g h \cos θ$

La energía E' tiene que ser mayor que el mínimo de la energía potencial efectiva. El mínimo se calcula derivado la energía potencial respecto de θ e igualando a cero

$\begin{array}{l} \frac{d E_{p} (θ)}{d θ} = 0 \\ (p_{φ} - p_{ψ} \cos θ) (p_{ψ} - p_{φ} \cos θ) - m g h I_{1} \sin^{4} θ = 0 \end{array}$

Cuando la energía E' es igual al mínimo de la energía potencial efectiva, el ángulo de inclinación del eje de la peonza es constante θ=θ₀

Representamos la energía potencial efectiva E_p(θ) con los siguientes datos

Constante, p_φ=2
Constante, p_ψ=2.5
Masa de la peonza, m=1
Momento de inercia, I₁=1
Distancia entre el c.m. y el punto de apoyo, h=1
Aceleración de la gravedad, g=1

El programa calcula, mediante la función fzero de MATLAB, el ángulo θ₀≈π/4, para el cual la energía potencial es mínima ≈0.76.

pg=2.5; %constantes del movimiento
pa=2;
m=1; %masa
g=1; %aceleración de la gravedad
I1=1; %momento de inercia
h=1; %distancia del c.m. al punto de apoyo
Ep=@(x) ((pa-pg*cos(x)).^2)./(2*I1*sin(x).^2)+m*g*h*cos(x);
f=@(x) (pa-pg*cos(x))*(pg-pa*cos(x))-I1*m*g*h*sin(x)^4;
min=fzero(f,[pi/6,pi/3]);

fplot(Ep,[pi/6,pi/3])
Em=Ep(min);
line([min,min],[0.75,Em],'lineStyle','--')
line([pi/6,min],[Em,Em],'lineStyle','--')
set(gca,'XTick',pi/6:pi/12:pi/3)
set(gca,'XTickLabel',{'\pi/6','\pi/4','\pi/3'})
grid on
xlabel('\beta')
ylabel('E_p')
title('Energía potencial efectiva')

Despejamos dθ/dt en la ecuación de la energía E

$\begin{array}{l} {(\frac{d θ}{d t})}^{2} \sin^{2} θ = \frac{2 E}{I_{1}} \sin^{2} θ - \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1}^{2}} - \frac{p_{ψ}^{2}}{I_{3} I_{1}} \sin^{2} θ - 2 \frac{m g h}{I_{1}} \cos θ \sin^{2} θ \\ {(\frac{d θ}{d t})}^{2} \sin^{2} θ = (\frac{2 E}{I_{1}} - \frac{p_{ψ}^{2}}{I_{3} I_{1}} - 2 \frac{m g h}{I_{1}} \cos θ) \sin^{2} θ - \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1}^{2}} \end{array}$

Llamamos u=cosθ y definimos los parámetros A, B, a y b

$\begin{array}{l} {(\frac{d θ}{d t})}^{2} \sin^{2} θ = (A - B \cos θ) \sin^{2} θ - {(b - a \cos θ)}^{2} \\ {(\frac{d u}{d t})}^{2} = (A - B u) (1 - u^{2}) - {(b - a u)}^{2} \\ {(\frac{d u}{d t})}^{2} = B u^{3} - (A + a^{2}) u^{2} + (2 a b - B) u + A - b^{2} \end{array}$

Integramos

$t = \int \frac{d u}{\sqrt{(A - B u) (1 - u^{2}) - {(b - a u)}^{2}}}$

que es una integral elíptica

Bajo la raíz tenemos un polinomio de tercer grado de coeficientes

$f (u) = u^{3} - \frac{A + a^{2}}{B} u^{2} + \frac{2 a b - B}{B} u + \frac{A - b^{2}}{B}$

El valor absoluto de u, |u|≤1. De las tres raíces, dos cumplirán

$- 1 \leq u_{1} \leq u_{2} \leq 1$

El ángulo del eje de la peonza con la vertical θ estará limitado a θ₁=arccos(u₁) y θ₂=arccos(u₂)

A=1.6;
B=2.0;
a=2.5;
b=1.7;
raices=roots([1,-(A+a^2)/B,(2*a*b-B)/B, (A-b^2)/B]);
f=@(u) (A-B*u).*(1-u.^2)-(b-a*u).^2;
hold on
fplot(f,[0,3])
for k=1:length(raices)
    plot(raices(k),0,'ro','markersize',3,'markerfacecolor','r')
end
hold off
grid on
xlabel('u')
ylabel('f(u)')
title('Inclinación máxima y mínima')

Los ángulos límite en el movimiento de precesión son: θ₁=40.8° y θ₂=72.7°

>> raices
raices =
    2.8714
    0.7568
    0.2968
>> acosd(raices(2))
ans =   40.8148
>> acosd(raices(3))
ans =   72.7340
>> acosd(b/a)
ans =   47.1564

El ángulo θ_φ=47.1° para el cual dφ/dt=0, está comprendido entre θ₁ y θ₂ (tercer caso)

Ecuación del movimiento en θ

Nos falta la ecuación del movimiento en θ

$\begin{array}{l} \frac{d}{d t} \frac{\partial L}{\partial \dot{θ}} - \frac{\partial L}{\partial θ} = 0 \\ I_{1} \frac{d^{2} θ}{d t^{2}} - I_{1} {(\frac{d φ}{d t})}^{2} \sin θ \cos θ + I_{3} (\frac{d φ}{d t} \cos θ + \frac{d ψ}{d t}) \frac{d φ}{d t} \sin θ - m g h \sin θ = 0 \\ I_{1} \frac{d^{2} θ}{d t^{2}} + (I_{3} - I_{1}) {(\frac{d φ}{d t})}^{2} \sin θ \cos θ + I_{3} \frac{d ψ}{d t} \frac{d φ}{d t} \sin θ - m g h \sin θ = 0 \end{array}$

Introducimos las constantes del movimiento

$\begin{array}{l} I_{1} \frac{d^{2} θ}{d t^{2}} + (I_{3} - I_{1}) {(\frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ})}^{2} \sin θ \cos θ + I_{3} (\frac{p_{γ}}{I_{3}} - \frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ} \cos θ) (\frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin^{2} θ}) \sin θ - m g h \sin θ = 0 \\ I_{1} \frac{d^{2} θ}{d t^{2}} - \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1} \sin^{3} θ} \cos θ + p_{γ} \frac{p_{φ} - p_{ψ} \cos θ}{I_{1} \sin θ} - m g h \sin θ = 0 \\ \frac{d^{2} θ}{d t^{2}} = - \frac{(p_{φ} - p_{ψ} \cos θ) (p_{ψ} - p_{φ} \cos θ)}{I_{1}^{2} \sin^{3} θ} + \frac{m g h}{I_{1}} \sin θ \\ \frac{d^{2} θ}{d t^{2}} = \frac{p_{φ}^{2} + p_{ψ}^{2}}{I_{1}^{2}} \frac{\cos θ}{\sin^{3} θ} - \frac{p_{φ} p_{ψ}}{I_{1}^{2}} \frac{3 + \cos (2 θ)}{2 \sin^{3} θ} + \frac{m g h}{I_{1}} \sin θ \end{array}$

Solución numérica de las ecuaciones del movimiento

Resolveremos mediante el procedimiento ode45 de MATLAB el sistema de tres ecuaciones diferenciales.

$\begin{array}{l} b = \frac{p_{φ}}{I_{1}}, a = \frac{p_{ψ}}{I_{1}}, B = 2 \frac{m g h}{I_{1}} \\ {\begin{cases} \frac{d^{2} θ}{d t^{2}} = (b^{2} + a^{2}) \frac{\cos θ}{\sin^{3} θ} - a b \frac{3 + \cos (2 θ)}{2 \sin^{3} θ} + \frac{B}{2} \sin θ \\ \frac{d φ}{d t} = \frac{b - a \cos θ}{\sin^{2} θ} \\ \frac{d ψ}{d t} = a \frac{I_{1}}{I_{3}} - \frac{b - a \cos θ}{\sin^{2} θ} \cos θ \end{cases} \end{array}$

con la condiciones iniciales

$t = 0, {\begin{cases} θ = θ_{0} \\ \frac{d θ}{d t} = 0 \\ φ = 0 \\ ψ = 0 \end{cases}$

Representaremos la trayectoria del extremo de la peonza z'=1; sobre una esfera de radio unidad. Para ello, es necesario hacer la transformación

$\begin{array}{l} (\begin{matrix} x \\ y \\ z \end{matrix}) = (\begin{matrix} \cos φ & - \sin φ & 0 \\ \sin φ & \cos φ & 0 \\ 0 & 0 & 1 \end{matrix}) (\begin{matrix} 1 & 0 & 0 \\ 0 & \cos θ & - \sin θ \\ 0 & \sin θ & \cos θ \end{matrix}) (\begin{matrix} \cos ψ & - \sin ψ & 0 \\ \sin ψ & \cos ψ & 0 \\ 0 & 0 & 1 \end{matrix}) (\begin{matrix} 0 \\ 0 \\ 1 \end{matrix}) \\ (\begin{matrix} x \\ y \\ z \end{matrix}) = (\begin{matrix} \cos φ & - \sin φ & 0 \\ \sin φ & \cos φ & 0 \\ 0 & 0 & 1 \end{matrix}) (\begin{matrix} 1 & 0 & 0 \\ 0 & \cos θ & - \sin θ \\ 0 & \sin θ & \cos θ \end{matrix}) (\begin{matrix} 0 \\ 0 \\ 1 \end{matrix}) \\ (\begin{matrix} x \\ y \\ z \end{matrix}) = (\begin{matrix} \cos φ & - \sin φ & 0 \\ \sin φ & \cos φ & 0 \\ 0 & 0 & 1 \end{matrix}) (\begin{matrix} 0 \\ - \sin θ \\ \cos θ \end{matrix}) \\ {\begin{cases} x = \sin φ \sin θ \\ y = - \cos φ \sin θ \\ z = \cos θ \end{cases} \end{array}$

Comprobaremos que la energía E se mentiene constante o bien, la cantidad e

$\begin{array}{l} E = \frac{1}{2} I_{1} {(\frac{d θ}{d t})}^{2} + \frac{1}{2} \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1} \sin^{2} θ} + \frac{p_{ψ}^{2}}{2 I_{3}} + m g h \cos θ \\ \frac{E}{I_{1}} - \frac{p_{ψ}^{2}}{2 I_{1} I_{3}} = \frac{1}{2} {(\frac{d θ}{d t})}^{2} + \frac{1}{2} \frac{{(p_{φ} - p_{ψ} \cos θ)}^{2}}{I_{1}^{2} \sin^{2} θ} + \frac{m g h}{I_{1}} \cos θ \\ e = 2 (\frac{E}{I_{1}} - \frac{a^{2} I_{1}}{2 I_{3}}) = {(\frac{d θ}{d t})}^{2} + \frac{{(b - a \cos θ)}^{2}}{\sin^{2} θ} + B \cos θ \end{array}$

Ejemplos

Ejemplo 1

Parámetros

A=1.6
B=2
a=2.5
b=1.7
I₁=2
I₃=1

Como hemos calculado anteriormente, las inclinaciones máxima y mínima del eje son θ₁=40.8° y θ₂=72.7°. El ángulo θ_φ=47.1° para el cual dφ/dt=0

Representamos la trayectoria del extremo del eje de la peonza para varios ángulos iniciales θ₀

Angulo inicial, θ₀=42°

B=2.0;
a=2.5;
b=1.7;
I1=2;
I3=1;
hold on
[X,Y,Z] = sphere;
h1=mesh(X,Y,Z);
set(h1,'EdgeColor',[0.6,0.6,0.6],'EdgeAlpha',0.5,'FaceAlpha',0.5)

% x(1) es beta, x(2) dbeta/dt, x(3) es alfa, x(4) gamma
fg=@(t,x)[x(2); (b^2+a^2)*cos(x(1))/sin(x(1))^3-a*b*(3+cos(2*x(1)))
/(2*sin(x(1))^3)+B*sin(x(1))/2; (b-a*cos(x(1)))/sin(x(1))^2; a*I1/I3
-(b-a*cos(x(1)))+cos(x(1))/sin(x(1))^2]; 
[~,x]=ode45(fg,[0,70],[42*pi/180, 0,0,0]); %modificar el ángulo inicial
%trayectoria
h1=line(sin(x(:,3)).*sin(x(:,1)),-cos(x(:,3)).*sin(x(:,1)),cos(x(:,1)));
set(h1,'Color',[.7,0,0],'LineWidth',1.5)
hold off

axis equal
grid on
xlabel('x')
ylabel('y')
zlabel('z')
title('Peonza')
view(12,16)
%energía
e=x(:,2).^2+((b-a*cos(x(:,1)))./sin(x(:,1))).^2+B*cos(x(:,1));

Comprobamos que la energía se mentiene constante

>> e(1:10)
ans =
    1.4926
    ...
    1.4926

>> e(200:220)
ans =
    1.4987
    1.4986
    ...
    1.5000
    ...
    1.4999

Angulo inicial, θ₀=45°

Angulo inicial, θ₀=47.16°

Angulo inicial, θ₀=50°

Angulo inicial, θ₀=55°

Angulo inicial, θ₀=65°

Angulo inicial, θ₀=70°

Ejemplo 2

Parámetros

A=1.6
B=2
a=2.5
b=2
I₁=2
I₃=1

Las inclinaciones máxima y mínima del eje son θ₁=36.9° y θ₂=53.7°. El ángulo θ_φ=36.9° para el cual dφ/dt=0

Ejemplo 3

Parámetros

A=2.8
B=2
a=1.75
b=2
I₁=2
I₃=1

Las inclinaciones máxima y mínima del eje son θ₁=24.0° y θ₂=65.9°. No hay ángulo θ_φ para el cual dφ/dt=0

Otros ejemplos, en el ártículo citado en segundo lugar, en las referencias

Referencias

Walter Greiner. Classical Mechanics, Systems of Particles and Hamiltonian Dynamics. Second Edition. Springer 2010 pp. 185-200, 238-241

Masatsugu Sei Suzuki. Chapter 13S Symmetric top . 2010