Ir al contenido
Background Image

Mecanismos en Ingeniería Mecánica con Python

1051 palabras·5 mins·
Métodos Numéricos - Este artículo es parte de una serie.
Parte 9: Este artículo

⚙️ Módulo 9: Aplicación Final - Análisis Cinemático de Mecanismos
#

Este módulo finaliza el curso aplicando todos los métodos numéricos aprendidos para analizar la posición, velocidad y aceleración de un mecanismo de cuatro barras, utilizando el potente enfoque de Coordenadas Cartesianas (Restricciones de Longitud al Cuadrado).


I. Fundamentos Teóricos: El Enfoque Cartesiano
#

En este análisis, las coordenadas \(x\) e \(y\) de los pivotes móviles se tratan como variables desconocidas. Las ecuaciones de restricción (o ecuaciones de cierre) garantizan que la longitud calculada entre dos pivotes sea siempre igual a la longitud fija del eslabón, formulado como: \(L^2 = \Delta x^2 + \Delta y^2\).

A. Nodos, Variables y Restricciones
#

Para el Mecanismo de Cuatro Barras con longitudes \(a_1, a_2, a_3, a_4\):

  • Pivotes Fijos (Bancada): \(A=(0, 0)\) y \(D=(a_1, 0)\).
  • Pivotes Móviles (Incógnitas): \(B=(x_B, y_B)\) y \(C=(x_C, y_C)\).
  • Vector Incógnita: \(\mathbf{x} = (x_B, y_B, x_C, y_C)^T\) (4 incógnitas).

Las 4 ecuaciones de restricción son:

$$\mathbf{F}(\mathbf{x}) = \begin{cases} f_1: (x_C - x_B)^2 + (y_C - y_B)^2 - a_3^2 = 0 & \text{(Restricción } \overline{BC}) \\ f_2: (a_1 - x_C)^2 + (0 - y_C)^2 - a_4^2 = 0 & \text{(Restricción } \overline{CD}) \\ f_3: x_B^2 + y_B^2 - a_2^2 = 0 & \text{(Restricción } \overline{AB}) \\ f_4: y_B - a_2 \sin(\theta_{in}) = 0 & \text{(Restricción de Entrada)} \end{cases}$$

II. Fase I: Cinemática de Posición (Problema No Lineal)
#

La posición del mecanismo se resuelve mediante el Método de Newton-Raphson Generalizado (Módulo 8).

A. La Matriz Jacobiana (\(\mathbf{J}\))
#

El Jacobiano es la matriz \(4 \times 4\) de derivadas parciales \(\partial f_i / \partial x_j\). La simplicidad de las ecuaciones cuadráticas facilita su cálculo analítico.

$$\mathbf{J} = \begin{pmatrix} 2(x_B - x_C) & 2(y_B - y_C) & 2(x_C - x_B) & 2(y_C - y_B) \\ 0 & 0 & 2(x_C - a_1) & 2y_C \\ 2x_B & 2y_B & 0 & 0 \\ 0 & 1 & 0 & 0 \end{pmatrix}$$

B. Solución Iterativa
#

El paso de corrección se obtiene de la solución del sistema lineal: \(\mathbf{J}(\mathbf{x}^k) \boldsymbol{\Delta}\mathbf{x} = -\mathbf{F}(\mathbf{x}^k)\).


III. Fase II: Cinemática de Velocidad (Problema Lineal)
#

Se obtiene al diferenciar \(\mathbf{F}(\mathbf{x}) = \mathbf{0}\) respecto al tiempo.

A. El Sistema de Velocidad
#

$$\mathbf{J} \mathbf{\dot{x}} = \mathbf{b}_{\text{vel}}$$

Donde \(\mathbf{\dot{x}} = (\dot{x}_B, \dot{y}_B, \dot{x}_C, \dot{y}_C)^T\) es el vector de velocidades nodales, y \(\mathbf{J}\) es la misma Matriz Jacobiana de la Fase I.

B. El Vector del Lado Derecho (\(\mathbf{b}_{\text{vel}}\))
#

