Sistemas de ecuaciones diferenciales

Sistemas de ecuaciones diferenciales acopladas, diagonalización y resolución de estados físicos múltiples.

Qué es

Un sistema de ecuaciones diferenciales ordinarias (EDO) de primer orden es un conjunto acoplado de $n$ ecuaciones dinámicas para $n$ funciones incógnitas $\mathbf{x}(t) = \left(x_1(t), x_2(t), \dots, x_n(t)\right)^T$: $\frac{d\mathbf{x}}{dt} = \mathbf{f}(t, \mathbf{x})$ Cualquier EDO escalar de orden $n$, $y^{(n)} = F(t, y, y', \dots, y^{(n-1)})$, se reduce de forma canónica a un sistema de primer orden de dimensión $n$ definiendo el vector de estado $x_1 = y, x_2 = y', \dots, x_n = y^{(n-1)}$.

Para el caso fundamental de sistemas lineales homogéneos con coeficientes constantes: $\frac{d\mathbf{x}}{dt} = \mathbf{A} \mathbf{x}, \quad \mathbf{A} \in \mathbb{R}^{n \times n}$ La solución general se construye mediante los autovalores $\lambda_i$ y autovectores $\mathbf{v}_i$ de la matriz del sistema $\mathbf{A}$, satisfaciendo $\mathbf{A} \mathbf{v}_i = \lambda_i \mathbf{v}_i$:

• Si $\mathbf{A}$ posee $n$ autovectores linealmente independientes: $\mathbf{x}(t) = \sum_{i=1}^n c_i e^{\lambda_i t} \mathbf{v}_i$
• Formalmente, la solución única sujeta a la condición inicial $\mathbf{x}(0) = \mathbf{x}_0$ se expresa de manera compacta mediante la matriz exponencial: $\mathbf{x}(t) = e^{\mathbf{A} t} \mathbf{x}_0, \quad e^{\mathbf{A} t} = \sum_{k=0}^\infty \frac{(\mathbf{A} t)^k}{k!} = \mathbf{P} e^{\mathbf{D} t} \mathbf{P}^{-1}$ donde $\mathbf{P}$ es la matriz de paso cuyas columnas son los autovectores y $\mathbf{D} = \operatorname{diag}(\lambda_1, \dots, \lambda_n)$.

Por qué aparece en física

En la naturaleza, los elementos rara vez se mueven en aislamiento; interactúan e intercambian energía continuamente:

• Osciladores acoplados: Dos masas unidas por resortes elásticos transfieren vibración entre sí ($m_1 \ddot{x}_1 = -k_1 x_1 + k_{12}(x_2 - x_1)$). La diagonalización del sistema desacopla el movimiento en modos normales independientes, cada uno vibrando a una frecuencia propia pura.
• Redes eléctricas multipolares: Circuitos con múltiples mallas y nodos regidos por las leyes de Kirchhoff para corrientes de malla y voltajes de nudo.
• Cinética de desintegración nuclear y cadenas de transporte: Cadenas radiactivas de Bateman ($A \to B \to C$), donde la población de un isótopo hijo depende de la tasa de desintegración del padre.
• Mecánica Hamiltoniana: La formulación más profunda de la física clásica reescribe la segunda ley de Newton como un sistema lineal o no lineal de $2n$ ecuaciones de primer orden: $\dot{q}_i = \frac{\partial H}{\partial p_i}, \quad \dot{p}_i = -\frac{\partial H}{\partial q_i}$

Idea visual o geométrica

El vector de estado $\mathbf{x}(t)$ describe una trayectoria continua en el espacio de estados o espacio de fases $\mathbb{R}^n$. En cada punto $\mathbf{x}$, la matriz $\mathbf{A}$ define un vector de velocidad instantánea $\mathbf{v} = \mathbf{A} \mathbf{x}$. Los autovectores representan direcciones invariantes: a lo largo de la recta generada por un autovector $\mathbf{v}_i$, el campo apunta exactamente en la misma dirección (o en sentido opuesto), de modo que si el sistema parte sobre dicho eje, nunca lo abandona, expandiéndose ($e^{\lambda t}$ con $\lambda > 0$) o colapsando ($e^{\lambda t}$ con $\lambda < 0$) como una línea recta. Cualquier otra trayectoria es una composición de estas direcciones fundamentales deformadas por el flujo lineal.

Lenguaje simbólico mínimo

• Matriz fundamental $\mathbf{\Phi}(t)$: Matriz de dimensión $n \times n$ cuyas columnas son $n$ soluciones linealmente independientes $\mathbf{x}^{(1)}(t), \dots, \mathbf{x}^{(n)}(t)$: $\frac{d\mathbf{\Phi}}{dt} = \mathbf{A} \mathbf{\Phi}(t), \quad W(t) = \det\mathbf{\Phi}(t) \neq 0$
• Ecuación de Liouville-Jacobi para el Wronskiano matricial: $W(t) = W(0) \exp(\operatorname{tr}(\mathbf{A}) t)$ donde $\operatorname{tr}(\mathbf{A}) = \sum \lambda_i$ mide la tasa de divergencia del volumen en el espacio de fases (teorema de Liouville).
• Sistemas inhomogéneos y fórmula de variación de constantes de Duhamel: $\mathbf{x}' = \mathbf{A} \mathbf{x} + \mathbf{g}(t) \implies \mathbf{x}(t) = e^{\mathbf{A} t} \mathbf{x}_0 + \int_0^t e^{\mathbf{A}(t - s)} \mathbf{g}(s)\,ds$

Qué significa el resultado

1. Desacoplamiento modal: Pasar de la base física $\mathbf{x}$ a la base modal $\mathbf{y} = \mathbf{P}^{-1} \mathbf{x}$ convierte un sistema fuertemente acoplado en $n$ ecuaciones escalares independientes $\dot{y}_i = \lambda_i y_i$, resolubles de forma aislada.
2. Conservación de volumen y divergencia: Si $\operatorname{tr}(\mathbf{A}) = 0$ (como en sistemas hamiltonianos conservativos sin fricción), el flujo en el espacio de fases preserva rígidamente el volumen, sin atractores ni colapsos. Si $\operatorname{tr}(\mathbf{A}) < 0$, el sistema es disipativo y los volúmenes colapsan hacia atractores.
3. Autovalores complejos: Un par $\lambda = \alpha \pm i\omega$ introduce rotación periódica con pulsación $\omega$ modulada por una envolvente exponencial $e^{\alpha t}$.

Cuándo no debes usarlo

• Matrices defectivas (no diagonalizables): Si un autovalor tiene multiplicidad algebraica mayor que su multiplicidad geométrica, no hay suficientes autovectores; es obligatorio emplear cadenas de autovectores generalizados y la forma canónica de Jordan.
• Sistemas no lineales lejos del equilibrio: La aproximación $\mathbf{x}' = \mathbf{A} \mathbf{x}$ solo es válida localmente en el entorno de un punto de equilibrio; fuera de esa vecindad, las no linealidades inducen fenómenos globales (bifurcaciones, caos, atractores extraños).
• Matrices dependientes del tiempo $\mathbf{A}(t)$ no conmutativas: Si $\mathbf{A}(t) \mathbf{A}(s) \neq \mathbf{A}(s) \mathbf{A}(t)$, la solución NO es $e^{\int \mathbf{A}(t) dt}$, requiriendo series de Dyson (o producto de Magnus).

Resumen operativo

1. Escribe el sistema en forma matricial canónica $\mathbf{x}' = \mathbf{A} \mathbf{x}$.
2. Calcula los autovalores resolviendo el polinomio característico $\det(\mathbf{A} - \lambda \mathbf{I}) = 0$.
3. Para cada autovalor $\lambda_i$, halla el autovector asociado resolviendo $(\mathbf{A} - \lambda_i \mathbf{I}) \mathbf{v}_i = \mathbf{0}$.
4. Si las raíces son complejas $\lambda = \alpha \pm i\omega$, toma la parte real e imaginaria de $e^{(\alpha + i\omega) t} \mathbf{v}$ para obtener soluciones reales ortogonales.
5. Construye la matriz fundamental $\mathbf{\Phi}(t)$ y aplica las condiciones iniciales $\mathbf{x}_0$ resolviendo $\mathbf{\Phi}(0) \mathbf{c} = \mathbf{x}_0$.

Procedimiento: Resolución de Sistemas Lineales de EDOs

Antes de calcular

1. Escribe el vector de estado: Define con claridad el vector $\mathbf{x}(t) = (x_1, x_2, \dots)^T$ y sus unidades.
2. Construye la matriz $\mathbf{A}$: Asegúrate de que las ecuaciones estén ordenadas de forma consistente con las filas correspondientes.
3. Calcula la traza y el determinante: $\operatorname{tr}(\mathbf{A}) = a_{11} + a_{22}$ y $\det(\mathbf{A}) = a_{11} a_{22} - a_{12} a_{21}$. Recuerda que el polinomio característico en $2 \times 2$ es siempre $\lambda^2 - \operatorname{tr}(\mathbf{A}) \lambda + \det(\mathbf{A}) = 0$.

Paso 1: identifica el objeto matemático

Escribe la matriz del sistema: $\mathbf{A} = \begin{pmatrix} a_{11} & a_{12} \\ a_{21} & a_{22} \end{pmatrix}$ Verifica si el sistema es homogéneo (sin vector término fuente $\mathbf{g}(t)$).

Paso 2: elige la operación

1. Resuelve la ecuación característica: $\det(\mathbf{A} - \lambda \mathbf{I}) = 0 \implies \lambda^2 - \operatorname{tr}(\mathbf{A}) \lambda + \det(\mathbf{A}) = 0$
2. Determina los autovectores para cada autovalor: $(\mathbf{A} - \lambda_k \mathbf{I}) \begin{pmatrix} v_{1} \\ v_{2} \end{pmatrix} = \begin{pmatrix} 0 \\ 0 \end{pmatrix}$
3. Tratamiento según el tipo de autovalores:
   • Reales distintos ($\lambda_1 \neq \lambda_2$): $\mathbf{x}(t) = c_1 e^{\lambda_1 t} \mathbf{v}_1 + c_2 e^{\lambda_2 t} \mathbf{v}_2$
   • Complejos conjugados ($\lambda = \alpha \pm i\beta$, autovector $\mathbf{v} = \mathbf{u} + i\mathbf{w}$): Una solución compleja es $\mathbf{z}(t) = e^{(\alpha + i\beta) t} (\mathbf{u} + i\mathbf{w})$. Separando en parte real e imaginaria mediante la fórmula de Euler: $\mathbf{x}_1(t) = \operatorname{Re}[\mathbf{z}(t)] = e^{\alpha t}\left(\mathbf{u} \cos(\beta t) - \mathbf{w} \sin(\beta t)\right)$ $\mathbf{x}_2(t) = \operatorname{Im}[\mathbf{z}(t)] = e^{\alpha t}\left(\mathbf{u} \sin(\beta t) + \mathbf{w} \cos(\beta t)\right)$ $\mathbf{x}(t) = c_1 \mathbf{x}_1(t) + c_2 \mathbf{x}_2(t)$
   • Raíz real doble repetida defectiva: Si solo hay un autovector independiente $\mathbf{v}$, halla el autovector generalizado $\mathbf{w}$ resolviendo $(\mathbf{A} - \lambda \mathbf{I}) \mathbf{w} = \mathbf{v}$. $\mathbf{x}(t) = c_1 e^{\lambda t} \mathbf{v} + c_2 e^{\lambda t} (t \mathbf{v} + \mathbf{w})$

Paso 3: conserva unidades y signos

• Los autovalores $\lambda_i$ tienen unidades de frecuencia o tasa inversa (ej. $\text{s}^{-1}$).
• Los autovectores $\mathbf{v}_i$ fijan la relación de amplitudes (adimensional o ratio de unidades físicas entre variables).

Paso 4: calcula

1. Agrupa las soluciones en la matriz fundamental $\mathbf{\Phi}(t) = [\mathbf{x}_1(t) \quad \mathbf{x}_2(t)]$.
2. Evalúa en $t = 0$: $\mathbf{\Phi}(0) \mathbf{c} = \mathbf{x}_0$.
3. Resuelve el sistema lineal algebraico para el vector de constantes $\mathbf{c} = \mathbf{\Phi}(0)^{-1} \mathbf{x}_0$.

Paso 5: interpreta el resultado

1. Si los autovalores tienen parte real negativa ($\operatorname{Re}(\lambda_i) < 0$), el origen es asintóticamente estable (atractor).
2. Si los autovalores son imaginarios puros ($\alpha = 0$), las trayectorias forman elipses concéntricas cerradas (centro conservativo, órbitas periódicas perpetuas).
3. Si los autovalores tienen signos opuestos ($\lambda_1 < 0 < \lambda_2$), el origen es un punto de ensilladura (inestable).

Checklist de autocontrol

• ¿Coincide la traza con la suma de autovalores $\lambda_1 + \lambda_2$?
• ¿Coincide el determinante con el producto de autovalores $\lambda_1 \lambda_2$?
• ¿He comprobado que $(\mathbf{A} - \lambda \mathbf{I}) \mathbf{v} = \mathbf{0}$ para cada autovector?
• ¿He separado correctamente las partes real e imaginaria sin introducir unidades imaginarias $i$ en la solución final física?
• ¿Al evaluar $\mathbf{x}(0)$ recupero exactamente el vector de condiciones iniciales $\mathbf{x}_0$?

Ejemplos resueltos: Sistemas de Ecuaciones Diferenciales

Ejemplo guiado paso a paso

Problema

Dos tanques de mezcla interconectados contienen salmuera. El tanque 1 contiene $V_1 = 100\,\text{L}$ y el tanque 2 contiene $V_2 = 100\,\text{L}$. Se bombea salmuera del tanque 1 al tanque 2 a razón de $r = 4\,\text{L}/\text{min}$, y del tanque 2 al tanque 1 a la misma tasa de $r = 4\,\text{L}/\text{min}$, manteniendo los volúmenes constantes. En $t = 0$, el tanque 1 contiene $x_1(0) = 20\,\text{kg}$ de sal disuelta y el tanque 2 contiene agua pura $x_2(0) = 0\,\text{kg}$.

1. Plantea el sistema matricial de EDOs para las masas de sal $x_1(t)$ y $x_2(t)$.
2. Calcula los autovalores y autovectores de la matriz del sistema.
3. Encuentra la solución analítica para $x_1(t)$ y $x_2(t)$.
4. Determina el estado límite de equilibrio cuando $t \to \infty$.

Solución detallada

Paso 1: Planteamiento del sistema matricial La tasa de variación de sal en cada tanque viene dada por el balance neto (caudal entrante $\times$ concentración entrante $-$ caudal saliente $\times$ concentración saliente): $\frac{dx_1}{dt} = -\frac{r}{V_1} x_1 + \frac{r}{V_2} x_2 = -\frac{4}{100} x_1 + \frac{4}{100} x_2 = -0.04 x_1 + 0.04 x_2$ $\frac{dx_2}{dt} = \frac{r}{V_1} x_1 - \frac{r}{V_2} x_2 = 0.04 x_1 - 0.04 x_2$ En forma matricial: $\frac{d}{dt} \begin{pmatrix} x_1 \\ x_2 \end{pmatrix} = \begin{pmatrix} -0.04 & 0.04 \\ 0.04 & -0.04 \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \end{pmatrix}$ La matriz es $\mathbf{A} = \begin{pmatrix} -0.04 & 0.04 \\ 0.04 & -0.04 \end{pmatrix}$.

Paso 2: Cálculo de autovalores y autovectores Polinomio característico: $\det(\mathbf{A} - \lambda \mathbf{I}) = \begin{vmatrix} -0.04 - \lambda & 0.04 \\ 0.04 & -0.04 - \lambda \end{vmatrix} = (-0.04 - \lambda)^2 - 0.04^2 = 0$ $\lambda^2 + 0.08 \lambda = 0 \implies \lambda(\lambda + 0.08) = 0$ Los autovalores son: $\lambda_1 = 0, \quad \lambda_2 = -0.08\,\text{min}^{-1}$

• Para $\lambda_1 = 0$: $(\mathbf{A} - 0 \mathbf{I}) \mathbf{v}_1 = \begin{pmatrix} -0.04 & 0.04 \\ 0.04 & -0.04 \end{pmatrix} \begin{pmatrix} v_{11} \\ v_{12} \end{pmatrix} = \begin{pmatrix} 0 \\ 0 \end{pmatrix} \implies -0.04 v_{11} + 0.04 v_{12} = 0 \implies v_{11} = v_{12}$ Elegimos el autovector: $\mathbf{v}_1 = \begin{pmatrix} 1 \\ 1 \end{pmatrix}$.
• Para $\lambda_2 = -0.08$: $(\mathbf{A} - (-0.08) \mathbf{I}) \mathbf{v}_2 = \begin{pmatrix} 0.04 & 0.04 \\ 0.04 & 0.04 \end{pmatrix} \begin{pmatrix} v_{21} \\ v_{22} \end{pmatrix} = \begin{pmatrix} 0 \\ 0 \end{pmatrix} \implies v_{21} + v_{22} = 0 \implies v_{21} = -v_{22}$ Elegimos el autovector: $\mathbf{v}_2 = \begin{pmatrix} 1 \\ -1 \end{pmatrix}$.

Paso 3: Solución general y aplicación de condiciones iniciales La solución general es: $\mathbf{x}(t) = c_1 e^{0 t} \begin{pmatrix} 1 \\ 1 \end{pmatrix} + c_2 e^{-0.08 t} \begin{pmatrix} 1 \\ -1 \end{pmatrix} = \begin{pmatrix} c_1 + c_2 e^{-0.08 t} \\ c_1 - c_2 e^{-0.08 t} \end{pmatrix}$ Aplicamos la condición inicial en $t = 0$: $\mathbf{x}(0) = \begin{pmatrix} 20 \\ 0 \end{pmatrix}$: $\begin{cases} c_1 + c_2 = 20 \\ c_1 - c_2 = 0 \end{cases} \implies 2 c_1 = 20 \implies c_1 = 10, \quad c_2 = 10$ Por tanto, la solución para cada tanque es: $x_1(t) = 10 + 10 e^{-0.08 t}\,\text{kg}$ $x_2(t) = 10 - 10 e^{-0.08 t}\,\text{kg}$

Paso 4: Estado asintótico y balance de masa

• Masa total en todo instante: $x_1(t) + x_2(t) = 20\,\text{kg}$ (conservación estricta de soluto).
• Cuando $t \to \infty$, el modo transitorio $e^{-0.08 t} \to 0$: $x_1(\infty) = 10\,\text{kg}, \quad x_2(\infty) = 10\,\text{kg}$ La concentración en ambos tanques se iguala a $0.1\,\text{kg}/\text{L}$, alcanzando el equilibrio térmico-químico completo.

Ejemplo inverso

Problema

En un sistema biológico bidimensional de interacción neuronal, las señales observadas en dos electrodos evolucionan según: $x_1(t) = 2 e^{-t} \cos(3 t) - e^{-t} \sin(3 t)$ $x_2(t) = e^{-t} \cos(3 t) + 2 e^{-t} \sin(3 t)$ Reconstruye la matriz del sistema $\mathbf{A} \in \mathbb{R}^{2 \times 2}$ que satisface $\frac{d\mathbf{x}}{dt} = \mathbf{A} \mathbf{x}$.

Solución detallada

Paso 1: Cálculo de las derivadas temporales Derivamos $x_1(t)$ y $x_2(t)$ respecto de $t$: $\dot{x}_1(t) = -e^{-t}(2 \cos 3t - \sin 3t) + e^{-t}(-6 \sin 3t - 3 \cos 3t) = e^{-t}\left((-2 - 3)\cos 3t + (1 - 6)\sin 3t\right) = e^{-t}(-5 \cos 3t - 5 \sin 3t)$ $\dot{x}_2(t) = -e^{-t}(\cos 3t + 2 \sin 3t) + e^{-t}(-3 \sin 3t + 6 \cos 3t) = e^{-t}\left((-1 + 6)\cos 3t + (-2 - 3)\sin 3t\right) = e^{-t}(5 \cos 3t - 5 \sin 3t)$

Paso 2: Relación con el vector de estado en $t = 0$ Evaluamos en $t = 0$: $\mathbf{x}(0) = \begin{pmatrix} 2 \\ 1 \end{pmatrix}, \quad \dot{\mathbf{x}}(0) = \begin{pmatrix} -5 \\ 5 \end{pmatrix}$ Como $\dot{\mathbf{x}}(0) = \mathbf{A} \mathbf{x}(0)$: $\mathbf{A} \begin{pmatrix} 2 \\ 1 \end{pmatrix} = \begin{pmatrix} -5 \\ 5 \end{pmatrix}$

Paso 3: Evaluación en $3 t = \frac{\pi}{2} \implies t_1 = \frac{\pi}{6}$ En $t_1$, $\cos(3 t_1) = 0$ y $\sin(3 t_1) = 1$: $\mathbf{x}(t_1) = e^{-\pi/6} \begin{pmatrix} -1 \\ 2 \end{pmatrix}, \quad \dot{\mathbf{x}}(t_1) = e^{-\pi/6} \begin{pmatrix} -5 \\ -5 \end{pmatrix}$ Cancelando el factor no nulo $e^{-\pi/6}$: $\mathbf{A} \begin{pmatrix} -1 \\ 2 \end{pmatrix} = \begin{pmatrix} -5 \\ -5 \end{pmatrix}$

Paso 4: Resolución matricial para $\mathbf{A}$ Agrupamos ambas relaciones vectoriales en forma matricial: $\mathbf{A} \begin{pmatrix} 2 & -1 \\ 1 & 2 \end{pmatrix} = \begin{pmatrix} -5 & -5 \\ 5 & -5 \end{pmatrix}$ La matriz invertida es $\begin{pmatrix} 2 & -1 \\ 1 & 2 \end{pmatrix}^{-1} = \frac{1}{5} \begin{pmatrix} 2 & 1 \\ -1 & 2 \end{pmatrix}$. Multiplicamos por la derecha: $\mathbf{A} = \begin{pmatrix} -5 & -5 \\ 5 & -5 \end{pmatrix} \frac{1}{5} \begin{pmatrix} 2 & 1 \\ -1 & 2 \end{pmatrix} = \begin{pmatrix} -1 & -1 \\ 1 & -1 \end{pmatrix} \begin{pmatrix} 2 & 1 \\ -1 & 2 \end{pmatrix} = \begin{pmatrix} -2 + 1 & -1 - 2 \\ 2 + 1 & 1 - 2 \end{pmatrix} = \begin{pmatrix} -1 & -3 \\ 3 & -1 \end{pmatrix}$ Comprobación de autovalores de $\mathbf{A}$: $\operatorname{tr}(\mathbf{A}) = -2, \quad \det(\mathbf{A}) = 1 - (-9) = 10 \implies \lambda^2 + 2 \lambda + 10 = 0 \implies \lambda = -1 \pm 3 i$ Coincide de forma exacta con la tasa de decaimiento $\alpha = -1$ y pulsación $\beta = 3$.

Errores comunes: Sistemas de Ecuaciones Diferenciales

Error 1: Emparejar erróneamente autovalores con autovectores de otro autovalor

Al formar la solución $\mathbf{x}(t) = c_1 e^{\lambda_1 t} \mathbf{v}_1 + c_2 e^{\lambda_2 t} \mathbf{v}_2$, multiplicar el factor temporal $e^{\lambda_1 t}$ por el autovector $\mathbf{v}_2$.

• Corrección: Cada modo dinámico $e^{\lambda_k t}$ está indisolublemente ligado a la dirección propia $\mathbf{v}_k$ generada por ese mismo autovalor. Mezclarlos produce un vector que no satisface la EDO.

Error 2: Dejar números imaginarios $i$ en la solución física real

Al resolver sistemas con autovalores complejos $\lambda = \alpha \pm i\beta$, escribir la solución general como combinación lineal con coeficientes complejos arbitrarios sin extraer la base ortogonal real.

• Corrección: Para observables físicos (voltajes, posiciones), la base debe descomponerse en parte real y parte imaginaria mediante la fórmula de Euler, obteniendo senos y cosenos reales ponderados por $e^{\alpha t}$.

Error 3: Asumir que toda matriz con raíz doble es diagonalizable

Si una matriz tiene un autovalor repetido $\lambda_1 = \lambda_2$, asumir sin comprobar que existen dos autovectores independientes y escribir $\mathbf{x}(t) = c_1 e^{\lambda t} \mathbf{v}_1 + c_2 e^{\lambda t} \mathbf{v}_2$.

• Consecuencia: Si el rango de $(\mathbf{A} - \lambda \mathbf{I})$ es 1 en una matriz $2 \times 2$, la multiplicidad geométrica es 1 (matriz defectiva). Se requiere un autovector generalizado $(t \mathbf{v} + \mathbf{w})$.

Error 4: Calcular $e^{\mathbf{A} t}$ aplicando la exponencial elemento a elemento

Creer erróneamente que $e^{\mathbf{A} t} = \begin{pmatrix} e^{a_{11} t} & e^{a_{12} t} \\ e^{a_{21} t} & e^{a_{22} t} \end{pmatrix}$.

• Corrección: La exponencial de matrices se define mediante la serie de Taylor $\sum \frac{(\mathbf{A} t)^k}{k!}$ o la diagonalización $\mathbf{P} e^{\mathbf{D} t} \mathbf{P}^{-1}$. La exponencial elemento a elemento carece de sentido físico y matemático.

Error 5: Olvidar el acoplamiento mutuo en las ecuaciones de balance

Al modelar tanques o circuitos, omitir el término que entra al tanque 2 procedente del tanque 1, tratando el sistema como si fueran subsistemas independientes desacoplados.

Regla rápida de control

Para cualquier sistema de EDOs lineales autónomo: calcula la traza de la matriz $\operatorname{tr}(\mathbf{A})$. La suma de los autovalores hallados DEBE coincidir con la traza: $\lambda_1 + \lambda_2 = a_{11} + a_{22}$, y el producto con el determinante: $\lambda_1 \lambda_2 = \det(\mathbf{A})$. Si no coinciden, repite el cálculo algebraico.

Aplicaciones físicas: Sistemas de Ecuaciones Diferenciales

1. Modos normales y espectroscopía vibracional molecular

En moléculas poliatómicas (como el agua $\text{H}_2\text{O}$ o el dióxido de carbono $\text{CO}_2$), los desplazamientos de los núcleos atómicos se rigen por: $\mathbf{M} \ddot{\mathbf{x}} + \mathbf{K} \mathbf{x} = \mathbf{0}$ La diagonalización simultánea de la matriz de masa $\mathbf{M}$ y la matriz de rigidez interatómica $\mathbf{K}$ determina los modos normales de vibración (estiramiento simétrico, asimétrico y flexión). Las frecuencias calculadas coinciden con los picos de absorción en espectroscopía infrarroja (FTIR) y Raman, permitiendo identificar firmas químicas en astronomía y farmacología.

2. Redes de distribución eléctrica y microredes inteligentes

El análisis de estabilidad transitoria en redes de transporte de energía eléctrica de alta tensión acopla generadores síncronos mediante las ecuaciones de oscilación: $M_i \ddot{\delta}_i + D_i \dot{\delta}_i = P_{m,i} - \sum_{j=1}^N V_i V_j Y_{ij} \sin(\delta_i - \delta_j)$ Linealizando este sistema alrededor del punto de operación en régimen permanente, se calculan los autovalores de la red para prevenir fenómenos de desconexión en cascada (blackouts) ante cortocircuitos o variaciones abruptas de carga.

3. Dinámica de poblaciones y competencia ecológica (Lotka-Volterra)

La interacción competitiva o depredador-presa entre dos especies biológicas con densidades $N_1$ y $N_2$ se describe por el sistema acoplado: $\frac{dN_1}{dt} = N_1 (r_1 - a_{11} N_1 - a_{12} N_2)$ $\frac{dN_2}{dt} = N_2 (r_2 - a_{21} N_1 - a_{22} N_2)$ La linealización matricial en torno a los puntos de equilibrio de coexistencia determina las condiciones de exclusión competitiva o coexistencia armónica estable en ecología teórica.

4. Transferencia de calor en intercambiadores de calor multietapa

En intercambiadores de calor industriales de carcasa y tubos a contracorriente, las temperaturas de las corrientes caliente ($T_h$) y fría ($T_c$) satisfacen el sistema de primer orden: $\dot{m}_h c_{p,h} \frac{dT_h}{dx} = -U P (T_h - T_c), \quad -\dot{m}_c c_{p,c} \frac{dT_c}{dx} = U P (T_h - T_c)$ La solución analítica del sistema matricial proporciona la diferencia de temperatura media logarítmica (LMTD) utilizada para dimensionar calderas, condensadores nucleares y plantas desalinizadoras.

5. Dinámica de guiado y control de satélites (Ecuaciones de Clohessy-Wiltshire)

El movimiento relativo entre dos naves espaciales en órbita circular cerrada (maniobras de aproximación y acoplamiento orbital con la Estación Espacial Internacional) se gobierna por las ecuaciones lineales de Clohessy-Wiltshire: $\ddot{x} - 2 \omega \dot{y} - 3 \omega^2 x = F_x / m$ $\ddot{y} + 2 \omega \dot{x} = F_y / m$ $\ddot{z} + \omega^2 z = F_z / m$ Al reescribirlo como un sistema de primer orden de dimensión 6, la matriz exponencial genera la matriz de transición de estado que guía los propulsores automáticos de acoplamiento espacial.

Conexión interdisciplinar

Los sistemas de EDOs constituyen el corazón de los modelos macroeconómicos de equilibrio dinámico estocástico general (DSGE), la transducción de señales en vías de señalización celular y la cinemática de brazos robóticos manipuladores con múltiples articulaciones.

Práctica guiada: Sistemas de Ecuaciones Diferenciales

Ejercicio 1: Modos normales de osciladores acoplados sin fricción

Dos masas idénticas $m = 1\,\text{kg}$ deslizan sobre una guía horizontal sin rozamiento. La masa 1 está unida a una pared fija por un resorte de constante $k = 3\,\text{N}/\text{m}$, ambas masas están conectadas entre sí por un resorte de acoplamiento $k_c = 2\,\text{N}/\text{m}$, y la masa 2 está libre en su extremo opuesto.

1. Plantea las ecuaciones de Newton y exprésalas en forma matricial para el vector de desplazamientos $\mathbf{x} = (x_1, x_2)^T$.
2. Halla las frecuencias propias cuadradas $\omega^2$ y los autovectores de modos normales.
3. Describe físicamente el movimiento en cada uno de los dos modos normales.

Solución

1. Ecuaciones dinámicas: $m \ddot{x}_1 = -k x_1 + k_c(x_2 - x_1) = -3 x_1 + 2(x_2 - x_1) = -5 x_1 + 2 x_2$ $m \ddot{x}_2 = -k_c(x_2 - x_1) = 2 x_1 - 2 x_2$ En forma de sistema de segundo orden $\ddot{\mathbf{x}} = \mathbf{K} \mathbf{x}$: $\begin{pmatrix} \ddot{x}_1 \\ \ddot{x}_2 \end{pmatrix} = \begin{pmatrix} -5 & 2 \\ 2 & -2 \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \end{pmatrix}$
2. Proponiendo $\mathbf{x}(t) = \mathbf{v} e^{i \omega t} \implies \ddot{\mathbf{x}} = -\omega^2 \mathbf{x}$, definimos $\lambda = -\omega^2$: $\det(\mathbf{K} - \lambda \mathbf{I}) = \begin{vmatrix} -5 - \lambda & 2 \\ 2 & -2 - \lambda \end{vmatrix} = (\lambda + 5)(\lambda + 2) - 4 = \lambda^2 + 7 \lambda + 6 = 0$ $(\lambda + 1)(\lambda + 6) = 0 \implies \lambda_1 = -1, \quad \lambda_2 = -6$
   • Frecuencia del modo 1: $\omega_1^2 = 1 \implies \omega_1 = 1\,\text{rad}/\text{s}$. $(\mathbf{K} - (-1)\mathbf{I})\mathbf{v}_1 = \begin{pmatrix} -4 & 2 \\ 2 & -1 \end{pmatrix} \begin{pmatrix} v_{11} \\ v_{12} \end{pmatrix} = \mathbf{0} \implies 2 v_{11} = v_{12} \implies \mathbf{v}_1 = \begin{pmatrix} 1 \\ 2 \end{pmatrix}$
   • Frecuencia del modo 2: $\omega_2^2 = 6 \implies \omega_2 = \sqrt{6} \approx 2.45\,\text{rad}/\text{s}$. $(\mathbf{K} - (-6)\mathbf{I})\mathbf{v}_2 = \begin{pmatrix} 1 & 2 \\ 2 & 4 \end{pmatrix} \begin{pmatrix} v_{21} \\ v_{22} \end{pmatrix} = \mathbf{0} \implies v_{21} = -2 v_{22} \implies \mathbf{v}_2 = \begin{pmatrix} -2 \\ 1 \end{pmatrix}$
3. Interpretación física:
   • Modo fundamental ($\omega_1 = 1$): Ambas masas oscilan en fase hacia la misma dirección, con la masa 2 teniendo el doble de amplitud que la masa 1.
   • Modo antisimétrico ($\omega_2 = \sqrt{6}$): Las masas oscilan en contrafase (sentidos opuestos), comprimiendo fuertemente el muelle de acoplamiento central.

Ejercicio 2: Cadena de desintegración nuclear radiactiva (Ecuaciones de Bateman)

Una sustancia radiactiva madre $A$ se desintegra a una tasa $\lambda_1 = 0.5\,\text{h}^{-1}$ generando un isótopo hijo inestable $B$, el cual a su vez se desintegra a tasa $\lambda_2 = 0.1\,\text{h}^{-1}$ formando un producto estable $C$: $\frac{dN_A}{dt} = -\lambda_1 N_A, \quad \frac{dN_B}{dt} = \lambda_1 N_A - \lambda_2 N_B$ En $t = 0$, $N_A(0) = N_0$ y $N_B(0) = 0$.

1. Escribe la matriz del sistema triangular $\mathbf{A}$.
2. Encuentra la función analítica $N_B(t)$.
3. Determina el instante de actividad máxima del isótopo intermedio $B$.

Solución

1. Sistema matricial: $\frac{d}{dt} \begin{pmatrix} N_A \\ N_B \end{pmatrix} = \begin{pmatrix} -0.5 & 0 \\ 0.5 & -0.1 \end{pmatrix} \begin{pmatrix} N_A \\ N_B \end{pmatrix}$ Como es triangular inferior, los autovalores son los elementos diagonales: $\lambda_1 = -0.5\,\text{h}^{-1}$ y $\lambda_2 = -0.1\,\text{h}^{-1}$.
2. Integración secuencial: $N_A(t) = N_0 e^{-\lambda_1 t} = N_0 e^{-0.5 t}$ Sustituyendo en la ecuación de $N_B$: $\frac{dN_B}{dt} + 0.1 N_B = 0.5 N_0 e^{-0.5 t}$ Factor integrante $\mu(t) = e^{0.1 t}$: $N_B(t) e^{0.1 t} = 0.5 N_0 \int e^{-0.4 t}\,dt = -\frac{0.5}{0.4} N_0 e^{-0.4 t} + C = -1.25 N_0 e^{-0.4 t} + C$ $N_B(t) = -1.25 N_0 e^{-0.5 t} + C e^{-0.1 t}$ Con $N_B(0) = 0 \implies C = 1.25 N_0$: $N_B(t) = 1.25 N_0 \left(e^{-0.1 t} - e^{-0.5 t}\right)$
3. Máximo de $N_B$: $\frac{dN_B}{dt} = 0 \implies -0.1 e^{-0.1 t} + 0.5 e^{-0.5 t} = 0$: $e^{-0.4 t} = \frac{0.1}{0.5} = 0.2 \implies -0.4 t = \ln(0.2) = -\ln 5$ $t_{\text{max}} = \frac{\ln 5}{0.4} = \frac{1.6094}{0.4} \approx 4.024\,\text{horas}$

Ejercicio 3: Exponencial matricial para un sistema 2x2

Dado el sistema $\mathbf{x}' = \mathbf{A} \mathbf{x}$ con $\mathbf{A} = \begin{pmatrix} 0 & 1 \\ -4 & 0 \end{pmatrix}$:

1. Calcula $\mathbf{A}^2, \mathbf{A}^3, \mathbf{A}^4$.
2. Obtén la matriz exponencial $e^{\mathbf{A} t}$ mediante la definición en serie de potencias.
3. Comprueba que coincide con la solución del oscilador armónico simple $\ddot{x} + 4 x = 0$.

Solución

1. Potencias sucesivas de $\mathbf{A}$: $\mathbf{A}^2 = \begin{pmatrix} 0 & 1 \\ -4 & 0 \end{pmatrix} \begin{pmatrix} 0 & 1 \\ -4 & 0 \end{pmatrix} = \begin{pmatrix} -4 & 0 \\ 0 & -4 \end{pmatrix} = -4 \mathbf{I} = -(2)^2 \mathbf{I}$ $\mathbf{A}^3 = \mathbf{A}^2 \mathbf{A} = -4 \mathbf{A}$ $\mathbf{A}^4 = (\mathbf{A}^2)^2 = (-4 \mathbf{I})^2 = 16 \mathbf{I} = (2)^4 \mathbf{I}$
2. Serie de la exponencial: $e^{\mathbf{A} t} = \sum_{k=0}^\infty \frac{(\mathbf{A} t)^k}{k!} = \mathbf{I} \sum_{m=0}^\infty \frac{(-1)^m (2 t)^{2m}}{(2m)!} + \frac{\mathbf{A}}{2} \sum_{m=0}^\infty \frac{(-1)^m (2 t)^{2m+1}}{(2m+1)!}$ Reconociendo las series de Taylor del coseno y del seno: $e^{\mathbf{A} t} = \cos(2 t) \mathbf{I} + \frac{1}{2} \sin(2 t) \mathbf{A} = \begin{pmatrix} \cos(2 t) & \frac{1}{2} \sin(2 t) \\ -2 \sin(2 t) & \cos(2 t) \end{pmatrix}$
3. Con $\mathbf{x}(0) = \begin{pmatrix} x_0 \\ v_0 \end{pmatrix}$: $\mathbf{x}(t) = e^{\mathbf{A} t} \mathbf{x}(0) = \begin{pmatrix} x_0 \cos(2 t) + \frac{v_0}{2} \sin(2 t) \\ -2 x_0 \sin(2 t) + v_0 \cos(2 t) \end{pmatrix}$ La primera componente es $x(t) = x_0 \cos(2 t) + \frac{v_0}{2} \sin(2 t)$, que es exactamente la solución general del oscilador $\ddot{x} + 4 x = 0$.