Integrador exponencial

Los integradores exponenciales son una clase de métodos numéricos para la solución de ecuaciones diferenciales ordinarias , específicamente problemas de valor inicial . Esta gran clase de métodos de análisis numérico se basa en la integración exacta de la parte lineal del problema de valor inicial. Debido a que la parte lineal se integra exactamente, esto puede ayudar a mitigar la rigidez de una ecuación diferencial. Los integradores exponenciales se pueden construir para que sean explícitos o implícitos para ecuaciones diferenciales ordinarias numéricas o para que sirvan como integrador temporal para ecuaciones diferenciales parciales numéricas .

Fondo

Estos métodos, que se remontan al menos a la década de 1960, fueron reconocidos por Certaine ^[1] y Pope. ^[2] Últimamente, los integradores exponenciales se han convertido en un área activa de investigación, consulte Hochbruck y Ostermann (2010). ^[3] Originalmente desarrollados para resolver ecuaciones diferenciales rígidas , los métodos se han utilizado para resolver ecuaciones diferenciales parciales, incluidos problemas hiperbólicos y parabólicos ^[4] como la ecuación del calor .

Introducción

Consideramos problemas de valor inicial de la forma,

${\displaystyle u'(t)=Lu(t)+N(u(t)),\qquad u(t_{0})=u_{0},\qquad \qquad (1)}$

donde se compone de términos lineales y se compone de términos no lineales . Estos problemas pueden surgir de un problema de valor inicial más típico ${\displaystyle L}$ ${\displaystyle N}$

${\displaystyle u'(t)=f(u(t)),\qquad u(t_{0})=u_{0},}$

después de linealizar localmente sobre un estado fijo o local : ${\displaystyle u^{*}}$

${\displaystyle L={\frac {\partial f}{\partial u}}(u^{*});\qquad N=f(u)-Lu.}$

Aquí, se refiere a la derivada parcial de con respecto a (el jacobiano de f). ${\displaystyle {\frac {\partial f}{\partial u}}}$ ${\displaystyle f}$ ${\displaystyle u}$

La integración exacta de este problema desde el tiempo 0 hasta un tiempo posterior se puede realizar utilizando matrices exponenciales para definir una ecuación integral para la solución exacta: ^[3] ${\displaystyle t}$

${\displaystyle u(t)=e^{Lt}u_{0}+\int _{0}^{t}e^{L(t-\tau )}N\left(u\left(\tau \right)\right)\,d\tau .\qquad (2)}$

Esto es similar a la integral exacta utilizada en el teorema de Picard-Lindelöf . En el caso de , esta formulación es la solución exacta de la ecuación diferencial lineal . ${\displaystyle N\equiv 0}$

Los métodos numéricos requieren una discretización de la ecuación (2). Pueden basarse en discretizaciones de Runge-Kutta , ^{[5] [6] [7]} métodos lineales de múltiples pasos o una variedad de otras opciones.

Métodos exponenciales de Rosenbrock

Se ha demostrado que los métodos exponenciales de Rosenbrock son muy eficientes para resolver grandes sistemas de ecuaciones diferenciales ordinarias rígidas, generalmente resultantes de la discretización espacial de ecuaciones diferenciales parciales (parabólicas) dependientes del tiempo. Estos integradores se construyen en base a una linealización continua de (1) a lo largo de la solución numérica. ${\displaystyle u_{n}}$

${\displaystyle u'(t)=L_{n}u(t)+N_{n}(u(t)),\qquad u(t_{0})=u_{0},\qquad (3)}$

donde Este procedimiento tiene la ventaja, en cada paso, de que Esto simplifica considerablemente la derivación de las condiciones de orden y mejora la estabilidad al integrar la no linealidad . Nuevamente, aplicando la fórmula de variación de constantes (2) se obtiene la solución exacta en el tiempo como ${\displaystyle L_{n}={\frac {\partial f}{\partial u}}(u_{n}),N_{n}(u)=f(u)-L_{n}u.}$ ${\displaystyle {\frac {\partial N_{n}}{\partial u}}(u_{n})=0.}$ ${\displaystyle N(u(t))}$ ${\displaystyle t_{n+1}}$

${\displaystyle u(t_{n+1})=e^{h_{n}L_{n}}u(t_{n})+\int _{0}^{h_{n}}e^{(h_{n}-\tau )L_{n}}N_{n}(u(t_{n}+\tau ))d\tau .\qquad (4)}$

La idea ahora es aproximar la integral en (4) mediante alguna regla de cuadratura con nodos y pesos ( ). Esto produce la siguiente clase de métodos exponenciales explícitos de Rosenbrock, véase Hochbruck y Ostermann (2006), Hochbruck, Ostermann y Schweitzer (2009): ${\displaystyle c_{i}}$ ${\displaystyle b_{i}(h_{n}L_{n})}$ ${\displaystyle 1\leq i\leq s}$ ${\displaystyle s-stage}$

${\displaystyle U_{ni}=e^{c_{i}h_{n}L_{n}}u_{n}+h_{n}\sum _{j=1}^{i-1}a_{ij}(h_{n}L_{n})N_{n}(U_{nj}),}$

${\displaystyle u_{n+1}=e^{h_{n}L_{n}}u_{n}+h_{n}\sum _{i=1}^{s}b_{i}(h_{n}L_{n})N_{n}(U_{ni})}$

con . Los coeficientes se eligen generalmente como combinaciones lineales de todas las funciones , respectivamente, donde ${\displaystyle u_{n}\approx u(t_{n}),U_{ni}\approx u(t_{n}+c_{i}h_{n}),h_{n}=t_{n+1}-t_{n}}$ ${\displaystyle a_{ij}(z),b_{i}(z)}$ ${\displaystyle \varphi _{k}(c_{i}z),\varphi _{k}(z)}$

${\displaystyle \varphi _{0}(z)=e^{z},\quad \varphi _{k}(z)=\int _{0}^{1}e^{(1-\theta )z}{\frac {\theta ^{k-1}}{(k-1)!}}d\theta ,\quad k\geq 1.}$

Estas funciones satisfacen la relación de recursión

${\displaystyle \varphi _{k+1}(z)={\frac {\varphi _{k}(z)-\varphi _{k}(0)}{z}},\ k\geq 0.}$

Al introducir la diferencia , se pueden reformular de una manera más eficiente para su implementación (ver también ^[3] ) como ${\displaystyle D_{ni}=N_{n}(U_{ni})-N_{n}(u_{n})}$

${\displaystyle U_{ni}=u_{n}+c_{i}h_{n}\varphi _{1}(c_{i}h_{n}L_{n})f(u_{n})+h_{n}\sum _{j=2}^{i-1}a_{ij}(h_{n}L_{n})D_{nj},}$

${\displaystyle u_{n+1}=u_{n}+h_{n}\varphi _{1}(h_{n}L_{n})f(u_{n})+h_{n}\sum _{i=2}^{s}b_{i}(h_{n}L_{n})D_{ni}.}$

Para implementar este esquema con tamaño de paso adaptativo, se pueden considerar, para fines de estimación de error local, los siguientes métodos integrados

${\displaystyle {\bar {u}}_{n+1}=u_{n}+h_{n}\varphi _{1}(h_{n}L_{n})f(u_{n})+h_{n}\sum _{i=2}^{s}{\bar {b}}_{i}(h_{n}L_{n})D_{ni},}$

que utilizan las mismas etapas pero con pesos . ${\displaystyle U_{ni}}$ ${\displaystyle {\bar {b}}_{i}}$

Para mayor comodidad, los coeficientes de los métodos exponenciales explícitos de Rosenbrock junto con sus métodos integrados se pueden representar utilizando la denominada tabla de Butcher reducida de la siguiente manera:

Condiciones de pedido estrictas

^{Además, en Luan y Ostermann (2014a) [8]} se muestra que el enfoque de reformulación ofrece una forma nueva y sencilla de analizar los errores locales y, por lo tanto, de derivar las condiciones de orden rígido para los métodos exponenciales de Rosenbrock hasta el orden 5. Con la ayuda de esta nueva técnica junto con una extensión del concepto de serie B, finalmente se ha dado una teoría para derivar las condiciones de orden rígido para los integradores exponenciales de Rosenbrock de orden arbitrario en Luan y Ostermann (2013). ^[9] Como ejemplo, en ese trabajo se han derivado las condiciones de orden rígido para los métodos exponenciales de Rosenbrock hasta el orden 6, que se indican en la siguiente tabla:

${\displaystyle {\begin{array}{|c|c|c| }\hline No.&{\text{Stiff order condition}}&{\text{Order}}\\\hline 1&\sum _{i=2}^{s}b_{i}(Z)c_{i}^{2}=2\varphi _{3}(Z)&3\\\hline 2&\sum _{i=2}^{s}b_{i}(Z)c_{i}^{3}=6\varphi _{4}(Z)&4\\\hline 3&\sum _{i=2}^{s}b_{i}(Z)c_{i}^{4}=24\varphi _{5}(Z)&5\\4&\sum _{i=2}^{s}b_{i}(Z)c_{i}K\psi _{3,i}(Z)=0&5\\\hline 5&\sum _{i=2}^{s}b_{i}(Z)c_{i}^{5}=120\varphi _{6}(Z)&6\\6&\sum _{i=2}^{s}b_{i}(Z)c_{i}^{2}M\psi _{3,i}(Z)=0&6\\7&\sum _{i=2}^{s}b_{i}(Z)c_{i}K\psi _{4,i}(Z)=0&6\\\hline \end{array}}}$

Aquí denotamos matrices cuadradas arbitrarias. ${\displaystyle Z,K,M.\ }$

Análisis de convergencia

Se demuestran los resultados de estabilidad y convergencia para los métodos exponenciales de Rosenbrock en el marco de semigrupos fuertemente continuos en algún espacio de Banach.

Ejemplos

Todos los esquemas presentados a continuación cumplen las condiciones de orden rígido y, por lo tanto, también son adecuados para resolver problemas rígidos.

Método de segundo orden

El método exponencial de Rosenbrock más simple es el esquema exponencial de Rosenbrock-Euler, que tiene orden 2, véase, por ejemplo, Hochbruck et al. (2009):

${\displaystyle u_{n+1}=u_{n}+h_{n}\ \varphi _{1}(h_{n}L_{n})f(u_{n}).}$

Métodos de tercer orden

En Hochbruck et al. (2009) se derivó una clase de métodos Rosenbrock exponenciales de tercer orden, denominada exprb32, que se presenta de la siguiente manera:

exprb32:

que se lee como

${\displaystyle U_{n2}=u_{n}+h_{n}\ \varphi _{1}(h_{n}L_{n})f(u_{n}),}$

${\displaystyle u_{n+1}=u_{n}+h_{n}\ \varphi _{1}(h_{n}L_{n})f(u_{n})+h_{n}\ 2\varphi _{3}(h_{n}L_{n})D_{n2},}$

dónde ${\displaystyle D_{n2}=N_{n}(U_{n2})-N_{n}(u_{n}).}$

Para una implementación de este esquema de tamaño de paso variable, se puede integrar con el Rosenbrock-Euler exponencial:

${\displaystyle {\hat {u}}_{n+1}=u_{n}+h_{n}\ \varphi _{1}(h_{n}L_{n})f(u_{n}).}$

Método ETDRK4 de cuarto orden de Cox y Matthews

Cox y Matthews ^[5] describen un método de diferenciación temporal exponencial (ETD) de cuarto orden que utilizaron Maple para derivar.

Usamos su notación y asumimos que la función desconocida es , y que tenemos una solución conocida en el tiempo . Además, haremos uso explícito de un lado derecho posiblemente dependiente del tiempo: . ${\displaystyle u}$ ${\displaystyle u_{n}}$ ${\displaystyle t_{n}}$ ${\displaystyle {\mathcal {N}}={\mathcal {N}}(t,u)}$

Primero se construyen tres valores de etapa:

${\displaystyle a_{n}=e^{Lh/2}u_{n}+L^{-1}\left(e^{Lh/2}-I\right){\mathcal {N}}(t_{n},u_{n})}$

${\displaystyle b_{n}=e^{Lh/2}u_{n}+L^{-1}\left(e^{Lh/2}-I\right){\mathcal {N}}(t_{n}+h/2,a_{n})}$

${\displaystyle c_{n}=e^{Lh/2}a_{n}+L^{-1}\left(e^{Lh/2}-I\right)\left(2{\mathcal {N}}(t_{n}+h/2,b_{n})-{\mathcal {N}}(t_{n},u_{n})\right)}$

La actualización final viene dada por,

${\displaystyle u_{n+1}=e^{Lh}u_{n}+h^{-2}L^{-3}\left\{\left[-4-Lh+e^{Lh}\left(4-3Lh+(Lh)^{2}\right)\right]{\mathcal {N}}(t_{n},u_{n})+2\left[2+Lh+e^{Lh}\left(-2+Lh\right)\right]\left({\mathcal {N}}(t_{n}+h/2,a_{n})+{\mathcal {N}}(t_{n}+h/2,b_{n})\right)+\left[-4-3Lh-(Lh)^{2}+e^{Lh}\left(4-Lh\right)\right]{\mathcal {N}}(t_{n}+h,c_{n})\right\}.}$

Si se implementa de manera ingenua, el algoritmo anterior sufre inestabilidades numéricas debido a errores de redondeo de punto flotante . ^[10] Para ver por qué, considere la primera función,

${\displaystyle \varphi _{1}(z)={\frac {e^{z}-1}{z}},}$

que está presente en el método de Euler de primer orden, así como en las tres etapas de ETDRK4. Para valores pequeños de , esta función sufre errores de cancelación numérica. Sin embargo, estos problemas numéricos se pueden evitar evaluando la función a través de un enfoque de integral de contorno ^[10] o mediante un aproximador de Padé . ^[11] ${\displaystyle z}$ ${\displaystyle \varphi _{1}}$

Aplicaciones

Los integradores exponenciales se utilizan para la simulación de escenarios rígidos en computación científica y visual , por ejemplo en dinámica molecular , ^[12] para simulación de circuitos VLSI , ^{[13] [14]} y en gráficos de computadora . ^[15] También se aplican en el contexto de métodos híbridos de Monte Carlo . ^[16] En estas aplicaciones, los integradores exponenciales muestran la ventaja de una gran capacidad de paso de tiempo y una alta precisión. Para acelerar la evaluación de funciones matriciales en escenarios tan complejos, los integradores exponenciales a menudo se combinan con métodos de proyección del subespacio de Krylov.

Véase también

• Métodos lineales generales

Notas

1. Cierta (1960)
2. Papa (1963)
3. Hochbruck y Ostermann (2010)
4. Hochbruck y Ostermann (2005a), Hochbruck y Ostermann (2005b)
5. Cox y Matthews (2002)
6. Tokman (2006)
7. Tokman (2011)
8. Luan y Ostermann (2014a)
9. Luan y Ostermann (2013)
10. Kassam y Trefethen (2005)
11. Berland, Skaflestad y Wright (2007)
12. Michels y Desbrun (2015)
13. Zhuang y otros (2014)
14. Weng, Chen y Cheng (2012)
15. Michels, Sobottka y Weber (2014)
16. Chao y otros (2015)

Referencias

• Berland, Havard; Owren, Brynjulf; Skaflestad, Bard (2005). "Series B y condiciones de orden para integradores exponenciales". Revista SIAM sobre análisis numérico . 43 (4): 1715–1727. CiteSeerX  10.1.1.216.5645 . doi :10.1137/040612683.
• Berland, Havard; Skaflestad, Bard; Wright, Will M. (2007). "EXPINT-A MATLAB Package for Exponential Integrators". ACM Transactions on Mathematical Software . 33 (1): 4–es. doi :10.1145/1206040.1206044. S2CID  1525599.
• Chao, Wei-Lun; Solomon, Justin; Michels, Dominik L.; Sha, Fei (2015). "Integración exponencial para el método Monte Carlo hamiltoniano". Actas de la 32.ª Conferencia internacional sobre aprendizaje automático (ICML-15) : 1142–1151.
• Certaine, John (1960). "La solución de ecuaciones diferenciales ordinarias con constantes de tiempo grandes". Métodos matemáticos para computadoras digitales . Wiley. págs. 128–132.
• Cox, SM; Matthews, PC (marzo de 2002). "Diferenciación temporal exponencial para sistemas rígidos". Journal of Computational Physics . 176 (2): 430–455. Bibcode :2002JCoPh.176..430C. doi :10.1006/jcph.2002.6995.
• Hochbruck, Marlis ; Ostermann, Alexander (mayo de 2010). "Integradores exponenciales". Acta Numerica . 19 : 209–286. Código Bibliográfico :2010AcNum..19..209H. CiteSeerX  10.1.1.187.6794 . doi :10.1017/S0962492910000048. S2CID  4841957.
• Hochbruck, Marlis ; Ostermann, Alexander (2005a). "Métodos de Runge-Kutta exponenciales explícitos para problemas parabólicos semilineales". Revista SIAM sobre análisis numérico . 43 (3): 1069–1090. CiteSeerX  10.1.1.561.5501 . doi :10.1137/040611434.
• Hochbruck, Marlis ; Ostermann, Alexander (mayo de 2005b). "Métodos exponenciales de Runge-Kutta para problemas parabólicos". Matemáticas numéricas aplicadas . 53 (2–4): 323–339. doi :10.1016/j.apnum.2004.08.005.
• Luan, Vu Thai; Ostermann, Alexander (2014a). "Métodos exponenciales de Rosenbrock de orden cinco: construcción, análisis y comparaciones numéricas". Revista de Matemática Computacional y Aplicada . 255 : 417–431. doi : 10.1016/j.cam.2013.04.041 .
• Luan, Vu Thai; Ostermann, Alexander (2014c). "Métodos de Runge-Kutta exponenciales explícitos de alto orden para problemas parabólicos". Revista de Matemática Computacional y Aplicada . 256 : 168–179. arXiv : 1307.0661 . doi :10.1016/j.cam.2013.07.027. S2CID  18448807.
• Luan, Vu Thai; Ostermann, Alexander (2013). "Serie B exponencial: el caso rígido". Revista SIAM sobre análisis numérico . 51 (6): 3431–3445. doi :10.1137/130920204.
• Luan, Vu Thai; Ostermann, Alexander (2014). "Condiciones de orden rígido para métodos exponenciales de Runge–Kutta de orden cinco". En Bock, Hans Georg; Hoang, Xuan Phu; Rannacher, Rolf; Schlöder, Johannes P. (eds.). Modelado, simulación y optimización de procesos complejos – HPSC 2012: Actas de la Quinta Conferencia Internacional sobre Computación Científica de Alto Rendimiento, 5–9 de marzo de 2012, Hanoi, Vietnam . Springer. págs. 133–143. doi :10.1007/978-3-319-09063-4_11. ISBN 978-3-319-09062-7.
• Luan, Vu Thai; Ostermann, Alexander (2016). "Métodos Rosenbrock exponenciales paralelos". Computadoras y matemáticas con aplicaciones . 71 (5): 1137–1150. doi : 10.1016/j.camwa.2016.01.020 .
• Michels, Dominik L.; Desbrun, Mathieu (2015). "Un enfoque semianalítico de la dinámica molecular". Journal of Computational Physics . 303 : 336–354. Bibcode :2015JCoPh.303..336M. doi : 10.1016/j.jcp.2015.10.009 .
• Michels, Dominik L.; Sobottka, Gerrit A.; Weber, Andreas G. (2014). "Integradores exponenciales para problemas elastodinámicos rígidos". ACM Transactions on Graphics . 33 : 7:1–7:20. doi :10.1145/2508462. S2CID  207207156.
• Pope, David A (1963). "Un método exponencial de integración numérica de ecuaciones diferenciales ordinarias". Comunicaciones de la ACM . 6 (8): 491–493. doi : 10.1145/366707.367592 . S2CID  18598461.
• Tokman, Mayya (octubre de 2011). "Una nueva clase de métodos iterativos de propagación exponencial de tipo Runge–Kutta (EPIRK)". Journal of Computational Physics . 230 (24): 8762–8778. Bibcode :2011JCoPh.230.8762T. doi :10.1016/j.jcp.2011.08.023.
• Tokman, Mayya (abril de 2006). "Integración eficiente de grandes sistemas rígidos de EDO con métodos iterativos de propagación exponencial (EPI)". Journal of Computational Physics . 213 (2): 748–776. Bibcode :2006JCoPh.213..748T. doi :10.1016/j.jcp.2005.08.032.
• Kassam, Aly-Khan; Trefethen, Lloyd N. (2005). "Pasos temporales de cuarto orden para ecuaciones diferenciales parciales rígidas". Revista SIAM de informática científica . 26 (4): 1214–1233. Bibcode :2005SJSC...26.1214K. CiteSeerX  10.1.1.15.6467 . doi :10.1137/S1064827502410633.
• Zhuang, Hao; Weng, Shih-Hung; Lin, Jeng-Hau; Cheng, Chung-Kuan (2014). "MATEX" (PDF) . Actas de la 51.ª Conferencia Anual de Automatización del Diseño sobre la Conferencia de Automatización del Diseño - DAC '14 . págs. 1–6. arXiv : 1511.04519 . doi :10.1145/2593069.2593160. ISBN 9781450327305.S2CID10585362  .​
• Weng, Shih-Hung; Chen, Quan; Cheng, Chung-Kuan (2012). "Análisis en el dominio del tiempo de circuitos a gran escala mediante el método exponencial matricial con control adaptativo". IEEE Transactions on Computer-Aided Design of Integrated Circuits and Systems . 32 (8): 1180–1193. doi :10.1109/TCAD.2012.2189396. S2CID  14977067.

Enlaces externos

• Integradores en GPGPU
• Código para un integrador exponencial sin malla