El vector \(\mathbf{b}_{\text{vel}}\) solo contiene términos dependientes de la velocidad de entrada \(\dot{\theta}_{in}\):

$$\mathbf{b}_{\text{vel}} = \begin{pmatrix} 0 \\ 0 \\ 0 \\ a_2 \cos(\theta_{in}) \dot{\theta}_{in} \end{pmatrix}$$

La solución se obtiene con numpy.linalg.solve (Módulo 1).


IV. Fase III: Cinemática de Aceleración (Segundo Problema Lineal)
#

Se obtiene al diferenciar las ecuaciones de velocidad respecto al tiempo.

A. El Sistema de Aceleración
#

$$\mathbf{J} \mathbf{\ddot{x}} = \mathbf{b}_{\text{acc}}$$

Donde \(\mathbf{\ddot{x}} = (\ddot{x}_B, \ddot{y}_B, \ddot{x}_C, \ddot{y}_C)^T\) es el vector de aceleraciones nodales.

B. El Vector del Lado Derecho (\(\mathbf{b}_{\text{acc}}\))
#

Este vector es el más complejo, ya que incorpora términos de aceleración centrípeta (velocidades al cuadrado):

$$\mathbf{b}_{\text{acc}} = \begin{pmatrix} -2 [ (\dot{x}_C - \dot{x}_B)^2 + (\dot{y}_C - \dot{y}_B)^2 ] \\ -2 [ \dot{x}_C^2 + \dot{y}_C^2 ] \\ -2 [ \dot{x}_B^2 + \dot{y}_B^2 ] \\ -a_2 \sin(\theta_{in}) (\dot{\theta}_{in})^2 \end{pmatrix}$$

(Se asume \(\ddot{\theta}_{in}=0\) para simplificar el análisis).


V. Integración, Barrido y Visualización (Implementación Final)
#

Implementaremos el bucle de \(360^\circ\), resolviendo Posición, Velocidad y Aceleración de forma secuencial en cada paso.

A. Código Base para la Solución Completa
#

import numpy as np
from scipy.optimize import root
import matplotlib.pyplot as plt

# --- 1. DATOS Y CONSTANTES ---
A1, A2, A3, A4 = 8, 2, 7, 6  # Longitudes de eslabones
XA, YA = 0.0, 0.0
XD, YD = A1, 0.0
DTHETA_IN = 10.0      # Velocidad angular de entrada (rad/s)
DDTHETA_IN = 0.0      # Aceleración angular de entrada (rad/s^2)

# 

def fun_posicion_cartesianas(x, theta_in):
    # (El código de fun_posicion_cartesianas y jac_posicion_cartesianas se define aquí
    # como en la respuesta anterior para F y J)
    xB, yB, xC, yC = x
    t_in_rad = np.deg2rad(theta_in)
    f1 = (xC - xB)**2 + (yC - yB)**2 - A3**2
    f2 = (XD - xC)**2 + (YD - yC)**2 - A4**2
    f3 = xB**2 + yB**2 - A2**2
    f4 = yB - A2 * np.sin(t_in_rad)
    return np.array([f1, f2, f3, f4])

def jac_posicion_cartesianas(x, theta_in):
    # (El código de jac_posicion_cartesianas se define aquí para J)
    xB, yB, xC, yC = x
    J = np.zeros((4, 4))
    J[0, 0], J[0, 1], J[0, 2], J[0, 3] = 2*(xB-xC), 2*(yB-yC), 2*(xC-xB), 2*(yC-yB)
    J[1, 2], J[1, 3] = 2*(xC-XD), 2*(yC-YD)
    J[2, 0], J[2, 1] = 2*xB, 2*yB
    J[3, 1] = 1.0
    return J

