You are currently browsing the tag archive for the ‘PDEs elípticas’ tag.

Una PDE cuasilineal de segundo orden en dos variables independientes x e y con función incógnita u(x,y) tiene la forma general:

a(x,y,u,u_x,u_y) u_{xx} + 2b(x,y,u,u_x,u_y) u_{xy}

+ c(x,y,u,u_x,u_y)u_{yy} + d(x,y,u,u_x,u_y) = 0,

donde a,b,c,d son funciones contínuas en un subconjunto abierto \mathcal{V} de \mathcal{U} \times \mathbb{R}^3 de las variables (x,y,u,u_x,u_y) donde \mathcal{U} es un abierto de \mathbb{R}^2.

Definimos el discriminante como

D(v^0) := a(v^0)c(v^0) - b^2(v^0),

con v^0 = (x^0,y^0,u^0,u_x^0,u_y^0) \in \mathcal{V}. Diremos que la ecuación anterior es:

  1. Elíptica en el punto v^0 si D(v^0) > 0,
  2. Parabólica en el punto v^0 si D(v^0) = 0,
  3. Hiperbólica en el punto v^0 si D(v^0) < 0.

Por tanto, el carácter elíptico, parabólico o hiperbólico depende no solo del punto (x^0,y^0) \in \mathcal{U} sino también del valor de una solución y sus derivadas parciales de primer orden en dicho punto. Además, en el caso de que la parte principal, los coeficientes que multiplican a las derivadas de segundo orden, sea de coeficientes constantes, el carácter se mantiene en todos los puntos donde esté definida la función d.

De esta manera, la ecuación de Laplace u_{tt} + u_{xx} = 0 es elíptica en todos los puntos; la ecuación del calor u_t - u_{xx} =0 es parabólica; y la ecuación de ondas u_{tt} - u_{xx} = 0 es hiperbólica.

En el caso de tener n variables independientes x_1, x_2, \ldots, x_n entonces la ecuación general tiene la forma:

a^{ij}(x_1,\ldots, x_n, u, u_{x_1},\ldots, u_{x_n}) u_{x_i, x_j} + \ldots =0,

donde a^{ij} es la parte principal y el resto son terminos de menor orden. En este caso, el carácter de la ecuación depende de la signatura de los valores propios de la matriz de coeficientes:

  1. Elíptica si los valores propios son todos positivos o todos negativos,
  2. Parabólica cuando todos los valores propios son positivos o negativos excepto uno que es zero,
  3. Hiperbólica si todos los valores propios son positivos excepto uno que es negativo o todos son negativos excepto uno que es positivo.

Finalmente, toda ecuación se puede reducir a una forma canónica, que corresponde a uno de los tres tipos clásicos: Laplace, calor u ondas.

En el libro de Análisis numérico de Burden y Faires aparecen una serie de ejemplos y ejercicios resueltos de PDEs elípticas en 2D. Vamos a intentar resolverlos utilizando Chombo…

La primera ecuación corresponde con un ejemplo y es

\frac{\partial^2}{\partial x^2}u(x,y) + \frac{\partial^2}{\partial y^2}u(x,y) = 0

en

\Omega = \{ (x,y): 0<x<0.5, 0<y<0.5\}

y cuyas condiciones en la frontera \partial \Omega son

u(x,0)=0, u(0,y)=0, u(x,0.5)=200x, u(0.5,y)=200y.

El resultado es:

ejeM1

Para el ejercicio 12.1.1

u_{xx} + u_{yy} = 4

en

\Omega = \{ (x,y): 0<x<1, 0<y<2\}

y cuyas condiciones en la frontera \partial \Omega son

u(x,0)=x^2, u(0,y)=y^2, u(x,2)=(x-2)^2, u(1,y)=(y-1)^2.

El resultado es:

ejeR1

A continuación, para el 12.1.3.a

u_{xx} + u_{yy} = 0

en

\Omega = \{ (x,y): 0<x<1, 0<y<1\}

y cuyas condiciones en la frontera \partial \Omega son

u(x,0)=0, u(0,y)=0, u(x,1)=x, u(1,y)=y.

El resultado es:

ejeR3a

En el 12.1.3.b encontramos

\Delta u = -(\cos(x+y)+\cos(x-y))

en

\Omega = \{ (x,y): 0<x<\pi, 0<y<\frac{\pi}{2}\}

y

