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]
Para una función diferenciable n veces, según el teorema de Taylor , la expansión en serie de Taylor viene dada como
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: Dividiendo por h se obtiene: Resolviendo para :
Suponiendo que sea suficientemente pequeña, la aproximación de la primera derivada de f es:
Esto es similar a la definición de derivada, que es: excepto por el límite hacia cero (el método lleva este nombre).
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).
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 . Utilizando la forma de Lagrange del resto del polinomio de Taylor , se puede descubrir cuál es el término dominante del error de truncamiento local. Por ejemplo, usando nuevamente la fórmula de diferencias directas para la primera derivada, sabiendo que , y con alguna manipulación algebraica, esto lleva a observar, además, que la cantidad de la izquierda es la aproximación del método de diferencias finitas y que la cantidad de la izquierda La derecha es la cantidad exacta de interés más un resto, claramente ese resto es el error de truncamiento local. Una expresión final de este ejemplo y su orden es:
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]
Por ejemplo, considere la ecuación diferencial ordinaria. El método de Euler para resolver esta ecuación utiliza el cociente de diferencias finitas para aproximar la ecuación diferencial sustituyéndola primero por u'(x) y luego aplicando un poco de álgebra (multiplicando ambos lados por h y luego sumando u(x) a ambos lados) para obtener La última ecuación es una ecuación en diferencias finitas y al resolverla se obtiene una solución aproximada a la ecuación diferencial.
Considere la ecuación de calor normalizada en una dimensión, con condiciones de contorno de Dirichlet homogéneas.
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
representará la aproximación numérica de
Usando una diferencia directa en el tiempo y una diferencia central de segundo orden para la derivada espacial en la posición ( FTCS ), se obtiene la ecuación de recurrencia:
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:
dónde
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.
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:
Usando la diferencia hacia atrás en el tiempo 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:
Este es un método implícito para resolver la ecuación del calor unidimensional .
Se puede obtener resolviendo un sistema de ecuaciones lineales:
El esquema es siempre 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:
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:
Esta fórmula se conoce como método de Crank-Nicolson .
Se puede obtener resolviendo un sistema de ecuaciones lineales:
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:
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.
con la condición de frontera
La solución exacta es
El operador de Laplace (continuo) en dimensiones viene dado por . El operador discreto de Laplace depende de la dimensión .
En 1D, el operador de Laplace se aproxima como Esta aproximación generalmente se expresa mediante la siguiente plantilla y representa una matriz tridiagonal simétrica. Para una cuadrícula equidistante se obtiene una 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 que generalmente se muestra en la siguiente plantilla
La coherencia de la aproximación mencionada anteriormente se puede demostrar para funciones muy regulares, como . La declaración es
Para demostrar esto, es necesario sustituir expansiones de la serie de Taylor hasta el orden 3 en el operador discreto de Laplace.
De manera similar a las funciones subarmónicas continuas, se pueden definir funciones subarmónicas para aproximaciones en diferencias finitas.
Se puede definir una plantilla general de tipo positivo mediante
Si es un subarmónico (discreto), entonces se cumple la siguiente propiedad del valor medio cuando la aproximación se evalúa en puntos de la cuadrícula y se supone que la plantilla es de tipo positivo.
Una propiedad similar del valor medio también se cumple para el caso continuo.
Para una función subarmónica (discreta) se cumple lo siguiente: ¿dónde están las discretizaciones del dominio continuo , respectivamente, el límite ?
Un principio máximo similar también se aplica al caso continuo.
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.