Hola a todos
Seguramente a muchos se nos viene a la cabeza preguntarnos que es una ecuación hiperbólica o parabólica. Voy a intentar detallar brevemente lo que son estos conceptos como los entiendo desde mi punto de vista.
Las ecuaciones en derivadas parciales (Partial differential equations en ingles, o PDE) se pueden clasificar atendiendo a varios criterios de clasificación, pero uno de los más utilizados divide las PDEs en 3 grandes grupos: parabólicas, elípticas e hiperbólicas. El criterio concreto para dividir los 3 grupos es la presencia o no en la ecuación de ciertos términos. Considerando que la ecuación se puede escribir de esta forma genérica:
$$A·u_{xx} + 2B·u_{xy} + C·u_{yy}+D·u_x+E·u_y+F = 0$$
Las ecuaciones se englobarán en una categoría u otra según los distintos coeficientes en letras mayúsculas en la ecuación sean cero o no. No obstante, mi criterio a la hora de recordar que tipo de ecuación tengo delante cuando la veo es guiarme por 3 criterios básicos y luego asociarlos.
Un método más intuitivo de aprender a indentificar cada tipo de ecuación a su forma matemática es tomando ejemplos clásicos de física, entendiendo lo que ocurre en la realidad y asociarle una forma matemática estándar que producirá dicho efecto numéricamente. Voy a poner un ejemplo para cada tipo:
· Ecuación parabólica: Un ejemplo clásico asociado a este tipo de ecuaciones es la ecuación del calor. Esta ecuación describe como se propaga el calor en un medio: si pensamos en un frigorífico cuyas paredes interiores están a una temperatura $T_{cold}$ e introducimos un vaso de agua caliente a temperatura $T_{hot}$, la intuición nos dice que el calor fluirá del vaso a través del interior del frigorífico hacia las paredes (porque están más frías) de forma gradual y sostenida en el tiempo, de tal modo que al principio el agua perderá calor rápidamente pero conforme esté más templada tardará más tiempo en seguir enfriándose hasta que finalmente el agua alcanza la misma temperatura que la nevera, y el proceso se detiene. Vemos que en este proceso intervienen 2 variables dimensionales importantes: el tiempo y la distancia entre el vaso y las paredes (no es lo mismo que estén muy cerca que estén muy separadas del vaso). Juntando estas, la ecuación se puede escribir como:
$$ \partial_t T = k \partial_{xx} T$$
con sus condiciones de contorno escritas de forma no ortodoxa (para que sea entendible):
$$T(x = pared )= T_{cold}$$
$$T(x = vaso) = T_{hot} $$
Lo que describe la ecuación es sencillamente como varía la temperatura en el interior de la nevera al pasar el tiempo en cada punto del interior de la nevera. Por lo tanto, si hablamos de ecuaciones parabólicas en procesos físicos donde las dimensiones espaciales y temporales están involucradas, entonces podemos afirmar que describen la evolución no estacionaria (hay un término con derivada temporal) de una magnitud (en este caso la temperatura) en un dominio dado (sería la zona del espacio donde quieres conocer la distribución de temperaturas). Esta ecuación da la temperatura para cada punto del espacio y para cada instante de tiempo.
· Ecuación elíptica: Un ejemplo para este tipo de ecuaciones es la ecuación de difusión. Consideremos una tubería rectangular llena de agua donde supongamos introducimos un tinte con una concentración C_{in} a través de uno de sus extremos. Supongamos también que por el otro lado de la tubería sale una concentración de tinte C_{out}. Si no hay corriente de agua alguna de líquido (el líquido esta en reposo), la pregunta que se nos plantea es como varíará la concentración de tinte a través de la tubería. La intuición nos vuelve a decir que el tinte "difundirá" progresivamente desde un lado de la tubería hasta el otro de tal manera que al observar la variación de la concentración a lo largo de la tubería, esta variará de forma progresiva.
La ecuación que describiría esta variación será:
$$ D \partial_{xx} C= 0 T$$
con sus condiciones de contorno:
$$C(x = entrada )= C_{in}$$
$$C(x = salida) = C_{out} $$
Observad que en este caso, la ecuación NO DICE NADA SOBRE LA EVOLUCIÓN TEMPORAL DE LA CONCENTRACIÓN. Sólo da la DISTRIBUCIÓN ESPACIAL DE CONCENTRACIONES. Por lo tanto, las ecuaciones elípticas en problemas de física con dimensiones espaciales están asociadas a distribuciones estacionarias de magnitudes. De hecho, si tomamos la ecuación parabólica del calor (que es no estacionaria), y hacemos que el tiempo tienda a infinito, se podrá observar que el término temporal tiende a anularse, y la solución tiende a ser la misma que la de una ecuación elíptica.
· Ecuación hiperbólica: El ejemplo clásico que se suele poner en los libros de texto es el de la ecuación de ondas. Sin embargo, esto a priori no debe decir mucho a alguien no familiarizado con el mundo de las ondas.
Veamos un ejemplo inverso para entenderlo. Primero asumiremos que un ejemplo de una ecuación de ondas puede ser el siguiente:
$H_t + c· H_x = 0$
con condición inicial
$H(x,0) = H_0(x)$
Esta ecuación en el ejemplo que tratamos referirá a como varía en el tiempo la altura de una cuerda. Si nos fijamos en la ecuación, hay 2 magnitudes independientes: tiempo y espacio. Luego la solución debe ser función de estas. Esta ecuación expresa que la variación en el tiempo la altura de la cuerda es igual a como varía en dirección longitudinal.
Sabemos que H(x,0) es la condición inicial del sistema, por lo que satisface la ecuación. Es fácil comprobar que una solución idéntica espacialmente a esta pero trasladada en el tiempo exactamente c veces un intervalo de tiempo t satisface también la ecuación, esto es: H(x-ct, t) también satisfará la ecuación. Esto se podría entender mejor si efectivamente asumimos que la solución posee una forma de onda: si a tiempo 0 la forma de la onda es mi solución, si yo me desplazo a una velocidad EXACTAMENTE IGUAL a la que viaja la onda, esta no ha variado su forma desde mi punto de vista, por lo que en cualquier instante de tiempo también será solución. Efectivamente si me quedo estático podré observar como la cuerda se elevará y ddescenderá físicamente formando una geometría igual a la de la onda original. La solución es una onda que se propaga físicamente por la cuerda, pero también en el tiempo, ya que a tiempo suficiente da igual en que punto de la cuerda este, veré llegar la onda y la veré marcharse a una velocidad c con exactamente la misma geometría que a instante t = 0.
Espero que os sirva de ayuda.
viernes, 20 de julio de 2012
lunes, 7 de mayo de 2012
Level - Set y método adjunto: comparativa
Estimados lectores
En este post voy a mostrar los resultados que se obtienen al utilizar dos métodos distintos para dibujar los resultados de un problema de optimización: en este caso se trata de la resolución de la "ecuación de radiación reducida", o también llamada ecuación de difusión (es el caso en el que el medio es ópticamente "grueso", esto es, que la absorción es despreciable frente al "scattering":
$$ \nabla \cdot (D(r) \nabla I(r)) + \mu_a (r)I(r) = 0 \text{ in } \Omega $$
$$ \mu_a (r)I(r) - \frac{2}{\pi} D(r) \nu(r_b) \nabla I(r)= f(r) \text{ in } \partial \Omega $$
He implementado dos soluciones: por un lado el método adjunto y por otro un método level-set. Las gráficas se anexan a continuación. Para el método level set se han utilizado dos pares de valores para el campo de valores de coeficiente de absorción del medio en la malla: por un lado un valor de 0.1 para el background y 0.03 para el objeto y para otro caso distinto un valor de 0.1 para el background y de 0.17 para el objeto. Esto se ha planteado así para analizar los resultados y ver las conclusiones.
En este post voy a mostrar los resultados que se obtienen al utilizar dos métodos distintos para dibujar los resultados de un problema de optimización: en este caso se trata de la resolución de la "ecuación de radiación reducida", o también llamada ecuación de difusión (es el caso en el que el medio es ópticamente "grueso", esto es, que la absorción es despreciable frente al "scattering":
$$ \nabla \cdot (D(r) \nabla I(r)) + \mu_a (r)I(r) = 0 \text{ in } \Omega $$
$$ \mu_a (r)I(r) - \frac{2}{\pi} D(r) \nu(r_b) \nabla I(r)= f(r) \text{ in } \partial \Omega $$
He implementado dos soluciones: por un lado el método adjunto y por otro un método level-set. Las gráficas se anexan a continuación. Para el método level set se han utilizado dos pares de valores para el campo de valores de coeficiente de absorción del medio en la malla: por un lado un valor de 0.1 para el background y 0.03 para el objeto y para otro caso distinto un valor de 0.1 para el background y de 0.17 para el objeto. Esto se ha planteado así para analizar los resultados y ver las conclusiones.
El resultado del método adjunto (resultado arriba ilustrado) muestra que en el dominio estudiado hay dos objetos
diagonalmente opuestos con unos coeficientes de difusión mayores y
menores que el valor de fondo del dominio, cuyo valor es de 0.1 . La
velocidad de convergencia fue muy lenta, del orden de horas.
Este resultado puede
confirmarse al desarrollar los métodos Level-Set:
- Al utilizar como
valores de referencia 0.1 para el fondo y 0.03 para el objeto, se
puede observar que el método level-set da como resultado un contorno
cercano al objeto con valor de coeficiente de absorción bajo en el
método adjunto, confirmando la ubicación predicha por el
primer método. El resultado se obtuvo en un tiempo de convergencia
mucho más bajo que en el caso del método adjunto. El resultado se ilustra a continuación:
- Al utilizar como
valores de referencia 0.1 para el fondo y 0.17 para el objeto se
visualiza que el método level-set da como resultado final el objeto
correspondiente a un coeficiente de absorción elevado obtenido en el
método adjunto. La velocidad de convergencia fue similar al
primer caso level-set.
Estos resultados permiten
concluir que:
- El método adjunto permite obtener en tu dominio todos los objetos que estén presentes
en el mismo, sea cual sea su coeficiente de absorción.
- El resultado de los
métodos level-set mostrarán un resultado u otro dependiendo del
valor escogido para el background y para el objeto. Si estos valores
no son cogidos de forma adecuada, puede que la solución obtenida sea
incompleta. Haría falta resolver el problema varias veces con
distintos valores para cerciorarse que el resultado obtenido es
correcto.
- El método adjunto es claramente más lento en términos de convergencia que los métodos
level-set, necesitando un número de iteraciones mucho mayor. Esto
concuerda con el hecho de que este método busca soluciones de coef.
de absorción en un campo de valores a priori infinito. Sin embargo,
el método level-set busca soluciones en un campo de coeficientes de
absorción limitado a dos valores: fondo u objeto. Por lo tanto se
espera que level- set sea mucho más rápido.
- El método adjunto da una solución difuminada en el espacio, esto es, la solución no
se sitúa en un punto sino en un área más o menos extensa.
- El método level-set
permite ubicar la posición del objeto de forma mucho más precisa,
ya que el campo de soluciones sólo permite dos valores para cada
punto de la malla.
- Como resumen de estas
conclusiones se puede comentar que el método level-set es un método
mucho más eficaz para si se conoce a priori el rango de valores que
puede adoptar el campo buscado. Si no se conoce información alguna,
el método adjunto puede resultar una opción acertada debido a
que permite averiguar todos los objetos con cualquier valor dentro
del dominio.
Si se deseara una buena
precisión, la opción preferida sería utilizar el método adjunto para acotar el valor que puede adoptar el campo estudiado y
posteriormente definir la posición de los objetos de forma mucho más
precisa con el método level-set.
Los resultados, como se pueden observar, no son exactos: como todo método de optimización, tu resultado final dependerá de lo preciso que sea el algoritmo de convergencia y de la estimación inicial que se realice. No obstante, este ejemplo permite visualizar la potencia de este tipo de técnicas.
Los resultados, como se pueden observar, no son exactos: como todo método de optimización, tu resultado final dependerá de lo preciso que sea el algoritmo de convergencia y de la estimación inicial que se realice. No obstante, este ejemplo permite visualizar la potencia de este tipo de técnicas.
Un saludo
lunes, 23 de abril de 2012
Métodos Level-Set
Estimados lectores
En esta ocasión voy a escribir sobre un tema sobre el cual estoy actualmente desarollando un trabajo y que me ha parecido ciertamente util: los métodos level-set.
Este tipo de métodos se engloban dentro de los métodos matemáticos de optimización, especialmente en aquellos destinados a la resolución de problemas inversos: ejemplos clásicos de problemas inversos pueden ser la reconstrucción de imágenes o la determinación de variables de campo implícitas mediante el uso de un número de datos obtenidos mediante la resolución del problema explícito (por ejemplo conocer el coeficiente de absorción de un dominio a partir de varios estímulos en la frontera con ciertas fuentes de luz de valor conocido).
Los métodos tipo level - set fueron planteados hace tiempo por Sethian y Onsager. Plantearon que resolviendo la siguiente ecuación, se podía conocer la evolución de la frontera de una función dada $\phi$ conocido el campo de velocidades $\vec{v}$ de cada uno de los puntos de la frontera:
$ \frac{\partial \phi}{\partial t} + \vec{v} \vec{\nabla} \phi = 0 $
con
$ \vec{\nabla} \phi = \vec{n} \mid \nabla \phi \mid $
El método en si es más sencillo de lo que parece. El punto clave consiste en que dada una función F en una dimensión arbitraria, por ejemplo dimensión 2, el problema se puede generalizar utilizando una función de una dimensión adicional $\phi$ (en el ejemplo dimensión 3) y buscar el valor de mi función cuando esta dimensión adicional adquiere un valor 0 "busco el nivel cero de la función $\phi$.
$$ F(x,y,t) \longrightarrow \phi(x,y,z,t) \longrightarrow z = \phi(x,y,t)$$
$$\Gamma (x,y,t) \equiv \phi(x,y,0,t) $$
En este caso, resulta que los valores de la función que hagan cero mi función son los valores de la frontera que estoy buscando.
Una vez conocida mi función $\phi$ de dimensión adicional, y conociendo el campo de velocidades $\vec{v}$, mediante la integración de la ecuación hiperbólica anteriormente descrita puedo hallar la evolución de mi frontera.
¿Como se traduce esto a efectos prácticos de resolver un problema concreto?
Es relativamente sencillo. Dada la función $F$ de tu problema dentro del dominio estudiado, se propone una función $\phi$ tal que $\Gamma = \phi = 0$
CUALQUIER FUNCIÓN PUEDE VALER. A priori vale cualquier forma siempre y cuando cumpla que en el instante inicial $\Gamma = \phi = 0$. Un ejemplo para el caso 3D puede ser un paraboloide de revolución clásico.
$$ z=(\frac{x}{a})^2 + (\frac{y}{b})^2 $$
Una vez se propone esta función, se comienza a integrar la PDE hiperbólica propuesta arriba utilizando el campo de velocidades del dominio. En cada paso "temporal", la función $\phi$ se irá deformando, pero sólo su proyección en el plano dimensional del problema adicional es nuestra solución ($\Gamma$). Así que representando al final de cada paso iterativo $\phi = 0$, obtendremos la deformación de la frontera.
Este tipo de métodos pueden ser combinado conjuntamente con otros métodos de optimización, como el método del gradiente. Utilizando la función de optimización respectiva y ajustando el tamaño del incremento mediante un parámetro escalar variable en magnitud en cada paso, se pueden resolver problemas inversos como se ha comentado anteriormente (como una reconstrucción de imágenes).
Un saludo
En esta ocasión voy a escribir sobre un tema sobre el cual estoy actualmente desarollando un trabajo y que me ha parecido ciertamente util: los métodos level-set.
Este tipo de métodos se engloban dentro de los métodos matemáticos de optimización, especialmente en aquellos destinados a la resolución de problemas inversos: ejemplos clásicos de problemas inversos pueden ser la reconstrucción de imágenes o la determinación de variables de campo implícitas mediante el uso de un número de datos obtenidos mediante la resolución del problema explícito (por ejemplo conocer el coeficiente de absorción de un dominio a partir de varios estímulos en la frontera con ciertas fuentes de luz de valor conocido).
Los métodos tipo level - set fueron planteados hace tiempo por Sethian y Onsager. Plantearon que resolviendo la siguiente ecuación, se podía conocer la evolución de la frontera de una función dada $\phi$ conocido el campo de velocidades $\vec{v}$ de cada uno de los puntos de la frontera:
$ \frac{\partial \phi}{\partial t} + \vec{v} \vec{\nabla} \phi = 0 $
con
$ \vec{\nabla} \phi = \vec{n} \mid \nabla \phi \mid $
El método en si es más sencillo de lo que parece. El punto clave consiste en que dada una función F en una dimensión arbitraria, por ejemplo dimensión 2, el problema se puede generalizar utilizando una función de una dimensión adicional $\phi$ (en el ejemplo dimensión 3) y buscar el valor de mi función cuando esta dimensión adicional adquiere un valor 0 "busco el nivel cero de la función $\phi$.
$$ F(x,y,t) \longrightarrow \phi(x,y,z,t) \longrightarrow z = \phi(x,y,t)$$
$$\Gamma (x,y,t) \equiv \phi(x,y,0,t) $$
En este caso, resulta que los valores de la función que hagan cero mi función son los valores de la frontera que estoy buscando.
Una vez conocida mi función $\phi$ de dimensión adicional, y conociendo el campo de velocidades $\vec{v}$, mediante la integración de la ecuación hiperbólica anteriormente descrita puedo hallar la evolución de mi frontera.
¿Como se traduce esto a efectos prácticos de resolver un problema concreto?
Es relativamente sencillo. Dada la función $F$ de tu problema dentro del dominio estudiado, se propone una función $\phi$ tal que $\Gamma = \phi = 0$
CUALQUIER FUNCIÓN PUEDE VALER. A priori vale cualquier forma siempre y cuando cumpla que en el instante inicial $\Gamma = \phi = 0$. Un ejemplo para el caso 3D puede ser un paraboloide de revolución clásico.
$$ z=(\frac{x}{a})^2 + (\frac{y}{b})^2 $$
Una vez se propone esta función, se comienza a integrar la PDE hiperbólica propuesta arriba utilizando el campo de velocidades del dominio. En cada paso "temporal", la función $\phi$ se irá deformando, pero sólo su proyección en el plano dimensional del problema adicional es nuestra solución ($\Gamma$). Así que representando al final de cada paso iterativo $\phi = 0$, obtendremos la deformación de la frontera.
Este tipo de métodos pueden ser combinado conjuntamente con otros métodos de optimización, como el método del gradiente. Utilizando la función de optimización respectiva y ajustando el tamaño del incremento mediante un parámetro escalar variable en magnitud en cada paso, se pueden resolver problemas inversos como se ha comentado anteriormente (como una reconstrucción de imágenes).
Un saludo
Disculpas a todos los lectores por la tardanza
Estimados lectores
Lamento profundamente no haber podido atender mi blog como es debido. Las obligaciones personales y con mi tesis doctoral no me han permitido actualizar los contenidos ni atender las preguntas a tiempo.
Espero en los próximos días ir actualizando con nuevo contenido mi blog.
Sin más, un saludo a todos
DAVID RODRIGUEZ
Lamento profundamente no haber podido atender mi blog como es debido. Las obligaciones personales y con mi tesis doctoral no me han permitido actualizar los contenidos ni atender las preguntas a tiempo.
Espero en los próximos días ir actualizando con nuevo contenido mi blog.
Sin más, un saludo a todos
DAVID RODRIGUEZ
lunes, 6 de febrero de 2012
Computación en paralelo I
Hola a todos
En esta ocasión voy a hablar muy escuetamente de lo que va a ser el futuro de la programación científica: la computación en paralelo.
Actualmente los procesadores tienden a tener un número de núcleos mayor de uno, tal y como era hace unos años: todos estamos familiarizados con los "Pentium core 2 duo" o los modernos "Pentium I3, I5, etc". Sin embargo, por lo que a mi respecta, la formación universitaria en ingeniería sólo cubre nociones básicas de programación en serie. Es obvio que, desde un punto de vista computacional, si se poseen varios núcleos de procesamiento, el realizar una programación capaz de utilizar simultáneamente la potencia de cálculo de esos cores va a ser mucho más ventajoso.
Sin embargo, hasta hace bien poco desconocía por completo un punto de partida que me pudiera orientar para empezar a caminar en este vasto campo. Pues bien, si queréis una piedra angular de referencia, vuestra palabra mágica es MPI.
El MPI se puede considerar como un anexo a los lenguajes de programación típicos en ciencias (Fortran, C) que permite procesar un código arrancando varios "threads", o hilos de computación. Realmente lo ideal es ejecutar tantos threads como núcleos tengas en tu cluster, pero no tiene por que ser así.
Para realizar un código en paralelo, hay que pasar de una mentalidad de cálculo en serie, y comenzar a pensar como se puede distribuir todas las tareas que tiene que hacer el código: habrá tareas, como una búsqueda de elementos en un fichero o matriz, que se puede "descomponer" en partes que pueden llevar a cabo cada núcleo. En cambio habrá otras que no se puedan descomponer de modo directo, como la escritura y lectura simultánea de un fichero.
Todas estas tareas deben ser coordinadas, esto es, normalmente se programa de forma que uno de los núcleos actúa como maestro y el resto como esclavos: el primero realizará tareas no paralelizables y se encargará de enviar y recibir información que los esclavos necesiten y/o proporcionen. Los cores esclavos se encargarán de realizar la computación respectiva, enviando o recibiendo información vital para el código.
MPI es una alternativa básica a la programación tradicional. Sin embargo, otras aplicaciones están desarrollándose rápidamente, con el afán de simplificar los procesos de comunicación entre cores y la distribución del trabajo: un ejemplo clásico es el OpenMP.
En próximos posts hablaré más en detalle de MPI y OpenMP.
Un saludo
En esta ocasión voy a hablar muy escuetamente de lo que va a ser el futuro de la programación científica: la computación en paralelo.
Actualmente los procesadores tienden a tener un número de núcleos mayor de uno, tal y como era hace unos años: todos estamos familiarizados con los "Pentium core 2 duo" o los modernos "Pentium I3, I5, etc". Sin embargo, por lo que a mi respecta, la formación universitaria en ingeniería sólo cubre nociones básicas de programación en serie. Es obvio que, desde un punto de vista computacional, si se poseen varios núcleos de procesamiento, el realizar una programación capaz de utilizar simultáneamente la potencia de cálculo de esos cores va a ser mucho más ventajoso.
Sin embargo, hasta hace bien poco desconocía por completo un punto de partida que me pudiera orientar para empezar a caminar en este vasto campo. Pues bien, si queréis una piedra angular de referencia, vuestra palabra mágica es MPI.
El MPI se puede considerar como un anexo a los lenguajes de programación típicos en ciencias (Fortran, C) que permite procesar un código arrancando varios "threads", o hilos de computación. Realmente lo ideal es ejecutar tantos threads como núcleos tengas en tu cluster, pero no tiene por que ser así.
Para realizar un código en paralelo, hay que pasar de una mentalidad de cálculo en serie, y comenzar a pensar como se puede distribuir todas las tareas que tiene que hacer el código: habrá tareas, como una búsqueda de elementos en un fichero o matriz, que se puede "descomponer" en partes que pueden llevar a cabo cada núcleo. En cambio habrá otras que no se puedan descomponer de modo directo, como la escritura y lectura simultánea de un fichero.
Todas estas tareas deben ser coordinadas, esto es, normalmente se programa de forma que uno de los núcleos actúa como maestro y el resto como esclavos: el primero realizará tareas no paralelizables y se encargará de enviar y recibir información que los esclavos necesiten y/o proporcionen. Los cores esclavos se encargarán de realizar la computación respectiva, enviando o recibiendo información vital para el código.
MPI es una alternativa básica a la programación tradicional. Sin embargo, otras aplicaciones están desarrollándose rápidamente, con el afán de simplificar los procesos de comunicación entre cores y la distribución del trabajo: un ejemplo clásico es el OpenMP.
En próximos posts hablaré más en detalle de MPI y OpenMP.
Un saludo
miércoles, 1 de febrero de 2012
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.
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.
domingo, 4 de diciembre de 2011
Cambio en server de LATEX
Hola
Debido a que ha habido cambios en el server de LATEX del cual suministré su enlace, he buscado otro server para poder mostrar LATEX en tu propia web sin cambios. El nuevo script que tenéis que introducir (sustituyendo al que ya tenías en caso de que siguierais las instrucciones tal y como expliqué en el anterior post) es el siguiente:
<script type="text/javascript" src="http://cdn.mathjax.org/mathjax/latest/MathJax.js">
MathJax.Hub.Config({
extensions: ["tex2jax.js","TeX/AMSmath.js","TeX/AMSsymbols.js"],
jax: ["input/TeX", "output/HTML-CSS"],
tex2jax: {
inlineMath: [ ['$','$'], ["\\(","\\)"] ],
displayMath: [ ['$$','$$'], ["\\[","\\]"] ],
},
"HTML-CSS": { availableFonts: ["TeX"] }
});
</script>
Con esto ya podréis escribir con LATEX.
Nota:
El script esta diseñado para que si se pone un dolar a cada lado de la escritura en latex ($ texto $), se mostrará el texto en línea, esto es, sin espacios. Si se escriben dobles dólares ($$ texto $$), se hace un punto y aparte y se separa del resto del cuerpo. Ejemplos:
Esto es una prueba de este texto $f=\sqrt{x^2+y^2}$ dentro de linea.
Esto es una prueba
$$ g=\frac{\partial}{\partial x}\left( z^2+6x \right) $$
de texto separado de cuerpo
Suscribirse a:
Entradas (Atom)