ECUACIONES DIFERENCIALES
PARCIALES
….1
El orden de una EDP está dado por la derivada parcial de mayor orden que
aparece en la ecuación. Por ejemplo, el orden de la ecuación (1) es k si F , como
función de alguna de las derivadas de orden k, no es constante.
8
A continuación, se presenta ejemplos de clasificación de EDPs:
Ejemplo 1
La
ecuación
es una EDP lineal no homogénea de primer orden.
Ejemplo 2
La
ecuación
conocida como la ecuación de Korteweg y de Vries, es semi-lineal de
tercer orden.
Ejemplo 3
La ecuación
conocida como la ecuación de Burger, con viscosidad v constante, es también
semi-lineal pero de segundo orden.
Ejemplo 4
La ecuación
denominada ecuación de Sine-Gordon, es otro ejemplo de semi-lineal de
segundo orden.
Ejemplo 5
La ecuación de Poisson
es una ecuación lineal de segundo orden no homogénea si la función h no
es idénticamente nula. En el caso en que h 0 la ecuación es homogénea, y
es llamada ecuación de Laplace.
Ejemplo 6
La ecuación de calor
Ejemplo 7
La ecuación de onda:
EDPs DE SEGUNDO ORDEN CON DOS VARIABLES
INDEPENDIENTES
10
10
11
(12)
(13)
(14)
(15)
(16)
Tomando en cuenta el criterio de clasificación expuesto, analicemos las
siguientes ecuaciones:
17
donde u y v son constantes reales;
18
15
Reordenando la ecuación (17) bajo la forma de (11), se tiene:
17
Sin embargo, tal como se muestra a seguir, la ecuación (17) tiene una mezcla
de propiedades parabólicas e hiperbólicas según si el término difusivo o el
advectivo sea dominante. Así, para un proceso puramente difusivo [término
difusivo dominante], resulta la ecuación de difusión pura (ecuación 18); y para
un proceso puramente advectivo [término advectivo dominante], resulta la
ecuación de advección pura (ecuación 15).
Entonces, reordenando (18) bajo la forma de (11), resulta:
19
15
20
21
21 11
Ejercicio 1
Clasifique la ecuación
Escribimos esta ecuación como
e identificamos de esta forma los coeficientes: A = 3 , B = 0 y C = 0 . En
vista de que B2 −4 AC = 0 , la ecuación es parabólica.
Ejercicio 2
Clasifique la ecuación
Arreglamos la ecuación bajo la forma
y vemos que: A = 1, B = 0 , C = −1 y B2 − 4AC =- 4(1)( -1) >
0. Por tanto, la ecuación es hiperbólica.
donde A, B y C son funciones de x y y, y D es una función de x, y, u, ∂u/∂x y
∂u/∂y. Dependiendo de los valores de los coeficientes de los términos de la
segunda derivada (A, B y C), la ecuación (7) se clasifica en una de tres
categorías (tabla 1). Esta clasificación, que se basa en el método de las
características (por ejemplo, véase Vichnevetsky, 1981, o Lapidus y Pinder,
1981), es útil debido a que cada categoría se relaciona con problemas de
ingeniería específicos y distintos, que demandan técnicas de solución
especiales.
4
5
4 5
7
Ecuaciones en Derivadas
Parciales
Diferencias finitas
Convergencia y estabilidad
Ecuaciones hiperbólicas: ecuación de
Ondas
Ecuaciones parabólicas: ecuación del
Calor
Ecuaciones elípticas: ecuación de
Laplace
Introducción
EDP de orden 2, lineales de coeficientes
constantes.
Auxx+Buxy+Cuyy+Dux+Euy+Fu=G
• Ecuación de Laplace uxx + uyy = 0
• Ecuación del Calor ut - cuxx = 0, c>0
• Ecuación de Ondas utt - c2uxx = 0
Condiciones iniciales y de contorno
Diferencias finitas
Discretización: EDP EDF
Métodos explícitos h→
yj+1
Sencillos k
Inestables yj
u i,j
Métodos implícitos yj - 1
Más complejos
x i-1 x i x i+1
Estables
Diferencias primeras
Hacia adelante
u( x i + h, y j ) - u( x i , y j ) u i +1, j - u i , j
u x (xi , y j ) =
h h
Error
u i +1, j - u i , j h
u x (xi , y j ) = - u xx ( x i + h, y j ), 0 1
h 2
Hacia atrás
u( x i , y j ) - u( x i - h, y j ) u i , j - u i -1, j
u x (xi , y j ) =
h h
Diferencias primeras (cont.)
Diferencias simétricas
u( x i + h, y j ) - u( x i - h, y j ) u i +1, j - u i -1, j
u x (xi , y j ) =
2h 2h
Error
u( x + h, y) - u ( x - h , y) h 2
u x ( x, y) = - u xxx (, y),
2h 6
]x - h, x + h[
Diferencias segundas
Diferencias simétricas
u x ( x i + h, y j ) - u x ( x i , y j ) u i +1, j - 2 u i , j + u i -1, j
u xx ( x i , y j )
h h2
Error
u( x - h, y) - 2 u( x, y) + u ( x + h, y) h 2
u xx ( x, y) = 2 - u xxxx (, y)
k 12
]x - h, x + h[
Convergencia y estabilidad
EDP F(x,y,u)=0 Solución: ~
u ( x, y)
EDF Gi,j(h,k,u)=0, para cada (i,j) u h ,k ( xi , y j )
Convergencia u h ,k ( xi , y j ) ⎯h⎯ ⎯0→ u( x i , y j )
,k→
Consistencia G i , j ( h, k , ~
u) ⎯h⎯ ⎯0→ 0
,k→
Estabilidad: Control del error de redondeo
Consistencia + Estabilidad Convergencia
Ecuaciones elípticas
Ecuación de Laplace
uxx + uyy = 0, 0 < x < a, 0 < y <b
Condiciones de contorno
u(x,0), u(x,b), u(0,y), u(a,y)
Discretización
u i+1, j - 2 u i , j + u i-1, j u i , j+1 - 2 u i , j + u i , j-1
2 + 2 =0
h k
Ecuación de Laplace
Ecuación en diferencias: a=k/h
a2(ui-1,j + ui+1,j) + ui,j-1 + ui,j+1 - 2(a2+1)ui,j = 0
Matriz del
sistema:
grande ,
dispersa
Caso h = k : ui-1,j + ui+1,j + ui,j-1 + ui,j+1 = 4ui,j
Ec. de Laplace: Métodos
iterativos
Método de Jacobi
( ( )
u (i ,kj) = a 2 u (i -k1-,1j) + u (i +k1-,1j) + u (i ,kj--11) + u (i ,kj+-11) ) (2 + 2a 2 )
Método de Gauss-Seidel
( ( )
u(i ,kj) = a 2 u(i -k1) , j + u(i +k1-,1j) + u(i ,kj)-1 + u(i ,kj-+11) ) (2 + 2a 2 )
Criterio de parada max u (i ,kj) - u (i ,kj-1)
i, j
Método de Sobrerrelajación
Idea:
ponderar el desplazamiento de Gauss-Seidel
Pasos
( ( )
u (i ,kj) = a 2 u (i -k1) , j + u (i +k1-,1j) + u(i ,kj)-1 + u(i ,kj-+11) ) ( 2 + 2a 2 )
( k - 1)
u (k)
i, j = (1 - w ) u i, j + wui , j
(k)
Si w = 1 coincide con Gauss-Seidel
Ecuación de Laplace.
Ejemplo
uxx+ uyy=0, 0 < x < 1
0 < y < 1,
n=4 m=4,
u(x, 0) = 0 u (x, 1) = 100x
u(0, y) = 0 u(1, y) = 100y
= 0.01
Ajuste de las condiciones de contorno:
ui(,0j ) = 25 i = 1,...,m-1 j= 1,...,n-1
Método de Jacobi.
Iteraciones: 8
Operaciones en coma flotante: 1142
y = 0.25 y = 0. 5 y = 0.75
x = 0.25 6.2561 12.5031 18.7500
x = 0.5 12.5031 25.0000 37.4969
x = 0.75 18.7500 37.4969 56.2439
Método de Gauss-Seidel.
Iteraciones: 11
Operaciones en coma flotante: 1378
y = 0.25 y = 0. 5 y = 0.75
x = 0.25 6.2447 12.4947 18.7473
x = 0.5 12.4947 24.9947 37.4973
x = 0.75 18.7473 37.4973 56.2487
Ecuaciones parabólicas
Ecuación ut = cuxx, 0 < x < L, t > 0
del Calor
Condición u(x, 0) = f(x)
inicial
Condiciones u(0, t) =T0 u(L, t) = TL
de contorno
Ecuación en diferencias
ui , j+1 - ui , j ui +1, j - 2ui , j + ui -1, j
=c
k h2
Ec. del Calor: Método explícito
Condición inicial
ui,0 = f(xi)
Condiciones de contorno
u0,j = T0 unx,t = TL para j>0
Pasos siguientes
ui,j+1 = a(ui+1,j+ui-1,j) +(1-2a)ui,j
Convergencia a 1/2 Óptimo a = 1/6
Ecuación del Calor. Método explícito.
Ejemplo
Hallar la temperatura para t = 0.3 de una barra
de 1m cuyos extremos se mantienen a 20ºC y a
40ºC. La temperatura inicial de la barra es de
100ºC y el coeficiente c = 0.1. Tomar Dx = 0.2
y Dt = 0.1. Justificar la aplicabilidad del método
explícito.
ck 0.1 0.1 1 1
a= 2= =
h 0.04 4 2
Ajuste de las condiciones iniciales y de
contorno:
u0,0 = 60, u1,0 = u2,0 = u3,0 = u4,0 = 100,
u5,0 = 70
Instante t = 0.1
u1,1 = (u0,0 + u2,0)/4 + u1,0/2
= (60+100)/4 + 100/2 = 90
u2,1 = u3,1 = 100
u4,1 = (u3,0 + u5,0)/4 + u4,0/2
= (100+70)/4 + 100/2 = 92.5
Instante t = 0.2 :
u1,2 = 75 u2,2 = 97.5
u3,2 = 98.125 u4,2 = 81.25
Instante t = 0.3:
u1,3 = 66.875 u2,3 = 92.0313
u3,3 = 93.75 u4,3 = 75.1563
Ec. del Calor: Método implícito
Idea: Diferencias hacia atrás
ui , j - ui , j-1 ui +1, j - 2ui , j + ui -1, j
=c 2
k h
Pasos
(1+2a)ui,j - a(ui-1,j + ui+1,j) = ui,j-1
Convergencia
para todo a
Ecuación del Calor. Método implícito
Se verifica la condición de convergencia :
a = 1/4 < 1/2
Diagonal principal: 1 + 2a = 3/2,
Diagonales contiguas -a = -1/4.
Para t = 0.1: 32 - 14 u1,1 u1, 0 + a u0, 0
- 14 3 2 - 14 u2,1 u2, 0
=
- 4
1 3 - 1 u u
2 4
3,1 3, 0
- 3 u u +a u
1
4 2 4, 0
4 ,1 5, 0
Valores obtenidos por este método:
x = 0.2 x = 0. 4 x = 0.6 x = 0.8
t = 0.1 86.2237 97.3423 97.8301 89.6384
t = 0.2 76.3776 93.3707 94.4771 82.1718
t = 0.3 69.0598 88.8487 90.5494 76.5394
Ecuaciones hiperbólicas
Ecuación de utt = c²uxx , 0 < x < L, t > 0
Ondas
Condiciones u(x, 0) = f(x) ut(x, 0) = g(x)
iniciales
Condiciones
u(0,t) = l(t) u(L,t) = r(t)
de contorno
Ecuación en diferencias finitas
ui , j+1 - 2 ui , j + ui , j-1 u
2 i + 1, j
- 2ui , j + ui -1, j
2 =c
k h2
Ec. de Ondas: Método explícito
Condiciones iniciales
ui,0 = fi y ui,1 - ui,-1 = 2kgi
Paso 1º
ui,1 = a2 (fi-1+fi+1)/2 + (1-a2)fi + kgi
Pasos siguientes
ui,j+1 = a2(ui+1,j + ui-1,j) +2(1 - a2)ui,j - ui,j-1
Convergencia a1
Ecuación de ondas. Método explícito.
Ejemplo utt = c²uxx , 0 < x < L, t > 0
c = 1, L=T=4, nx=4, nt=8,
u(x, 0) = 2-|x-2| ut(x, 0) = 0
u(0,t) = 0 u(L,t) = 0
Condición de convergencia : a = ck 1 0 .5 1
= = 1
h 1 2
Instante t = 0:
u0,0 = f(x0) = 2 - |x0 - 2| = 2 - |0 - 2| = 0 = f(x4)
u1,0 = f(x1) = 2 - |x1 - 2| = 2 - |1 - 2| = 1 = f(x3)
u2,0 = f(x2) = 2 - |x2 - 2| = 2 - |2 - 2| = 2
Instante t=1:
ui,1 = a2·(ui-1,0+ui+1,0)/2 + (1 - a2)·ui,0 + k·g(xi)
donde a2 = 1/4, 1 - a2 = 3/4:
u1,1 = (1/4)(u0,0 + u2,0)/2 + (3/4)u1,0 =
(1/4)(0 + 2)/2 + (3/4)1 = 1 = u3,1
u2,1 = (1/4)(u1,0 + u3,0)/2 + (3/4)u2,0 =
(1/4)(1 + 1)/2 + (3/4)2 = 7/4
Aplicando la fórmula genérica
ui,j+1 = a2·(ui-1,j + ui+1,j) + 2·(1 - a2)·ui,j - ui,j-1
con lo que, para t = 1 obtenemos:
u1,2 = (1/4)(u0,1 + u2,1) + (3/2)u1,1 - u1,0
= (1/4)(0 + 7/4) + (3/2)1 - 1 = 15/16
= u3,2
u2,2 = (1/4)(u1,1 + u3,1) + (3/2)u2,1 - u2,0
= (1/4)(1 + 1) + (3/2)(7/4) - 2 = 9/8
Procediendo análogamente
x=0 x=1 x=2 x=3 x=4
t=0 0 1.0000 2.0000 1.0000 0
t = 0.5 0 1.0000 1.7500 1.0000 0
t=1 0 0.9375 1.1250 0.9375 0
t = 1.5 0 0.6875 0.4063 0.6875 0
t=2 0 0.1953 -0.1719 0.1953 0
t = 2.5 0 -0.4375 -0.5664 -0.4375 0
t=3 0 -0.9932 -0.8965 -0.9932 0
t = 3.5 0 -1.2764 -1.2749 -1.2764 0
t=4 0 -1.2401 -1.6541 -1.2401 0
Ec. de Ondas: Método implícito
Idea
ui,j+1 - 2ui,j + ui,j-1 = a2 [(ui+1,j+1 - 2ui,j+1 + ui-1,j+1)
+ (ui+1,j-1 - 2ui,j-1 + ui-1,j-1)]/2
Pasos
(1+a2)ui,j+1 - a2(ui+1,j+1 + ui-1,j+1)/2 =
2ui,j + a2(ui+1,j-1 + ui-1,j-1)/2 - (1+a2)ui,j-1
Convergencia para todo a
Algoritmo del método implícito
Truco ecuación implícita
- a2( ui-1,j-1 + ui-1,j+1)/4
+ (1 + a2)(ui,j-1 + ui,j+1)/2
- a2(ui+1,j-1 + ui+1,j+1)/4 = ui,j .
Sistema Aw = v, v = (u1,j,u2,j,...,unx-1,j)'
tridiagonal ui,j+1 = wi - ui,j-1
Factorización LU Lz = v
Uw = z
Método implícito. 58 -1 0 x1 1
16
-116 16 x1 = 4
5 -1 7
8
Resolución del sistema 0 -1 5 x 1
16 8 1
Sustitución
u1, 2 x1 u1, 0 1.9184 1 0.9184
u2, 2 = x1 - u2, 0 = 3.1837 - 2 = 1.1837
u x u 1.9184 1 0.9184
3, 2 1 3, 0
Factorización LU
0.625 0 0 1 - 0.1 0
L = - 0.0625 0.61875 0 U = 0 1 - 0.10
- 0
0 0 . 0625 0 . 6186 0 1
x=0 x=1 x=2 x=3 x=4
t=0 0 1.0000 2.0000 1.0000 0
t = 0.5 0 1.0000 1.7500 1.0000 0
t=1 0 0.9184 1.1837 0.9184 0
t = 1.5 0 0.6926 0.4824 0.6926 0
t=2 0 0.2912 -0.1699 0.2912 0
t = 2.5 0 -0.2449 -0.6647 -0.2449 0
t=3 0 -0.7996 -0.9953 -0.7996 0
t = 3.5 0 -1.2231 -1.2214 -1.2231 0
t=4 0 -1.3966 -1.3981 -1.3966 0