def b_aceleracion(x, dx, theta_in, dtheta_in):
    """Calcula el vector b_acc de 4x1 para el sistema de aceleración."""
    xB, yB, xC, yC = x
    dxB, dyB, dxC, dyC = dx
    t_in_rad = np.deg2rad(theta_in)
    
    b1 = -2 * ((dxC - dxB)**2 + (dyC - dyB)**2)
    b2 = -2 * (dxC**2 + dyC**2)
    b3 = -2 * (dxB**2 + dyB**2)
    b4 = -A2 * np.sin(t_in_rad) * (dtheta_in**2) # + A2 * cos(t_in_rad) * DDTHETA_IN
    
    return np.array([b1, b2, b3, b4])

# --- 2. BUCLE DE INTEGRACIÓN (360 GRADOS) ---
THETA_IN_RANGE = np.linspace(0, 360, 100)
posiciones_C, velocidades_C, aceleraciones_C = [], [], []

x_guess_k = np.array([1.0, 1.732, 6.0, 1.0]) # Suposición inicial (~60 deg)

for theta_in_k in THETA_IN_RANGE:
    # FASE I: POSICIÓN (Newton-Raphson)
    sol_pos = root(fun_posicion_cartesianas, x_guess_k, args=(theta_in_k,), jac=jac_posicion_cartesianas, method='hybr')
    if not sol_pos.success: break
    x_k = sol_pos.x
    posiciones_C.append(x_k[2:4])
    x_guess_k = x_k 
    
    J_k = jac_posicion_cartesianas(x_k, theta_in_k)
    t_in_rad = np.deg2rad(theta_in_k)
    
    # FASE II: VELOCIDAD (J * dx = b_vel)
    b_vel = np.zeros(4)
    b_vel[3] = A2 * np.cos(t_in_rad) * DTHETA_IN
    try:
        dx_k = np.linalg.solve(J_k, b_vel)
        velocidades_C.append(dx_k[2:4])
    except np.linalg.LinAlgError:
        dx_k = np.array([np.nan] * 4)
        velocidades_C.append(dx_k[2:4]) 
        
    # FASE III: ACELERACIÓN (J * ddx = b_acc)
    if np.any(np.isnan(dx_k)):
         aceleraciones_C.append([np.nan, np.nan])
         continue
         
    b_acc = b_aceleracion(x_k, dx_k, theta_in_k, DTHETA_IN)
    try:
        ddx_k = np.linalg.solve(J_k, b_acc)
        aceleraciones_C.append(ddx_k[2:4])
    except np.linalg.LinAlgError:
        aceleraciones_C.append([np.nan, np.nan])

posiciones_C = np.array(posiciones_C)
velocidades_C = np.array(velocidades_C)
aceleraciones_C = np.array(aceleraciones_C)

# --- 3. VISUALIZACIÓN ---
fig, axes = plt.subplots(1, 3, figsize=(18, 5))

# Trayectoria
axes[0].plot(posiciones_C[:, 0], posiciones_C[:, 1], 'r-')
axes[0].set_title('Trayectoria de Posición C')
axes[0].set_xlabel('\(x\) (m)'); axes[0].set_ylabel('\(y\) (m)'); axes[0].axis('equal')

# Velocidad
velocidad_resultante_C = np.linalg.norm(velocidades_C, axis=1)
axes[1].plot(THETA_IN_RANGE, velocidad_resultante_C, 'g-')
axes[1].set_title('Magnitud de Velocidad \(|\\mathbf{v}_C|\)')
axes[1].set_xlabel('Ángulo de Entrada \(\\theta_{in}\) (grados)'); axes[1].set_ylabel('\(|\\mathbf{v}_C|\) (m/s)')

# Aceleración
aceleracion_resultante_C = np.linalg.norm(aceleraciones_C, axis=1)
axes[2].plot(THETA_IN_RANGE, aceleracion_resultante_C, 'b-')
axes[2].set_title('Magnitud de Aceleración \(|\\mathbf{a}_C|\)')
axes[2].set_xlabel('Ángulo de Entrada \(\\theta_{in}\) (grados)'); axes[2].set_ylabel('\(|\\mathbf{a}_C|\) (m/s\(^2\))')

for ax in axes: ax.grid(True, linestyle='--')
plt.tight_layout()
plt.show()
Métodos Numéricos - Este artículo es parte de una serie.
Parte 9: Este artículo