u(x,0)=\cos x , u(0,y)=\cos y, u(x,\frac{\pi}{2})=0, u(\pi,y)=-\cos y,

obteniendo:

 ejeR3b ejeR3b2

En 1D tenemos:

u_{xx} = f(x),

con 0 < x < 1 y u(0)= u(1)=0. Particionamos el dominio en n subintervalos introduciendo los puntos x_i = ih con h=\frac{1}{n}. De esta manera, para cada uno de los puntos interiores reemplazamos las derivadas por su aproximación de segundo orden en diferencias finitas (a partir de Taylor alrededor de x_i):

\frac{\partial^2}{\partial x^2} u(x_i) = \frac{u_{i-1} -2 u_i + u_{i+1}}{h^2} - \frac{h^2}{12} \frac{\partial^4}{\partial x^4}u(\xi_i)

con \xi_i \in (x_{i-1}, x_{i+1}) por lo que:

\frac{v_{i-1}-2v_i+v_{i+1}}{h^2} = f_i.

con un error de truncamiento local del orden de O(h^2) y donde utilizamos la v para hablar de la aproximación a la u. Si lo escribimos de manera matricial, tenemos:

\frac{1}{h^2} \begin{bmatrix} -2 & 1 & & & & & \\ 1 & -2 & 1 & & & & \\ & 1 & -2 & 1 & & & \\ & & \cdot & \cdot & \cdot & & \\ & & & \cdot & \cdot & \cdot & \\ & & & & 1 & -2 & 1 \\ & & & & & 1 & -2\end{bmatrix} \begin{bmatrix} v_1 \\ v_2 \\ v_3 \\ \cdot \\ \cdot \\ v_{n-2} \\ v_{n-1} \end{bmatrix} = \begin{bmatrix} f_1 \\ f_2 \\ f_3 \\ \cdot \\ \cdot \\ f_{n-2} \\ f_{n-1} \end{bmatrix}

Existen otros esquemas de mayor orden para la aproximación, por ejemplo el esquema central de cuarto orden:

\frac{\partial^2}{\partial x^2} u(x_i) = \frac{-u_{i-2}+16u_{i-1}-30u_{i}+16u_{i+1}-u_{i+2}}{12 h^2} + \frac{h^4}{90} \frac{\partial^6}{\partial x^6}u(\xi_i)

con \xi_i \in (x_{i-2},x_{i+2}).

En 2D tenemos:

u_{xx}+u_{yy} = f(x,y),

con 0 < x < 1, 0 < y < 1 y u=0 en el cuadrado unidad. Particionamos el dominio en n_x \times n_y subintervalos introduciendo los puntos (x_i,y_j) = (ih_x, jh_y) con h_x=\frac{1}{n_x} y h_y = \frac{1}{n_y}. De esta manera, para cada uno de los puntos interiores reemplazamos las derivadas por su aproximación de segundo orden en diferencias finitas:

\frac{v_{i-1\,j}-2v_{i\,j}+v_{i+1\,j}}{h_x^2}+\frac{v_{i\,j-1}-2v_{i\,j}+v_{i\,j+1}}{h_y^2} = f_{i\,j}

con v_{0\,j} = v_{n_x\,j} = v_{i\,0} = u_{i\,n_y} = 0. Tenemos un error local de truncamiento del orden de O(h_x^2+h_y^2).

Si hacemos h_x = h_y = h nos queda una plantilla, que representa un operador que interactua localmente en la malla, como esta:

A = \frac{1}{h^2} \begin{bmatrix} & 1 & \\ 1 & -4 & 1 \\ & 1 & \end{bmatrix}

En 3D tenemos:

u_{xx} + u_{yy} + u_{zz} = f(x,y,z),

con 0 < x, y, z < 1 y u=0 en la frontera del cubo unitario. Procediento de la misma manera que antes tenemos:

\frac{v_{i-1\,j\,k}-2v_{i\,j\,k}+v_{i+1\,j\,k}}{h_x^2}+\frac{v_{i\,j-1\,k}-2v_{i\,j\,k}+v_{i\,j+1\,k}}{h_y^2} + \frac{v_{i\,j\,k-1}-2v_{i\,j\,k}+v_{i\,j\,k+1}}{h_z^2}= f_{i\,j\,k}

En este caso, la manera mas sencilla de representar el operador local es mediante el tensor:

