Introducción sobre elementos finitos II

Hola a todos

Tras un largo parón navideño (de casi 2 meses), quería volver a retomar mi actividad hablando un poco más de los elementos finitos.

Como recordatorio breve, comentar que el método de los elementos finitos es una técnica matemática que permite obtener el cálculo aproximado de la solución de una ecuación determinada sin llegar a resolverla. Este método permite la implementación de la formulación débil de las ecuaciones diferenciales presentes en los modelos físicos: lo que quiero decir es que este método satisface forma integral de nuestro problema de contorno, y no el problema de contorno en si.

El método se basa, esencialmente, en aproximar nuestra solución por una combinación lineal de ciertas funciones que nosotros propondremos y que se denominan "funciones de base". Esto es, si nuestro problema es el siguiente:

$$ L(u)=f $$

donde u es nuestra variable y L(u) es un operador aplicado a u, aproximamos la solución de u por $\hat{u}$:

$$ \hat{u} = \Sigma a_i · N_i $$

donde $N_i$ son las funciones de forma, que a priori son conocidas.

Existen una variedad muy grande de funciones de forma: normalmente, estas misteriosas funciones que aparecen en todos los libros como la panacea para resolver todos nuestros problemas por arte de magia no son más que el fruto del estudio de lo que se denomina teoría de espacios funcionales.
Este campo se dedica a hallar funciones, a priori desconocidas, que cumplan una serie de propiedades que puedan ser útiles para la resolución de problemas matemáticos: cumplir con condiciones de contorno, con condiciones de integrabilidad, topología de dimensión determinada.

En los libros suelen aparecer familias de funciones más o menos complejas: normamente poseen una forma de spline (funciones polinomiales con n puntos conocidos), si bien todas ellas suelen ser de cuadrado integrable: esta condición asegura que la integral de dicha función no va a ser infinito, lo cual es condición necesaria en la resolución de cualquier problema físico real. Ejemplos clásicos pueden ser los polinomios de Lagrange o los polinomios de Hermite. Sin embargo, otras funciones de base pueden ser utilizados en ciertos problemas, como las series de Fourier.

Por lo tanto, conocida la forma de dichas funciones, el objetivo del método, esencialmente, consistirá en buscar aquellos valores que cumplan que:

$L(\hat{u}) -   f  \approx 0$

Sabemos que, al no ser u solución exacta, no lo cumple. Como ya dijimos, entonces debe cumplir:

$L(\hat{u}) -   f  = r $

Con r, función residuo no nulo de la ecuación. La formulación integral del problema entra aquí. Nosotros vamos a imponer que el valor del residuo en promedio a lo largo de toda la función va a ser cero. Para ello, es necesario integrar el valor de r a lo largo del dominio.

$\int w(u) r (u) d\Omega = 0$

Esta condición se denomina en matemáticas "Método de los residuos ponderados".

La función w(u) se denomina función peso. La pregunta que surge es obvia: ¿Por que utilizamos esta fórmula con una función misteriosa w(u)? ¿Que motivo real existe para utilizar dicha expresión?. Bien, en el fondo la razón que subyace es la de que queremos que nuestra solución sea una buena solución. El método de elementos finitos busca una SOLUCIÓN APROXIMADA a la ecuación diferencial, cuyas soluciones existen en un subespacio concreto. Nosotros no sabemos cual es en concreto (sólo podemos conocer el ERROR que cometemos entre nuestra aproximación y el valor real),  pero si podemos establecer un CRITERIO DE APROXIMACIÓN, esto es, establecer que aproximación puede ser válida para nosotros.
 Esto se realiza matemáticamente definiendo la función residuo r(u) y la función peso w(u) de tal modo que cumpla nuestro criterio de aproximación: por ejemplo, si quiero seguir un criterio de mínima distancia geométrica entre mi solución real y la aproximada, definiré mi función residuo como el cuadrado de la diferencia (utilizaría el método de mínimos cuadrados), o si quiero que mi función valga exactamente cero en una serie de puntos, utilizaré unas funciones peso tipo delta de Dirac en esos puntos (método de colocación). Al establecer este criterio de aproximación estamos acotando matemáticamente un subespacio de soluciones válido dentro del subespacio que forman mis funciones de forma. El método se encargará de "ajustar" tus resultados a las condiciones de tu problema de la forma que hayas elegido.

En el caso del método de Galerkin (uno de los métodos más usados), las funciones peso que se escogen son las propias funciones de forma del subespacio que forma mi solución. Esta condición que se establece con el método de residuos ponderados en Galerkin es totalmente lógica: yo impongo que la mejor solución aproximada de la ecuación diferencial a resolver es su PROYECCIÓN ORTOGONAL sobre mi espacio de soluciones que me he definido con las funciones de forma que he escogido (esto matemáticamente se expresa como que su producto escalar debe ser nulo) . De este modo las funciones de forma se escogen para que sean ortogonales al propio error.

Todo este proceso, que en apariencia es relativamente sencillo, a la hora de realizar una implementación numérica (me refiero a utilizar un ordenador para realizar los cálculos) resulta hartamente más complejo. La razón es la propia naturaleza de los ordenadores: no están diseñados para trabajar con variables continuas. Por ello, lo que para nosotros en la teoría es una función continua, al ordenador hay que introducirlo como una serie discreta de puntos que representen dichas funciones. Esto obliga en la práctica a utilizar una notación matricial para resolver el problema.

En el próximo artículo voy a detallar un poco más del proceso de implementación de elementos finitos en un ordenador.