Breathers en la ecuación de sine-Gordon con laplaciano fraccional

Universidad de Sevilla Doble Grado en Física y Matemáticas Trabajo de Fin de Grado Breathers en la ecuación de sine-Gordon con laplaciano fraccional Autor: Jorge Catarecha Otero Saavedra Tutores: Jesús Cuevas Maraver Ana María Gómez Ramírez Departamento de Física Atómica, Molecular y Nuclear Facultad de Física 23 de junio de 2022 i Agradecimientos Con estas líneas me gustaría agradecer a todas las personas que, directa o indirectamente, han contribuido a este trabajo de fin de grado. En primer lugar, a Jesús Cuevas Maraver quien con su tiempo, sabiduría y comprensión no sólo me ha enseñado sobre la física no lineal sino también acerca de cómo es iniciarse en el mundo de la investigación científica. De igual forma, a Ana María Gómez Ramírez por su ayuda para poder presentar este trabajo en el departamento de Físicca Atómica, Molecular y Nuclear. Querría agradecer asimismo al Instituto de Matemáticas de la Universidad de Sevilla, ya que con su programa de iniciación a la investigación, PI3, me puso en contacto con el proyecto que ofertaba Jesús Cuevas sobre la ecuación de sine-Gordon, y que fue el germen de este trabajo. Debo también agradecer igualmente al Ministerio de Educación y Formación Profesional por concederme la beca de colaboración en el departamento de Física Aplicada I de la Universidad de Sevilla, permitiéndome desarrollar este trabajo. También quisiera agradecer a mis compañeros de clase; quienes, aún no conociendo el tema con anterioridad, me han ayudado de buena fe siempre que han podido. Por último, no puedo pasar sin agradecer especialmente a mi familia y resto de personas allegadas por su gran apoyo siempre que me ha hecho falta. Sin ellos, el camino recorrido para llevar a buen puerto este trabajo se me habría hecho mucho más arduo. ii Resumen El objetivo de este proyecto es estudiar las soluciones tipo breather en la ecuación de sine-Gordon clásica y la ecuación con laplaciano de orden fraccionario. Hemos dividido el trabajo en cuatro capítulos. En el capítulo 1 presentamos el interés de este trabajo e introducimos los conceptos básicos con los que trabajamos: ecuación de sine-Gordon, breathers y derivadas de orden fraccionario. El segundo capítulo lo dedicamos a analizar los métodos numéricos que implementaremos en MATLAB para calcular breathers: el método de disparo y un método basado en el desarrollo de Fourier de la solución. Un concepto de gran importancia es la estabilidad de las soluciones a la que le dedicamos el capítulo 3. En dicho capítulo estudiamos los modos lineales de vibración (una de las causas principales de la aparición de inestabilidades) y exponemos los resultados principales de la teoría de Floquet, teoría que nos permite dar una condición necesaria y suficiente para conocer la estabilidad lineal de los breathers. Por último, reservamos el cuarto capítulo para exponer los resultados principales de este trabajo: primeramente, reproducimos los resultados clásicos ya conocidos para breathers en la ecuación de sine-Gordon y, en segundo lugar, presentamos los resultados novedosos acerca de la dependencia de las propiedades de los breathers con el orden fraccionario del laplaciano considerado. iii Abstract The aim of this project is to study breather-like solutions in the sine-Gordon equation with the classical Laplacian and the fractional one. We have divided it into four chapters. Chapter 1 is devoted to motivating the interest of this work as well as to introducing the main concepts we are working with: the sine-Gordon equation, breathers and fractional derivatives. In chapter 2, we analyze the numerical methods we will implement in MATLAB to calculate breathers: shooting method and a method based on expressing the solution as Fourier series. Another main concept is stability of solutions, thus we will consecrate chapter 3 to it. In the aforementioned chapter, we will study linear vibration modes (which are one of the main causes of instability) and we will present as well the main results of Floquet theory. This will allow us to state a necessary and sufficient condition for breathers lineal stability. Lastly, chapter 4 will be dedicated to presenting our results: firstly, we will reproduce classical, well-known results about the sine-Gordon equation; and, secondly, we will show new results about the relation between breathers properties and the fractional-order Laplacian. Índice general 1. Introducción 1 1.1. Ecuación de sine-Gordon . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2. Breathers .................................... 2 1.3. Derivadas de orden fraccionario . . . . . . . . . . . . . . . . . . . . . . . . 5 2. Resolución numérica 9 2.1. Métododedisparo ............................... 11 2.2. Método de la serie de Fourier . . . . . . . . . . . . . . . . . . . . . . . . . 12 3. Estabilidad de soluciones 15 3.1. ModosNormales ................................ 15 3.2. AnálisisdeFloquet............................... 18 4. Resultados 25 4.1. ResultadosClásicos............................... 25 4.2. Resultados Fraccionarios . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 5. Conclusiones 48 iv Capítulo 1 Introducción 1.1. Ecuación de sine-Gordon Una de las ecuaciones más conocidas en el mundo de la física para modelar la propagación de ondas dispersivas es la ecuación de Klein-Gordon: d2u dt2−∆u+u= 0,(1.1) donde u(x, t)es la perturbación como función de la posición y el tiempo, y ∆el operador laplaciano. La razón de que sea una ecuación dispersiva es que si proponemos soluciones del tipo u(x, t) = exp(−i(ωt − kx)), la condición necesaria para que sea solución de (1.1) es la relación de dispersión: ω2= 1 + k2.(1.2) Esta relación es, por tanto, distinta a la propia de ondas no dispersivas: ω=c·k. El hecho de ser lineal es una de las razones para que esta ecuación sea tan utilizada, especialmente, en los ámbitos de la mecánica cuántica relativista y la teoría cuántica de campos, donde se utiliza para modelar los campos mesónicos (véase el capítulo 4 de [16]). Sin embargo, también existen generalizaciones no lineales con propiedades muy interesantes desde el punto de vista físico-matemático. En este trabajo nosotros nos centramos en generalización no lineal conocida como ecuación de sine-Gordon: d2u dt2−∆u+ sinu= 0.(1.3) 1 CAPÍTULO 1. INTRODUCCIÓN 2 Como la versión linealizada de esta ecuación (que recordemos que es la que siguen las soluciones de pequeña amplitud) coincide con (1.1), la relación de dispersión para los modos normales sigue siendo (1.2). Esta es una ecuación completamente integrable, permitiendo encontrar amplias familias de soluciones analíticas y cantidades conservadas. Algunas herramientas matemáticas que se han usado para ello son las transformaciones Bäcklund [14], que datan del siglo XIX y se basan en un estudio geométrico-diferencial, y las transformadas de scattering inverso (inverse scattering transforms), llevadas a cabo en la década de 1970 [23]. Históricamente, esta ecuación la encontramos por primera vez en un trabajo de Jacques Edmond Émile Bour de 1862, en el que estudiaba las ecuaciones de Gauss-Codazzi para superficies bidimensionales de curvatura gaussiana negativa inmersas en el espacio euclídeo tridimensional [5]. Posteriormente, en 1939, volvió a ganar popularidad gracias al modelo de Frenkel-Kontorova para explicar dislocaciones en sólidos. En este modelo se usaba una versión unidimensional y discretizada de (1.3): d2un dt2−C(un+1 +un−1−2un) + sinun= 0 (1.4) donde unes la perturbación en la posición del átomo n-ésimo de la red y C=1 h2, con hla distancia entre átomos de la red. Esta versión discretizada será con la que trabajaremos nosotros [6]. Algunos otros ejemplos donde se han encontrado aplicaciones para esta ecuación son el modelo de Skyrme de nucleones en un fluido piónico (1958), las uniones de Josephson en superconductores (1967), la transparencia autoinducida (1967) [13] o el más reciente modelado de estrellas de axiones [15]. El hecho de que esta ecuación surja en ámbitos de la física de tan diversa índole no es más que una clara comprobación de la utilidad y flexibilidad de la misma. Por ello, en este trabajo estudiaremos la versión discreta de dicha ecuación y su versión de orden fraccional. Este concepto nunca antes se ha investigado aplicado a la ecuación de sine-Gordon discreta. 1.2. Breathers Como ya hemos mencionado, la total integrabilidad de la ecuación de sine-Gordon permite la obtención analítica de una amplia familia de soluciones. Algunos ejemplos CAPÍTULO 1. INTRODUCCIÓN 3 de estas son solitones topológicos (kinks yanti-kinks) y breathers. De entre todas ellas nosotros nos centraremos en el estudio de los breathers. Los breathers son soluciones localizadas exponencialmente en el espacio y periódicas en el tiempo. En la literatura también se denominan breathers móviles a soluciones localizadas que se desplazan en el espacio; sin embargo, en este trabajo (salvo que se indique lo contrario) siempre nos referiremos a breathers estacionarios. En la figura 1.1a, vemos el perfil de un breather permitiéndonos observar cómo la perturbación está localizada (el mapa de color de 1.1c nos muestra la evolución temporal). Para que una ecuación diferencial tenga soluciones tipo breather es necesario que sea no lineal y dispersiva. Esta condición diferencia estas soluciones localizadas de las que surgen, por ejemplo, debido a impurezas en una cadena monoatómica de osciladores lineales. Algunas ecuaciones diferenciales en derivadas parciales que presentan soluciones de este tipo son la ecuación de sine-Gordon y la ecuación de Schrödinger no lineal [11]. En concreto, el breather de la ecuación de sine-Gordon toma la forma: u= 4arctan  √1−ω2cos(ωt) ωcosh√1−ω2x ,(1.5) donde ωes la frecuencia de oscilación del breather y debe cumplir ω < 1. Que una ecuación diferencial tenga soluciones de este tipo no es habitual; sin embargo, en redes discretas la situación es completamente antagónica. Este hecho es fundamental tenerlo en cuenta ya que, cuando se resuelven numéricamente ecuaciones diferenciales en derivadas parciales, es necesario discretizar el dominio en nodos acoplados entre sí por el esquema numérico considerado, siendo en esta red discreta donde sí pueden aparecer breathers en condiciones bastante generales. Este resultado, debido a MacKay y Aubry (1994) [18], enuncia que en redes discretas de osciladores existen breathers cuando el acoplamiento entre nodos es suficientemente débil. La razón para que el acoplamiento sea suficientemente débil subyace en la necesidad de que ninguno de los armónicos del breather resuene con la banda de fonones (modos lineales). Es decir, ningún múltiplo de la frecuencia del breather debe coincidir con una frecuencia de oscilación de un modo lineal. Además, el extremo superior de la banda de fonones aumenta conforme lo hace el acoplamiento, lo que justifica por qué es necesario que el acoplamiento sea débil. CAPÍTULO 1. INTRODUCCIÓN 4 0 5 10 15 20 25 30 35 40 n -0.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 un max Perfil en máxima elongación, C = 1.26 (a) Perfil de un breather. 0 5 10 15 20 25 30 35 40 n -0.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 un max Perfil en máxima elongación, C = 1.27 (b) Perfil de un nanóptero. Desplazamiento de nodos, C = 1.26 5 10 15 20 25 30 35 40 Nodos 0 5 10 15 Tiempo -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 (c) Evolución temporal del breather. Desplazamiento de nodos, C = 1.27 5 10 15 20 25 30 35 40 Nodos 0 5 10 15 Tiempo -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 (d) Evolución temporal del nanóptero. Se aprecia cómo la perturbación afecta toda la red. Figura 1.1: Comparación entre un breather y un nanóptero. CAPÍTULO 2. RESOLUCIÓN NUMÉRICA 11 2.1. Método de disparo Nuestro sistema es un problema de contorno en el tiempo, con condiciones de contorno periódicas. El método de disparo trata de simplificar la dificultad que surge de estas condiciones de contorno. La idea del método es transformar el problema de contorno en el tiempo en un problema de valores iniciales y la minimización de una función. El esquema del método es el siguiente: Seleccionamos un vector inicial o semilla U(0) = U0(recordemos que U(t)es un vector de 2Ncomponentes: las Nposiciones y las Nvelocidades de los nodos). Resolvemos numéricamente la ecuación diferencial de (2.1) con la condición inicial U(0) = U0. Evaluamos la solución obtenida en t=T= 2π/ω, y estimamos la diferencia con la condición inicial: ∆ = U0−U(T, U0). Ajustamos las condiciones iniciales para minimizar ∆. Es decir, podemos pensar nuestro problema como la búsqueda del vector U0tal que ∆(U0) = U0−U(T;U0) = 0, siendo ahora un problema de encontrar los ceros de dicha función. Hay que destacar que la solución nula es solución del nuestro problema; sin embargo, no es la que deseamos encontrar. Por ello, y debido a la convergencia local de este método, a la hora de implementarlo numéricamente es necesario partir de una solución suficientemente cerca de la deseada. Este es el motivo de realizar la continuación de la solución previamente explicada ya que, al variar ligeramente la constante de acoplamiento, esperamos que la solución buscada haya variado poco, siendo la solución con el valor anterior de la constante un buen punto de partida. Este método, sin embargo, presenta algún inconveniente desde el punto de vista computacional. Por una parte, en la búsqueda de la solución se obtendrán velocidades iniciales no nulas, lo que posteriormente aumentaría la complejidad del cálculo de la energía de las soluciones al deber tener en cuenta este término cinético. Además, el método de minimización de la función ∆se basa en el método de Newton-Raphson a la hora de elegir la nueva condición inicial en cada iteración. Este método es un método iterativo de segundo orden para encontrar los ceros de un sistema de ecuaciones  F(x). En dicho método, a partir de un CAPÍTULO 2. RESOLUCIÓN NUMÉRICA 12 vector inicial x0, las sucesivas iteraciones se calcucan como xn+1 =xn−JF(xn)−1F(xn), donde JFes el jacobiano de la función F. Por tanto, utilizar este método conlleva la construcción de la matriz jacobiana de la función ∆, lo cual presenta un coste computacional importante (será necesario calcular 2N derivadas numéricas de la función ∆). Por ello, otra opción es considerar una aproximación finita de la matriz jacobiana, requiriendo trabajar con precisiones numéricas muy altas a fin de que la aproximación sea buena y no se vea perjudicada la velocidad de convergencia del método. 2.2. Método de la serie de Fourier El segundo método que usaremos se desarrolla en torno a la idea de imponer la periodicidad de la solución directamente expresándola como una serie de Fourier en cada nodo. Este método, que resumimos a continuación, se puede consultar en detalle en [12]. De esta manera, la perturbación en cada nodo toma la forma: un(t) = zo n+ 2 km X k=1 zk ncos(kωbt),con n= 1, . . . , N (2.5) Aunque teóricamente la serie debe ser infinita, a la hora de calcularla numéricamente será necesario truncarla. Este es el motivo de que el sumatorio acabe en el término km. Esto será una fuente de error numérico; sin embargo, es un parámetro fácil de controlar, lo que nos permitirá en la práctica reducir el error por debajo de otros errores de integración. Otra ventaja es que al imponer que la serie de Fourier solo tenga cosenos nos aseguramos que las velocidades sean nulas en el instante inicial. De esta forma, los breathers calculados estarán en la posición de elongación máxima para t= 0. Sabemos que cada nodo cumple la ecuación: d2un dt2+ sinun−CX i=n Kα(|n−i|)(ui−un) = 0,(2.6) donde para el caso clásico basta tomar como valores de Kαlos valores límites calculados en (1.11). Introduciendo ahora la expresión (2.5), multiplicando por cos(kωbt)e integrando entre 0 y T, nos permite obtener el siguiente sistema de ecuaciones para los coeficientes CAPÍTULO 2. RESOLUCIÓN NUMÉRICA 13 de cada armónico: Fk n:=−k2ω2 bzk n+V′k n−CX i=n Kα(|n−i|)zk i−zk n= 0,con k= 0, . . . , km.(2.7) V′k nno es más que el k-ésimo coeficiente coseno de Fourier de la función V′(t) = sin(un(t)) (que presenta la misma periodicidad en el tiempo que u(t)). El objetivo de este método es encontrar los coeficientes de la serie de Fourier como solución a un sistema de ecuaciones algebraicas. Para ello, lo siguiente que necesitamos es relacionar V′k ncon los coeficientes del desarrollo en serie de u(t). Necesitaremos definir las siguientes cantidades: un(ti) = zo n+ 2 km X k=1 zk ncos(kωbti), ti=2πi ωb(2km+ 1),con i= 0, . . . , km.(2.8) Con estas evaluaciones equiespaciadas de u(t)podemos recurrir ahora a la herramienta matemática conocida como transformada discreta de Fourier (DFT), que permite recuperar los primeros coeficientes del desarrollo de Fourier a partir de únicamente la evaluación en una cantidad finita de puntos. Esto es suficiente, pues recordemos que no nos interesan todos los coeficientes al estar trabajando con series truncadas. De esta manera valen: V′k n=1 2km+ 1  sin(un(0)) + 2 km X i=1 sin(un(ti))cos(kωbti) ,(2.9) dependientes implícitamente de los coeficientes zk n. Con esto hemos logrado que, si consideramos las ecuaciones (2.7) para todos los nodos, tengamos un sistema de ecuaciones implícitas para todos los coeficientes. Este conjunto de ecuaciones se puede resolver de nuevo por el método de Newton-Raphson, ya que es equivalente a encontrar los ceros de la función  F(zk n)definida en (2.7). La principal ventaja de este método es que, al contrario que en el método de disparo, la matriz jacobiana Jque se utiliza en el método de Newton sí se puede calcular de forma analítica. Por tanto, se le puede proporcionar al algoritmo numérico para reducir drásticamente el tiempo de cómputo, así como los errores provocados por la aproximación numérica de J. Procedamos a calcular las componentes de J, que las podemos escribir CAPÍTULO 2. RESOLUCIÓN NUMÉRICA 14 como: Jk,k′ n,n′=∂Fk n ∂zk′ n′ ,(2.10) Gracias a cómo aparecen los coeficientes en las ecuaciones, podemos separar en una parte diagonal y otra no diagonal: J=JD+JND,(2.11) siendo cada una de estas:            JD=(" −k2ω2 b+CP i=n Kα(|n−i|)!δk,k′+∂V ′k n ∂zk′ n#δn,n′) JND =(−C P i=n Kα(|n−i|)δi,n′!δk,k′).(2.12) En esta ecuación la derivada ∂V ′k n ∂zk′ nse calcula derivando (2.9), obteniendo: ∂V ′k n ∂zk′ n =1 1 + δk′,0V′′|k−k′| n+V′′(k+k′) n,(2.13) donde los coeficientes V′′ kson los coeficientes de la transformada de Fourier discreta de cos(un(t)), calculados de forma análoga a (2.9). Para esto es importante que, debido al término V′′ (k+k′), sea necesario considerar coeficientes hasta 2km. En la práctica, estos valores se calcularán considerando 0 los valores de zk ncon k > km, estando considerado este error de truncamiento en el error de truncamiento de la serie original de u(t). Con esto, ya tendríamos deducido analíticamente el jacobiano necesario para el solucionador numérico. Sin embargo, este método presenta también el inconveniente de que hemos aumentado la dimensión del sistema a N·(km+ 1) (km+ 1 coeficientes por cada nodo). Por tanto, la inversión de la matriz jacobiana para el método de Newton puede llegar a consumir mucho tiempo en sistemas con una gran cantidad de nodos, a pesar de que la matriz jacobiana en sí puede ser calculada de forma analítica. En estas condiciones será preferible utilizar el método de disparo. Capítulo 3 Estabilidad de soluciones 3.1. Modos Normales Como ya adelantamos en la introducción, una de las razones de que los breathers discretos se desestabilicen es la resonancia con los modos lineales de vibración de la red. Según el teorema de McKay-Aubry [18], la inestabilidad puede ocurrir cuando la frecuencia de alguno de los modos normales de vibración es múltiplo de la frecuencia de vibración del breather. Debido, por tanto, al interés que presenta, dedicaremos esta sección al cálculo de estos modos normales y cuándo pueden ocurrir estas resonancias, que dan lugar a hibridaciones breather-fonón. Los modos normales cumplen la versión linealizada de la ecuación (2.6): d2un dt2+un−CX i=n Kα(|n−i|)(ui−un)=0.(3.1) Y se calculan imponiendo soluciones armónicas un(t) = e−i(ωt−kn): −ω2·e−i(ωt−kn)=−e−i(ωt−kn)+CX j=n Kα(|n−j|)e−i(ωt−kj)−e−i(ωt−kn), ω2= 1 −CX j=n Kα(|n−j|)eik(j−n)−1, ω2= 1 −C∞ X m=1 Kα(m)eikm −2 + e−ikm, ω2= 1 + 4C∞ X m=1 Kα(m)sin2 mk 2!.(3.2) Esta última fórmula es la relación de dispersión de la red, que la hemos supuesto infinita 15 CAPÍTULO 3. ESTABILIDAD DE SOLUCIONES 16 para poder llegar a expresiones analíticas. En el caso clásico esto no supone ninguna diferencia ya que la interacción es de corto alcance, y en el caso fraccionario debido al decaimiento del núcleo, si la red es suficientemente larga estaremos despreciando la cola de la serie convergente (que se puede hacer suficientemente despreciable). La relación de dispersión nos permite calcular el modo normal asociado a cada número de onda k. Al ser periódica, nos podemos restringir a k∈[0, π], donde alcanza el máximo en k=πy el mínimo ω= 1 en k= 0. Si queremos una fórmula cerrada analítica, podemos reescribir la relación de dispersión como: ω2=1 + 16Γ(α+ (1/2)) √πΓ(1 + α)1−exp(−ik)αΓ(1 + α) hR(1,1−α, 2 + α;exp(−ik)) + exp(2ik)R(1,1−α, 2 + α;exp(ik))i, donde R(a, b, c;z) = 2F1(a, b, c;z)/Γ(c)es la función hipergeométrica regularizada [21]. En general, nos interesará el valor máximo que toman las frecuencias de los modos normales ya que indicará cuándo alguno de los armónicos del breather entra en la banda de fonones (como la frecuencia inferior de la banda es 1, nosotros trabajaremos con breathers de frecuencias menores). El valor máximo depende del orden fraccionario considerado; nosotros trabajaremos entre α= 0,5(ya que para valores menores el acoplamiento es de tan largo alcance que se vuelve inviable para simular numéricamente) y α= 1, el caso clásico. En estos dos casos la expresión del extremo superior (k=π) de la banda de frecuencias toma una expresión particularmente sencilla: ωClasica max =√1+4C, (3.3) ωα=0,5 max =√1+2C. (3.4) Si consideramos órdenes fraccionarios intermedios, el extremo superior de la banda se situará entre estos dos valores. Una consecuencia directa de estas expresiones es que en el caso fraccionario la banda de fonones es más estrecha que en el caso clásico. De hecho, para α= 0 la banda es infinitamente estrecha formándose lo que se denomina una ‘flat band’. De aquí deducimos que necesitaremos acoplamientos mayores (Cmás grande) para llegar al caso en el que el extremo superior de la banda alcance un armónico del breather; y, por tanto, las inestabilidades causadas por las hibridaciones se retrasarán. CAPÍTULO 3. ESTABILIDAD DE SOLUCIONES 17 Este enfoque nos ha permitido obtener la relación de dispersión; sin embargo, cuando queremos obtener las frecuencias de los modos normales de nuestra red de forma numérica es más recomendable actuar de otra forma. El principal motivo es que en el cálculo previo consideramos una red infinita para poder llegar a expresiones analíticas cerradas, aunque en la práctica nuestra red siempre será finita. Por tanto, debemos proponer soluciones de la forma un=An·e−iωt (igual que antes pero incluyendo la fase relativa que introduce la posición dentro de la amplitud An). Con un desarrollo análogo llegamos a que se deben cumplir el siguiente conjunto de ecuaciones: An 1 + CX i=n Kα(|n−i|) −CX i=n Kα(|n−i|)Ai=ω2An,(3.5) para todo ndesde 1 hasta N. Este es un conjunto de ecuaciones lineales en los coeficientes An; de hecho, viendo los coeficientes como un vector  Ade Ncomponentes, hemos llegado a una ecuación de autovalores. Por tanto, las frecuencias al cuadrado de la banda de fonones no son más que los autovalores de la matriz de coeficientes del término de la izquierda de la ecuación (3.5). Es fácil ver que las frecuencias estarán bien definidas, ya que la matriz es hermítica y diagonal dominante (i.e. la suma de los elementos de fuera de la diagonal en valor absoluto de cualquier fila es menor que el elemento diagonal de dicha fila). La primera propiedad nos asegura que los autovalores son reales y la segunda que es una matriz semidefinida positiva y la parte real de sus autovalores será no negativa (como corolario del teorema de los círculos de Gershgorin). Por tanto, todos los autovalores son números reales no negativos y las frecuencias que obtendremos serán también reales no negativas. En la figura 3.1a hemos representado la banda de frecuencias de modos lineales para C= 1,5calculadas como autovalores del sistema (3.5). Se observa cómo la frecuencia aumenta conforme lo hace k(kpuede ser recuperado a partir de la periodicidad espacial de los autovalores del sistema (3.5)). La analogía con la física del estado sólido es habernos restringido a la mitad de la primera zona de Brilloin. En la figura 3.1b, hemos representado la dependencia del valor máximo de la banda de frecuencias con la constante de acoplamiento para distintos valores de α. Como habíamos predicho, la banda se ensancha conforme aumenta el valor del orden fraccionario. Además, hemos representado las expresiones analíticas que habíamos deducido para el caso clásico y el fraccionario (3.3) y (3.4), CAPÍTULO 3. ESTABILIDAD DE SOLUCIONES 18 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 k/ 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 Frecuencias Modos Normales C = 1.5 (a) Banda de frecuencias de los modos lineales de una red con laplaciano clásico y 41 nodos. 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 Constante de acoplamiento, C 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 3.0 Frecuencia Extremo superior de la banda de modos lineales Clásico (b) Extremo superior de la banda de frecuencias de una red de 41 nodos en función del acoplamiento para distintos valores de α. También se incluyen las expresiones analíticas (3.3) y (3.4). Figura 3.1 respectivamente. Se observa cómo ajustan de forma casi perfecta los valores numéricos obtenidos; las ligeras diferencias ocurren debido a que en la expresión numérica hemos usado redes de 41 nodos y las expresiones analíticas se dedujeron para redes infinitas. Las mayores diferencias se alcanzan en el caso de α= 0,5ya que la interacción es de más largo alcance y se aprecia más el hecho de que la red no sea infinita. 3.2. Análisis de Floquet Ya hemos presentado una de las causas por las que las soluciones pueden desestabilizarse. En esta sección presentaremos de forma sucinta el marco teórico que usaremos para determinar la estabilidad de las soluciones. Al ser nuestras soluciones periódicas usaremos la teoría de Floquet. Una solución se considera linealmente estable si, al perturbarla ligeramente, la perturbación crece a lo sumo linealmente con el tiempo. El elemento central del análisis es el operador de Floquet que se define como:    ξ(T) π(T)  =F   ξ(0) π(0)   ,(3.6) siendo ξla perturbación u→u+ξyπla derivada temporal de la perturbación (cada CAPÍTULO 3. ESTABILIDAD DE SOLUCIONES 19 uno de ellos un vector de Ncomponentes). Es decir, el operador de Floquet es aquel que a partir de la perturbación inicial te devuelve cómo ha evolucionado un periodo más tarde. En [12] se encuentran demostradas numerosas propiedades de este operador; nosotros nos restringiremos a mencionar y comentar las más importantes. En primer lugar, estudiaremos la ecuación diferencial que cumple la perturbación. Para ello, basta con sustituir en (2.6) la solución periódica (breather en nuestro caso) upor u+ξy quedarnos en primer orden de la perturbación: 0 =d2(un+ξn) dt2+ sin(un+ξn)−CX i=n Kα(|n−i|)(ui−un+ξi−ξn)≈ ≈d2(un+ξn) dt2+ sinun+ cos(un)·ξn−CX i=n Kα(|n−i|)(ui−un+ξi−ξn) = =d2ξn dt2+ cos(un)·ξn−CX i=n Kα(|n−i|)(ξi−ξn)≡ Nξ, (3.7) donde hemos aprovechado para definir operador de Newton. Por cómo la hemos obtenido, esta ecuación es lineal en ξy de esta manera el operador de Floquet es un operador lineal; su representación matricial se conoce con el nombre de matriz de monodromía (dimensiones 2N×2N). El principal resultado que usaremos de la teoría de Floquet es que nuestra solución será linealmente estable si ninguno de los autovalores del operador de Floquet tiene módulo mayor a 1. No entraremos en la demostración de este resultado, pero la idea es que si existiera un autovalor λcon módulo estrictamente mayor que uno, la perturbación asociada a dicho autovalor a lo largo de n periodos sería: ξ(n·T) = Fnξ(0) = λnξ. (3.8) Y, por tanto, el módulo de esta perturbación crecería de forma exponencial, provocando que nuestra solución no sea linealmente estable. En la sección 4 de [2] se encuentra recogida la propiedad de que el operador de Floquet en esta ecuación es simpléctico, lo que conlleva que si λes autovalor entonces λ−1tam- bién lo será. Además, el operador es real (pues en su definición no aparecen operaciones complejas); por tanto, si λes autovalor entonces λ∗también lo será. En resumen, es un operador cuyos autovalores aparecen agrupados por cuadrupletes {λ, λ−1, λ∗, λ−1∗}. Con esto y el teorema principal de la teoría de Floquet podemos enunciar el siguiente resultado CAPÍTULO 3. ESTABILIDAD DE SOLUCIONES 20 para los breathers: Un breather será linealmente estable si y solo si todos los autovalores de su operador de Floquet asociado tienen módulo igual a 1 (i.e. son de la forma λ=eiθ con θreal, donde θrecibe el nombre de argumento de Floquet). Por tanto, cuando un breather pasa de estable a inestable corresponde a que algún autovalor se ha salido de la circunferencia unidad. Pero, debido a la estructura que tiene el espectro no es posible que únicamente salga un autovalor, sino que lo harán los cuatro del mismo cuadruplete a la vez (o una pareja si lo hacen por θ=0óπ). Para ello, es necesario que dos de ellos colisionen en el mismo punto de la circunferencia unidad antes de salirse. Un concepto útil para entender cuándo pueden ocurrir estas colisiones es la signatura de Krein. A cada autovalor con argumento de Floquet θle podemos asignar su signatura de Krein como: κ(θ) = sgn([Re(Ω(t)),Im(Ω(t))]) = sgn"iX n (ξn(t)π∗ n(t)−ξ∗ n(t)πn(t))#,(3.9) siendo Ωel autovector correspondiente al autovalor asociado a dicho argumento de Floquet, y ξyπlas perturbaciones en las posiciones y velocidades, respectivamente. Gracias a la propiedad de que el operador de Floquet es simpléctico, se puede demostrar que la signatura de Krein está bien definida y no depende del instante de tiempo en que se evalúe. Con esta definición podemos enunciar el criterio de Krein que afirma que para que dos autovalores colisionen en la circunferencia unidad y se salgan de ella creando una inestabilidad es necesario que tengan signaturas de Krein opuestas (en el caso de que ambos sean reales, tendrán signatura nula los dos). A la hora de implementar este análisis de estabilidad de forma numérica en MATLAB, el procedimiento será el siguiente: Consideraremos las 2Nperturbaciones correspondientes los vectores con todas las componentes nulas menos una que será igual a uno. Para cada una de ellas, resolveremos numéricamente la ecuación (3.7) y evaluamos ξyπen t=T: Con estas evaluaciones por columnas construimos la matriz de monodromía a la que le calculamos los autovalores. CAPÍTULO 4. RESULTADOS 27 0 5 10 15 20 25 30 35 40 n -1.0 -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 1.0 un(0) Perfil inicial, C = 0.5 (a) Perfil inicial del breather obtenido con el método de disparo. 0 5 10 15 20 25 30 35 40 n 0.0 0.5 1.0 1.5 2.0 2.5 3.0 un(0) Perfil inicial, C = 0.5 (b) Perfil inicial del breather obtenido con el método de la serie de Fourier. Desplazamiento de nodos, C = 0.5 5 10 15 20 25 30 35 40 Nodos 0 5 10 15 Tiempo -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 (c) Desplazamiento de cada nodo a lo largo de dos periodos para el breather calculado con el método de disparo. Desplazamiento de nodos, C = 0.5 -20 -15 -10 -5 0 5 10 15 20 Nodos 0.0 5.0 10.0 15.0 Tiempo -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 (d) Desplazamiento de cada nodo a lo largo de dos periodos para el breather calculado con el método de la serie de Fourier. 0 5 10 15 20 25 30 35 40 n 0.0 0.5 1.0 1.5 2.0 2.5 3.0 un max Perfil en máxima elongación, C = 0.5 (e) Perfil de máxima elongación del breather calculado con el método de disparo Figura 4.1: Comparación de las soluciones obtenidas por los dos métodos numéricos presentados. CAPÍTULO 4. RESULTADOS 28 no lineal, el segundo al acoplamiento entre primeros vecinos del laplaciano clásico y el último al acoplamiento debido a las condiciones de contorno periódicas entre el primer y último nodo. En la figura 4.2 hemos representado la energía de cada solución frente a su constante de acoplamiento. Primeramente observamos una tendencia creciente. Este comportamiento es el esperado pues a mayor constante de acoplamiento más nodos son afectados por la perturbación y más energía tienen los breathers. Además observamos que hay algunas soluciones que tienen energía superior a la de sus adyacentes, saliéndose de la tendencia general. Estas soluciones son debidas a bifurcaciones; es decir, ocurren bifurcaciones en la rama de soluciones que estamos obteniendo y nos pasamos a otra rama. En la figura 4.3 observamos cómo sería la reconstrucción de todas las ramas, pudiendo identificar las bifurcaciones como bifurcaciones tipo pitchfork rotas (ya que las ramas no se cruzan). El hecho de que regresemos a la rama inferior se explica por la gran pendiente de la rama superior, lo que provoca que el algoritmo rápidamente encuentre soluciones por debajo de la rama intermedia (breather inestable) convergiendo a la inferior. Estas bifurcaciones se deben a la coincidencia de frecuencia de armónicos del breather con modos lineales, en nuestro caso con el tercer armónico. La razón es que el primero se encuentra por debajo del límite inferior de la banda (que está fijo), y debido a las simetrías de nuestro potencial (es espacialmente par i.e. V(−u) = V(u)) nuestra solución no presenta términos pares de la serie de Fourier (gracias a ser una función impar en el intervalo [0, T/2]). Si calculamos para qué constante de acoplamiento el extremo superior de la banda de fonones alcanza el tercer armónico de nuestro breather (que tenía frecuencia ωb= 0,8), obtenemos C0= 1,20 lo que es consistente con que la primera hibridación haya aparecido en C= 1,27. Para poder entender mejor en qué consisten estas bifurcaciones debidas a hibridaciones con los modos lineales hemos realizado la comparación en la figura 4.4. En las figuras 4.4a y 4.4b, observamos el perfil de dos breathers con constantes de acoplamiento parecidas (C= 4,80 yC= 4,82). El primero coincide con una de las soluciones de la figura 4.2 que presenta más energías que sus adyacentes. Esto se entiende al ver que el breather está hibridado con un modo lineal de vibración (i.e. superpuesto con una componente cosenoidal). Este tipo de soluciones reciban el nombre de nanóptero o breather fantasma. CAPÍTULO 4. RESULTADOS 29 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 Constante de acoplamiento, C 0 5 10 15 20 25 Energía Energía en función del acoplamiento Figura 4.2: Energía de cada solución respecto a la constante de acoplamiento. 0.80 0.85 0.90 0.95 1.00 1.05 1.10 1.15 Constante de acoplamiento, C 9 10 11 12 13 14 15 16 17 18 Energía Figura 4.3: Energía en función de la constante de acoplamiento para breathers con frecuencia ω= 0,7en una red de 30 nodos. Imagen generada a partir de datos suministrados por el tutor de este trabajo, Jesús Cuevas Maraver. CAPÍTULO 4. RESULTADOS 30 0 5 10 15 20 25 30 35 40 n -0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 un max Perfil en máxima elongación, C = 4.8 (a) Perfil de máxima elongación de breather hibridado, C = 4,80. 0 5 10 15 20 25 30 35 40 n 0.0 0.5 1.0 1.5 2.0 2.5 3.0 un max Perfil en máxima elongación, C = 4.82 (b) Perfil de máxima elongación de breather sin hibridar, C = 4,82. Desplazamiento de nodos, C = 4.8 5 10 15 20 25 30 35 40 Nodos 0 5 10 15 Tiempo -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 (c) Mapa de color del breather hibridado. Desplazamiento de nodos, C = 4.82 5 10 15 20 25 30 35 40 Nodos 0 5 10 15 Tiempo -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 (d) Mapa de color de un breather sin hibridar. Figura 4.4: Comparación entre breather hibridado (nanóptero) y breather normal. En el mapa de color 4.4c se observa además cómo la hibridación no está localizada sino que afecta a toda la red y tiene frecuencia distinta a la del breather. De hecho, estas pequeñas ‘alas’ (de ahí el nombre de nanóptero), vibran con frecuencia 3ωbal haber ocurrido la hibridación con el tercer armónico del breather. En oposición, en el mapa de color 4.4d de la solución sin hibridar se observa que el breather sí está localizado. El hecho de que al variar tan poco la constante de acoplamiento regresemos a soluciones sin hibridar se debe a que hemos saltado a una rama cercana. Habiendo presentado estas dos soluciones, veamos a continuación las diferencias que presentan los espectros de sus operadores de Floquet. Estos espectros los hemos representado en la figura 4.5, siendo las figuras 4.5a y 4.5b los espectros asociados a la solución hibridada y sin hibridar, respectivamente. La principal diferencia es que en el espectro asociado al nanóptero hay dos autovalores que se han salido de la circunferencia unidad, CAPÍTULO 4. RESULTADOS 31 -1.0 -0.5 0.0 0.5 1.0 Parte Real -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 Parte Imaginaria Autovalores Operador de Floquet (a) Autovalores del operador de Floquet para el breather con C = 4,80. Los asteriscos azules, cruces negras y círculos rojos denotan signatura de Krein 1, 0 y -1 respectivamente. -1.0 -0.5 0.0 0.5 1.0 Parte Real -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 Parte Imaginaria Autovalores Operador de Floquet (b) Autovalores del operador de Floquet para el breather con C = 4,82. Los asteriscos azules, cruces negras y círculos rojos denotan signatura de Krein 1, 0 y -1 respectivamente. Figura 4.5: Comparación del espectro de autovalores del operador de Floquet para un breather hibridado y otro sin hibridar. lo que nos indica que dicha solución no es linealmente estable. Hay que destacar que lo que introduce y presenta la inestabilidad es el modo lineal con el que se ha hibridado, no el breather en sí. Queremos mencionar que, como demostramos en el fundamento teórico, hemos obtenido el 1 como autovalor asociado a los modos de fase y crecimiento. De hecho, podemos hacer el estudio del espectro de autovalores del operador de Floquet para todas las soluciones en función de la constante de acoplamiento; lo mostramos en las figuras 4.6 y 4.8. La primera corresponde al módulo de los autovalores y la segunda a su argumento. Recordemos que por cada solución (i.e. valor de la constante de acoplamiento) tenemos 82 autovalores. Respecto al módulo, lo primero que observamos es que, como ya enunciamos, si λes autovalor, también lo es λ−1; esta es la causa de la cuasi-simetría vertical de la gráfica. En segundo lugar, cada vez que hay autovalores fuera de la línea de módulo unidad nos indica que la solución es inestable. Para valores altos del acoplamiento hay soluciones inestables correspondientes a hibridaciones con fonones de la red que concuerdan con una bifurcación o salto de rama en la gráfica 4.2. Sin embargo, aparecen también inestabilidades que no se corresponden con las hibridaciones comentadas. El primer tipo de estas otras inestabilidades corresponden a autovalores espurios que no siguen un patrón definido y que se desplazan de la posición de módulo 1 ligeramente (en torno a 1,006). Estas inestabilidades no son propias del breather en sí, sino de las colas, CAPÍTULO 4. RESULTADOS 32 0.0 0.5 1.0 1.5 2.0 2.5 3.0 Constante de acoplamiento, C 0.85 0.90 0.95 1.00 1.05 1.10 1.15 | | Evolución Módulo Autovalores Clásico Figura 4.6: Módulo de los autovalores del operador de Floquet en función de la constante de acoplamiento. y se deben al tamaño finito de la red. Cuando aumenta el tamaño de la red, el numero de inestabilidades de este tipo aumenta aunque su amplitud disminuye. En el límite de considerar la red infinita, estas inestabilidades tendrían amplitud cero y, por tanto, no existirían. Para más información sobre las inestabilidades provocadas por el tamaño finito de la red véase [19]. Estas inestabilidades, en la representación de autovalores de Floquet, se reflejan la mayoría de las veces como pequeñas bifurcaciones de Hopf (los autovalores no se separan por el eje real). Dichas inestabilidades comienzan a ocurrir cuando las bandas de autovalores se cruzan en θ=π, ya que antes el breather es estable, como mencionamos en el fundamento teórico. Debemos destacar que dichas inestabilidades son tan pequeñas que, para apreciar la cuasiperiodicidad que introducen, haría falta realizar simulaciones numéricas a lo largo de muchos periodos. Esto en nuestro método de cálculo no provoca ningún inconveniente ya que, al estar usando métodos de punto fijo, nos estamos restringiendo siempre a soluciones periódicas de frecuencia ωb= 0,8. Por otra parte, encontramos una bifurcación que ocurre en C= 0,81 con un patrón definido en el que conforme aumenta C, aumenta el módulo de los autovalores hasta un punto en el que comienza a decrecer. Esta inestabilidad está relacionada con la ruptura CAPÍTULO 4. RESULTADOS 33 Figura 4.7: Imagen de una bifurcación de intercambio de estabilidad en una red discreta con acoplamiento a primeros vecinos y un potencial tipo Morse: V(r) = De1−e−a(r−r0)2. Se muestran también los pérfiles de las soluciones de cada rama. Extraída de la sección 8 de [12]. Creative Commons Attribution License (CC BY). https://creativecommons.org/licenses/ by/4.0/ de la simetría en las soluciones. Nosotros trabajamos siempre con soluciones simétricas. Sin embargo, en este punto la bifurcación ocurre porque aparecen además dos soluciones inestables con simetría especular una respecto a la otra (en las que ya no hay simetría respecto al nodo central). Un ejemplo de las distintas ramas que aparecen en esta bifurcación tipo pitchfork se encuentra en la figura 4.7. En ella observamos cómo el breather simétrico pierde estabilidad, dando lugar además a dos ramas con breathers intermedios que, finalmente, se estabilizan en breathers con dos nodos inicialmente excitados (en nuestro caso inicial sólo es uno, el central); por esta razón también se denomina bifurcación de intercambio de estabilidad. Si perturbamos nuestras soluciones en la dirección del autovector asociado al autovalor con módulo mayor conseguiremos que el breather se mueva por la red obteniendo un breather móvil [3]. El cálculo y la existencia de breather móviles en distintos tipos de redes son objeto de estudio desde la década de los 90 del siglo pasado. Para más información sobre ello se puede consultar el capítulo 8 de [12]. Es importante mencionar que el autovalor que provoca la inestabilidad decrece en módulo hasta llegar de nuevo a 1, donde vuelve a estabilizarse la rama del breather simétrico. Sin embargo, no observamos el decrecimiento de forma continua hasta uno en la figura 4.6, ya que en esos CAPÍTULO 4. RESULTADOS 34 0.0 0.5 1.0 1.5 2.0 2.5 3.0 Constante de acoplamiento, C -3 -2 -1 0 1 2 3 Evolución Argumento Autovalores Clásico Figura 4.8: Argumento de los autovalores del operador de Floquet en función de la constante de acoplamiento. valores del acoplamiento estamos teniendo también bifurcaciones debidas a hibridaciones conllevando los consecuentes saltos de rama. Estos saltos de rama provocan que el autovalor no se mueva de forma continua, y, por ende, no observemos claramente su regreso a 1. Otro aspecto interesante de esta bifurcación es que no se produce por la colisión de dos autovalores de la banda de autovalores (es decir, los que se mueven en conjunto por la circunferencia unidad), sino por dos autovalores o modos localizados separados del comportamiento de la banda. Esto se observa claramente en la figura 4.8, donde los autovalores de la banda se corresponden con la evolución en línea rectas y los modos localizados con un comportamiento distinto. Observamos en C= 0,30 dos autovalores que se separan de la banda y comienzan a aproximarse al 1 (argumento 0) al que llegan en C= 0,80 (corresponden al comportamiento curvo que se aprecia). Estos al encontrarse en 1 colisionan y dan lugar a la bifurcación de ruptura de simetría antes explicada (nótese que en C= 0,80 llegan a 1 y la bifurcación aparece en C= 0,81, justo cuando se salen de la circunferencia unidad). Aprovechamos la figura 4.8, para fijarnos que en el límite anticontinuo los autovalores tienen argumento 0 (correspondiente a los modos de fase y crecimiento) o ±π 2, concordante con que los autovalores de la monodromía son ±icomo CAPÍTULO 4. RESULTADOS 35 vimos en el fundamento teórico. Además se aprecia bien cómo las bandas de autovalores surgen de estos últimos valores, mientras que los del 1 se mantienen fijos. 4.2. Resultados Fraccionarios Tras haber estudiado distintas propiedades de los breathers soluciones de la ecuación clásica de sine-Gordon, en esta sección procederemos, usando las mismas técnicas, al estudio de las soluciones en la ecuación con laplaciano de orden fraccionario. De esta manera podremos estudiar si existen y comprobar cómo el orden fraccionario ha afectado a sus características. En primer lugar, recordemos que la discretización del laplaciano de orden fraccionario se puede entender como una interacción de largo alcance, en la que el decaimiento del kernel es algebraico. Esto provoca que sea necesario trabajar con redes más largas para que en los extremos de la red la oscilación localizada haya decaído hasta la precisión de máquina. En nuestro caso, hemos pasado de los 41 nodos que utilizábamos hasta a redes de 750 nodos. Esto corresponde a ecuaciones diferenciales de 1500 variables (posiciones y velocidades) donde todos los nodos están acoplados entre sí. Desde el punto de numérico esto conlleva un gran aumento del coste computacional, por lo que sólo consideraremos órdenes fraccionarios mayores o iguales a 0,5(ya que de usar valores menores habría que aumentar aún más la longitud de la red). En este caso sí se aprecia que el método de cálculo de los breathers a partir de su serie de Fourier es mucho más rápido, por ende todas las soluciones que presentaremos en esta sección han sido calculadas de dicha forma. Lo primero que nos planteamos es que, para calcular un breather con constante de acoplamiento Cy orden fraccionario α, ahora podemos actuar de dos formas distintas. Por una parte, podemos fijar αy realizar la continuación en el acoplamiento como antes: desde C= 0 hasta el deseado. Sin embargo, también podemos partir de la solución ya calculada en el caso clásico con el Cdeseado y α= 1, y realizar la continuación en el orden fraccionario, es decir, ir variándolo lentamente usando cada solución como semilla para la solución con α−∆α(es negativo porque la solución de la que partimos tiene αmáximo e igual a 1). En la figura 4.9 hemos comparado las soluciones para el caso C= 1,5yα= 0,5siguiendo ambas rutas. La figura 4.9a corresponde al caso en el que hemos realizado la continuación en αy la figura 4.9b al caso de haber realizado CAPÍTULO 4. RESULTADOS 36 0 10 20 30 40 50 60 70 80 90 100 n 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 un max Perfil en máxima elongación, C = 1.5 (a) Breather con C= 1,5yα= 0,5, calculado a partir del breather clásico con C= 1,5y α= 1. 0 10 20 30 40 50 60 70 80 90 100 n 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 un(0) Perfil breather, C = 1.5 (b) Breather con C= 1,5yα= 0,5, calculado a partir del breather con C= 0 yα= 0,5. Figura 4.9: Comparación de un mismo breather calculado desde dos límites distintos. la continuación en el acoplamiento. Además, el primero de ellos está calculado usando el método de disparo mientras que el segundo usando el método de la serie de Fourier. Como observamos, ambos breathers son idénticos salvo que el de la figura 4.9a ha sufrido un desplazamiento temporal en el cálculo y no está representado en el instante inicial sino en el de máxima elongación (para poder compararlo correctamente). Por tanto, con estas imágenes dejamos patente la existencia de breathers en la ecuación de sine-Gordon fraccional (pues hemos sido capaces de calcularlos). También apreciamos un efecto de la fraccionalidad en las soluciones y es que la oscilación localizada afecta a un mayor número de nodos que en el caso clásico; basta comparar las imágenes de las figuras 4.9 y 4.4b. Esta es la causa de que en este caso nos sea necesario considerar redes más largas. Para apreciar más claramente este efecto, presentamos la figura 4.10 en la que hemos representado dos breathers para una misma constante de acoplamiento pero dos órdenes fraccionarios distintos. Dicha figura nos confirma nuestra intuición previa: a un orden fraccionario menor le corresponde una interacción de alcance más largo, provocando que una mayor cantidad de nodos se vean afectados. Para esa figura, hemos considerado una red de 200 nodos, adecuada para el breather con α= 0,95 pero observamos que no para el caso α= 0,5ya que si quisiéramos aumentar la constante de acoplamiento la solución no podríamos considerarla localizada pues no decaería lo suficiente antes de llegar a los extremos de la red. Por esto, trabajaremos a partir de ahora con una red de 750 nodos, siempre que necesitemos llegar a constantes de acoplamiento más altas (y lo CAPÍTULO 4. RESULTADOS 43 -1.0 -0.5 0.0 0.5 1.0 Parte Real -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 Parte Imaginaria Autovalores Operador de Floquet Figura 4.18: Autovalores del operador de Floquet para C= 1,58, uno de los breathers que presenta inestabilidades debido al tamaño finito de la red. Los asteriscos azules, cruces negras y círculos rojos denotan signatura de Krein 1, 0 y -1 respectivamente. diferencia es que las inestabilidades debidas a hibridaciones con los modos lineales (recordemos que sólo pueden aparecer por encima de C= 2,39) tienen mucha menor amplitud que antes (debida al tamaño del modo lineal hibridado), llegando a confundirse con las inestabilidades por el tamaño finito de la red. Estas últimas podríamos pensar que no debieran aparecer al estar considerando una red más larga (200 nodos en dicha gráfica); sin embargo, que el acoplamiento sea de mayor alcance ha provocado finalmente que aumenten en tamaño (alcanzando valores cercanos a módulo igual a 1,02, mayores que en el caso clásico). Sería necesario aumentar aún más el tamaño de la red para disminuir su amplitud. En la figura 4.18, hemos representado el espectro del operador de Floquet para el breather α= 0,5yC= 1,58, uno de los que en la figura 4.17 presenta este tipo de inestabilidades. Podemos observar cómo, además de la pareja que se han salido por el eje real (debido a la bifurcación por intercambio de estabilidad), también hay autovalores cercanos al uno fuera de la circunferencia. Estas son las bifurcaciones de Hopf que hemos comentado que introducen cuasiperiodicidades debidas al tamaño finito de la red. CAPÍTULO 4. RESULTADOS 44 Figura 4.19: Evolución del argumento de los autovalores del operador de Floquet para α= 0,5en función del acoplamiento. Se observa cómo dos autovalores localizados se acercan a argumento igual a 0 para colisionar en C= 1,10. Queremos destacar que en este caso sí se aprecia perfectamente el punto en el que la rama simétrica de la bifurcación por intercambio de estabilidad vuelve a estabilizarse, en C= 1,82. El hecho de que podamos apreciarlo se debe a que, en este caso, no estamos experimentando saltos de rama en el cálculo de nuestras soluciones (gracias a que las hibridaciones ocurren por encima de C= 2,39). Además, también podemos entender qué ocurre antes y después de esta bifurcación si observamos la evolución del argumento de los autovalores. En la figura 4.19, observamos cómo dos autovalores localizados convergen hasta colisionar en C= 1,10 en el valor 1 (argumento igual a 0). Estos autovalores al colisionar son los que provocan que surja la bifurcación. Además, en la figura 4.20, podemos observar que justo en C= 1,82, hay dos autovalores que se localizan desde el argumento igual a 0; estos son los que habían provocado la bifurcación por intercambio de estabilidad que al regresar a la circunferencia unidad por el eje real se localizan. De nuevo, podemos apreciar tan claramente este hecho gracias a que no estamos sufriendo saltos de rama de soluciones. En la figura 4.21 hemos realizado el mismo tipo de representación, pero en este caso con α= 0,95. De nuevo, observamos la bifurcación de intercambio de estabilidad, en este caso en C= 0,84 (más cercano al valor clásico, como esperábamos). De nuevo, al igual que CAPÍTULO 4. RESULTADOS 45 Figura 4.20: Evolución del argumento de los autovalores del operador de Floquet para α= 0,5 en función del acoplamiento. Se observa la localización de dos autovalores en C= 1,82 a partir de argumento igual a 0. en el caso clásico, no es posible seguir la trayectoria continua de los autovalores debido a los saltos de rama provocados por hibridaciones. Dichas hibridaciones en este caso sí son mayores, ya que para este orden fraccionario (y, en general, para órdenes fraccionarios cercanos a uno) los modos lineales con los que aparecen las hibridaciones tienen mayor amplitud (como ocurre en el caso clásico). Por último, podemos comentar que al seguir trabajando con una red de 200 nodos y una interacción de menor alcance que en el caso α= 0,5, las inestabilidades debidas al tamaño finito de la red han disminuido en amplitud (medidas en módulo de los autovalores que las provocan) desde valores cercanos a 1,02 (caso α= 0,5) hasta valores alrededor a 1,001 (caso α= 0,95), siendo hasta menores que en el caso clásico. De nuevo, en la figura 4.22 hemos hecho la representación del argumento de los autovalores, donde (al igual que en casos anteriores) podemos ver cómo la bifurcación se debe a dos autovalores que se localizan de la banda y colisionan en el 1. Este último caso, gracias a sus similitudes con el caso clásico nos permite establecer el acuerdo entre las soluciones con el laplaciano clásico y las soluciones asociadas a laplacianos de orden fraccionario pero cercano a uno. Para terminar de observar cómo recuperamos el caso clásico al acercarnos a α= 1, hemos representado en la figura 4.23 el acoplamiento para el cual ocurre la bifurcación relacionada con los breathers móviles en CAPÍTULO 4. RESULTADOS 46 Figura 4.21: Módulo de los autovalores del operador de Floquet en función de la constante de acoplamiento para α= 0,95. Figura 4.22: Argumento de los autovalores del operador de Floquet en función de la constante de acoplamiento para α= 0,95. Se observa cómo dos modos se localizan para colisionar en C= 0,83 CAPÍTULO 4. RESULTADOS 47 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00 Orden fraccionario, 0.80 0.85 0.90 0.95 1.00 1.05 1.10 1.15 Constante de acoplamiento, C Acoplamiento en el que ocurre la bifuración Figura 4.23: Acoplamiento en el que ocurre la bifurcación por aparición de soluciones no simétricas frente al orden fraccionario del laplaciano. función del orden fraccionario del laplaciano. Hemos tomado como criterio para considerar que ha ocurrido la bifurcación que un autovalor tenga módulo mayor a 1,05. Observamos una tendencia general de que a menor orden fraccionario más tarde ocurre la bifurcación (esto puede ocurrir debido a razones análogas a las que provocan que la hibridación con los modos lineales ocurra más tarde). Además, el último punto de la gráfica, con α= 1, ha sido calculado de forma clásica, lo que no es sino otra prueba más de que el paso del caso fraccionario al paso clásico es una evolución continua. Capítulo 5 Conclusiones Con este trabajo hemos logrado el objetivo principal de nuestra investigación que consistía en aprender sobre los breathers en la ecuación de sine-Gordon clásica y estudiar qué ocurría al considerar la ecuación de orden fraccionario. Llevarlo a cabo no ha sido un trabajo directo pues ha sido necesario formarse de distintos campos para poder aplicar este conocimiento en conjunto. Por una parte, hemos aprendido sobre distintos métodos numéricos que pueden ser implementados en MATLAB para resolver ecuaciones diferenciales, pudiendo compararlos entre ellos y decidir cuál era el mejor atendiendo al caso en el que nos encontráramos. Por otra parte, también ha sido necesario instruirse en el estudio de estabilidad de soluciones periódicas a ecuaciones diferenciales mediante la teoría de Floquet (que luego ha sido necesario implementar de forma numérica). Con estos ingredientes hemos podido abordar el estudio de los breathers en la ecuación de sine-Gordon clásica, siendo capaces de reproducir los resultados conocidos en la literatura: obtener los breathers y estudiar las distintas bifurcaciones que ocurren desde el punto de vista de los autovalores del operador de Floquet. Tras esto, hemos podido aplicar y adaptar un procedimiento análogo para estudiar el caso de la ecuación de orden fraccionario y obtener los resultados principales de este trabajo. Hemos sido capaces de constatar la existencia de breathers en la ecuación de sine- Gordon de orden fraccionario, comprendiendo el papel que juega la fraccionalidad como una interacción entre nodos de largo alcance. También se ha logrado estudiar el efecto que tiene en las inestabilidades, observando cómo siguen ocurriendo hibridaciones pero de menor amplitud y que sigue mantiéndose la bifurcación por intercambio de estabilidad re- 48 CAPÍTULO 5. CONCLUSIONES 49 lacionada con los breathers móviles. Por último, hemos conseguido calcular la dependencia con la fraccionalidad del acoplamiento en el que ocurre esta última bifurcación. En resumen, hemos logrado los objetivos que nos habíamos propuesto, sin embargo, este es un campo muy amplio donde aún quedan muchos frentes abiertos. Algunos de los posibles proyectos a realizar en el futuro podrían ser los siguientes: Estudio teórico de cómo afecta la fraccionalidad a las hibridaciones, demostrando que efectivamente se ven reducidas en amplitud como hemos observado en este trabajo. Cálculo de breathers móviles en la red con acoplamiento de largo alcance, investigando cómo dependen sus propiedades de la fraccionalidad. Al ocurrir menos hibridaciones, esperamos que las soluciones que se hallen de forma numérica radien menos energía y, por tanto, perduren más tiempo. Al igual que hemos considerado la fraccionalidad en el laplaciano, dando lugar a interacciones de largo alcance, también es posible considerar la derivada temporal de orden fraccionario. La interpretación sería la de un material con efecto de memoria, de forma parecida a los materiales que presentan histéresis. Bibliografía [1] Abdon Atangana y Dumitru Baleanu. «New fractional derivatives with nonlocal and non-singular kernel: theory and application to heat transfer model». En: arXiv preprint arXiv:1602.03408 (2016). [2] Serge Aubry. «Breathers in nonlinear lattices: Existence, linear stability and quantization». En: Physica D: Nonlinear Phenomena 103.1-4 (1997), págs. 201-250. [3] Serge Aubry y Thierry Cretegny. «Mobility and reactivity of discrete breathers». En: Physica D: Nonlinear Phenomena 119.1-2 (1998), págs. 34-46. [4] C Baesens, S Kim y RS MacKay. «Localised modes on localised equilibria». En: Physica D: Nonlinear Phenomena 113.2-4 (1998), págs. 242-247. [5] E Bour. Journ. de l’École Imperial Polytechnique 19. 1862. [6] Oleg M Braun y Yuri S Kivshar. «Nonlinear dynamics of the Frenkel–Kontorova model». En: Physics Reports 306.1-2 (1998), págs. 1-108. [7] Roland Bulirsch, Josef Stoer y J Stoer. Introduction to numerical analysis. Vol. 3. Springer, 2002. [8] Michele Caputo. «Linear models of dissipation whose Q is almost frequency independent—II». En: Geophysical Journal International 13.5 (1967), págs. 529-539. [9] Michele Caputo y José M Carcione. «Hysteresis cycles and fatigue criteria using anelastic models based on fractional derivatives». En: Rheologica acta 50.2 (2011), págs. 107-115. [10] Michele Caputo y Mauro Fabrizio. «A new definition of fractional derivative without singular kernel». En: Progress in Fractional Differentiation & Applications 1.2 (2015), págs. 73-85. 50 BIBLIOGRAFÍA 51 [11] Robert Conte. «Explicit breather solution of the nonlinear Schrödinger equation». En: Theoretical and Mathematical Physics 209.1 (2021), págs. 1357-1366. [12] Jesús Cuevas Maraver. «Localización y Transferencia de Energía en Redes Anarmónicas No Homogéneas». Tesis doct. 2003. [13] Jesús Cuevas-Maraver, Panayotis G Kevrekidis y Floyd Williams. «The sine-Gordon model and its applications». En: Nonlinear systems and complexity 10 (2014). [14] RK Dodd y RK Bullough. «Bäcklund transformations for the sine–Gordon equations». En: Proceedings of the Royal Society of London. A. Mathematical and Physical Sciences 351.1667 (1976), págs. 499-523. [15] Joshua Eby et al. «Global view of QCD axion stars». En: Physical Review D 100.6 (2019), pág. 063002. [16] Franz Gross. Relativistic quantum mechanics and field theory. John Wiley & Sons, 2008. [17] R Hilfer. «Fractional diffusion based on Riemann-Liouville fractional derivatives». En: The Journal of Physical Chemistry B 104.16 (2000), págs. 3914-3917. [18] RS MacKay y S Aubry. «Proof of existence of breathers for time-reversible or Hamiltonian networks of weakly coupled oscillators». En: Nonlinearity 7.6 (1994), pág. 1623. [19] JL Marín y S Aubry. «Finite size effects on instabilities of discrete breathers». En: Physica D: Nonlinear Phenomena 119.1-2 (1998), págs. 163-174. [20] Mario I. Molina. «The fractional discrete nonlinear Schrödinger equation». En: Physics Letters A 384.8 (2020), pág. 126180. [21] Mario I. Molina. «The two-dimensional fractional discrete nonlinear Schrödinger equation». En: Physics Letters A 384.33 (2020), pág. 126835. [22] Anna Maria Morgante et al. «Breather–phonon resonances in finite-size lattices:‘phantom breathers’?» En: Journal of Physics A: Mathematical and General 35.24 (2002), pág. 4999. [23] L Pilloni y D Levi. «The inverse scattering transform for solving the discrete sine- Gordon equation». En: Physics Letters A 92.1 (1982), págs. 5-8. BIBLIOGRAFÍA 52 [24] Harvey Segur y Martin D Kruskal. «Nonexistence of small-amplitude breather solutions in phi 4 theory». En: Physical review letters 58.8 (1987), pág. 747. [25] Igor M Sokolov, Joseph Klafter y Alexander Blumen. «Fractional kinetics». En: Physics Today 55.11 (2002), págs. 48-54. [26] Vasily E Tarasov. «Non-linear macroeconomic models of growth with memory». En: Mathematics 8.11 (2020), pág. 2078.