A = \frac{1}{h^2} a_{ijk} con a_{0jk}=\begin{bmatrix} & & \\ & 1 & \\ & & \end{bmatrix}, a_{1jk}=\begin{bmatrix} & 1 & \\ 1 & -6 & 1 \\ & 1 & \end{bmatrix} y a_{2jk} = \begin{bmatrix} & & \\ & 1 & \\ & & \end{bmatrix}

si h_x = h_y = h_z = h.

Aquí está el artículo donde aparece el nuevo esquema en el que el sistema se desacopla de manera jerárquica:

(1) Conocidas las cantidades hidrodinámicas conservadas, resolver:

\Delta X^i + \frac{1}{3} \mathcal{D}^i \mathcal{D}_j X^j = 8 \pi f^{ij} S_j^*

para encontrar

\hat{A}^{ij} \approx (LX)^{ij} = \mathcal{D}^i X^j + \mathcal{D}^j X^i - \frac{2}{3} \mathcal{D}_k X^k f^{ij}.

(2) Resolver la ecuación:

\Delta \psi = -2 \pi \psi^{-1} E^{*} - \psi^{-7} \frac{ f_{il} f_{jm} \hat{A}^{lm} \hat{A}^{ij} }{8}

para encontrar \psi, donde la unicidad local ahora esta garantizada. Podemos calcular S^* de manear consistente.

(3) Resolver la ecuación:

\Delta(\psi N) = 2 \pi N \psi^{-1} (E^* + 2 S^*) + N \psi^{-7} \frac{7 f_{il} f_{jm} \hat{A}^{lm} \hat{A}^{ij} }{8}

para N \psi, una ecuación lineal donde podemos aplicar el principio del máximo con lo que, con las codiciones de contorno apropiadas, se sigue la unicidad y existencia.

(4) Finalmente, resolver:

\Delta \beta^i + \frac{1}{3} \mathcal{D}^i ( \mathcal{D}_j \beta^j ) = D_j( 2 N \psi^{-6} \hat{A}^{ij} )

para encontrar \beta^i.

Además, en este otro artículo, presentan una manera de reducir una ecuación elíptica vectorial, un complicado sistema acoplado de PDEs, a un conjunto de ecuaciones Poisson escalares desacopladas. Para el caso del shift, la \beta anterior, por ejemplo, en coordenadas esféricas, tendríamos:

(1) Resolver ecuación:

\Delta \mu = \mu_S

que corresponde a la parte toroidal, para la resolución de la parte angular se introducen un potencial poloidal \eta y  un potencial toroidal \mu de manera que \boldsymbol{\beta} = , y está desacoplada del resto para obtener \mu.

(2) Resolver la también desacoplada ecuación para la divergencia (de \boldsymbol{\beta} respecto de la conexión plana \mathcal{D}):

\Delta \Theta = \frac{3}{4} \mathcal{D}_{\hat{k}} S(\boldsymbol{\beta}^{\hat{k}}).

(3) Obtener \beta^r a partir de una de las siguiente ecuaciones:

(i) \frac{\partial^2 \beta^r}{\partial r^2} + \frac{4}{r} \frac{\partial \beta^r}{\partial r} + \frac{2 \beta^r}{r^2} + \frac{1}{r^2}\Delta_{\theta \varphi} \beta^r = S(\boldsymbol{\beta})^r - \frac{1}{3} \frac{\partial \Theta}{\partial r} + \frac{2}{r} \Theta

(ii) \Delta \chi = r S(\boldsymbol{\beta})^r - \frac{r}{3} \frac{\partial \Theta}{\partial r} + 2 \Theta , donde \chi = r \beta^r

(4) Deducir \eta de una de las siguientes ecuaciones:

(i) \Delta_{\theta \varphi} \eta = r \Theta - r \frac{\partial \beta^r}{\delta r} - 2 \beta^r, que tiene la ventaja de que solo requiere una división por -l (l+1) de los coeficientes de la expansión por armónicos esféricos pero la desventaja de que utiliza la derivada radial de \beta^r que puede tener problemas con el orden.

(ii) \Delta \eta = \eta_S - \frac{2 \beta^r}{r^2} - \frac{1}{3} \frac{\Theta}{r}, que requiere la resolución de otra ecuación de Poisson adicional.

agosto 2017
L M X J V S D
« Jul    
 123456
78910111213
14151617181920
21222324252627
28293031