método de diferencias finitas

Clase de técnicas numéricas.

En análisis numérico , los métodos de diferencias finitas ( FDM ) son una clase de técnicas numéricas para resolver ecuaciones diferenciales mediante la aproximación de derivadas con diferencias finitas . Tanto el dominio espacial como el dominio del tiempo (si corresponde) están discretizados o divididos en un número finito de intervalos, y los valores de la solución en los puntos finales de los intervalos se aproximan resolviendo ecuaciones algebraicas que contienen diferencias finitas y valores de puntos cercanos. .

Los métodos de diferencias finitas convierten ecuaciones diferenciales ordinarias (ODE) o ecuaciones diferenciales parciales (PDE), que pueden ser no lineales , en un sistema de ecuaciones lineales que pueden resolverse mediante técnicas de álgebra matricial . Las computadoras modernas pueden realizar estos cálculos de álgebra lineal de manera eficiente, lo que, junto con su relativa facilidad de implementación, ha llevado al uso generalizado de FDM en el análisis numérico moderno. ^[1] Hoy en día, los FDM son uno de los enfoques más comunes para la solución numérica de PDE, junto con los métodos de elementos finitos . ^[1]

Deducir el cociente de diferencias del polinomio de Taylor

Para una función diferenciable n veces, según el teorema de Taylor, la expansión en serie de Taylor viene dada como

