21.4: IVPs - Sistema de n ODEs Lineales

Consideramos aquí por simplicidad una familia particular de problemas: los osciladores\(n/2\) acoplados. Esta familia de sistemas puede ser descrita por el conjunto de ecuaciones.

\ begin {reunió}
\ frac {d^ {2} u^ {(1)}} {d t^ {2}} =h^ {(1)}\ left (\ frac {d u^ {(j)}} {d t}, u^ {(j)}, 1\ leq j\ leq n/2\ derecha) +f^ {(1)} (t)\
\ frac {d^ {2} u^ {(2)}} {d t^ {2}} =h^ {(2)}\ left (\ frac {d u^ {(j)}} {d t}, u^ {(j)}, 1\ leq j\ leq n/2\ derecha) +f^ {(2)} (t)\
\ vdots\\
\ frac {d^ {2} u^ {(n/2)}} {d t^ {2}} =h^ {(n/2)}}\ izquierda (\ frac {d u^ {(j)}} {d t}, u^ {(3)}, 1\ leq j\ leq n/2\ derecha) +f^ {(n/2)} (t)
\ fin {reunidos}

donde\(h^{(k)}\) se supone que es una función lineal de todos sus argumentos.

Figura 21.21: Comparación de esquemas de integración numérica para un sistema amortiguador-masa-amortiguador sobreamortiguado con\(ω_n = 1.0\) y\(ζ = 50\). Obsérvese que los pasos de tiempo utilizados para los esquemas explícitos son diferentes de los de los esquemas implícitos.

Primero convertimos este sistema de ecuaciones en forma de espacio de estado. Identificamos

\(w_{1}=u^{(1)}, \quad w_{2}=\frac{d u^{(1)}}{d t}\)

\(w_{3}=u^{(2)}, \quad w_{4}=\frac{d u^{(2)}}{d t}\)

\(w_{n-1}=u^{(n / 2)}, \quad w_{n}=\frac{d u^{(n / 2)}}{d t}\)

Entonces podemos escribir nuestro sistema -usando el hecho de que\(h\) es lineal en sus argumentos- como

\[\dfrac{d w}{d t}=A w+F \tag{21.25}\]

\[w(0)=w_{0} \tag{21.25}\]

donde\(h\) determina\(A, F\) es dado por\(\left(\begin{array}{llll}0 & f^{(1)}(t) & 0 & f^{(2)}(t) \cdots 0 f^{(n/2)}(t) \end{array}\right)^T.\)

\(w_{0}=\left(u^{(1)}(0) \frac{d u^{(1)}}{d t}(0) \quad u^{(2)}(0) \frac{d u^{(2)}}{d t}(0) \quad \ldots \quad u^{(n / 2)}(0) \frac{d u^{(n / 2)}}{d t}(0)\right)^{\mathrm{T}}\)

Ahora hemos reducido nuestro problema a una forma abstracta idéntica a (21.20) y por lo tanto podemos aplicar cualquier esquema\(\mathbb{S}\) a (21.25) de la misma manera que a (21.20).

Por ejemplo, Euler Adelante, Euler Hacia Atrás, Crank-Nicolson, y AB2 aplicados a (21.25) tomarán la misma forma (21.21), (21.22), (21.23), (21.24), respectivamente, excepto que ahora\(w \in \mathbb{R}^{n}\),\(A \in \mathbb{R}^{n \times n}, F \in \mathbb{R}^{n}, w_{0} \in \mathbb{R}^{n}\) se dan en (21.25)\(n / 2\), donde, el número de osciladores (o masas) en nuestro sistema, ya no está restringido a\(n / 2=1\) (es decir,\(n=2\)). De manera similar podemos aplicar AB3 o BD1 o RK4.

Nuestro criterio de estabilidad también se extiende fácilmente. Primero notamos que ahora\(A\) tendrá en general\(n\) valores propios,\(\lambda_{1}, \lambda_{2}, \ldots, \lambda_{n}\). (En ciertos casos múltiples valores propios pueden crear dificultades; no consideramos estos casos típicamente bastante raros aquí). Nuestra condición de estabilidad se afirma entonces simplemente: un paso de tiempo\(\Delta t\) conducirá a un comportamiento estable si y solo si\(\lambda_{i} \Delta t\) está en\(\mathcal{R}_{\mathbb{S}}\) para todos\(i, 1 \leq i \leq n\). Si no se cumple esta condición entonces habrá uno (o más) modos que explotarán, llevándose consigo (o ellos) toda la solución. (Para ciertas condiciones iniciales muy especiales -en las que el\(w_{0}\) se elige de tal manera que todos los modos peligrosos son inicialmente exactamente cero- esta explosión podría evitarse con infinita precisión; pero en precisión finita todavía estaríamos condenados). Para esquemas explícitos,\(\Delta t_{\mathrm{cr}}\) es el mayor paso de tiempo tal que todos los rayos\(\left[0, \lambda_{i} \Delta t\right], 1 \leq i \leq n\), se encuentran dentro\(\mathcal{R}_{\mathbb{S}}\).

Ciertamente hay dificultades computacionales que surgen para grandes\(n\) que no son un problema para\(n=2\) (o pequeñas\(n\)). En primer lugar, para esquemas implícitos, la división necesaria -solución en lugar de evaluación de ecuaciones matriz-vector- se volverá considerablemente más costosa. Segundo, para esquemas explícitos, determinación de\(\Delta t_{\mathrm{cr}}\), o un atado\(\Delta t_{\mathrm{cr}}^{\text {conservative }}\) tal que\(\Delta t_{\mathrm{cr}}^{\text {conservative }} \approx t_{\mathrm{cr}}\) y\(\Delta t_{\mathrm{cr}}^{\text {conservative }} \leq \Delta t_{\mathrm{cr}}\), puede ser difícil. Como ya se mencionó, la descomposición modal completa puede ser costosa. Afortunadamente, para determinar\(\Delta t_{\mathrm{cr}}\), a menudo solo necesitamos como estimación para decir el valor propio real más negativo, o el valor propio imaginario más grande (en magnitud); estos valores propios extremos a menudo se pueden estimar de manera relativamente eficiente.

Finalmente, observamos que en la práctica a menudo se utilizan esquemas adaptativos en los que la estabilidad y la precisión son monitoreadas y\(∆t\) modificadas adecuadamente. Estos métodos también pueden abordar problemas no lineales —en los que ya\(h\) no depende linealmente de sus argumentos—.