Sea un robot de dos grados de libertad de accionamiento directo con articulaciones RR como el de la figura. Supóngase que las piezas tienen una densidad constante.
Las piezas A y B son iguales en longitud (\(L_1=L_2=2L\)) con masa \({^AM} = 1.0\text{ Kg}\), tienen la misma densidad y momento de inercia que la B. Para dicho robot, y en función de los parámetros del mismo, se pide:
Calculo_NE.m como ayuda.robot_RR.mdl. Se pide programar la función ModeloDinamico_R.m y crear un simulador del robot propuesto utilizando las ecuaciones obtenidas del apartado anterior.graficas, siempre y cuando se respeten los nombres de las variables del espacio de trabajo de Matlab.
El manipulador es un robot espacial de 2 GDL con cinemática RR (Rotacional-Rotacional) de accionamiento directo (direct-drive, relación de transmisión \(r = 1\)):
| Elemento | Parámetro Físico | Símbolo | Valor Numérico | Unidades |
|---|---|---|---|---|
| Eslabones A y B | Longitud total (\(2L\)) | \(L_1, L_2\) | \(1.0\) | \(\text{m}\) |
| Eslabones A y B | Masa total | \(m_1, m_2\) | \(1.0\) | \(\text{kg}\) |
| Eslabones A y B | Distancia al centro de masas (\(L\)) | \(l_{c1}, l_{c2}\) | \(0.5\) | \(\text{m}\) |
| Pieza B (y A) | Inercia transversal baricéntrica | \(I_{\perp} = {^B I_{11}} = {^B I_{33}}\) | \(0.5\) | \(\text{kg}\cdot\text{m}^2\) |
| Pieza B (y A) | Inercia longitudinal baricéntrica | \(I_{\parallel} = {^B I_{22}}\) | \(0.00021\) | \(\text{kg}\cdot\text{m}^2\) |
| Motores Direct-Drive | Inercia de rotor | \(J_{m1}, J_{m2}\) | \(0.025\) | \(\text{kg}\cdot\text{m}^2\) |
| Motores Direct-Drive | Fricción viscosa de rotor | \(B_{m1}, B_{m2}\) | \(2.0 \times 10^{-5}\) | \(\text{N}\cdot\text{m}/(\text{rad/s})\) |
| Motores Direct-Drive | Constante de par electromagnético | \(K_{t1}, K_{t2}\) | \(10.0\) | \(\text{N}\cdot\text{m/A}\) |
1Eslabón 1 (Pieza A):
Al estar en el plano horizontal a altura \(z = 0\), su centro de masas se ubica en \(\mathbf{r}_{c1} = [l_{c1}\cos q_1, l_{c1}\sin q_1, 0]^T\). Su velocidad lineal es \(v_{c1}^2 = l_{c1}^2 \dot{q}_1^2\). Su velocidad angular es \(\boldsymbol{\omega}_1 = [0, 0, \dot{q}_1]^T\). El eje vertical es transversal a la barra, por lo que su momento de inercia respecto al eje de giro es \(I_{\perp} = 0.5\):
\[ K_1 = \frac{1}{2} m_1 l_{c1}^2 \dot{q}_1^2 + \frac{1}{2} I_{\perp} \dot{q}_1^2 + \frac{1}{2} J_{m1} \dot{q}_1^2 = \frac{1}{2} \left( m_1 l_{c1}^2 + I_{\perp} + J_{m1} \right) \dot{q}_1^2 = \frac{1}{2} (0.775) \dot{q}_1^2 \]2Eslabón 2 (Pieza B):
La posición de su centro de masas en el espacio tridimensional es:
\[ \mathbf{r}_{c2} = \begin{bmatrix} (L_1 + l_{c2}\cos q_2)\cos q_1 \\ (L_1 + l_{c2}\cos q_2)\sin q_1 \\ l_{c2}\sin q_2 \end{bmatrix} \implies v_{c2}^2 = (L_1 + l_{c2}\cos q_2)^2 \dot{q}_1^2 + l_{c2}^2 \dot{q}_2^2 \]La velocidad angular total de la barra 2 es \(\boldsymbol{\omega}_2 = \dot{q}_1 \hat{\mathbf{k}} + \dot{q}_2 \hat{\mathbf{e}}_2\). Proyectando sobre sus ejes baricéntricos principales:
Como \(I_{11} = I_{33} = I_{\perp}\), la energía cinética rotacional de la barra 2 es:
\[ K_{2,\text{rot}} = \frac{1}{2} I_{\parallel} (\dot{q}_1 \sin q_2)^2 + \frac{1}{2} I_{\perp} (\dot{q}_1^2 \cos^2 q_2 + \dot{q}_2^2) \]Añadiendo la inercia del rotor del motor 2 montado en la articulación 2: \(\frac{1}{2} J_{m2} \dot{q}_2^2\).
Sumando todas las componentes, la energía cinética total no contiene términos cruzados \(\dot{q}_1 \dot{q}_2\), resultando una matriz de inercia diagonal:
\[ M(q) = \begin{bmatrix} M_{11}(q_2) & 0 \\ 0 & M_{22} \end{bmatrix} \]Donde:
\[ \begin{aligned} M_{11}(q_2) &= m_1 l_{c1}^2 + I_{\perp} + J_{m1} + m_2 (L_1 + l_{c2}\cos q_2)^2 + I_{\perp}\cos^2 q_2 + I_{\parallel}\sin^2 q_2 \\ &= 0.775 + (1.0 + \cos q_2 + 0.25\cos^2 q_2) + 0.5\cos^2 q_2 + 0.00021\sin^2 q_2 \\ &= 1.775 + \cos q_2 + 0.75\cos^2 q_2 + 0.00021\sin^2 q_2 \quad [\text{kg}\cdot\text{m}^2] \\[6pt] M_{22} &= m_2 l_{c2}^2 + I_{\perp} + J_{m2} = 1.0(0.25) + 0.5 + 0.025 = 0.775 \quad [\text{kg}\cdot\text{m}^2] \end{aligned} \]Derivando \(M_{11}(q_2)\) respecto a \(q_2\):
\[ \frac{\partial M_{11}}{\partial q_2} = -\sin q_2 - 1.5\cos q_2 \sin q_2 + 0.00042\sin q_2 \cos q_2 \approx -(1 + 1.49958\cos q_2)\sin q_2 \]Aplicando las ecuaciones de Euler-Lagrange:
\[ c_1(q_2, \dot{q}_1, \dot{q}_2) = \frac{\partial M_{11}}{\partial q_2}\dot{q}_1 \dot{q}_2, \qquad c_2(q_2, \dot{q}_1) = -\frac{1}{2}\frac{\partial M_{11}}{\partial q_2}\dot{q}_1^2 \]La gravedad actúa verticalmente (\(-\hat{\mathbf{k}}\)). El eslabón 1 se mueve horizontalmente (\(g_1 = 0\)). El eslabón 2 se eleva con altura de CoM \(z = l_{c2}\sin q_2\):
\[ U = m_2 g l_{c2}\sin q_2 \implies g_2(q_2) = \frac{\partial U}{\partial q_2} = m_2 g l_{c2}\cos q_2 = (1.0)(9.81)(0.5)\cos q_2 = 4.905\cos q_2 \quad [\text{N}\cdot\text{m}] \]En accionamiento directo (\(r=1\)), los pares entregados a las articulaciones son proporcionales a las corrientes eléctricas de control:
\[ \tau_1 = K_{t1} i_1 = 10 i_1, \qquad \tau_2 = K_{t2} i_2 = 10 i_2 \]Incorporando la fricción viscosa \(B_{m1}\dot{q}_1\) y \(B_{m2}\dot{q}_2\), el modelo dinámico directo resultante es:
Para su inclusión en el bloque de MATLAB Function del archivo robot_RR.mdl, se ha creado el archivo oficial Matlab/03_Dinamica_y_Modelado/ModeloDinamico_R.m:
Dado que la matriz de inercia es diagonal (\(M_{12} = M_{21} = 0\)), el sistema dinámico se encuentra naturalmente desacoplado a bajas velocidades. Para cada articulación \(k\), la planta linealizada se modela como un doble integrador accionado en corriente:
\[ J_{k,\text{nom}} \ddot{q}_k = K_{tk} i_k \implies \frac{q_k(s)}{i_k(s)} = \frac{K_{tk}}{J_{k,\text{nom}} s^2} \]Se implementa una ley de control PID para cada eje en términos de la corriente del motor:
\[ i_k(t) = K_{pk} e_k(t) + K_{ik} \int_0^t e_k(\tau) d\tau - K_{dk} \dot{q}_k(t) \quad [+ \text{feedforward de gravedad para eje 2}] \]La ecuación característica de lazo cerrado deseada para un comportamiento críticamente amortiguado (\(\xi = 1.0\)) con frecuencia natural \(\omega_n\) y polo rápido auxiliar \(p_o = 4\omega_n\) es:
\[ s^3 + \frac{K_{tk} K_{dk}}{J_k} s^2 + \frac{K_{tk} K_{pk}}{J_k} s + \frac{K_{tk} K_{ik}}{J_k} = (s^2 + 2\xi\omega_n s + \omega_n^2)(s + p_o) = s^3 + 6\omega_n s^2 + 9\omega_n^2 s + 4\omega_n^3 \]Igualando coeficientes, las ganancias de control en corriente son:
\[ K_{pk} = \frac{J_{k,\text{nom}}}{K_{tk}} (9\omega_n^2), \qquad K_{dk} = \frac{J_{k,\text{nom}}}{K_{tk}} (6\omega_n), \qquad K_{ik} = \frac{J_{k,\text{nom}}}{K_{tk}} (4\omega_n^3) \]| Eje Articular | Inercia Nominal \(J_{k,\text{nom}}\) | Frecuencia \(\omega_{nk}\) | \(K_{pk}\) (A/rad) | \(K_{dk}\) (A\(\cdot\)s/rad) | \(K_{ik}\) (A/(rad\(\cdot\)s)) | Compensación Gravitatoria |
|---|---|---|---|---|---|---|
| Eje 1 (Giro Base) | \(2.525\text{ kg}\cdot\text{m}^2\) | \(3.0\text{ rad/s}\) | \(20.4525\) | \(4.5450\) | \(27.2700\) | Ninguna (\(g_1 = 0\)) |
| Eje 2 (Elevación) | \(0.775\text{ kg}\cdot\text{m}^2\) | \(4.0\text{ rad/s}\) | \(11.1600\) | \(1.8600\) | \(19.8400\) | \(i_{2,\text{grav}} = \frac{4.905\cos q_2}{10.0} = 0.4905\cos q_2\) A |
Simulando la respuesta ante escalones de referencia \([q_{1,\text{ref}}, q_{2,\text{ref}}] = [45^\circ, 30^\circ]\) partiendo del reposo: