Resolución numérica del problema de Bratu: diferencias finitas y métodos de continuación

Traballo Fin de Grao Resolución numérica del problema de Bratu: diferencias nitas y métodos de continuación Mateo Pérez Rodríguez 2019/2020 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA GRAO DE MATEMÁTICAS Traballo Fin de Grao Resolución numérica del problema de Bratu: diferencias nitas y métodos de continuación Mateo Pérez Rodríguez Junio, 2020 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Trabajo propuesto Área de Conocimiento: Matemática Aplicada Título: Resolución numérica del problema de Bratu: diferencias - nitas y métodos de continuación Breve descripción del contenido El problema de Bratu es una ecuación diferencial elíptica con una no linealidad de forma exponencial que aparece, entre otros, en modelos de ignición. Dependiendo de los valores del coeciente del término no lineal, la ecuación puede tener dos soluciones, una o ninguna. Por tanto, la completa resolución del problema conlleva al cálculo de las soluciones dependientes de dicho parámetro. La resolución numérica se abordará mediante la discretización con diferencias nitas del operador diferencial combinada con un método de continuación, para lo cual será necesario introducir una coordenada curvilínea, denida mediante una ecuación escalar, que se añadirá al sistema discretizado del modelo. Recomendaciones Haber cursado las materias de Métodos Numéricos en Optimización y Ecuaciones Diferenciales y estar cursando Análisis Numérico de Ecuaciones en Derivadas Parciales iii Índice general Resumen viii Introducción xi 1. Discretización numérica del problema de Bratu 1 1.1. Resolución del problema de Bratu mediante el método de Newton . . . . . . 3 2. Parametrización de las curvas de soluciones de problemas no lineales dependientes de un parámetro 5 2.1. Parametrización natural . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.2. Parametrización en un entorno de un punto singular simple . . . . . . . . . 8 2.3. Parametrización en un entorno de un punto de retorno . . . . . . . . . . . . 11 2.4. Parametrización por longitud de arco . . . . . . . . . . . . . . . . . . . . . . 14 3. Aproximación numérica de las curvas de soluciones 17 3.1. Aproximación numérica en un entorno de un punto regular . . . . . . . . . . 17 3.2. Aproximación numérica en un entorno de un punto singular simple . . . . . 21 3.3. Aproximación numérica en un entorno de un punto de retorno . . . . . . . . 23 3.3.1. Puntos de retorno no degenerados . . . . . . . . . . . . . . . . . . . 25 4. Detección de puntos de retorno 29 4.1. Funcionestest................................... 29 4.1.1. Determinantes .............................. 30 4.1.2. Matrices ampliadas . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 5. Resultados numéricos 41 5.1. Solución analítica del problema de Bratu . . . . . . . . . . . . . . . . . . . . 41 5.2. Solución numérica del problema de Bratu mediante diferencias nitas . . . . 43 5.3. Tablasdeerror .................................. 44 v 5.3.1. Error en norma L∞(0,1) ......................... 45 5.3.2. Error en norma L2(0,1) ......................... 46 5.3.3. Error en norma H1 0 ............................ 47 5.4. Aproximación numérica de la autofunción en un entorno del punto de retorno 48 Conclusiones y nuevas líneas de investigación 49 5.5. Conclusiones de la memoria . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 5.6. Problema de Bratu en dimensión superior . . . . . . . . . . . . . . . . . . . 52 5.7. Discretizaciones alternativas del problema de Bratu . . . . . . . . . . . . . . 52 5.8. Aproximación numérica en un entorno de un punto de bifurcación simple . . 53 Bibliografía 55 vi xiv Capítulo 1 Discretización numérica del problema de Bratu Consideremos el problema de contorno no lineal (−u00(x) = λeu(x)∀x∈[0,1], u(0) = 0 = u(1). (1.1) Con objetivo de discretizar el dominio, denamos los nodos x0, . . . , xn+1 de manera que 0 = x0< x1< ...... < xn< xn+1 = 1, que además supondremos equiespaciados, esto es xi=i h , i = 0, ..., n + 1, siendo h=1 n+ 1 . Sea u∈ C4([0,1]) . Empleando la fórmula de derivación numérica centrada u00(xi) = u(xi+h)−2u(xi) + u(xi−h) h2+O(h2), se obtiene la siguiente aproximación de la derivada segunda, de orden dos, u00(xi)≈u(xi+h)−2u(xi) + u(xi−h) h2. Esto sugiere la siguiente aproximación de la ecuación diferencial (1.1) en el punto xi −ui+1 + 2ui−ui−1 h2−λeui= 0, 1 donde ui+1 , ui y ui−1 representan las aproximaciones de u(xi+h) , u(xi) y u(xi−h) , respectivamente. Por tanto, como se tienen las condiciones de contorno homogéneas u(0) = u0= 0 = un+1 =u(1) , tomamos uh= (u1, . . . , un)T , de manera que el problema original se puede discretizar y escribir en notación matricial como A0huh−λeuh=O, (1.2) donde A0h=1 h2          2−1· · · 0 0 −1 2 · · · 0 0 . . .. . ..... . . 0 0 · · · 2−1 0 0 · · · − 1 2          ∈ Mn×n (1.3) y euh= (eu1, eu2, . . . , eun−1, eun)T∈Rn. De esta forma, si se dene la función Fh:Rn−→ Rn Z= (z1,· · · , zn)T7−→ Fh(Z) = A0hZ−λeZ, es claro que uh es solución del sistema n -dimensional no lineal Fh(uh) = O. (1.4) Nótese que el orden de aproximación del esquema numérico construido es dos, como se recoge en el siguiente resultado: Teorema 1.1. Sea g una función de clase C2([0,1]) y supongamos que la solución u de la ecuación lineal −u00 =g con condiciones de contorno u(0) = α y u(1) = β es una función de clase C4([0,1]) . Entonces kuexa −uhk∞=m´ax i=1,...,n |u(xi)−ui| ≤ K sup x∈[0,1] u(4)(x)h2, siendo uexa = (u(x1), . . . , u(xn)) ∈Rn y K una constante independiente de h . Demostración. Véase [CMT89]. 2 1.1. Resolución del problema de Bratu mediante el método de Newton Una vez discretizada la ecuación diferencial (1.1), el siguiente paso será resolver la ecuación (1.4) haciendo uso de un método iterativo, como es el método de Newton vectorial. Formalmente, se construye el siguiente proceso iterativo (u(0) hdado u(k) h=u(k−1) h−DFh(u(k−1) h)−1Fh(u(k−1) h), k = 1,2,3,· · · donde, en nuestro caso DFh(u(k) h) = A0h−λ       eu(k) 10· · · 0 0 0eu(k) 2· · · 0 0 . . .. . ..... . . 0 0 · · · 0eu(k) n        . Sin embargo, a la hora de implementar este método en un ordenador, es preferible, en términos de eciencia, evitar el cálculo de la inversa de DFh(u(k−1) h) resolviendo un único sistema, dado por (u(0) hdado DFh(u(k−1) h)∆u(k) h=−Fh(u(k−1) h), k = 1,2,3,· · · donde ∆u(k) h=u(k) h−u(k−1) h . Como test de parada se puede emplear ||∆u(k) h||∞≤ε con ε > 0 previamente jado. 3 4 Capítulo 2 Parametrización de las curvas de soluciones de problemas no lineales dependientes de un parámetro Sea F:R×Rn−→ Rn una aplicación de clase C∞ . Se considera la ecuación F(λ, u) = O, (2.1) donde u∈Rn es un vector que describe el estado de un sistema denido por el parámetro λ∈R . Con objeto de detallar las distintas situaciones que se analizan y se resuelven numéricamente en este trabajo, introducimos las siguientes deniciones (véase [MEI13], [BRR80] y [BRR81a]). Denición 2.1. Sea F:Rn+1 −→ Rn una aplicación de clase C∞ . Un punto g∈Rn se dice valor regular de F si rang (DF(λ, u))) = n , ∀(λ, u)∈F−1(g), siendo DF (λ, u) : Rn+1 −→ Rn la aplicación diferencial de F en (λ, u) . Un punto g∈Rn que no es valor regular de F se dice valor singular . Observación 2.2 . Nótese que en la denición anterior se entiende que rang(DF(λ, u)) = dim(Im(DF(λ, u))) . Por tanto, un punto g∈Rn se dirá valor regular de F si la aplicación DF (λ, u) : Rn+1 −→ Rn es sobreyectiva, para todo (λ, u)∈F−1(g) . 5 Denición 2.3. Sea F:Rn+1 −→ Rn una aplicación de clase C∞ . Diremos que el punto (λ0, u0)∈Rn+1 es un punto regular de F si la aplicación DuF0:= DuF(λ0, u0) : Rn−→ Rn es un isomorsmo. Denición 2.4. Un punto (λ0, u0)∈Rn+1 se dice punto singular de F si rang(DF0)≤n−1. En particular, diremos que (λ0, u0) es un punto singular simple de F si: 1. F(λ0, u0) = O , 2. DuF0∈L(Rn,Rn) es singular, con rang(DuF0) = n−1 . (λ0, u0) es un punto de retorno (o punto límite simple ) de F si: 1. F(λ0, u0) = O , 2. DuF0∈L(Rn,Rn) es singular, con rang(DuF0) = n−1 , 3. DλF0:= DλF(λ0, u0)/∈Im(DuF0) . Observación 2.5 . Como veremos más adelante, los puntos regulares verican las hipótesis del teorema de la función implícita, de manera que podremos obtener en estos puntos u como función explícita de λ , a diferencia de lo que sucede con los puntos singulares simples y los puntos de retorno. Nótese además que tanto los puntos regulares como los puntos de retorno pueden aplicarse en valores regulares de F . Observación 2.6 . Las anteriores deniciones, así como todos los resultados que se presentan a lo largo de este Capítulo, se aplican de manera inmediata al problema de Bratu sin más que tomar F(λ, u) = −u00 −λeu. En este trabajo consideraremos únicamente el caso en el que O ∈ Rn es un valor regular de F . El objetivo será estudiar procedimientos numéricos para el cálculo aproximado de las familias de soluciones de (2.1). 6 Denotemos por M0:= (λ, u)∈Rn+1 :F(λ, u) = O, rang(DF(λ, u)) = n (2.2) al conjunto de soluciones de (2.1) para las cuales la matriz jacobiana DF (λ, u) tiene rango máximo. En este capítulo, discutiremos cómo parametrizar el conjunto M0 . 2.1. Parametrización natural Consideremos el problema uniparamétrico con variables (λ, u)∈Rn+1 , donde u representa una variable de estado y λ su parámetro de control. Dado un punto regular (λ0, u0) solución de (2.1), el teorema de la función implícita (que se puede consultar en [ROD03]) establece las condiciones bajo las cuales λ puede ser utilizado como parámetro de la curva de soluciones que pasa por (λ0, u0) en un entorno de dicho punto. En otras palabras, establece cuándo el conjunto de soluciones M0 denido en (2.2) puede ser parametrizado localmente respecto al parámetro λ . Lema 2.7. Sea (λ0, u0)∈Rn+1 tal que F(λ0, u0) = O y rang(DuF0) = n . Entonces existe una constante maximal δ > 0 y una única función u:λ0−δ, λ0+δ−→ Rn tal que (λ, u (λ)) : λ−λ0< δ , u λ0=u0⊂ M0. (2.3) Observación 2.8 . El lema anterior extiende la parametrización a ambos lados de λ0 hasta que DF (λ, u(λ)) deje de ser inversible. Notemos además que, como hemos tomado δ de forma que sea maximal, el rango de DuF(λ+δ, u(λ+δ)) o DuF(λ−δ, u(λ−δ)) es menor que n . Nótese que, si O ∈ Rn es un valor regular de F y en (λ0, u0) se produce una deciencia del rango de DuF0 , entonces DuF0:Rn−→ Rn deja de ser un isomorsmo y el lema anterior ya no se aplica. En este caso, basta cambiar el parámetro λ con otra componente de u para conseguir una parametrización de la curva de soluciones en un entorno de (λ0, u0) . En los dos apartados siguientes veremos cómo realizar una parametrización global de la curva de soluciones de (2.1) cuando rang(DuF0)< n . 7 2.2. Parametrización en un entorno de un punto singular simple Dada una aplicación F:Rn+1 −→ Rn de clase C∞ , el objetivo de este apartado es estudiar las soluciones de la ecuación F(λ, u) = O (2.4) en un entorno de un punto singular simple (λ0, u0) , esto es, de un punto (λ0, u0)∈Rn+1 que verique las dos siguientes condiciones: 1. F(λ0, u0) = O , 2. DuF0 es singular, con rang(DuF0) = n−1 (o, equivalentemente, singular y tal que 0 es un autovalor simple de DuF0 ). En esta sección seguiremos el formalismo de [BRR81a], adaptado al caso de dimensión nita. Denotemos por k·k a la norma euclídea en Rn y consideremos el producto de dualidad h·,·i :Rn×Rn→R . A continuación presentamos un resultado que será de gran utilidad a lo largo de todo el capítulo. Lema 2.9. Sea F:R×Rn−→ Rn una aplicación de clase C∞ tal que rang(DuF0) = n−1 . Entonces existen ϕ0, ϕ0∗∈Rn tales que DuF0ϕ0= 0 ,kϕ0k= 1, (DuF0)∗ϕ0∗= 0 ,ϕ0, ϕ0∗= 1, (2.5) siendo (DuF0)∗≡DuF0T la matriz adjunta de DuF0 . Además, si denotamos por V1≡Ker(DuF0) = ϕ0, V2≡Im(DuF0) = w∈Rn:w, ϕ0∗= 0, (2.6) entonces se verica que Rn=V1⊕V2 (2.7) y DuF0V2 es un isomorsmo de V2 . Observación 2.10 . Nótese que en el problema de Bratu el operador DuF0 es autoadjunto (esto es, verica que DuF0= (DuF0)∗ ). Por tanto, en el teorema anterior se tiene que ϕ0=ϕ0∗ y, además, la suma (2.7) es ortogonal. 8 Denotemos por L=DuF0V2−1 a la inversa del isomorsmo DuF0V2:V2−→ V2 . Si denimos el operador proyección Q:Rn−→ V2 como Q(v) = v−v, ϕ0∗ϕ0, v ∈Rn, (2.8) entonces la ecuación (2.1) es equivalente al sistema (Q F(λ, u) = O ∈ V2, (I−Q)F(λ, u) = O ∈ V1, (2.9) donde (I−Q) : Rn−→ V1 . De lo anterior se sigue que, dado u∈Rn , si (λ0, u0) es un punto singular simple, entonces existe una única descomposición de la forma u=u0+α ϕ0+v , α ∈R, v ∈V2, ya que Rn=ϕ0⊕V2 . Por tanto, si escribimos λ=λ0+ξ , ξ ∈R, entonces el sistema (2.9) se puede escribir como (F(ξ, α, v) = O ∈ V2, (I−Q)F(λ0+ξ, u0+α ϕ0+v) = O ∈ V1 o, equivalentemente, (F(ξ, α, v) = O ∈ V2, F(λ0+ξ, u0+αϕ0+v), ϕ0∗= 0 ∈R, (2.10) siendo F:R2×V2−→ V2 (ξ, α, v)7−→ F(ξ, α, v) = Q F(λ0+ξ, u0+α ϕ0+v) una función que verica: 1. F∈ C∞(R2×V2) , 2. F(0,0,O) = Q F(λ0, u0) = Q(O) = O , 3. DvF(0,0,O) = DuF0V2∈L(V2, V2) es un isomorsmo (por el lema 2.9). Por tanto, del teorema de la función implícita (véase [ROD03]) se sigue el siguiente lema. 9 respecto de s , obtenemos DλF(λ(s), u(s)) ˙ λ(s) + DuF(λ(s), u(s)) ˙u(s) = DF (λ(s), u(s)) ˙ λ(s),˙u(s)=O, esto es, la derivada de la parametrización longitud de arco es el vector tangente orientado, ˙ λ(s),˙u(s) . Además, por construcción en (2.22), se verica que k˙u(s)k2+|˙ λ(s)|2= 1. En las secciones anteriores se ha mostrado que un parámetro longitud de arco puede ser tomado en la dirección de ξ (o, equivalentemente, de λ ) en el caso de un punto regular; y en la dirección de α en el caso de un punto de retorno no degenerado. 16 Capítulo 3 Aproximación numérica de las curvas de soluciones En el capítulo anterior se ha estudiado la parametrización de las curvas de soluciones del problema continuo de dimensión nita F(λ, u) = O, λ ∈R, u ∈Rn, siendo F:R×Rn−→ Rn una aplicación de clase C∞ , en un entorno de soluciones regulares, de puntos singulares simples y de puntos de retorno no degenerados. En particular, estos resultados se aplican también a las aproximaciones nito-dimensionales construidas con esquemas de diferencias nitas de ecuaciones diferenciales. En los artículos ya clásicos de Brezzi, Rappaz y Raviart ([BRR80], [BRR81a] y [BRR81b]), se analiza el orden de aproximación que se obtiene al aproximar una curva de soluciones del problema continuo por la curva de soluciones del problema discreto. En este capítulo se presentan dichos resultados de aproximación numérica de problemas no lineales adaptados al marco funcional adecuado para el análisis numérico del problema de Bratu. 3.1. Aproximación numérica en un entorno de un punto regular Denamos el conjunto H1 0(0,1) := g: (0,1) −→ R:g∈H1(0,1) y g(0) = 0 = g(1), donde H1(0,1) = g∈L2(0,1) : g0∈L2(0,1) . 17 Se tiene que H1 0(0,1) es un espacio de Hilbert (y, por tanto, de Banach) de dimensión innita, dotado del producto escalar hu, vi=Z1 0 u(x)v(x)dx +Z1 0 u0(x)v0(x)dx, donde u0 y v0 denotan las derivadas débiles 1 de u y v , respectivamente. Nótese que dicho producto escalar induce la norma kukH1 0(0,1) =Z1 0 u(x)2dx +Z1 0 u0(x)2dx1/2 . A lo largo de esta sección emplearemos la notación H1 0≡H1 0(0,1) . Sean entonces G:R×H1 0−→ H1 0 una aplicación de clase C∞ y T∈L(H1 0, H1 0) un operador lineal tales que F(λ, u) = u+TG(λ, u)∈L(R×H1 0, H1 0), λ ∈R, u ∈H1 0 y consideremos el problema no lineal F(λ, u) = O. (3.1) Denotemos por DlF(λ, u)∈Ll(R×H1 0, H1 0) , l≥2 , a la derivada l− ésima de F . Consideremos una malla x0, x1, . . . , xn+1 vericando 0 = x0< x1< . . . < xn+1 = 1, que, además, supondremos uniforme (esto es, tal que xi−xi−1=h , para i= 1,...n+ 1 ). Para cada valor del parámetro h > 0 , que tenderá a cero, consideremos el subespacio de dimensión nita Vh⊂H1 0 , formado por todas las funciones anes a trozos en [0,1] Vh=nvh∈ C0([0,1]) : vh|[xi−1,xi]∈P1, i = 1, . . . , n + 1 ; vh(0) = 0 = vh(1)o donde por P1 denotamos al espacio de polinomios de grado menor o igual que uno. Nótese que cada elemento de Vh queda unívocamente determinado por los valores en los nodos xi . 1 El concepto de derivada débil o derivada en el sentido de las distribuciones constituye una generalización del concepto de derivada, en el sentido de que es aplicable a funciones no diferenciables pero sí integrables (véase [BRE10]). Más concretamente, dada una función u∈L1(0,1) , se dice que u0∈L1(0,1) es una derivada débil de u si y solo si R1 0u(t)ψ0(t)dt =−R1 0u0(t)ψ(t)dt , para todo ψ∈ C∞ . 18 Sea Th∈L(H1 0, Vh) un operador lineal y denamos la aplicación Fh(λ, u) = u+ThG(λ, uh), λ ∈R, u ∈H1 0. El problema aproximado consiste en resolver la ecuación Fh(λ, uh) = O, (3.2) esto es, en encontrar los pares (λ, uh)∈R×Vh soluciones de (3.2). Nótese que esto equivale a resolver la ecuación (3.2) en R×H1 0 . A continuación enunciamos un resultado previo que será de gran utilidad a lo largo de todo este capítulo. Teorema 3.1. Supongamos que u:λ∈R−→ u(λ)∈H1 0 es una función acotada de clase Cr tal que, para cada λ∈R , verica: 1. F(λ, u(λ)) = O , 2. DuF(λ, u(λ)) : H1 0−→ H1 0 es un isomorsmo tal que  DuF(λ, u(λ))−1 L(H1 0,H1 0)≤c. Para cada valor del parámetro h , sea Fh:R×H1 0−→ Vh⊂H1 0 una aplicación de clase Cr tal que, para todo subconjunto acotado B⊂R×H1 0 , verique: 1. l´ım h→0"sup (λ,u)∈B  DlF(λ, u)−DlFh(λ, u)  Ll(R×H1 0,H1 0)#= 0,0≤l < r, 2. sup (λ,u)∈B kDrFhkLr(R×H1 0,H1 0)≤˜c, con ˜c independiente de h. Entonces existen dos constantes a , h0∈R , h0>0 , y, para h≤h0 sucientemente pequeño, una única aplicación uh:R−→ Vh⊂H1 0 de clase Cr tal que, para todo λ∈R se tiene 1. Fh(λ, uh(λ)) = O, 2. kuh(λ)−u(λ)kH1 0≤a. Además, para todo λ, λ∗∈R y todo m∈Z tal que 0≤m<r se tiene la siguiente cota de error 19 1. kDmuh(λ∗)−Dmu(λ)kLm(R,H1 0)≤ ≤K(|λ∗−λ|R+ m X l=0     dl dλl(F(λ, u(λ)) −Fh(λ, u(λ)))   Ll(R,H1 0)), 2. sup λ∈R kDruhkLr(R,H1 0)≤K, donde Dmuh y Dmu denotan las derivadas m -ésimas de uh y u , respectivamente, y K > 0 es una constante independiente de h . Este resultado se aplica directamente a las parametrizaciones de las soluciones en un entorno de un punto regular del problema de Bratu, sin más que tomar G(λ, u) = λeu, T ∈L(H1 0, H1 0), siendo T tal que Tg =w es solución de de la ecuación w” = g . En efecto, se tiene que −u00 =λeu=⇒ −u=Tλeu=⇒u+Tλeu=O, (3.3) luego, efectivamente, basta tomar F(λ, u) = u+Tλeu. Además, de (3.3) se sigue que la solución de problema de Bratu es u=−Tλeu. En cuanto al problema aproximado, de manera análoga se tiene Fh(λ, uh) = uh+Thλeuh. En este caso, el operador Th es, precisamente, la solución del problema aproximado con diferencias nitas, uh=−Thλeuh. Observación 3.2 . Nótese que, así denido, T es un operador compacto. En efecto, la aplicación T:L2(0,1) −→ H1 0 es lineal y continua, y la inclusión H1 0,→L2(0,1) es compacta (véase [BRE10]), por lo que la composición es un operador lineal compacto. 20 3.2. Aproximación numérica en un entorno de un punto singular simple En esta sección estudiaremos la aproximación numérica de la ecuación F(λ, u) = O (3.4) en un entorno de un punto singular simple (λ0, u0) . Sea entonces (λ0, u0) un punto singular simple de F , esto es, tal que 1. F(λ0, u0) = O, 2. DuF0=I+TDuG0∈L(H1 0, H1 0) tiene a 0 como autovalor simple (o, equivalentemente, −1 es autovalor simple del operador compacto TDuG0 ), siendo ϕ0 la autofunción asociada. Observación 3.3 . Nótese que, en el problema de Bratu, el operador DuF0=I+TDuG0 es autoadjunto. En efecto, consideremos la notación Tf =v y Tw =u , donde v” = f y u” = w (ya que, como hemos dicho anteriormente, así es como se dene el operador T en el problema de Bratu). Entonces hTw, fi=hu, v”i=Z1 0 u(x)v”(x)dx. Integrando por partes dos veces y aplicando condiciones de contorno Dirichlet homogéneas (pues u, v ∈H1 0 ) resulta Z1 0 u”(x)v(x)dx =hu”, vi=hw, T ∗fi. Basta recordar, nalmente, que la denición de operador adjunto es hTw, fi=hw, T∗fi , lo cual completa la demostración. Por tanto, ϕ0=ϕ0∗ es también la autofunción asociada de DuF0∗ . Observación 3.4 . Nótese que podemos escribir la siguiente descomposición del espacio H1 0 H1 0=ϕ0⊕ϕ0⊥. 21 Al igual que sucedía en el caso continuo, el problema aproximado Fh(λ, uh) = O, (3.5) es equivalente a sistema (Q Fh(λ, uh) = Fh(λ, uh)−Fh(λ, uh), ϕ0∗ϕ0=O ∈ ϕ0⊥, (I−Q)Fh(λ, u) = Fh(λ, uh), ϕ0∗ϕ0=O ∈ ϕ0, (3.6) siendo Q:Rn−→ ϕ0⊥ el operador proyección denido en (2.8). Si escribimos λ=λ0+ξ , ξ ∈R, uh=u0+α ϕ0+vh, α ∈R, vh∈V2:= ϕ0⊥, entonces el sistema (3.6) se puede escribir como (Fh(ξ, α, vh) = O (I−Q)Fh(λ0+ξ, u0+α ϕ0+vh) = O ∈ ϕ0, (3.7) siendo Fh la función de clase C∞ denida como Fh:R2×V2−→ V2 (ξ, α, v)7−→ Fh(ξ, α, v) = Q Fh(λ0+ξ, u0+α ϕ0+v). Consideremos, nalmente, la aplicación J∈ C∞ denida por J: (−ξ0, ξ0)×(−α0, α0)−→ H1 0 (ξ, α)7−→ J(ξ, α) = G(λ0+ξ, u0+α ϕ0+v(ξ, α)), donde v:R2−→ V2 denota la función proporcionada por el lema (2.11) del capítulo anterior. El siguiente resultado de [BRR81a] establece la existencia de una familia de soluciones del problema discretizado (3.5). Teorema 3.5. Sea (λ0, u0) un punto singular simple de F y supongamos que G es una aplicación de clase Cp (con p > 1 ) tal que DpG está acotada. Si se verica que l´ım h→0kT−ThkL(H1 0,H1 0)= 0, entonces existen tres constantes positivas ξ0 , α0 , a∈R y, para h≤h0 sucientemente pequeño, una única aplicación vh: (−ξ0, ξ0)×(−α0, α0)−→ V2 de clase Cp tal que 1. Fh(ξ, α, vh(ξ, α)) = O, 2. kvh(ξ, α)−v(ξ, α)kH1 0≤a, |ξ|< ξ0,|α|< α0. 22 Además, existe una constante K > 0 independiente de h tal que, para todo ξ , ξ∗∈(−ξ0, ξ0) , todo α , α∗∈(−α0, α0) y todo m∈Z tal que 0≤m<p , se tiene la siguiente cota de error 1. kDmvh(ξ∗, α∗)−Dmv(ξ, α)kLm(R2,H1 0)≤ ≤K(|ξ∗−ξ|+|α∗−α|+ m X l=0   (T−Th)DlJ(ξ, α)  Ll(R2,H1 0)), 2. kDpvh(ξ∗, α∗)kLp(R2,H1 0)≤K. Observación 3.6 . Nótese que en el problema de Bratu se tiene que G:R×H1 0−→ H1 0 es de clase C∞ y, por tanto, el teorema anterior proporciona una aproximación de las derivadas parciales de v de todos los órdenes. Gracias a este teorema queda probado que, al igual que sucedía en el caso continuo, resolver el problema aproximado (3.5) en un entorno de un punto singular simple (λ0, u0) es equivalente a resolver la ecuación de bifurcación aproximada fh(ξ, α) = Fh(λ0+ξ, u0+αϕ0+vh(α, ξ)), ϕ0∗= 0 en un entorno del origen. 3.3. Aproximación numérica en un entorno de un punto de retorno Para nalizar este capítulo, estudiaremos la aproximación numérica del problema continuo F(λ, u) = O (3.8) en un entorno de un punto de retorno (λ0, u0) . Sea entonces (λ0, u0) un punto de retorno de F , esto es, un punto vericando 1. F(λ0, u0) = O, 2. DuF0=I+TDuG0∈L(H1 0, H1 0) tiene a 0 como autovalor simple (o, equivalentemente, −1 es autovalor simple del operador compacto TDuG0 ), siendo ϕ0=ϕ0∗ la autofunción asociada (pues, como ya se vió en la observación 3.3, el operador DuF del problema de Bratu es autoadjunto). 3. DλF0/∈Im(DuF0) o, equivalentemente, DλF0, ϕ0∗6= 0 . 23 Se ha visto en el capítulo anterior que ∂f ∂ξ (0,0) 6= 0 , lo cual permite construir una curva de soluciones (ξ(α), α) , obtenida mediante el teorema de la función implícita, sin más que resolver la ecuación de bifurcación f(ξ(α), α) = F(λ0+ξ(α), u0+αϕ0+v(ξ(α), α)), ϕ0∗= 0. A continuación enunciamos un resultado que establece la existencia de una rama de soluciones del problema aproximado Fh(λ, uh) = O (3.9) en un entorno de la rama de soluciones (λ(α), u(α)) : |α| ≤ α0 del problema continuo (3.8), para un h≤h0 sucientemente pequeño. Lema 3.7. Sea (λ0, u0) un punto de retorno de F . Entonces existen dos constantes positivas α0 , b∈R y, para h≤h0 sucientemente pequeño, una única aplicación ξh:α∈ (−α0, α0)−→ ξh(α)∈R de clase Cp tal que, para |α|< α0 , se tiene 1. fh(ξh(α), α)=0 , 2. |ξh(α)−ξ(α)| ≤ b. Además, existe una constante K > 0 independiente de h tal que, para todo α∈(−α0, α0) y todo entero m∈Z tal que 0≤m<p , se tiene la siguiente estimación del error 1.  dm dαm(ξh(α)−ξ(α)) ≤K m X l=0     (T−Th)dl dαlG(λ(α), u(α))   H1 0 , 2.  dp dαpξh(α) ≤K. Demostración. Véase [BRR81a]. Por tanto, si denimos el par de funciones de clase C∞ , λh:α∈(−α0, α0)−→ λh(α)∈R, uh:α∈(−α0, α0)−→ uh(α)∈H1 0 como λh(α) = λ0+ξh(α) uh(α) = u0+α ϕ0+vh(ξh(α), α), 24 entonces resolver el problema aproximado (3.9) en un entorno de un punto de retorno equivale a resolver la ecuación de bifurcación aproximada fh(ξh(α), α) = Fh(λ0+ξh(α), u0+αϕ0+vh(α, ξh(α))), ϕ0∗= 0, (3.10) de manera que (λh(α), uh(α)) : |α|< α0 es una rama de soluciones de (3.9). Denotemos por λ(m)(α) , u(m)(α) y λ(m) h(α) , u(m) h(α) a las derivadas m -ésimas de las funciones λ(α) , u(α) y λh(α) , uh(α) , respectivamente. El siguiente resultado, que constituye el resultado principal de este apartado, justica la aproximación de la curva de soluciones del problema continuo (3.8), obtenida en el Lema 2.13 del capítulo anterior, (λ(α), u(α)) , por la del problema discreto (3.9), (λh(α), uh(α)) , en un entorno del punto de retorno. Teorema 3.8. Sea (λ0, u0) un punto de retorno de F . Entonces, el problema aproximado Fh(λ, uh) = O tiene una única rama de soluciones (λh(α), uh(α)) : |α|< α0 en un entorno de la rama de soluciones (λ(α), u(α)) : |α|< α0 del problema continuo F(λ, u) = O. Además, estas ramas de soluciones son de clase C∞ y, para todo α∈(−α0, α0) y todo m∈Z tal que m≥0 , se tiene la siguiente estimación del error λ(m) h(α)−λ(m)(α)+  u(m) h(α)−u(m)(α)  H1 0 ≤K m X l=0     (T−Th)dl dαlG(λ(α), u(α))   H1 0 . 3.3.1. Puntos de retorno no degenerados En este apartado probaremos que si la curva de soluciones del problema continuo (3.8), (λ(α), u(α)) : |α|< α0 , presenta un punto de retorno no degenerado en α= 0 , entonces la curva de soluciones del problema aproximado (3.9), (λh(α), uh(α)) : |α|< α0 , también presenta un punto de retorno no degenerado en un punto λh(α0 h), uh(α0 h) . Supongamos entonces que (λ0, u0) es un punto de retorno no degenerado de F , esto es, un punto de retorno que, además, verica la condición DuuF0ϕ0ϕ0, ϕ0∗6= 0. 25 Como por hipótesis A es no singular, entonces el sistema (4.5) tiene solución para cualquier (y, η) y, por tanto, también la tiene la ecuación (4.7). Tomando y= 0 en dicha ecuación, tenemos (C−DB−1A)ζ=η, lo cual, tomando η6= 0 , implica que C−DB−1A6= 0 . 2. Caso 2: B singular y dim(Ker(B)) = 1 . Supongamos ciertas las igualdades (4.4) y veamos que entonces A es no singular. Nótese que si se satisfaciese el sistema (4.6), entonces A ζ =O=B x, pues, en caso contrario, si A ζ 6=O entonces B x 6=O , lo cual contradiría el hecho de que Im(A)∩Im(B) = {O}. Por tanto, como dim(Im(A)) = 1 y ζ∈R , entonces ζ= 0 . Sustituyendo esto en la segunda ecuación de (4.6) obtenemos D x =O, y como Ker(B)∩Ker(D) = {O} , entonces x=O , quedando probada la existencia y unicidad de soluciones de (4.6). Supongamos ahora A no singular. Entonces usando la primera ecuación de (4.5) se tiene que, dado y∈Rn , existe ζ∈R tal que y−A ζ =B x ∈Im(B). (4.8) Como dim(Ker(B)) = 1 , entonces dim(Im(B)) = n−1 y, por tanto, de (4.8) se deduce que dim(Im(A)) = 1 y, además, Im(A)∩Im(B) = {O}. Supongamos ahora, por reducción al absurdo, que Ker(B)∩Ker(D)6={O} . Entonces, para algún x06=O se tendría que B x0=O=D x0. Sustituyendo esto en (4.6), obtenemos (A ζ =O C ζ = 0 32 y, por tanto, (x, ζ) = (x0,0) sería una solución no trivial del sistema (4.6), contradiciendo la unicidad de soluciones de (4.6) (pues estamos suponiendo A no singular). Por tanto Ker(B)∩Ker(D) = {O} Por otro lado, sea Ker(B) = hφi. Notemos que A 0 φ!= O η!, donde η=D φ. Como A es no singular, entonces η6= 0 , de donde se deduce que dim(Im(D)) = 1. 3. Caso 3: B singular y dim(Ker(B)) >1 . Para probar que A es singular, basta construir una solución no trivial des sistema homogéneo. En este caso, existen φ1, . . . , φp∈Rn , p > 1 , vectores linealmente independientes tales que B φi=O,∀i∈1, . . . , p. Pero como el conjunto {D φi:i∈1, . . . , p} ⊂ R debe ser linealmente dependiente, entonces han de existir escalares αi∈R tales que x=Pp i=1 αiφi6=O satisfaga el sistema (B x =O D x =O y, por tanto, A es singular. Observación 4.4 . Para los dos primeros casos anteriormente estudiados podemos obtener explícitamente las soluciones del sistema (4.5). 33 1. Caso 1: B no singular. Supongamos B no singular y obtengamos la expresión explícita de la única solución (η, x) del sistema (4.5). Despejando x en la primera ecuación del sistema (4.5): x=B−1(y−Aζ), y sustituyendo en la segunda ecuación se obtiene C ζ +DB−1(y−A ζ) = η⇐⇒ (C−DB−1A)ζ=η−DB−1y, de donde ζ= (C−DB−1A)−1(η−DB−1y). 2. Caso 2: B singular y dim(Ker(B)) = 1 . Supongamos ahora B singular y tal que dim(Ker(B)) = 1 y obtengamos la expresión explícita de la única solución (η, x) del sistema (4.5). Como dim(Ker(B)) = 1 , entonces Rn=R⊕Im(B). Ahora bien, como dim(Im(A)) = 1 y Im(A)∩Im(B) = {O} , se tiene que R≈Im(A) y así Rn=Im(A)⊕Im(B). De esta forma, dado y∈Rn , existen x0∈Rn y ζ∈R tales que y=A ζ +B x0⇐⇒ B x0=y−A ζ. Calculemos, para este y , la solución general de la primera ecuación de (4.5) A ζ +B x =A ζ +B x0=⇒B x =B x0=⇒B(x−x0) = 0 =⇒x−x0∈Ker(B), esto es x=x0+z0, z0∈Ker(B). Para resolver la segunda ecuación del sistema, basta encontrar un z0∈Ker(B) tal que C ζ +D(x0+z0) = η⇐⇒ Dz0=η−C ζ −D x0. Nótese que, dado η∈R , la existencia de un z0 en tales condiciones queda garantizada por ser Ker(B)∩Ker(D) = {O} y dim(Im(D)) = 1 .  34 Si se aplica 1. del lema 4.3 a la matriz r DuF(λ, u) 0l∗! (4.9) obtenemos que, si la matriz DuF(λ, u) es no singular, entonces la matriz ampliada (4.9) es no singular si y solo si DuF−1(λ, u)r, l∗6= 0. (4.10) Si DuF(λ, u) , por el contrario, es singular y dim(Ker(DuF)) = 1 , entonces aplicando 2. del lema 4.3 a la matriz (4.9), resulta que la matriz ampliada (4.9) es no singular si y solo si r /∈Im(DuF)⇐⇒ r, ϕ0∗6= 0, l∗/∈Im(DuF∗)⇐⇒ ϕ0, l∗6= 0, (4.11) donde ϕ0 y ϕ0∗ son las denidas en (4.2). Observación 4.5 . En este capítulo se emplean las siguientes notaciones equivalentes en Rn v∗v=hv, v∗i, identicando así el producto escalar con el producto de dualidad. Teorema 4.6. La matriz ampliada r DuF(λ, u) 0l∗!∈ M(n+1)×(n+1) (4.12) es no singular si y solo si la matriz (I−rr∗)DuF+rl∗∈ Mn×n es no singular. Además, la solución del sistema implícito (4.3) viene dada por v(λ, u) = [(I−rr∗)DuF(λ, u) + rl∗]−1r∈Rn, σ(λ, u) = − hDuF(λ, u)v(λ, u), r∗i ∈ R. (4.13) De esta forma, (λ0, u0) es un punto de retorno si y solo si σ(λ0, u0) = 0, y se tiene que d dsσ(λ(s), u(s))s=0 =KDuuF0ϕ0ϕ0, ϕ0∗r, ϕ0∗−1, K ∈R\ {0}. 35 Demostración. Supongamos que la matriz ampliada (4.12) es no singular y consideremos la ecuación homogénea [(I−rr∗)DuF+rl∗]x=DuF x −rr∗DuF x +rl∗x =DuF x −r(r∗DuF x −l∗x) = O. (4.14) Multiplicando la ecuación (4.14) por r∗ y teniendo en cuenta que krk2=r∗r= 1 , resulta hDuF x −r(r∗DuF x −l∗x), r∗i=hDuF x , r∗i−krk2(r∗DuF x −l∗x) =hDuF x , r∗i−hDuF x , r∗i+l∗x=l∗x= 0, esto es, l∗x= 0 . Por tanto, la ecuación (4.14) se puede escribir como DuF x −rr∗DuF x =O, (4.15) de donde se deduce inmediatamente que r DuF(λ, u) 0l∗! −r∗DuF x x!= DuF x −rr∗DuF x l∗x!=O. (4.16) Como por hipótesis la matriz de coecientes es regular, entonces la igualdad (4.16) implica que x=O . Acabamos de probar entonces que [(I−rr∗)DuF+rl∗]x=O=⇒x=O, lo cual implica que la matriz (I−rr∗)DuF+rl∗ es no singular. Recíprocamente, supongamos que la matriz (I−rr∗)DuF+rl∗ es no singular y consideremos el sistema homogéneo de dos incógnitas, (x, α) , dado por (DuF x +rα =O, l∗x= 0. (4.17) Multiplicando la primera ecuación del sistema (4.17) por r∗ y teniendo en cuenta que krk2=r∗r=hr, r∗i= 1 , se deduce que hDuF x +rα , r∗i=hDuF x , r∗i+krk2α= 0 =⇒ hDuF x , r∗i+α=0=⇒ =⇒α=− hDuF x , r∗i. Sustituyendo esto en la primera ecuación de (4.17) y sumándole el vector cero r(l∗x) obtenemos DuF x −rhDuF x , r∗i+rl∗x= [(I−rr∗)DuF+rl∗]x=O. (4.18) 36 Como, por hipótesis, la matriz de coecientes (I−rr∗)DuF+rl∗ es no singular, entonces el sistema (4.18) tiene como única solución la solución trivial, x=O . Sustituyendo esto en (4.17), obtenemos que α= 0 . Por tanto, (x, α) = (O,0) es la única solución del sistema (4.17) , de donde se deduce inmediatamente que la matriz ampliada (4.12) es no singular. Por otro lado, es fácil comprobar que (4.13) es la única solución del sistema implícito r DuF(λ, u) 0l∗! σ v!= O 1!. (4.19) En efecto, multiplicando la primera ecuación del sistema (4.19) por r∗ y teniendo en cuenta que krk2=r∗r=hr, r∗i= 1 , resulta hrσ +DuF(λ, u)v, r∗i=krk2σ+hDuF v, r∗i=σ+hDuF v, r∗i=O, de donde σ=− hDuF(λ, u)v, r∗i. (4.20) Si ahora tomamos v:= [(I−rr∗)DuF(λ, u) + rl∗]−1r, operando se tiene que [(I−rr∗)DuF(λ, u) + rl∗]v=r=⇒DuF(λ, u)v−rhDuF(λ, u)v, r∗i+rl∗v=r. Sustituyendo (4.20) en la anterior expresión, resulta DuF(λ, u)v+rσ +rl∗v=r y, como DuF(λ, u)v+rσ =O (primera ecuación del sistema (4.19)), entonces rl∗v=r=⇒l∗v= 1, vericándose así la segunda ecuación del sistema (4.19). A continuación estudiamos el comportamiento de la función test σ en un entorno de un punto de retorno cuadrático. Sea (λ0, u0) un punto de retorno no degenerado de F y consideremos la primera ecuación del sistema (4.19) que dene de manera implícita a la función test σ≡σ(λ(s), u(s)) a lo largo de la curva de soluciones (λ(s), u(s)) , suponiendo que (λ(0), u(0)) = (λ0, u0) : r σ(λ(s), u(s)) + DuF(λ(s), u(s)) v(λ(s), u(s)) = O. (4.21) 37 Evaluando (4.21) en s= 0 resulta r σ(λ0, u0) + DuF0v(λ0, u0) = O. (4.22) Como en un punto de retorno se verica rang(DuF0) = n−1 , entonces, por (4.11), se tiene que r /∈Im(DuF) , luego de (4.22) se sigue σ(λ0, u0)=0 (4.23) de donde v(λ0, u0)∈Ker(DuF0) = ϕ0=⇒v(λ0, u0) = K ϕ0, K ∈R\ {0}. (4.24) Por otro lado, si derivamos (4.21) respecto de s rd dsσ(λ(s), u(s)) + DuλF(λ(s), u(s)) ˙ λ(s) + DuuF(λ(s), u(s)) ˙u(s)v(λ(s), u(s))+ +DuF(λ(s), u(s)) d dsv(λ(s), u(s)) = O. (4.25) Evaluando (4.25) en s= 0 rd dsσ(λ(s), u(s))s=0 +DuλF0˙ λ(0) + DuuF0˙u(0)v(λ0, u0)+ +DuF0d dsv(λ(s), u(s))s=0 =O. (4.26) Calculemos, en primer lugar, las expresiones de ˙ λ(0) y de ˙u(0) . Para ello, basta usar la regla de la cadena en la ecuación F(λ(s), u(s)) = O en el punto s= 0 , obteniendo d dsF(λ(s), u(s))s=0 =DλF0˙ λ(0) + DuF0˙u(0) = O. (4.27) Como en un punto de retorno se tiene que DλF0/∈Im(DuF0) , entonces de (4.27) se deduce ˙ λ(0) = 0, (4.28) de donde ˙u(0) ∈Ker(DuF0) = ϕ0=⇒˙u(0) = C ϕ0, C ∈R\ {0}. (4.29) 38 Sustituyendo estos valores en (4.26) rd dsσ(λ(s), u(s))s=0 +DuuF0C ϕ0v(λ0, u0) + DuF0d dsv(λ(s), u(s))s=0 =O y multiplicando ambos lados de la igualdad por ϕ0∗ y usando (4.24), resulta rd dsσ(λ(s), u(s))s=0 , ϕ0∗+DuuF0CK ϕ0ϕ0, ϕ0∗+ +DuF0d dsv(λ(s), u(s))s=0 , ϕ0∗=O. (4.30) Ahora bien, como Im(DuF0) = ϕ0∗ , entonces (4.30) se puede escribir como rd dsσ(λ(s), u(s))s=0 , ϕ0∗+DuuF0CK ϕ0ϕ0, ϕ0∗=O. (4.31) Como r, ϕ0∗6= 0 por ser r /∈Im(DuF0) , se puede despejar de (4.31), obteniendo d dsσ(λ(s), u(s))s=0 =K2DuuF0ϕ0ϕ0, ϕ0∗r, ϕ0∗−1, K2:= −CK ∈R\ {0}, de donde ˙σ(λ0, u0)≡d dsσ(λ(s), u(s))s=0 6= 0, (4.32) ya que en un punto de retorno no degenerado se verica que DuuF0ϕ0ϕ0, ϕ0∗6= 0 . Por tanto, de (4.23) y de (4.32) se sigue inmediatamente que, en un punto de retorno no degenerado, la función test σ(λ, u) cambia estrictamente de signo. 39 40 Capítulo 5 Resultados numéricos 5.1. Solución analítica del problema de Bratu Como ya se ha mencionado a lo largo de este trabajo, el problema de Bratu se dene como el problema de contorno no lineal (−u00(x) = λeu(x)∀x∈[0,1], u(0) = 0 = u(1). (5.1) Dado λ∈R , la expresión analítica de las soluciones de (5.1) es conocida y se dene en función de un parámetro auxiliar, µ , de la siguiente manera: Si λ≥0 , entonces la solución analítica viene dada por u(x) = ln cosh2µ 2·cosh−2µx−1 2, (5.2) donde el parámetro µ≥0 es solución de la ecuación trascendente λ= 2µ2cosh−2µ 2. (5.3) Si λ < 0 , entonces u(x) = ln cos2µ 2·cos−2µx−1 2, siendo µ∈[0, π) solución de λ=−2µ2cos−2µ 2. (5.4) Observación 5.1 . Las ecuaciones no lineales (5.3) y (5.4) se pueden resolver fácilmente haciendo uso de un método iterativo como, por ejemplo, el método de Newton. 41 Para λ=−2,5 : Paso Error ndim h= 0,1 2,0662 ·10−2 9 h= 10−22,8543 ·10−3 99 h= 10−32,9389 ·10−4 999 Tabla 5.5: Tabla de error calculada con k·kH1 0 , para λ=−2,5 . Para λ= 2,5 : Paso Error ndim h= 0,1 3,8522 ·10−2 9 h= 10−25,0750 ·10−3 99 h= 10−35,2000 ·10−4 999 Tabla 5.6: Tabla de error calculada con k·kH1 0 , para λ= 2,5 . Como vemos, el empleo de una fórmula de derivación numérica de orden uno provoca que el orden del método descienda. 5.4. Aproximación numérica de la autofunción en un entorno del punto de retorno Finalmente, en el punto de retorno (λ0, u0) también se puede obtener numéricamente la autofunción ϕ0 denida en el Capítulo 2. Para ello, basta tener en cuenta que ϕ0 es solución del problema de contorno lineal (−ϕ000 (x) = λ0eu0(x)ϕ0(x)∀x∈[0,1], ϕ0(0) = 0 = ϕ0(1). De esta manera, ϕ0 se puede obtener numéricamente sin más que resolver el siguiente sistema lineal con matriz singular (obtenido tras discretizar el operador derivada segunda con diferencias nitas) (A0h−D)vh=O, (5.8) 48 siendo D∈ Mn×n la matriz diagonal D=λ0       eu0 10· · · 0 0 0eu0 2· · · 0 0 . . .. . ..... . . 0 0 · · · 0eu0 n        . En la Figura 5.3 se representa la autofunción aproximada, vh , obtenida tras resolver el sistema (5.8) por remonte, para lo cual se ha empleado la factorización de Gauss de la matriz A0h−D . 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Nodos (x) 0 0.05 0.1 0.15 vh Autofuncion 0 aproximada Figura 5.3: Aproximación numérica de la autofunción ϕ0 de DuF(λ0, u0) . 49 50 Conclusiones y nuevas líneas de investigación 5.5. Conclusiones de la memoria En este trabajo se ha visto que el problema de Bratu unidimensional (−u00(x) = λeu(x)∀x∈[0,1], u(0) = 0 = u(1). (5.9) puede presentar dos soluciones, una o ninguna, dependiendo del valor del parámetro de control λ∈R . Además, deniendo la función F:R×V−→ V de manera adecuada, se ha probado que la curva de soluciones de (5.9) en función del parámetro λ se puede obtener de manera equivalente resolviendo el problema no lineal F(λ, u) = O. (5.10) Sin embargo, la curva de soluciones de (5.10) presenta un punto de retorno no degenerado, lo cual diculta la parametrización de la misma. Por ello, mientras que en un entorno de un punto regular se ha podido emplear la parametrización natural (esto es, expresar u en función de λ ) sin más que usar el teorema de la función implícita, en un entorno del punto de retorno ha sido necesario emplear una parametrización respecto al parámetro longitud de arco. También se ha probado que el problema discretizado obtenido mediante un esquema de diferencias nitas de orden dos presenta un comportamiento similar al del problema continuo, de manera que la curva de soluciones del problema aproximado también posee un punto de retorno no degenerado. Además, se han ofrecido las correspondientes cotas de error de la aproximación de la curva de soluciones del problema continuo (5.10) por las del problema discreto correspondiente. 51 Finalmente, se han presentado algunos tipos de funciones escalares (llamadas funciones test ) que permiten detectar, de manera sencilla, la presencia de dichos puntos de retorno. 5.6. Problema de Bratu en dimensión superior Aunque en este trabajo se ha estudiado únicamente el caso unidimensional del problema de Bratu, también es de gran importancia su versión en dimensión superior, denido como la siguiente ecuación diferencial no lineal, elíptica, con dominio acotado Ω⊂Rn , con n≤3 , y condiciones de contorno Dirichlet homogéneas en la frontera ∂Ω (−∆u(x) = λeu(x)∀x∈Ω, u(x) = O ∀x∈∂Ω, donde ∆u≡ ∇2u denota el laplaciano de u . De nuevo, este problema aparece en numerosos ámbitos de la ciencia, como son la teoría de Frank-Kamenetskii sobre la combusión térmica de una mezcla homogénea de reactivos mantenida dentro de un recipiente cerrado con paredes de temperatura constante; o la ecuación de Emden-Chandrasekhar que, además de describir la expansión de la nebulosa planetaria, también permite describir la carga espacial de electricidad alrededor de un cable incandescente. 5.7. Discretizaciones alternativas del problema de Bratu Aunque en este trabajo la discretización del problema unidimensional de Bratu se ha llevado a cabo mediante diferencias nitas de orden dos, también se podría haber abordado mediante un esquema de elementos nitos Lagrange P1 . Sin embargo, todos los resultados que se han presentado a lo largo de los Capítulos 2 y 3 son también válidos para problemas variacionales, de manera que se aplicarían de manera inmediata. Por otro lado, en la discretización de la derivada segunda se podría haber optado por un esquema de diferencias nitas de mayor orden, como el método de Numerov. Este método numérico, empleado para resolver problemas de contorno como el que se expone en este trabajo, es un método de diferencias nitas de cuarto orden utilizado generalmente en problemas de bifurcación, por lo que su uso habría aumentado el orden de la discretización considerablemente. 52 Por último, cabe mencionar que existen formas de obtener soluciones aproximadas del problema de Bratu sin necesidad de ningún tipo de discretización, como es el caso del VIM ( variational iteration method ). Este método, basado en los multiplicadores de Lagrange, es el que se expone con detalle en [BAT10]. 5.8. Aproximación numérica en un entorno de un punto de bifurcación simple Dada una función arbitraria F:R×V−→ V de clase C∞ , el problema no lineal F(λ, u) = O (5.11) puede presentar otros tipos de singularidades que provoquen una multiplicidad de soluciones. Este es el caso de, por ejemplo, los llamados puntos de bifurcación simples. Denición 5.5. Diremos que (λ0, u0)∈Rn+1 es un punto de bifurcación simple de F si: 1. F(λ0, u0) = O , 2. DuF0∈L(Rn,Rn) es singular, con rang(DuF0) = n−1 , 3. DλF0:= DλF(λ0, u0)∈Im(DuF0) . Observación 5.6 . Nótese que un punto de bifurcación simple y un punto de retorno se diferencian, únicamente, en que en este último se tiene que DλF0/∈Im(DuF0) . Por este motivo, en el caso de un punto de bifurcación simple, DF 0 no es de rango máximo, a diferencia de lo que sucede en un punto de retorno. Si además se impone una condición de intersección transversal sobre la derivada segunda de F , se puede probar que el problema (5.11) presenta dos ramas de soluciones distintas que se intersecan transversalmente en (λ0, u0) . En [CR71] se puede encontrar un mayor análisis de este tipo de puntos y en [BRR81b] se estudia la aproximación numérica de la curva de soluciones del problema continuo (5.11) por la del correspondiente problema aproximado en un entorno de un punto de bifurcación simple. 53 54 Bibliografía [BAT10] B. BATIHA. Numerical solution of bratu-type equations by the variational iteration method. Hacettepe Journal of Mathematics and Statistics , 39(1):23 29, 2010. [BRE10] H. BREZIS. Functional analysis, Sobolev spaces and partial dierential equations . Springer Science & Business Media, 2010. [BRR80] F. BREZZI, J. RAPPAZ, and P.A. RAVIART. Finite dimensional approximation of nonlinear problems. Part I: Branches of nonsingular solutions. Numer. Math. , (36):pp. 125, 1980. [BRR81a] F. BREZZI, J. RAPPAZ, and P.A. RAVIART. Finite dimensional approximation of nonlinear problems. Part II: Limit points. Numer. Math. , (37):pp. 128, 1981. [BRR81b] F. BREZZI, J. RAPPAZ, and P.A. RAVIART. Finite dimensional approximation of nonlinear problems. Part III: Simple Bifurcation points. Numer. Math. , (38):pp. 130, 1981. [CCÖ + 09] H. CAGLAR, N. CAGLAR, M. ÖZER, A. VALARISTOS, A. N. MILIOU, and A. N. ANAGNOSTOPOULOS. Dynamics of the solution of bratu's equation. Nonlinear Analysis: Theory, Methods & Applications , 71(12):e672e678, 2009. [CMT89] P. G. CIARLET, B. MIARA, and J. M. THOMAS. Introduction to numerical linear algebra and optimisation . Cambridge University Press, 1989. [CR71] M. G. CRANDALL and P. H. RABINOWITZ. Bifurcation from simple eigenvalues. Journal of Functional Analysis , 8(2):321340, 1971. [CR90] M. CROUZEIX and J. RAPPAZ. On numerical approximation in bifurcation theory . Springer Berlin, 1990. 55 [DK80] D. W. DECKER and H. B. KELLER. Multiple limit point bifurcation. Journal of mathematical analysis and applications , 75(2):417430, 1980. [GR84] A. GRIEWANK and G. W. REDDIEN. Characterization and computation of generalized turning points. SIAM journal on numerical analysis , 21(1):176185, 1984. [MEI13] Z. MEI. Numerical bifurcation analysis for reaction-diusion equations , volume 28. Springer Science & Business Media, 2013. [ROD03] G. RODRÍGUEZ. Diferenciación de funciones de varias variables reales , volume 4. Univ Santiago de Compostela, 2003. 56