$${\displaystyle f(x_{0}+h)=f(x_{0})+{\frac {f'(x_{0})}{1!}}h+{\frac {f^{(2)}(x_{0})}{2!}}h^{2}+\cdots +{\frac {f^{(n)}(x_{0})}{n!}}h^{n}+R_{n}(x),}$$

donde n ! denota el factorial de n , y R _n ( x ) es un término restante, que denota la diferencia entre el polinomio de Taylor de grado n y la función original.

El siguiente es el proceso para derivar una aproximación para la primera derivada de la función f truncando primero el polinomio de Taylor más el resto:

$${\displaystyle f(x_{0}+h)=f(x_{0})+f'(x_{0})h+R_{1}(x).}$$

h

$${\displaystyle {f(x_{0}+h) \over h}={f(x_{0}) \over h}+f'(x_{0})+{R_{1}(x) \over h}}$$

${\displaystyle f'(x_{0})}$

$${\displaystyle f'(x_{0})={f(x_{0}+h)-f(x_{0}) \over h}-{R_{1}(x) \over h}.}$$

Suponiendo que sea suficientemente pequeña, la aproximación de la primera derivada de f es: ${\displaystyle R_{1}(x)}$

$${\displaystyle f'(x_{0})\approx {f(x_{0}+h)-f(x_{0}) \over h}.}$$

Esto es similar a la definición de derivada, que es:

$${\displaystyle f'(x_{0})=\lim _{h\to 0}{\frac {f(x_{0}+h)-f(x_{0})}{h}}.}$$

Precisión y orden

El error en la solución de un método se define como la diferencia entre la aproximación y la solución analítica exacta. Las dos fuentes de error en los métodos de diferencias finitas son el error de redondeo , la pérdida de precisión debido al redondeo de cantidades decimales por computadora, y el error de truncamiento o error de discretización , la diferencia entre la solución exacta de la ecuación diferencial original y la cantidad exacta suponiendo aritmética perfecta (sin redondeo).

El método de diferencias finitas se basa en discretizar una función en una cuadrícula.

Para utilizar un método de diferencias finitas para aproximar la solución a un problema, primero se debe discretizar el dominio del problema. Esto generalmente se hace dividiendo el dominio en una cuadrícula uniforme (ver imagen). Esto significa que los métodos de diferencias finitas producen conjuntos de aproximaciones numéricas discretas a la derivada, a menudo en forma de "pasos de tiempo".

Una expresión de interés general es el error de truncamiento local de un método. Generalmente expresado usando la notación Big-O , el error de truncamiento local se refiere al error de una sola aplicación de un método. Es decir, es la cantidad si se refiere al valor exacto y a la aproximación numérica. El término restante del polinomio de Taylor se puede utilizar para analizar el error de truncamiento local . Usando la forma de Lagrange del resto del polinomio de Taylor para , que es ${\displaystyle f'(x_{i})-f'_{i}}$ ${\displaystyle f'(x_{i})}$ ${\displaystyle f'_{i}}$ ${\displaystyle f(x_{0}+h)}$

$${\displaystyle R_{n}(x_{0}+h)={\frac {f^{(n+1)}(\xi )}{(n+1)!}}(h)^{n+1}\,,\quad x_{0}<\xi <x_{0}+h,}$$

${\displaystyle f(x_{i})=f(x_{0}+ih)}$

$${\displaystyle f(x_{0}+ih)=f(x_{0})+f'(x_{0})ih+{\frac {f''(\xi )}{2!}}(ih)^{2},}$$

$${\displaystyle {\frac {f(x_{0}+ih)-f(x_{0})}{ih}}=f'(x_{0})+{\frac {f''(\xi )}{2!}}ih,}$$

$${\displaystyle {\frac {f(x_{0}+ih)-f(x_{0})}{ih}}=f'(x_{0})+O(h).}$$

En este caso, el error de truncamiento local es proporcional a los tamaños de paso. La calidad y duración de la solución FDM simulada depende de la selección de la ecuación de discretización y los tamaños de paso (pasos de tiempo y espacio). La calidad de los datos y la duración de la simulación aumentan significativamente con un tamaño de paso más pequeño. ^[2] Por lo tanto, para un uso práctico es necesario un equilibrio razonable entre la calidad de los datos y la duración de la simulación. Los pasos de tiempo grandes son útiles para aumentar la velocidad de la simulación en la práctica. Sin embargo, los intervalos de tiempo demasiado grandes pueden crear inestabilidades y afectar la calidad de los datos. ^{[3] [4]}

Los criterios de von Neumann y Courant-Friedrichs-Lewy a menudo se evalúan para determinar la estabilidad del modelo numérico. ^{[3] [4] [5] [6]}

Ejemplo: ecuación diferencial ordinaria

Por ejemplo, considere la ecuación diferencial ordinaria

$${\displaystyle u'(x)=3u(x)+2.}$$

método de Euler

$${\displaystyle {\frac {u(x+h)-u(x)}{h}}\approx u'(x)}$$

$${\displaystyle u(x+h)\approx u(x)+h(3u(x)+2).}$$

Ejemplo: la ecuación del calor.

Considere la ecuación de calor normalizada en una dimensión, con condiciones de contorno de Dirichlet homogéneas.

$${\displaystyle {\begin{cases}U_{t}=U_{xx}\\U(0,t)=U(1,t)=0&{\text{(boundary condition)}}\\U(x,0)=U_{0}(x)&{\text{(initial condition)}}\end{cases}}}$$

Una forma de resolver numéricamente esta ecuación es aproximar todas las derivadas por diferencias finitas. Primero divida el dominio en el espacio usando una malla y en el tiempo usando una malla . Supongamos una partición uniforme tanto en el espacio como en el tiempo, por lo que la diferencia entre dos puntos espaciales consecutivos será hy entre dos puntos temporales consecutivos será k . Los puntos ${\displaystyle x_{0},\dots ,x_{J}}$ ${\displaystyle t_{0},\dots ,t_{N}}$

$${\displaystyle u(x_{j},t_{n})=u_{j}^{n}}$$

representará la aproximación numérica de ${\displaystyle u(x_{j},t_{n}).}$

método explícito

La plantilla para el método explícito más común para la ecuación del calor.

Usando una diferencia directa en el tiempo y una ${\displaystyle t_{n}}$diferencia central de segundo orden para la derivada espacial en la posición ( FTCS ), se obtiene la ecuación de recurrencia: ${\displaystyle x_{j}}$

$${\displaystyle {\frac {u_{j}^{n+1}-u_{j}^{n}}{k}}={\frac {u_{j+1}^{n}-2u_{j}^{n}+u_{j-1}^{n}}{h^{2}}}.}$$

Este es un método explícito para resolver la ecuación del calor unidimensional .

Se pueden obtener de los demás valores de esta manera: ${\displaystyle u_{j}^{n+1}}$

$${\displaystyle u_{j}^{n+1}=(1-2r)u_{j}^{n}+ru_{j-1}^{n}+ru_{j+1}^{n}}$$

dónde ${\displaystyle r=k/h^{2}.}$

Entonces, con esta relación de recurrencia, y conociendo los valores en el tiempo n , se pueden obtener los valores correspondientes en el tiempo n +1. y debe ser reemplazado por las condiciones de contorno; en este ejemplo, ambas son 0. ${\displaystyle u_{0}^{n}}$ ${\displaystyle u_{J}^{n}}$

Se sabe que este método explícito es numéricamente estable y convergente siempre que . ^[7] Los errores numéricos son proporcionales al paso de tiempo y al cuadrado del paso de espacio: ${\displaystyle r\leq 1/2}$

$${\displaystyle \Delta u=O(k)+O(h^{2})}$$

Método implícito

La plantilla del método implícito.

Usando la diferencia hacia atrás en el tiempo ${\displaystyle t_{n+1}}$ y una diferencia central de segundo orden para la derivada espacial en la posición (El método del espacio centrado en el tiempo hacia atrás "BTCS") se obtiene la ecuación de recurrencia: ${\displaystyle x_{j}}$

$${\displaystyle {\frac {u_{j}^{n+1}-u_{j}^{n}}{k}}={\frac {u_{j+1}^{n+1}-2u_{j}^{n+1}+u_{j-1}^{n+1}}{h^{2}}}.}$$

Este es un método implícito para resolver la ecuación del calor unidimensional .

Se puede obtener resolviendo un sistema de ecuaciones lineales: ${\displaystyle u_{j}^{n+1}}$

$${\displaystyle (1+2r)u_{j}^{n+1}-ru_{j-1}^{n+1}-ru_{j+1}^{n+1}=u_{j}^{n}}$$

El esquema siempre es numéricamente estable y convergente, pero normalmente es más intensivo numéricamente que el método explícito, ya que requiere resolver un sistema de ecuaciones numéricas en cada paso de tiempo. Los errores son lineales en el paso del tiempo y cuadráticos en el paso del espacio:

$${\displaystyle \Delta u=O(k)+O(h^{2}).}$$

Método Crank-Nicolson

Finalmente, usando la diferencia central en el tiempo y una diferencia central de segundo orden para la derivada espacial en la posición ("CTCS") se obtiene la ecuación de recurrencia: ${\displaystyle t_{n+1/2}}$ ${\displaystyle x_{j}}$

$${\displaystyle {\frac {u_{j}^{n+1}-u_{j}^{n}}{k}}={\frac {1}{2}}\left({\frac {u_{j+1}^{n+1}-2u_{j}^{n+1}+u_{j-1}^{n+1}}{h^{2}}}+{\frac {u_{j+1}^{n}-2u_{j}^{n}+u_{j-1}^{n}}{h^{2}}}\right).}$$

Esta fórmula se conoce como método de Crank-Nicolson .

La plantilla de Crank-Nicolson.

Se puede obtener resolviendo un sistema de ecuaciones lineales: ${\displaystyle u_{j}^{n+1}}$

$${\displaystyle (2+2r)u_{j}^{n+1}-ru_{j-1}^{n+1}-ru_{j+1}^{n+1}=(2-2r)u_{j}^{n}+ru_{j-1}^{n}+ru_{j+1}^{n}}$$

El esquema es siempre numéricamente estable y convergente, pero normalmente es más intensivo numéricamente ya que requiere resolver un sistema de ecuaciones numéricas en cada paso de tiempo. Los errores son cuadráticos tanto en el paso de tiempo como en el paso de espacio:

$${\displaystyle \Delta u=O(k^{2})+O(h^{2}).}$$

Comparación

En resumen, normalmente el esquema de Crank-Nicolson es el esquema más preciso para pasos de tiempo pequeños. Para intervalos de tiempo mayores, el esquema implícito funciona mejor ya que es menos exigente desde el punto de vista computacional. El esquema explícito es el menos preciso y puede ser inestable, pero también es el más fácil de implementar y el menos intensivo numéricamente.

Aquí hay un ejemplo. Las siguientes figuras presentan las soluciones dadas por los métodos anteriores para aproximar la ecuación del calor.

$${\displaystyle U_{t}=\alpha U_{xx},\quad \alpha ={\frac {1}{\pi ^{2}}},}$$

con la condición de frontera

$${\displaystyle U(0,t)=U(1,t)=0.}$$

La solución exacta es

$${\displaystyle U(x,t)={\frac {1}{\pi ^{2}}}e^{-t}\sin(\pi x).}$$

Comparación de métodos de diferencias finitas

Ejemplo: el operador de Laplace

El operador de Laplace (continuo) en dimensiones viene dado por . El operador discreto de Laplace depende de la dimensión . ${\displaystyle n}$ ${\displaystyle \Delta u(x)=\sum _{i=1}^{n}\partial _{i}^{2}u(x)}$ ${\displaystyle \Delta _{h}u}$ ${\displaystyle n}$

En 1D, el operador de Laplace se aproxima como

$${\displaystyle \Delta u(x)=u''(x)\approx {\frac {u(x-h)-2u(x)+u(x+h)}{h^{2}}}=:\Delta _{h}u(x)\,.}$$

plantilla

$${\displaystyle \Delta _{h}={\frac {1}{h^{2}}}{\begin{bmatrix}1&-2&1\end{bmatrix}}}$$

matriz de Toeplitz

El caso 2D muestra todas las características del caso n-dimensional más general. Cada segunda derivada parcial debe aproximarse de manera similar al caso 1D

$${\displaystyle {\begin{aligned}\Delta u(x,y)&=u_{xx}(x,y)+u_{yy}(x,y)\\&\approx {\frac {u(x-h,y)-2u(x,y)+u(x+h,y)}{h^{2}}}+{\frac {u(x,y-h)-2u(x,y)+u(x,y+h)}{h^{2}}}\\&={\frac {u(x-h,y)+u(x+h,y)-4u(x,y)+u(x,y-h)+u(x,y+h)}{h^{2}}}\\&=:\Delta _{h}u(x,y)\,,\end{aligned}}}$$

plantilla

$${\displaystyle \Delta _{h}={\frac {1}{h^{2}}}{\begin{bmatrix}&1\\1&-4&1\\&1\end{bmatrix}}\,.}$$

Consistencia

La coherencia de la aproximación mencionada anteriormente se puede demostrar para funciones muy regulares, como . La declaración es ${\displaystyle u\in C^{4}(\Omega )}$

$${\displaystyle \Delta u-\Delta _{h}u={\mathcal {O}}(h^{2})\,.}$$

Para demostrar esto, es necesario sustituir expansiones de la serie de Taylor hasta el orden 3 en el operador discreto de Laplace.

Propiedades

subarmónico

De manera similar a las funciones subarmónicas continuas , se pueden definir funciones subarmónicas para aproximaciones en diferencias finitas. ${\displaystyle u_{h}}$

$${\displaystyle -\Delta _{h}u_{h}\leq 0\,.}$$

Valor medio

Se puede definir una plantilla general de tipo positivo mediante

$${\displaystyle {\begin{bmatrix}&\alpha _{N}\\\alpha _{W}&-\alpha _{C}&\alpha _{E}\\&\alpha _{S}\end{bmatrix}}\,,\quad \alpha _{i}>0\,,\quad \alpha _{C}=\sum _{i\in \{N,E,S,W\}}\alpha _{i}\,.}$$

Si es un subarmónico (discreto), entonces se cumple la siguiente propiedad del valor medio ${\displaystyle u_{h}}$

$${\displaystyle u_{h}(x_{C})\leq {\frac {\sum _{i\in \{N,E,S,W\}}\alpha _{i}u_{h}(x_{i})}{\sum _{i\in \{N,E,S,W\}}\alpha _{i}}}\,,}$$

Una propiedad similar del valor medio también se cumple para el caso continuo.

Principio máximo

Para una función subarmónica (discreta) se cumple lo siguiente ${\displaystyle u_{h}}$

$${\displaystyle \max _{\Omega _{h}}u_{h}\leq \max _{\partial \Omega _{h}}u_{h}\,,}$$

${\displaystyle \Omega _{h},\partial \Omega _{h}}$ ${\displaystyle \Omega }$ ${\displaystyle \partial \Omega }$

Un principio de máximo similar también se aplica al caso continuo.

El método SBP-SAT

El método SBP-SAT ( suma por partes - término de aproximación simultánea ) es una técnica estable y precisa para discretizar e imponer condiciones de contorno de una ecuación diferencial parcial bien planteada utilizando diferencias finitas de alto orden. ^{[8] [9]}

El método se basa en diferencias finitas donde los operadores de diferenciación exhiben propiedades de suma por partes . Por lo general, estos operadores consisten en matrices de diferenciación con plantillas de diferencia central en el interior con plantillas de límites unilaterales cuidadosamente elegidas diseñadas para imitar la integración por partes en el entorno discreto. Utilizando la técnica SAT, las condiciones límite del PDE se imponen débilmente, donde los valores límite son "tirados" hacia las condiciones deseadas en lugar de cumplirse exactamente. Si los parámetros de sintonización (inherentes a la técnica SAT) se eligen adecuadamente, el sistema resultante de ODE exhibirá un comportamiento energético similar al del PDE continuo, es decir, el sistema no tendrá crecimiento de energía no físico. Esto garantiza la estabilidad si se utiliza un esquema de integración con una región de estabilidad que incluya partes del eje imaginario, como el método de Runge-Kutta de cuarto orden . Esto hace que la técnica SAT sea un método atractivo para imponer condiciones de contorno para métodos de diferencias finitas de orden superior, en contraste, por ejemplo, con el método de inyección, que normalmente no será estable si se utilizan operadores de diferenciación de orden superior.

Ver también

• Método de elementos finitos
• Diferencia finita
• Dominio del tiempo de diferencias finitas
• método de diferencias infinitas
• Plantilla (análisis numérico)
• Coeficientes de diferencias finitas
• Plantilla de cinco puntos
• Teorema de Lax-Richtmyer
• Métodos de diferencias finitas para la valoración de opciones
• Esquema de diferenciación contra el viento para convección.
• Esquema de diferenciación central
• Ecuación de Poisson discreta
• Operador discreto de Laplace

Referencias

1. Christian Grossmann; Hans-G. Roos; Martín Stynes ​​(2007). Tratamiento Numérico de Ecuaciones Diferenciales Parciales . Medios de ciencia y negocios de Springer. pag. 23.ISBN _ 978-3-540-71584-9.
2. Arieh Iserles (2008). Un primer curso de análisis numérico de ecuaciones diferenciales . Prensa de la Universidad de Cambridge. pag. 23.ISBN _ 9780521734905.
3. Hoffman JD; Frankel S (2001). Métodos numéricos para ingenieros y científicos . Prensa CRC, Boca Ratón.
4. Jaluria Y; Atluri S (1994). "Transferencia de calor computacional". Mecánica Computacional . 14 (5): 385–386. Código Bib : 1994CompM..14..385J. doi :10.1007/BF00377593. S2CID  119502676.
5. Majumdar P (2005). Métodos computacionales para la transferencia de calor y masa (1ª ed.). Taylor y Francis, Nueva York.
6. Smith GD (1985). Solución numérica de ecuaciones diferenciales parciales: métodos de diferencias finitas (3ª ed.). Prensa de la Universidad de Oxford.
7. Crank, J. Las matemáticas de la difusión . Segunda edición, Oxford, 1975, pág. 143.
8. Bo Strand (1994). "Suma por partes para aproximaciones en diferencias finitas para d/dx". Revista de Física Computacional . 110 (1): 47–67. Código bibliográfico : 1994JCoPh.110...47S. doi :10.1006/jcph.1994.1005.
9. Mark H. Carpintero; David I. Gottlieb; Saúl S. Abarbanel (1994). "Condiciones de contorno estables en el tiempo para esquemas de diferencias finitas que resuelven sistemas hiperbólicos: metodología y aplicación a esquemas compactos de alto orden". Revista de Física Computacional . 111 (2): 220–236. Código Bib : 1994JCoPh.111..220C. doi :10.1006/jcph.1994.1057. hdl : 2060/19930013937 .

Otras lecturas

• KW Morton y DF Mayers, Solución numérica de ecuaciones diferenciales parciales, Introducción . Prensa de la Universidad de Cambridge, 2005.
• Autar Kaw y E. Eric Kalu, Métodos numéricos con aplicaciones , (2008) [1]. Contiene una breve introducción orientada a la ingeniería a FDM (para ODE) en el Capítulo 08.07.
• John Strikwerda (2004). Esquemas en diferencias finitas y ecuaciones diferenciales parciales (2ª ed.). SIAM. ISBN 978-0-89871-639-9.
• Smith, GD (1985), Solución numérica de ecuaciones diferenciales parciales: métodos de diferencias finitas, 3ª ed. , Prensa de la Universidad de Oxford
• Peter Olver (2013). Introducción a las ecuaciones diferenciales parciales. Saltador. Capítulo 5: Diferencias finitas. ISBN 978-3-319-02099-0..
• Randall J. LeVeque , Métodos de diferencias finitas para ecuaciones diferenciales ordinarias y parciales , SIAM, 2007.
• Sergey Lemeshevsky, Piotr Matus, Dmitriy Poliakov (Eds): "Esquemas exactos de diferencias finitas", De Gruyter (2016). DOI: https://doi.org/10.1515/9783110491326.