Sistemas dinámicos discretos unidimensionais: estabilidade, bifurcacións e aplicacións en bioloxía

Traballo Fin de Grao Sistemas dinámicos discretos unidimensionais: estabilidade, bifurcacións e aplicacións en bioloxía Borja Cabadas Dosantos Setembro, 2022 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA GRAO EN MATEMÁTICAS Traballo Fin de Grao Sistemas dinámicos discretos unidimensionais: estabilidade, bifurcacións e aplicacións en bioloxía Borja Cabadas Dosantos Setembro, 2022 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Agradecementos Quero agradecer a todos os profesores que tiven, tanto na universidade como fóra dela, por ter contribuído, aínda que uns máis que outros, ao desenvolvemento da miña curiosidade. Tamén quero agradecer ás miñas titoras deste Traballo Fin de Grao, Cristina Lois Prados e Rosana Rodríguez López, por guiarme nesta experiencia facendo un traballo que nunca será o sucientemente valorado. Por último, quero dar as grazas á miña familia e ás persoas máis próximas durante esta etapa educativa, polo apoio que me deron, e que me xeron ser o que son hoxe. iii Traballo proposto Área de Coñecemento: Análise Matemática Título: Sistemas dinámicos discretos unidimensionais: estabilidade, bifurcacións e aplicacións en bioloxía Breve descrición do contido As ecuacións diferenciais ordinarias serven para modelar o comportamento dunha especie cando esta evoluciona no tempo de maneira continua. Sen embargo, cando o crecemento da poboación é un proceso estacional, resulta máis axeitado traballar con sistemas dinámicos discretos unidimensionais. O obxectivo principal deste traballo é introducir ao alumno no estudo cualitativo dos sistemas dinámicos discretos unidimensionais e a súa aplicación. Un esquema aproximado dos contidos podería ser o seguinte. O primeiro paso sería establecer un marco teórico para levar a cabo o estudo da estabilidade local [7, Capítulo 1], estabilidade global [2, 5, 15] e as bifurcacións [16, Sección 3.2]. En segundo lugar, aplicaríase a teoría desenvolvida a modelos de poboación clásicos, coma o de Beverton-Holt ou o de Ricker. Os resultados teóricos acompañaríanse con grácas e simulacións numéricas realizadas co programa MATLAB. Para rematar, faríase unha interpretación dos resultados obtidos en termos biolóxicos. v Índice Agradecementos iii Resumo ix Introdución xi 1. Resultados preliminares 1 2. Conceptos e resultados sobre sistemas dinámicos discretos 5 2.1. Conceptos básicos e estabilidade local sobre puntos xos . . . . . . . . . . . . . . 6 2.2. Estabilidadeglobal................................... 9 2.3. Bifurcacións....................................... 22 3. Aplicacións a modelos biolóxicos 31 3.1. Descrición dalgúns modelos clásicos . . . . . . . . . . . . . . . . . . . . . . . . . . 31 3.1.1. Modelo de Beverton-Holt . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 3.1.2. ModelodeRicker................................ 33 3.2. Existenciadepuntosxos............................... 34 3.2.1. Modelo de Beverton-Holt . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 3.2.2. ModelodeRicker................................ 35 3.3. Análise da estabilidade dos puntos de equilibrio . . . . . . . . . . . . . . . . . . . 35 3.3.1. Modelo de Beverton-Holt . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 vii 2 1. Resultados preliminares Denición 1.4. Sexa I intervalo aberto e f:I−→ R unha función. Diremos que f é de clase r en I (denotado por f∈ Cr(I) ) se existe f(r)(x) e é continua para todo x∈I . O seguinte resultado, a regra da cadea, úsase para demostrar a fórmula da derivada Schwarziana dunha composición de funcións. Teorema 1.5 (Regra da cadea) . Sexan I , J intervalos en R e sexan g:I−→ R e f:J−→ R funcións tales que f(J)⊆I . Se f é derivable en J e g é derivable en f(J) , entón a función composta g◦f é derivable en J e (g◦f)′= (g′◦f)·g′. Corolario 1.6. Sexan I⊂R un intervalo e f:I−→ I . Se f é derivable en a∈I , entón, para cada n∈N , fn+1 =f◦n+1) ··· ◦f é derivable en a con (fn+1)′(a) = f′(fn(a)) ···f′(f2(a)) ·f′(f(a)) ·f′(a). Ademais, se f∈ Cr(I) , entón fn∈ Cr(I) . Demostración. Probamos o resultado facendo uso do método de indución. Comezamos co caso n= 1 , (f2)′(a) = f′(f(a)) ·f′(a) (regra da cadea para f=g ) . Supoñamos que, para n≥2 , (fn)′(a) = f′(fn−1(a)) ···f′(f2(a)) ·f′(f(a)) ·f′(a) e, a partir da fórmula para n e da regra da cadea, demostramos para n+ 1 , (fn+1)′(a) = (fn◦f)′(a) = ((fn)′◦f)(a)·f′(a) =f′(fn−1(f(a))) ···f′(f2(f(a))) ·f′(f(f(a))) ·f′(f(a)) ·f′(a) =f′(fn(a)) ···f′(f2(a)) ·f′(f(a)) ·f′(a). Isto proba que a propiedade é certa para todo n∈N . O seguinte resultado dinos que, para unha función continua nun intervalo pechado e derivable no correspondente intervalo aberto, existe un punto do intervalo aberto no que a recta tanxente é paralela a recta que pasa polos puntos da gráca correspondentes aos extremos do intervalo pechado, o que quere dicir que as rectas teñen a mesma pendente. Este teorema úsase para demostrar resultados de estabilidade local, estabilidade global e inestabilidade. 3 Teorema 1.7 (Teorema do valor medio) . Sexa f: [a, b]⊂R−→ R tal que f é continua en [a, b] e derivable en (a, b) , entón existe c∈(a, b) tal que f(b)−f(a) = f′(c)(b−a). Outro concepto que necesitamos denir é o de punto crítico, que se corresponde cun punto da curva no que a pendente é paralela ao eixe das abscisas. Con él veremos os conceptos de monotonía e extremos relativos, que serán relevantes para o estudo das funcións. Denición 1.8. Sexan I intervalo e f:I−→ R . Un punto crítico x é aquel que verica f′(x) = 0 . Denición 1.9. Sexan I un intervalo e f:I−→ R unha función. De xeito análogo, dise que f é ( estritamente ) crecente en I se, sempre que x1, x2∈I e x1(<)≤x2 , entón f(x1)(<)≤f(x2) . Dise que f é ( estritamente ) decrecente en I se, sempre que x1, x2∈I e x1(<)≤x2 , entón f(x1) (>)≥f(x2) . Se unha función é (estritamente) crecente ou (estritamente) decrecente en I , dicimos que é ( estritamente ) monótona en I . Denición 1.10. Sexan I un intervalo e f:I−→ R unha función. Dicimos que f ten un máximo ( mínimo ) relativo en c∈I se existe unha veciñanza V de c tal que f(x)≤(≥)f(c) para todo x∈V∩I . Dicimos que hai un extremo relativo en c∈I se f ten un máximo relativo ou un mínimo relativo en c . O seguinte resultado será de utilidade para estudar a monotonía e os extremos a través da primeira e segunda derivada. Proposición 1.11. Sexan I un intervalo e f:I−→ R unha función. Se f é derivable en I , f é (estritamente) crecente en I se, e só se, f′(x) (>)≥0 para todo x∈I . f é (estritamente) decrecente en I se, e só se, f′(x) (<)≤0 para todo x∈I . Se c∈˚ I , f é derivable en c e f ten un extremo relativo en c , entón f′(c) = 0 . Se f é continua en I= [a, b] , c∈˚ I , e f é derivable en (a, c)∪(c, b) , f ten un máximo (mínimo) relativo en c , se existe unha veciñanza (c−δ, c +δ)⊆I tal que f′(x)≥(≤) 0 para todo x∈(c−δ, c) , e f′(x)≤(≥) 0 para todo x∈(c, c +δ) . Alternativamente, se existe e é continua f′ nunha veciñanza de c∈˚ I tal que f′(c) = 0 , f ten un máximo (mínimo) relativo en c se f′′(c)<(>) 0 . 4 1. Resultados preliminares O último resultado que imos enunciar neste capítulo é o Teorema da función implícita, que se usará para deducir as condicións que caracterizan as bifurcacións. Teorema 1.12 (Teorema da función implícita) . Sexa un aberto Ω⊂R2 , (x0, y0)∈Ω e sexa unha función F: Ω ⊂R2−→ R tal que: 1. F∈ C1(Ω) , 2. F(x0, y0)=0 , 3. ∂F ∂y (x0, y0)= 0 . Baixo estas condicións, existen α, β > 0 e φ: [x0−α, x0+α]−→ [y0−β, y0+β] tales que: (i) φ(x0) = y0 , (ii) F(x, φ(x)) = 0 , para todo x∈[x0−α, x0+α] , (iii) se, para x∈[x0−α, x0+α] , existe y∈[y0−β, y0+β] tal que F(x, y) = 0 , entón y=φ(x) , (iv) φ é continua en [x0−α, x0+α] , φ∈ C1((x0−α, x0+α)) e ademais φ′(x) = −∂F ∂y (x, φ(x))−1∂F ∂x (x, φ(x)) para todo x∈(x0−α, x0+α). Capítulo 2 Conceptos e resultados sobre sistemas dinámicos discretos No que vén a continuación, traballaremos cun intervalo I⊂R (limitado ou non) e unha aplicación f:I−→ I continua, que nos dene o seguinte sistema dinámico discreto unidimensional xn+1 =f(xn), (2.1) onde n∈N∪{0} e x0∈I é coñecido. Como exemplo, na Figura 2.1, ilustramos o modelo de Beverton-Holt (3.8) e o modelo de Ricker (3.14) nos que afondaremos no seguinte capítulo. 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 (a) Modelo de Beverton-Holt dado por f(x) = rx 1+x , r= 2 . 0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5 (b) Modelo de Ricker, dado por f(x) = xer(1−x) , r= 1,5 . Figura 2.1: Exemplos de sistemas dinámicos discretos. 5 6 2. Conceptos e resultados sobre sistemas dinámicos discretos 2.1. Conceptos básicos e estabilidade local sobre puntos xos Os contidos desta sección pódense consultar en [7, Capítulo 1 ] e [8, Sección 2 ]. En primeiro lugar, presentamos os conceptos de punto xo e punto periódico, falamos tamén da órbita dun punto e exemplicamos un tipo de órbita particular, a órbita periódica. Denición 2.1. Dise que xf∈I é un punto xo de f se cumpre que f(xf) = xf . Dise que p∈I é un punto periódico de período m se cumpre que fm(x) = x . Dirase de período primo m cando m é o menor natural que cumpre a igualdade anterior. Denición 2.2. A órbita dun punto x0∈I é aquela formada polo conxunto de puntos fn(x0) , con n∈Z , sempre e cando exista a inversa de f e esta sexa continua. A nós interesaranos a semiórbita positiva ( O+ ), formada polos puntos fn(x0) tales que n∈N∪{0} . Exemplo 2.3. En particular, a semiórbita positiva dun punto xo é o propio punto e a dun punto periódico p está formada polas iteracións de p , referímonos a ela como órbita periódica . Para representar as órbitas dun sistema dinámico, faremos uso dos diagramas Cobweb . Para construír este tipo de diagramas, teremos que debuxar as curvas y=x e y=f(x) , cuxa intersección nos dá os puntos xos do sistema dinámico asociado. A continuación, dado un punto x0 , proxéctase unha recta vertical dende o punto (x0, x0) da diagonal ata a función para acadar a primeira iteración f(x0) . Agora, proxectamos unha recta horizontal dende o punto (x0, f(x0)) cara a diagonal, obtendo o punto (f(x0), f(x0)) . Continuando este proceso, teremos unha representación das iteracións de x0 representadas na diagonal y=x . Respresentamos estes pasos na Figura 2.2. Continuamos con algúns conceptos básicos sobre estabilidade, en particular, deniremos puntos estables, inestables, atractores locais e globais. Tamén veremos que é un punto hiperbólico, así como dous resultados relevantes que nos serven para deducir a estabilidade asintótica local ou a inestabilidade, para puntos xos hiperbólicos, en función do valor absoluto da derivada. Denición 2.4. Diremos que un punto xo xf∈I é: Un punto estable se, para cada veciñanza V de xf en I , existe unha veciñanza U de xf en I tal que fn(x)∈V , para todo x∈U, n ≥1 . Un punto inestable se xf non é estable. Un atractor local se existe un entorno V⊂I de xf tal que l´ım n→∞ fn(x) = xf , para todo x∈V . Un atractor global se l´ım n→∞ fn(x) = xf , para todo x∈I . 2.1. Conceptos básicos e estabilidade local sobre puntos xos 7 Un punto localmente asintoticamente estable (L.A.S.) se é estable e, ademais, atractor local. Un punto globalmente asintoticamente estable (G.A.S.) se é estable e, ademais, atractor global. 0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5 (a) Comezamos a primeira iteración coa recta vertical dende o punto (x0, x0) da identidade ata a función f(x) . 0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5 (b) Completamos a primeira iteración, engadindo a recta horizontal dende o punto (x0, f(x0)) da función f ata a identidade. 0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5 (c) A segunda iteración comeza cunha recta vertical dende o punto (f(x0), f(x0)) da identidade ata a función f(x) . 0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5 (d) Completamos a segunda iteración cunha recta horizontal dende o punto (f(x0), f2(x0)) da función f ata a identidade. Figura 2.2: Diagrama Cobweb paso a paso co punto inicial x0= 2 marcado cun asterisco vermello, a recta y=x representada en verde, a función de Ricker f(x) = xe1−x 2 para x∈[0,2,5] representada en azul e as iteracións representadas en vermello. Observación 2.5 . Os conceptos denidos para un punto xo xf pódense adaptar a un punto periódico p de período primo m , xa que este é un punto xo da aplicación fm . Ademais, un punto periódico de período primo m= 1 é un punto xo. 8 2. Conceptos e resultados sobre sistemas dinámicos discretos Denición 2.6. Sexa p un punto periódico de período primo m∈N . O punto p é hiperbólico se |(fm)′(p)| = 1 . Se |(fm)′(p)|= 1 , entón p dise non hiperbólico . As aplicacións con puntos periódicos hiperbólicos son máis doadas de analizar, veremos a continuación os dous resultados que nos proporcionan información sobre a dinámica local para este tipo de puntos. Teorema 2.7. Supoñamos que f∈ C1(I) e que xf é un punto xo hiperbólico tal que |f′(xf)|<1 . Entón, existe un intervalo aberto J⊂I tal que xf∈J e l´ım n→∞ fn(x) = xf,∀x∈J, (2.2) é dicir, xf é un atractor local. Demostración. Como f∈ C1(I) , xf∈I con |f′(xf)|<1 e a aplicación |·|:R−→ [0,∞) é continua, entón existen ε > 0 e A∈(0,1) tales que |f′(x)|<A<1 para todo x∈(xf−ε, xf+ε) =: J. Dado x∈J con x=xf , por ser A∈(0,1) , temos polo Teorema 1.7 (Teorema do valor medio) que: |f(x)−xf|=|f(x)−f(xf)|=|f′(cx,xf)||x−xf| ≤ A|x−xf|<|x−xf|< ε. Polo tanto, f(x)∈J e, de feito, está máis cerca de xf que o propio x . Polo mesmo argumento, para n∈N, n ≥2 : |fn(x)−xf|=|fn(x)−f(xf)|=|f(fn−1(x)) −f(xf)|=|f′(cfn−1(x),xf)| |fn−1(x)−xf| ≤A|fn−1(x)−xf| ≤ ··· ≤ AAn−1|x−xf|=An|x−xf|<|x−xf|< ε. Polo tanto, |fn(x)−xf| ≤ An|x−xf| para calquera n∈N , sendo x=xf e A∈(0,1) , e tamén se cumpre para x=xf (pois fn(xf) = xf para todo n∈N ), logo séguese que l´ım n→∞ An|x−xf|=|x−xf|l´ım n→∞ An= 0 e, entón, existe o límite l´ım n→∞ fn(x) = xf para todo x∈J , tal e como queríamos demostrar. O seguinte corolario é consecuencia directa do anterior teorema dado que un punto periódico de período m é un punto xo para fm . Corolario 2.8. Sexa f∈ C1(I) e p un punto periódico hiperbólico de período m con |(fm)′(p)|<1 . Entón existe un intervalo aberto J⊂I tal que p∈J e l´ım n→∞ (fm)n(x) = p , para todo x∈J . 2.2. Estabilidade global 9 A continuación, veremos que o comportamento da dinámica preto dos puntos xos nos que |f′(xf)|>1 é bastante distinto ao amosado no resultado anterior. Teorema 2.9. Sexa xf un punto xo hiperbólico con |f′(xf)|>1 . Entón existe un intervalo aberto J⊂I tal que xf∈J e de xeito que, para todo x∈J con x=xf , existe k > 0 tal que fk(x)/∈J. (2.3) Demostración. Como f∈ C1(I) , xf∈I con |f′(xf)|>1 e a aplicación |·|:R−→ [0,∞) é continua, entón existen ε > 0 e A > 1 tales que |f′(x)|>A>1 para todo x∈(xf−ε, xf+ε) =: J. Facemos a demostración por redución ao absurdo. Dado x∈J , x=xf , supoñamos que para todo k∈N se cumpre fk(x)∈J , entón, polo Teorema 1.7 (Teorema do valor medio), para cada n∈N , temos: |fn(x)−xf|=|fn(x)−f(xf)|=|f(fn−1(x)) −f(xf)|=|f′(cfn−1(x),xf)| |fn−1(x)−xf| > A |fn−1(x)−xf|>··· > AAn−1|x−xf|=An|x−xf|>|x−xf|. Como temos que |fn(x)−xf|> An|x−xf| para calquera n∈N , sendo x=xf e A > 1 , cúmprese: l´ım n→∞ |fn(x)−xf| ≥ l´ım n→∞ An|x−xf|=|x−xf|l´ım n→∞ An=∞. Polo tanto, existe algún k∈N tal que fk(x)/∈J= (xf−ε, xf+ε) . 2.2. Estabilidade global Esta sección recolle diversos resultados de estabilidade global para aplicacións cun único punto xo (positivo). Comezamos establecendo un par de resultados cuxa hipótese principal é a ben coñecida ausencia de órbitas de período 2 . O primeiro deles, debido a Coppel [4], é válido cando I é un intervalo pechado e limitado. O segundo, recollido en [8], é aplicable para calquera intervalo de R , tamén para os non limitados, pero necesita engadir unha hipótese adicional. Teorema 2.10 ([4]) . Sexan I= [a, b] e f: [a, b]−→ [a, b] unha función continua e limitada cun único punto xo xf no intervalo [a, b] . Temos que xf é un atractor global en [a, b] se, e só se, f non ten órbitas de período 2 en [a, b]\{xf} . Teorema 2.11 ([8]) . Sexa I⊂R un intervalo (non necesariamente limitado) e f:I−→ I unha función continua cun punto xo xf . Tense que xf é un atractor global se, e só se, f2(x)=x para todo x∈I\{xf} (é dicir, f non ten órbitas de período 2) e existe k≥1 tal que xf é estable para fk . 10 2. Conceptos e resultados sobre sistemas dinámicos discretos Observación 2.12 . Non incluímos aquí a demostración completa dos Teoremas 2.10 e 2.11. Porén, xusticamos que a non existencia de órbitas de período 2 é condición necesaria para que xf sexa un atractor global. En efecto, sexa xf∈I . Como f(xf) = xf , dedúcese que f2(xf) = f(f(xf)) = f(xf) = xf . Supoñamos que a ecuación f2(x) = f(f(x)) ten unha solución p∈I que non é solución de f(x) = x . Se tomamos x0=p , temos que x1=f(x0) = f(p)=p , pero x2=f(x1) = f(f(p)) = p , e os seguintes termos da sucesión {xn} repiten esta secuencia. Temos entón que a sucesión {xn} oscila e non converxe a xf , polo tanto, xf non pode ser un atractor global. A condición establecida nos resultados previos sobre a ausencia de órbitas de período 2 non é sinxela de comprobar. Establecemos a continuación outros resultados que serán de maior utilidade na práctica. Será na demostración dalgún destes novos resultados onde necesitaremos os Teoremas 2.10 e 2.11. O primeiro resultado, que se pode consultar en [2], require unhas condicións moi sinxelas, que aparecen ilustradas na Figura 2.3 e fan referencia á posición da función f con respecto á recta y=x e ao punto xo xf . Teorema 2.13 ([2]) . Sexan I= [0,∞) e f: [0,∞)−→ [0,∞) unha función continua tal que f(0) = 0 . Se xf é o único punto xo positivo de f e se cumpre que x < f(x)< xf en (0, xf) e x>f(x)> xf en (xf,∞) , entón xf é un atractor global en (0,∞) . Demostración. Sexa x0∈(0,∞) , probaremos que xn=fn(x0) converxe a xf cando n tende a ∞ . Distinguiremos 2 casos: Se, para algún j∈N , temos 0< xj≤xf , entón, por hipótese, para todo n≥j tamén se satisfai 0< xn< xf . Ademais, xa que f(x)> x para 0< x < xf , dedúcese que {xn}∞ n=j é crecente e acotada superiormente, polo tanto, ten un límite d , que cumpre a propiedade 0< d ≤f(d)≤xf . Pola Proposición 1.2, tomando os límites en ambos lados en xn+1 =f(xn) , obtemos d=f(d) , xa que f é continua. Por hipótese, xf é o único punto xo positivo de f , entón l´ım n→∞ xn=d=xf . Se xn> xf para calquera n∈N , entón xn é decrecente e acotada inferiormente pois x > f(x)> xf en (xf,∞) , polo tanto, ten un límite d que cumpre d≥f(d)≥xf . De novo, pola Proposición 1.2 e a continuidade de f , temos que f(d) = d ( d é un punto xo) e como, por hipótese, só hai un punto xo positivo, xa temos que l´ım n→∞ xn=d=xf . 2.2. Estabilidade global 11 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 2.5 3 Figura 2.3: Ilustración das hipóteses do Teorema 2.13 mediante a función f(x) = xe(1−x) 2 , x∈[0,∞) , en cor azul, así como do comportamento asintótico de dúas condicións iniciais que converxen ao punto xo positivo xf= 1 . En cor verde está representada a identidade, en vermello as iteracións do diagrama de Cobweb e os asteriscos vermellos indican as condicións iniciais. O seguinte teorema é unha pequena modicación dun resultado que se recolle en [5]. Neste, engádese unha hipótese que nos facilita os argumentos da demostración. Esta inclúese para poder utilizar o Teorema 2.11, en lugar do Teorema 2.10, pois este último é válido para intervalos pechados e limitados, mentres o novo resultado permite traballar en calquera intervalo e, en particular, en [0,∞) . As hipóteses que se requiren neste resultado tamén teñen que ver coa posición da función f con respecto á recta y=x , así como coa existencia dunha función adicional que envolve a función f . O concepto de función envolvente defínese a continuación e ilústrase na Figura 2.4. Denición 2.14. Sexan I= [0,∞) , ϕ: [0,∞)−→ [0,∞) e f: [0,∞)−→ [0,∞) tal que xf é o seu único punto xo positivo. Dicimos que a función ϕ é unha función envolvente de f ou que ϕ envolve á función f se (i.1) ϕ(x)> f(x) , para todo x∈(0, xf) , (i.2) ϕ(x)< f(x) , para todo x>xf tal que ϕ(x)>0 e f(x)>0 . Teorema 2.15. Sexan I= [0,∞) e f: [0,∞)−→ [0,∞) unha función continua que satisfai f(0) = 0 , xf é o seu único punto xo positivo, existe n≥1 tal que xf é estable para fn e, ademais, (a) f(x)> x , para todo x∈(0, xf) , (b) f(x)< x , para todo x∈(xf,∞) , 18 2. Conceptos e resultados sobre sistemas dinámicos discretos Por contrarrecíproco, supoñemos que f(x)< x , para todo x∈(a, xf) e f(x)> x , para todo x∈(xf, b) . Como tamén se cumpre que f(xf) = xf e f′(xf)=1 , necesariamente existe ε > 0 tal que f′(x)>1 para todo x∈(xf−ε, xf+ε)\{xf} . Polo tanto, f′ ten un mínimo local positivo en xf , o que supón unha contradición coa conclusión da Proposición 2.19. Supoñamos que se cumpre algunha das dúas condicións recollidas en (2.7): Sexa a∈[0,+∞) o valor máis pequeno para o que se cumpre a desigualdade f(x)> x , para todo x∈(a, xf) . Por como o denimos, f(a) = a e, usando o Teorema 1.7 (Teorema do valor medio) en [a, xf] , existe ca,xf∈(a, xf) tal que f(xf)−f(a) = f′(ca,xf)(xf−a) . Como a e xf son puntos xos, f′(ca,xf) = 1 . Entón, como f(x)> x en (a, xf), f′(a)≥1 para que f(x) estea por enriba da recta y=x para todo x∈(a, xf) . Usando agora que f′(ca,xf) = f′(xf)=1 e f′(a)≥1 , ao ser f′ continua e f non ter mínimos locais positivos, necesariamente existe xc∈(a, xf) tal que f′(xc)=0 . Tomamos o punto crítico xc∈(a, xf) máis próximo a xf e vemos que l´ım n→∞ fn(x) = xf , para todo x∈[xc, xf) . Como x<f(x) para todo x∈[xc, xf)⊂(a, xf) e f′(x)>0 para todo x∈(xc, xf) (pois f′(xc) = 0 , f′(xf) = 1 e xc é o punto crítico máis cercano a xf ), temos que x<f(x)< xf , para todo x∈[xc, xf) . Para cada x∈[xc, xf) , consideramos a sucesión {fn(x)}n∈N , que é monótona crecente e limitada superiormente por xf , polo tanto, é converxente e existe L∈(xc, xf] tal que l´ımn→∞ fn(x) = L . Agora ben, como f é continua, facendo uso da Proposición 1.2, temos que L= l´ım n→∞ fn(x) = l´ım n→∞ f(fn−1(x)) = fl´ım n→∞ fn−1(x)=f(L). En consecuencia, f(L) = L , é dicir, L é punto xo, pero o único punto xo en (xc, xf] é xf , entón L=xf . Consideramos agora b∈[0,+∞) o valor máis grande para o que se cumpre f(x)< x , para todo x∈(xf, b] . Por como o denimos, f(b) = b e, usando o Teorema 1.7 (Teorema do valor medio) en [xf, b] , existe cxf,b ∈(xf, b) tal que f(b)−f(xf) = f′(cxf,b)(b−xf) . Como b e xf son puntos xos, f′(cxf,b) = 1 . Ademais, como f(x)< x en (xf, b), f′(b)≥1 para que f(x) estea por debaixo da recta y=x para todo x∈(xf, b) e que f(b) = b . Como f′(cxf,b) = f′(xf) = 1 e f′(b)≥1 ao ser f′ continua e f non ter mínimos locais positivos, existe xc∈(xf, b) tal que f′(xc)=0 . Tomamos o punto crítico xc∈(xf, b) máis próximo a xf , vexamos que l´ım n→∞ fn(x) = xf , para todo x∈(xf, xc] . Como x>f(x) para todo x∈(xf, xc]⊂(xf, b) e f′(x)>0 para todo x∈(xf, xc) (pois f′(xc) = 0 , f′(xf) = 1 e xc é o punto crítico máis próximo a xf ), temos que x>f(x)> xf , para todo x∈(xf, xc] . Para cada x∈(xf, xc] , consideramos a sucesión {fn(x)}n∈N , que é monótona decrecente e limitada inferiormente por xf , polo tanto, é converxente e existe L∈[xf, xc) tal que 2.2. Estabilidade global 19 l´ım n→∞ fn(x) = L . Agora ben, como f é continua, volvendo a utilizar a Proposición 1.2, temos que L= l´ım n→∞ fn(x) = l´ım n→∞ f(fn−1(x)) = fl´ım n→∞ fn−1(x)=f(L). En consecuencia, f(L) = L , é dicir, L é punto xo, pero o único punto xo en [xf, xc) é xf , entón L=xf . Como consecuencia dos Lemas 2.21 e 2.23, temos o resultado que buscabamos, que permitirá abordar o estudo da estabilidade global en modelos con propiedades similares ás do modelo de Ricker. Teorema 2.24. Sexa f: [0,+∞)−→ [0,+∞), f ∈ C3([0,+∞)) , unha función cun único punto crítico xc e tal que: (i.1) f(0) = 0 , f(x)>0 en (0,∞) , (i.2) l´ım x→+∞f(x) = 0 , (i.3) f′(x)>0 en (0, xc) , (i.4) f′(x)<0 en (xc,∞) , (i.5) f′′(x)<0 en (0, xc) , (i.6) Sf(x)<0 en [0,+∞)\{xc} . Se f ten un único punto xo xf∈(0,+∞) e se cumpre que |f′(xf)| ≤ 1 , entón xf é un atractor global en (0,+∞) . Demostración. Consideramos por separado os casos |f′(xf)|<1 e o caso |f′(xf)|= 1 , para os que usaremos, respectivamente, o Lema 2.21 e o Lema 2.23. Supoñamos que |f′(xf)|<1 , polo Teorema 2.7, temos que xf é un atractor local. Como, respectivamente, por (i.3) e (i.4) , f′(x)>0 con x∈(0, xc) e f′(x)<0 con x∈(xc,∞) , entón xc é un máximo. En consecuencia, para cada par de valores a , b∈[0,+∞) tales que 0< a < m´ın{xc, xf} ≤ m´ax{xc, xf}< b , tendo en conta a condición (i.6) , f|[a,b] cumpre as condicións do Lema 2.21, así xf é un atractor global en [a, b] . Agora ben, para probar que xf é un atractor global en (0,+∞) , facemos uso das propiedades (i.1)  (i.3) e (i.5) . 20 2. Conceptos e resultados sobre sistemas dinámicos discretos Usando as condicións (i.1) , (i.3) e (i.5) , como xf é punto xo, existe ε > 0 tal que f(x)> x para todo x∈(0, ε) . Polo tanto, tomando a=ε 2 , temos asegurado que para todo x0∈(0, ε) existe k∈N tal que fk(x0)∈[a, xf]⊂[a, b] . Así, a sucesión {fn(x0)}n∈N converxe a xf . Usando agora a condición (i.2) , l´ım x→∞ f(x) = 0 , temos que para todo ˜ε, existe M˜ε tal que, se x>M˜ε, entón |f(x)|<˜ε. Tomando o valor de ε xado no punto previo, podemos considerar b=Mε , dado que, se x0∈(Mε,+∞) , entón f(x0)∈(0, ε) e estamos nas condicións do punto anterior. Supoñamos agora que |f′(xf)|= 1 , f cumpre as hipóteses do Lema 2.23 e temos: (a) se xc< xf , l´ım n→∞ fn(x) = xf , para todo x∈[xc, xf) , (b) se xc> xf , l´ım n→∞ fn(x) = xf , para todo x∈(xf, xc] . Sexa x0∈(0,+∞) arbitrariamente xado, probemos que {fn(x0)}n∈N converxe a xf en cada un destes 2 casos: se xc< xf , polas condicións (i.1)  (i.5) , temos unha función cunha gráca similar á da Figura 2.5. 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 Figura 2.5: Ilustración do diagrama de Cobweb para o modelo de Ricker con r= 1,8 , xc< xf , onde a recta verde representa á identidade, a recta azul á función de Ricker, o asterisco vermello é o punto inicial x0= 0,5 e as liñas vermellas son as iteracións. Consideramos o punto f−1 +(xc) = m´ax{x∈[0,+∞) : f−1(x) = xc} . Por (a), podemos armar que xf é un atractor global en [xc, f−1 +(xc)] . 2.2. Estabilidade global 21 Sexa agora x0∈(0, xc) , pola forma que ten a gráca de f e a denición de f−1 +(xc) , existirá k > 0 tal que fk(x0)∈[xc, f−1 +(xc)] (tal e como se ilustra na Figura 2.5). En consecuencia, {fn(x0)}n∈N converxe a xf . Por último, se x0∈(f−1 +(xc),+∞) , temos que f(x0)∈(0, xc) e, utilizando o caso anterior, concluímos que {fn(x0)}n∈N converxe a xf . Como (0,+∞) = (0, xc)∪[xc, f−1 +(xc)] ∪(f−1 +(xc),+∞) , podemos concluír que xf é un atractor global en (0,+∞) . se xc> xf , polas condicións (i.1)  (i.5) , temos unha función cuxa representación gráca é como a da Figura 2.6. 012345678910 0 0.5 1 1.5 Figura 2.6: Ilustración do diagrama de Cobweb para o modelo de Ricker con r= 0,7 , xc> xf , onde a recta verde representa á identidade, a recta azul á función de Ricker, o asterisco vermello é o punto inicial x0= 0,5 e as liñas vermellas son as iteracións. Por (b), podemos armar que xf é un atractor global en [xf, xc] . Sexa agora x0∈(0, xf) , como xc> xf , pola forma que ten a gráca de f (ilustrada na Figura 2.6, tense que {fn(x0)}n∈N é monótona crecente e limitada superiormente por xf . Seguindo argumentos similares aos de demostracións previas, tense que {fn(x0)}n∈N converxe a xf . Por último, se x0∈(xc,+∞) , temos que f(x0)∈(0, xc) e, utilizando os casos anteriores, concluímos que {fn(x0)}n∈N converxe a xf . Como (0,+∞) = (0, xf)∪[xf, xc]∪(xc,+∞) , podemos concluír que xf é un atractor global en (0,+∞) . 22 2. Conceptos e resultados sobre sistemas dinámicos discretos 2.3. Bifurcacións Os contidos recollidos nesta sección foron consultados en [16]. Supoñamos que a función f:I−→ I ( I⊂R intervalo) que dene o sistema dinámico xn+1 =f(xn) depende dun parámetro r∈R . Expresemos dita dependencia escribindo xn+1 =f(xn) = f(xn, r) ou f(x) = f(x, r) . Denotamos un punto xo de f por (xf, rf) , que cumpre que f(xf, rf) = xf . Para os distintos valores de r∈R , xorden naturalmente as seguintes cuestións: 1. Os puntos xos de f(x, r) son L.A.S. ou inestables? 2. Como se ve afectada a estabilidade cando r varía? Para comezar a responder estas preguntas, podemos realizar un estudo da primeira derivada da función no punto xo, pois, na Sección 2.1, vimos que o comportamento local do sistema dinámico preto dun punto xo (xf, rf) queda totalmente determinado se este é hiperbólico, é dicir, se |f′(xf, rf)| = 1 . En concreto, se |f′(xf, rf)|<1 , entón (xf, rf) é un atractor local e, se |f′(xf, rf)|>1 , entón (xf, rf) é inestable. En consecuencia, a estabilidade dun punto xo (xf, rf) (ou dun punto periódico) adoita verse afectada nos puntos non hiperbólicos, é dicir, cando |f′(xf, rf)|= 1 . No que vén a continuación, imos describir os cambios no comportamento a longo prazo das solucións da ecuación xn+1 =f(xn, r) nunha veciñanza dun punto xo non hiperbólico. Distinguimos dous casos: f′(xf, rf) = 1 e f′(xf, rf) = −1 . A estes cambios na dinámica asintótica referirémonos como bifurcacións, sendo puntos de bifurcación aqueles nos que se producen. Para cada tipo de cambio ou bifurcación, comezamos proporcionando a descrición e continuamos cun exemplo de función que a cumpre. Ademais, no caso da bifurcación de tipo sela-nodo, deduciremos as condicións analíticas necesarias para que se produza. Nas restantes situacións poderían facerse razoamentos análogos, facendo uso do Teorema 1.12 (Teorema da función implícita). Comezamos co caso f′(xf, rf) = ∂f ∂x (xf, rf) = 1 , no que poden aparecer 3 tipos de bifurcacións diferentes: 1. Bifurcación de tipo sela-nodo Unha bifurcación de tipo sela-nodo ocorre cando dous puntos xos coexistentes (un deles localmente asintoticamente estable e outro inestable) coliden nun único punto xo xf tal que f′(xf)=1 e despois desaparecen. 2.3. Bifurcacións 23 Exemplo 2.25. Considerando a función (x, r)−→ f(x, r) = x+r−x2 , temos que ∂f ∂x (x, r)=1−2x , polo que (xf, rf) = (0,0) é un punto xo non hiperbólico xa que f(0,0) = 0 e ∂f ∂x (0,0) = 1 . Estamos interesados nos puntos xos nunha veciñanza de (xf, rf) = (0,0) . Sabemos que os puntos xos son da forma f(x, r)−x=r−x2= 0 , é dicir, x1 f= +√r, x2 f=−√r . Vemos que os puntos xos existen se, e só se, r≥0 , coincidindo se r= 0 . Ademais, substituíndo en ∂f ∂x (x, r) , temos que ∂f ∂x(x1 f, r)=1−2√r < 1, para todo r > 0, ∂f ∂x(x2 f, r)=1−2(−√r) = 1 + 2√r > 1, para todo r > 0. En consecuencia, para r > 0 sucientemente pequeno, x1 f é un atractor local e x2 f é inestable. Deste xeito, acabamos de comprobar que existe unha bifurcación de tipo sela-nodo en (xf, rf) = (0,0) , como se representa na Figura 2.7. 5 10 15 20 r -4 -3 -2 -1 1 2 3 4 x Figura 2.7: Ilustración da bifurcación de tipo sela-nodo para f(x, r) = x+r−x2 , en vermello represéntanse os puntos xos inestables ( xf=−√r para r > 0 ) e, en azul, os localmente asintoticamente estables ( xf= +√r para r > 0 ). Para rematar, vexamos as condicións analíticas para que se dea unha bifurcación de tipo sela-nodo. Teorema 2.26. Sexa f:R2−→ R , (x, r)7−→ f(x, r) unha función en C2(R2) . Existe unha bifurcación de tipo sela-nodo no punto xo (xf, rf) = (0,0) se f cumpre as seguintes condicións: f(0,0) = 0,∂f ∂x(0,0) = 1, (2.8) ∂f ∂r (0,0) = 0,∂2f ∂x2(0,0) = 0. (2.9) 24 2. Conceptos e resultados sobre sistemas dinámicos discretos Coas condicións de (2.9) , pode ocorrer que −∂2f ∂x2(0,0) ∂f ∂r (0,0) >0 (2.10) ou que −∂2f ∂x2(0,0) ∂f ∂r (0,0) <0. (2.11) En cada un destes casos debe cumprirse algunha das seguintes condicións adicionais para r > 0 ou r < 0 , respectivamente:  ∂f ∂x(0,0)<1, x > 0 e  ∂f ∂x(0,0)>1, x < 0 (2.12) ou  ∂f ∂x(0,0)>1, x > 0 e  ∂f ∂x(0,0)<1, x < 0. (2.13) Demostración. Consideramos a función auxiliar h:R2−→ R dada por h(x, r) = f(x, r)−x . Os puntos xos da función f cumpren h(x, r)=0 , polo tanto, h(0,0) = 0 . Por outra banda, da primeira condición en (2.9), temos que ∂h ∂r (0,0) = ∂f ∂r (0,0) = 0 e o Teorema 1.12 (Teorema da función implícita) implica que unha única curva de puntos xos de f pasa por (xf, rf) = (0,0) ; ademais, para x sucientemente pequeno, esta curva de puntos xos pode representarse como un gráco sobre a variable x , é dicir, existe unha función única r(x)∈ C2(−ε, ε) , con ε sucientemente pequeno, de xeito que h(x, r(x)) = f(x, r(x)) −x= 0. (2.14) As condicións (2.8) e (2.9) permiten xusticar que ∂r ∂x(0) = 0 e ∂2r ∂x2(0) = 0, (2.15) o que asegura a existencia de dúas curvas de puntos xos como as da Figura 2.8. Facemos entón uso das condicións (2.8) e (2.9) para probar que se cumpre (2.15). En primeiro lugar, derivando implicitamente con respecto a x en (2.14), obtemos: ∂f ∂x(x, r(x)) −1 + ∂f ∂r (x, r(x)) ·∂r ∂x(x) = 0, (2.16) que, avaliando en x= 0 e despexando, implica que: ∂r ∂x(0) = −(∂f ∂x (0,0) −1) ∂f ∂r (0,0) = 0. 2.3. Bifurcacións 25 5 10 15 20 r -4 -3 -2 -1 1 2 3 4 x (a) ∂2r ∂x2(0) >0 . -20 -15 -10 -5 r -4 -3 -2 -1 1 2 3 4 x (b) ∂2r ∂x2(0) <0 . Figura 2.8: Ilustración das curvas de puntos xos r(x) , para x sucientemente pequeno, que verican: r(0) = 0,∂r ∂x (0) = 0 e ∂2r ∂x2(0) = 0 . En segundo lugar, derivando implicitamente con respecto a x en (2.16) tense: ∂2f ∂x2(x, r(x)) + ∂2f ∂x∂r(x, r(x)) ·∂r ∂x(x) + ∂2f ∂r∂x(x, r(x)) ·∂r ∂x(x) +∂2f ∂r2(x, r(x)) ·∂r ∂x(x)2 +∂f ∂r (x, r(x)) ·∂2r ∂x2(x) = 0, que, avaliando en x= 0 e usando que ∂r ∂x (0) = 0 , permite deducir: ∂2f ∂x2(0,0) + ∂f ∂r (0,0) ·∂2r ∂x2(0) = 0 e, despexando, obtense: ∂2r ∂x2(0) = −∂2f ∂x2(0,0) ∂f ∂r (0,0) = 0. Unha vez probado que as curvas de puntos xos son da forma representada na Figura 2.8, é claro que as condicións (2.12) ou (2.13) implican os requirimentos de estabilidade da bifurcación de tipo sela-nodo. 2. Bifurcación transcrítica Unha bifurcación transcrítica ocorre cando dous puntos xos coexistentes (un deles localmente asintoticamente estable e outro inestable) coliden nun único punto xo xf tal que f′(xf)=1 e despois intercambian a súa estabilidade. Exemplo 2.27. Considerando a función (x, r)−→ f(x, r) = x+rx −x2 , temos que ∂f ∂x (x, r) = 1 + r−2x , polo que é fácil vericar que (xf, rf) = (0,0) é un punto xo non hiperbólico xa que f(0,0) = 0 e ∂f ∂x (0,0) = 1 . Estamos interesados nos puntos xos nunha veciñanza de (xf, rf) = (0,0) . Sabemos que 26 2. Conceptos e resultados sobre sistemas dinámicos discretos os puntos xos son da forma f(x, r)−x=rx −x2= 0 , é dicir, x1 f= 0, x2 f=r . Vemos que os puntos xos existen para todo r∈R , coincidindo se r= 0 . Ademais, substituíndo en ∂f ∂x (x, r) , temos que ∂f ∂x(x1 f, r) = 1 + r   <1 se r < 0, >1 se r > 0, ∂f ∂x(x2 f, r)=1−r   <1 se r > 0, >1 se r < 0. En consecuencia, para r > 0 sucientemente pequeno, x1 f é inestable e x2 f é un atractor local. Porén, para r < 0 sucientemente pequeno, x1 f é un atractor local e x2 f é inestable. Así, acabamos de comprobar que existe unha bifurcación transcrítica en (xf, rf) = (0,0) , representada na Figura 2.9. -10 -8 -6 -4 -2 0 2 4 6 8 10 -8 -6 -4 -2 0 2 4 6 8 x r Figura 2.9: Ilustración da bifurcación transcrítica para f(x, r) = x+rx −x2 . En vermello, represéntase os puntos xos inestables ( xf=r para r < 0 e xf= 0 para r > 0 ) e, en azul, os localmente asintoticamente estables ( xf= 0 para r < 0 e xf=r para r > 0 ). 3. Bifurcación de tipo pitchfork Unha bifurcación de tipo pitchfork ocorre cando tres ramas de puntos xos coliden nun punto xo xf tal que f′(xf)=1 e entón só queda unha. No caso supercrítico, dous puntos xos L.A.S. chocan cun punto xo inestable e o punto xo resultante é L.A.S. No caso subcrítico, dous puntos xos inestables chocan cun punto xo L.A.S. e o punto xo resultante é inestable. Exemplo 2.28. Considerando a función (x, r)−→ f(x, r) = x+rx −x3 , temos que ∂f ∂x (x, r) = 1 + r−3x2 , polo que é fácil vericar que (xf, rf) = (0,0) é un punto xo non hiperbólico xa que f(0,0) = 0 e ∂f ∂x (0,0) = 1 . 2.3. Bifurcacións 27 Estamos interesados nos puntos xos nunha veciñanza de (xf, rf) = (0,0) . Sabemos que os puntos xos son da forma f(x, r)−x=rx −x3= 0 , é dicir, x1 f= 0, x2 f= +√r, x3 f=−√r . Vemos que os puntos xos existen se, e só se, r∈R para x1 f e r > 0 para x2 f e x3 f , coincidindo se r= 0 . Ademais, substituíndo en ∂f ∂x (x, r) , temos que ∂f ∂x(x1 f, r) = 1 + r   <1 se r < 0, >1 se r > 0, ∂f ∂x(x2 f, r) = 1 + r−3(+√r)2= 1 −2r < 1, para todo r > 0, ∂f ∂x(x3 f, r) = 1 + r−3(−√r)2= 1 + r−3r= 1 −2r < 1, para todo r > 0, En consecuencia, para r > 0 sucientemente pequeno, x1 f é inestable e x2 f, x3 f son atractores locais. Sen embargo, para r < 0 sucientemente pequeno, x1 f é un atractor local, como se representa na Figura 2.10. -10 -5 0 5 10 15 -8 -6 -4 -2 0 2 4 6 8 x r Figura 2.10: Ilustración da bifurcación pitchfork (caso supercrítico) para f(x, r) = x+xr−x3 . En vermello, represéntanse os puntos xos inestables ( xf= 0 para r > 0 ) e, en azul, os localmente asintoticamente estables ( xf= 0 para r < 0 e xf=±√r para r > 0 ). Consideramos agora o caso f′(xf, rf) = ∂f ∂x (xf, rf) = −1 no que só pode aparecer un tipo de bifurcación: 1. Bifurcación de dobramento de período Unha bifurcación de dobramento de período ocorre cando un punto xo xf , cumprindo que f′(xf) = −1 , cambia a súa estabilidade e aparece unha órbita de período 2 . No caso supercrítico, un punto xo L.A.S. pasa a ser inestable e aparece unha órbita de período 2 L.A.S. No caso subcrítico, un punto xo inestable pasa a ser L.A.S. e aparece unha órbita de período 2 inestable. 34 3. Aplicacións a modelos biolóxicos Facendo uso de que N0=ˆ fS e que T é a idade de recrutamento, tense: R≡N(T)=(ˆ fe−a2T)Se−b2T S ou R=αSe−βS, (3.12) sendo α=ˆ fe−a2T o parámetro de produtividade e β=b2T o parámetro que controla o nivel de dependencia da densidade. Se consideramos a función ˜ f: [0,∞)−→ [0,∞) , denida mediante x7−→ ˜ f(x) = αxe−βx , temos o sistema dinámico discreto unidimensional de Ricker , dado por: wn+1 =˜ f(wn) = αwne−βwn. (3.13) Esta relación ten forma de cúpula en función dos reprodutores, cun máximo en (S, R)=(α βe ,1 β) . Tamén sucede que, cando S tende a ∞ , R se aproxima a 0 . Polo tanto, o efecto dunha inuencia negativa dos reprodutores nas primeiras etapas de vida, mediante canibalismo, é producir un efecto cada vez máis prexudicial sobre o recrutamento a medida que aumenta a poboación reprodutora. De xeito similar ao mencionado para o modelo de Beverton-Holt, o modelo de Ricker segue sendo válido se o efecto de mortandade non afecta a todas as etapas temperás da historia da vida ata o recrutamento e a relación entre as restantes etapas non depende da densidade. Como se xo co modelo de Beverton-Holt, o modelo de Ricker establecido en (3.13) pode ser reescrito en función dun só parámetro r > 0 . Utilizando o cambio de variable zn=βwn , obtemos zn+1 =αzne−zn . Agora, utilizando α=er con r > 0 , quédanos zn+1 =zner−zn . Se, nalmente, empregamos o cambio de variable xn=zn r , chegamos a xn+1 =yner(1−xn) . Reescribimos o modelo como segue: xn+1 =f(xn) = xner(1−xn). (3.14) Polo tanto, pódese considerar f(x) = xer(1−x) . Esta función está representada para r= 1,5 na Figura 2.1b. 3.2. Existencia de puntos xos Imos analizar a existencia de puntos xos non negativos para as funcións que denen ambos modelos. 3.3. Análise da estabilidade dos puntos de equilibrio 35 3.2.1. Modelo de Beverton-Holt Consideramos o modelo de Beverton-Holt (3.8), dado pola función f(x) = r x 1+x , x∈[0,∞) , onde r > 0 . Calculamos os puntos xos a través da denición, é dicir, os puntos x∈[0,∞) tales que f(x) = x : rx 1 + x=x⇔rx = (1 + x)x⇔   x= 0, r= 1 + x, (3.15) quedando x= 0 ou x=r−1 . Polo tanto, para que existan puntos xos non negativos, necesariamente debemos considerar r > 1 . 3.2.2. Modelo de Ricker Consideramos o modelo de Ricker (3.14), dado pola función f(x) = xer(1−x) , x∈[0,∞) , onde r > 0 . Calculamos os puntos xos de f , é dicir, os puntos x∈[0,∞) tales que f(x) = x : xer(1−x)=x⇔   x= 0 er(1−x)= 1 ⇔x= 1, (3.16) isto é, x= 0 ou x= 1 . Neste caso, para todo r > 0 existe un único punto xo positivo x= 1 que non depende do parámetro. 3.3. Análise da estabilidade dos puntos de equilibrio Veremos a estabilidade asintótica local e global dos puntos xos non negativos para as funcións que denen os modelos anteriores. 3.3.1. Modelo de Beverton-Holt Para o modelo de Beverton-Holt, estudamos a estabilidade local e global: Estabilidade local: Como f′(x) = r (1+x)2 , x∈[0,∞) e r > 1 , |f′(0)|=r > 1 , polo tanto, x= 0 é inestable en virtude do Teorema 2.9 e, como |f′(r−1)|=1 r∈(0,1) , entón x=r−1 é asintoticamente estable polo Teorema 2.7. Estabilidade global: Como f: [0,∞)−→ [0,∞) , x7−→ f(x) = rx 1+x é unha función continua tal que f(0) = 0 , 36 3. Aplicacións a modelos biolóxicos x=r−1 é o único punto xo positivo e cumpre que x < f(x)< r −1 en (0, r −1) e x>f(x)> r −1 en (r−1,∞) , entón r−1 é un atractor global en (0,∞) , en virtude do Teorema 2.13. Este resultado de estabilidade global aparece ilustrado na Figura 3.1. En efecto, para x > 0 , comprobamos as condicións descritas anteriormente. Por unha banda, consideramos o caso x∈(0, r −1) . En primeiro lugar, f(x) = rx 1+x> x equivale a x<r−1 . En segundo lugar, como f é estritamente crecente, xa que f′(x)>0 , e r−1 é punto xo, temos que f(x)< f(r−1) = r−1 . Polo tanto, x<f(x)< r −1 para todo x∈(0, r −1) . Por outra banda, consideramos o caso x∈(r−1,∞) . En primeiro lugar, f(x) = rx 1+x< x equivale a x>r−1 . En segundo lugar, como f é estritamente crecente, xa que f′(x)>0 , e r−1 é punto xo, temos que r−1 = f(r−1) < f(x) . Polo tanto, r−1< f(x)< x para todo x∈(r−1,∞) . 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 (a) Converxencia pola esquerda, x0= 0,2∈(0,1) . 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 (b) Converxencia pola dereita, x0= 1,8∈(1,∞) . Figura 3.1: Ilustración da estabilidade global en (0,∞) do punto xo positivo x= 1 do modelo de Beverton-Holt dado por f(x) = rx 1+x , x∈[0,2] , con r= 2 (curva azul). A recta verde representa á función identidade e, en vermello, ilústranse as iteracións do diagrama Cobweb para dúas condicións iniciais distintas x0∈(0,∞) . 3.3.2. Modelo de Ricker Analogamente, para o modelo de Ricker, estudamos a estabilidade local e global: Estabilidade local: Da mesma forma que xemos para o modelo de Beverton-Holt, utilizaremos o Teorema 2.7 e o Teorema 2.9. 3.3. Análise da estabilidade dos puntos de equilibrio 37 Como f′(x) = (1 −rx)er(1−x), x ∈[0,∞) , r > 0 , entón |f′(0)|=er>1 , polo que x= 0 é inestable e |f′(1)|=|1−r| . Se |f′(1)|>1 , é dicir, se r > 2 , temos que x= 1 é inestable. Por outra banda, se |f′(1)|<1 , é dicir, se r∈(0,2) , entón x= 1 é asintoticamente estable. No caso |f′(1)|= 1 , que se corresponde con r= 2 , o criterio non decide. Este resultado aparece ilustrado na Figura 3.2, onde a gráca 3.2a representa o caso asintoticamente estable e a gráca 3.2b o caso inestable. Estabilidade global: Estudemos a estabilidade global para o modelo de Ricker de dúas formas distintas, nunha utilizaremos o Teorema 2.15, onde non se incluiría r= 2 e, para a outra forma, usaremos o Teorema 2.24, onde si se incluiría r= 2 . Método 1 : Xusticación mediante o Teorema 2.15 ou método envolvente. A función f: [0,∞)−→ [0,∞) , x7−→ f(x) = xer(1−x) con r > 0 , é unha función continua e x= 1 é o único punto xo positivo de f , con f(x)> x en (0,1) , f(x)< x en (1,∞) e f(x)>0 en (0,∞) . En efecto, para x > 0 , f(x) = xer(1−x)< x se, e só se, x > 1 . De xeito similar, tense que para x > 0 , f(x) = xer(1−x)> x se, e só se, x < 1 . Ademais, para x > 0 , como ey>0 para todo y∈R , tense que f(x)>0 . Entón pódese tomar ϕ(x) = 2 −x , que é monótona decrecente (xa que ϕ′(x) = −1 ), con ϕ2(x) = x e ϕ positiva no intervalo (0, x−) = (0,2) . Ademais, verica que ϕ(x)> f(x) en (0,1) e que ϕ(x)< f(x) en (1,∞) . Polo tanto, quedaría probado que xf= 1 é un atractor global en (0, ∞ ), en virtude do Teorema 2.15. Comprobemos estas dúas últimas desigualdades entre f e ϕ . En primeiro lugar, temos que f e ϕ son continuas, ϕ(0) = 2 >0 = f(0) , ϕ(1) = f(1) = 1 , ϕ′(x) = −1 e, por último, f′(x) = (1 −rx)er(1−x) para todo x∈[0,∞) . Interésanos estudar a pendente da función f , para isto calculamos f′′(x) = −2rer(1−x)+r2xer(1−x)= 0 ⇔x=2 r. Así, a mínima pendente de f acádase para x=2 r e vale f′(2 r) = −er−2 . Polo tanto, temos que f′2 r>−1⇔er−2<1⇔r−2<log(1) = 0 ⇔r < 2, ademais, f′(2 r) = −1 se, e só se, r= 2 . En consecuencia, podemos armar que, se 0< r < 2 , entón f(x)< ϕ(x) para todo x∈(0,1) e f(x)< ϕ(x) para todo x∈(1,∞) . Método 2 : Estudo baseado no Teorema 2.24 para funcións con derivada Schwarziana negativa. Comezamos estudando a existencia de puntos críticos de f . Como f′(x) = (1 −rx)(er(1−x))=0⇔1−rx = 0 ⇔x=1 r, 38 3. Aplicacións a modelos biolóxicos vemos que f só ten un único punto crítico xc=1 r . Ademais, a función de Ricker cumpre que f(0) = 0 , f(x)>0 en (0,∞) e l´ım x→∞ f(x)=0 . Comprobemos esta última propiedade tendo en conta que r > 0 : l´ım x→∞ xer(1−x)= l´ım x→∞ x 1 er(1−x) L'Hôpital = l´ım x→∞ 1 rer(x−1) = 0. Traballamos agora coas condicións sobre a primeira derivada, é dicir, f′(x)>0 en (0, xc) e f′(x)<0 en (xc,∞) : f′(x) = (1 −rx)er(1−x) | {z } >0 >0⇔1−rx > 0⇔x < 1 r, polo tanto, f′(x)>0 en (0, xc) = (0,1 r) e, ademais, f′(x) = (1 −rx)er(1−x) | {z } >0 <0⇔1−rx < 0⇔x > 1 r, polo que f′(x)<0 en (xc,∞)=(1 r,∞) . Respecto á segunda derivada, temos que f′′(x)<0 en (0, xc) , f′′(x)=(r2x−2r)er(1−x) | {z } >0 <0⇔(r2x−2r)<0⇔ ⇔r(rx −2) <0⇔(rx −2) <0⇔x < 2 r. Como 1 r<2 r , temos que f′′(x)<0 en (0, xc) = (0,1 r) . Por último, faltaría ver a hipótese da derivada Schwarziana, comprobemos que Sf(x)<0 en [0,+∞)\{xc} . Primeiro, necesitamos calcular a derivada terceira, f′′′(x) = r2er(1−x)+ (r2x−2r)(−rer(1−x)) = (−r3x+ 3r2)er(1−x). Estamos en condicións de calcular a derivada Schwarziana mediante a expresión (2.4): Sf(x) = f′′′(x) f′(x)−3 2f′′(x) f′(x)2 =(−r3x+ 3r2)er(1−x) (1 −rx)(er(1−x))−3 2 (r2x−2r)er(1−x) (1 −rx)(er(1−x))!2 =−r3x+ 3r2 1−rx −3 2 (r2x−2r)2 (1 −rx)2=−r3x+ 3r2 1−rx −3 2r4x2−4r3x+ 4r2 (1 −rx)2 =−r3x+ 3r2 1−rx −3r4x2−12r3x+ 12r2 2(1 −rx)2 =2(1 −rx)(−r3x+ 3r2)−3r4x2+ 12r3x−12r2 2(1 −rx)2 =−2r3x+ 6r2+ 2r4x2−6r3x−3r4x2+ 12r3x−12r2 2(1 −rx)2 = <0 z}| { −6r2−r4x2+ 4r3x 2(1 −rx)2 | {z } >0 <0. 3.4. Diagramas de bifurcación 39 Dado que r > 0 , como 6r2+r4x2>4r3x , temos que −6r2−r4x2+ 4r3x < 0 e, polo tanto, Sf(x)<0 en [0,+∞)\{1 r} . Así, estamos en condicións de poder aplicar o Teorema 2.24 e, como f(x) ten un único punto xo en (0,∞) , tal e como demostramos no método anterior, entón xf= 1 é un atractor global en (0,∞) para r∈(0,2] . 0 0.5 1 1.5 2 2.5 3 3.5 4 0 0.5 1 1.5 2 2.5 3 3.5 4 (a) Globalmente asintoticamente estable: r= 1,5≤ 2 , x0= 3 . 0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5 (b) Inestable: r= 3 >2 , x0= 1,4 . Figura 3.2: Ilustración da estabilidade en (0,∞) do punto xo positivo x= 1 do modelo de Ricker dado por f(x) = f(x) = xer(1−x) , x∈[0,4] , para distintos valores de r (curva azul). A recta verde representa á función identidade e, en vermello, ilústranse as iteracións do diagrama Cobweb para dúas condicións iniciais distintas x0∈(0,∞) . 3.4. Diagramas de bifurcación Un diagrama de bifurcación contén información para predicir a dinámica asintótica dun sistema, pois neste represéntase o comportamento a longo prazo das solucións. Dita dinámica pode cambiar de forma bastante drástica co cambio dos parámetros do sistema e é importante saber que estes cambios non son só cuantitativos, como, por exemplo, o cambio na localización dun punto xo, senón tamén cualitativos, é dicir, os puntos xos pódense crear ou destruír, a súa estabilidade pode cambiar, o comportamento do sistema pode cambiar de regular (estabilidade dun punto xo ou órbita periódica) a irregular. A estes cambios cualitativos na dinámica do sistema denotámolos bifurcacións e os valores dos parámetros nos que se producen son os puntos de bifurcación. Nos diagramas de bifurcación que aparecen a continuación, representamos no eixe vertical varios valores de xn , para n∈N grande, fronte a varios valores do parámetro r ao longo do eixe horizontal. Unha mostra dos programas que elaboramos en MATLAB para obter ditos grácos pode consultarse nos Programas A.13 e A.14, nestes, para cada valor de r xado, comezamos 40 3. Aplicacións a modelos biolóxicos realizando Npre iteracións a partir dunha condición inicial x0= 0,5 , gardamos a última iteración e realizamos, a partir dela, outras N iteracións adicionais, que serán as que utilicemos como mostra do comportamento a longo prazo. Estes valores quedan almacenados nun vector x que se debuxa fronte ao valor correspondente de r . Traballaremos cos modelos de Beverton-Holt e Ricker. 3.4.1. Modelo de Beverton-Holt A continuación, na Figura 3.3, móstrase a representación realizada co programa A.13 que constrúe un diagrama de bifurcación para o modelo de Beverton-Holt (3.8) co parámetro r no intervalo (0,5] . Utilizamos Npre = 200 e N= 100 . Como vemos na gura, cando r∈(0,1] , o Figura 3.3: Ilustración do diagrama de bifurcación para o modelo de Beverton-Holt f(x) = r x 1+x , x∈[0,4] , con r∈(0,5] . A liña de puntos en vermello representa o punto xo inestable x1 f= 0 , mentres que os puntos de cor azul representan a dinámica asintótica. único punto xo x1 f= 0 é localmente asintoticamente estable, e cando r > 1 , temos dous puntos xos, x1 f= 0 , que é inestable e x2 f=r−1 , que é globalmente asintoticamente estable. O estudo da estabilidade no caso r > 1 xa foi realizado na Sección 3.3. Vexamos a estabilidade global para o punto xo x1 f= 0 cando r∈(0,1] . Como f(0) = 0 e x>f(x)>0 para todo x∈(0,∞) , temos que a sucesión {fn(x)} é monótona decrecente e acotada inferiormente, polo que usando a Proposición 1.2 podemos concluír que x1 f= 0 é atractor global en [0,∞) . Probemos entón as desigualdades x>f(x)>0 en (0,∞) : por ser f′(x) = r (1+x)2 para x∈[0,∞) e r > 0 , f é unha función estritamente crecente e, posto que f(0) = 0 , temos que f(x)>0 . Como f′(x)<1 para todo r∈(0,1] , x∈(0,∞) e f(0) = 0 , f está por debaixo da recta y=x para r∈(0,1] , entón, f(x)< x en (0,1] . De todo isto, concluímos que, en (x1 f, r1 f) = (0,1) , temos unha bifurcación de tipo pitchfork 3.4. Diagramas de bifurcación 41 (caso subcrítico), xa que dous puntos xos coexistentes (un L.A.S. e outro inestable) coliden nun único punto xo (x1 f, r1 f) = (0,1) tal que f′(x1 f, r1 f) = f′(0,1) = 1 e despois só se mantén o que era L.A.S., que pasa a ser inestable. Cabe mencionar que a bifurcación non é exactamente da maneira que a describimos na Sección 2.3, pois aquí coliden dúas ramas de puntos xos en lugar de tres. 3.4.2. Modelo de Ricker Recollemos agora, na Figura 3.4, a representación do programa A.14, que constrúe un diagrama de bifurcación para o modelo de Ricker (3.14) co parámetro r no intervalo (0,4] . Neste caso, tamén utilizamos Npre = 200 e N= 100 . Para a interpretación da gura, cómpre recordar que 0 0.5 1 1.5 2 2.5 3 3.5 4 r 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 xn Figura 3.4: Ilustración do diagrama de bifurcación para o modelo de Ricker f(x) = xer(1−x) , x∈[0,5] , tomando r∈(0,4] . As liñas de puntos en vermello representan puntos xos inestables, mentres que os puntos azuis representan o comportamento a longo prazo das órbitas. cada liña vertical contén N puntos, entón, se vemos só un pequeno número de puntos aparecendo para un valor de r xado, signica que estes puntos son visitados repetidamente pola órbita e, polo tanto, haberá una órbita periódica L.A.S. Con respecto á estabilidade de puntos xos, temos que, para r > 0 , o punto xo x1 f= 0 é inestable, mentres que o punto xo x2 f= 1 é localmente asintoticamente estable para r∈(0,2) , pasando a ser inestable cando r > 2 . Ademais, como vimos con anterioridade, o punto xo x2 f= 1 é atractor global para r∈(0,2] . A primeira bifucación ocorre en (x2 f, r2 f) = (1,2) , trátase dunha bifurcación de dobramento de período, xa que o punto xo x2 f= 1 , tal que f′(x2 f, r) = f′(1,2) = −1 ( f′(x) = (1 −rx)er(1−x) ), cambia a súa estabilidade e aparece unha órbita de período 2 L.A.S. Como o punto xo pasa de ser L.A.S. a ser inestable e a órbita de período 2 que aparece é L.A.S., estamos no caso 42 3. Aplicacións a modelos biolóxicos supercrítico. En segundo lugar, ilustramos e comentamos outros detalles do diagrama de bifurcación que se poden consultar en [6]. Nos diagramas que se presentan nas Figuras 3.4 e 3.5, podemos cotexar que, para 2< r < 2,525 (aprox.), as órbitas do sistema converxen a unha órbita de período 2 . Aproximadamente en r= 2,525 , hai outra bifurcación de dobramento de período, onde unha órbita de período 2 L.A.S. pasa a ser inestable e aparece unha órbita de período 4 L.A.S. De feito, unha mirada máis atenta revela que o aumento de r dá como resultado toda unha fervenza de bifurcacións que duplican os períodos que se producen cada vez máis preto unhas das outras e producen órbitas dos períodos 8 , 16 , 32 , 64 , etc., ilústramolo nas Figuras 3.5a e 3.5b. A maiores, (a) Representación das órbitas periódicas de período 4 e 8 para r∈[2,5,3] . (b) Representación das órbitas periódicas de período 8 e 16 para r∈[2,65,2,7] . (c) Representación dun diagrama para r∈[3,6,3,64] , idéntico ao diagrama completo. Figura 3.5: Ilustración de varios diagramas de bifurcación para o modelo Ricker. As liñas de puntos en vermello representan puntos xos inestables, mentres que os puntos azuis ilustran a dinámica asintótica. na Figura 3.5c, vemos unha copia case exacta de todo o diagrama. En realidade, o diagrama 3.5. Interpretación biolóxica 43 contén innitas copias deste tipo de si mesmo, pódese ver que esta autosemellanza se repite con resolucións cada vez máis nas. Tal comportamento é característico das entidades xeométricas chamadas fractais, por iso se di que o diagrama de bifurcación do modelo Ricker é un fractal. 3.5. Interpretación biolóxica Nesta sección, imos facer uso dos diagramas de bifurcación para interpretar, en termos biolóxicos, o comportamento asintótico do número de individuos da especie nos modelos de BevertonHolt e Ricker prestando especial atención á inuencia dos parámetros. Antes de comezar, recordemos o que modelan cada un destes sistemas dinámicos discretos unidimensionais de tipo Reprodutores-Recrutas, pódense consultar os detalles na Sección 3.1. En ambos modelos, o número de ovos é proporcional ao número de reprodutores e, no transcurso de ovo a recruta, temos efectos independentes e dependentes da densidade que provocarán a mortandade dos individuos nas fases temperás evitando que aumente sen límite a poboación recruta. No modelo de Beverton-Holt, estes efectos dependentes da densidade serán causados pola competencia de comida e espazo, mentres que no modelo de Ricker se deben ao canibalismo. Vexamos que consecuencias teñen estas relacións entre as diferentes etapas de evolución da especie na dinámica asintótica. 3.5.1. Modelo de Beverton-Holt Na sección 3.1, para o modelo de Beverton-Holt, chegamos a unha primeira formulación dada pola función ˜ f: [0,∞)−→ [0,∞) , denida mediante x7−→ ˜ f(x) = αx 1+βx , cuxo sistema dinámico discreto unidimensional asociado dado por (3.7) é: wn+1 =˜ f(wn) = αwn 1 + βwn , no que recordamos que α é o parámetro de produtividade (número de recrutas aos que dá lugar cada reprodutor) e β é o parámetro que controla o nivel de dependencia da densidade de poboación. Para simplicar o estudo do modelo, xemos unha reescritura deste en función dun único parámetro r=α , mediante o cambio de variable xn=βwn , obtendo a formulación (3.8): xn+1 =f(xn) = r xn 1 + xn . A Figura 3.3 ilustra o comportamento a longo prazo deste último modelo, do que se interpreta que, para r∈(0,1] , tendo calquera número positivo de reprodutores S , como 0 é atractor global, a poboación de individuos que chegan a recrutas R tende a cero e, polo tanto, a especie extínguese. Para r > 1 e calquera número positivo de reprodutores S , como r−1 é atractor global, o número de recrutas estabilízase en r−1 individuos, que aumenta con r . 50 A. Código de MATLAB 12 % debuxamos a grafica de y con respecto a x en azul , a recta dende xmin a xmax en y=0 en negro , a recta en x=0, dende min(y) a max(y) e a funcion identidade en verde 13 14 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 15 16 Y=iter ates ( f , x0 ,N) ; % chama a funcion iterates , que aplica repetidamente a funcion f a partir do valor i n i c i a l x0 e realiza N iteracions 17 18 YY(1)=Y(1) ; % almacena o primeiro valor da iteracion nunha nova variable YY=Y(1) 19 20 for i =1:N % crea un bucle dende 1 ata o numero de iteracions N 21 22 XX(2 * i = 1)=Y( i ) ; % almacena nas entradas impares o noso novo valor de x , que denotamos por XX, cos valores do vector Y das iteracions 23 24 XX(2 * i )=Y( i ) ; % almacena as entradas pares , XX=(Y(1) ,Y(1) ,Y(2) ,Y(2) , . . . ,Y(N) ,Y(N) ) 25 26 YY(2 * i )=Y( i+1) ; % almacena nas entradas pares o noso novo valor de Y , que denotamos por YY, cos valores do vector Y das iteracions 27 28 YY(2 * i+1)=Y( i+1) ; % almacena as entradas impares menos a primeira entrada , YY=(Y(1) ,Y(2) ,Y(2) ,Y(3) ,Y(3) , . . . ,Y(N+1) ,Y(N+1)) 29 end 30 31 XX(2 * N+1)=Y(N+1) ; % falta ria por completar un valor no vector XX para que tenan a mesma lonxitude XX e YY, en particular , fa ltari a a entrada N+1, co ultimo valor da iteracion , XX=(Y(1) ,Y(1) ,Y(2) ,Y (2) , . . . ,Y(N) ,Y(N) ,Y(N+1)) 32 33 plot (XX,YY, ' r ' ,x0 ,0 , ' r * ')% representamos XX sobre YY en cor vermello e no punto (x0 ,0) facemos unha marca cun asterisco de cor vermello ( rectas verticais e horizantais entre f e a funcion A.2. Diagrama Cobweb 51 identidade ) . Programa A.6: Código representación dos modelos 2 (chamadas.m) 1 % Representacions graficas 2 clc , clear ; 3 4 % Grafica paso a paso do diagrama de Cobweb de Ricker ( r =1.5) 5 Ricker1=@(x) x. * exp(1.5 * (1. = x) ) ; 6 % Paso 1 7 representarci (Ricker1 ,0 ,2.5 , 'b ' ) ; 8 plot ([2 2] ,[2 Ricker1 (2) ] , ' r ' ,2 ,0 , ' r * ',[2.5 0] ,[2.5 2.5] , 'k ' ,[2.5 2.5] ,[0 2.5] , 'k ' ) 9 % Paso 2 10 cobweb( Ricker1 ,2 ,1 ,0 ,2.5) ; 11 % Paso 3 12 cobweb( Ricker1 ,2 ,1 ,0 ,2.5) ; 13 plot ([0.44626 0.44626] , [0.44626 Ricker1 (0.44626) ] , 'r') 14 % Paso 4 15 cobweb( Ricker1 ,2 ,2 ,0 ,2.5) ; 16 17 % Ilustracion Teorema 2.13 , Ricker ( r =0.5) con diagrama Cobweb con 3 iteracions en [0 ,3] cos puntos i n i c i a i s x01=0.3 e x02=1.5 18 Ricker5=@(x) x. * exp(0.5 * (1. = x) ) ; 19 cobweb( Ricker5 , 0.3 , 8 , 0 , 3) ; 20 cobweb( Ricker5 , 1.5 , 3 , 0 , 3) ; 21 22 % Representacions funcion envolvente e funcion de Ricker ( r =1.5) coa identidade 23 Envoltura=@(x) 2 = x ; 24 Ricker1=@(x) x. * exp(1.5 * (1. = x) ) ; 25 representarsi ( Envoltura ,0 ,2 , ' = r ' ) ; 26 representarci (Ricker1 ,0 ,2 , ' = b ' ) ; 27 28 % Ilustracion Teorema 2.24 , Ricker ( r =1.8) , xc < xf con diagrama Cobweb con 3 iteracions en [0 ,3] , punto i n i c i a l x0=0.5 29 Ricker8=@(x) x. * exp (1.8. * (1. = x) ) ; 30 representar2 ( Ricker8 ,0 ,3 , ' = b ' ) ; 31 ylim ([0 1.5]) ; 52 A. Código de MATLAB 32 cobweb( Ricker8 ,0.5 ,3 ,0 ,3) ; 33 ylim ([0 1.5]) ; 34 % Ilustracion Teorema 2.24 , Ricker ( r =0.7) , xc > xf con diagrama Cobweb con 3 iteracions en [0 ,10] , punto i n i c i a l x0=0.5 35 Ricker7=@(x) x. * exp (0.7. * (1. = x) ) ; 36 representarci ( Ricker7 ,0 ,10 , ' = b ' ) ; 37 cobweb( Ricker7 ,0.5 ,3 ,0 ,10) ; 38 ylim ([0 1.5]) ; A.3. Representación de bifurcacións Programa A.7: Código das curvas de puntos xos da bifurcación de tipo sela-nodo na Figura 2.8 (CurvaPtsFixos.m) 1 % Curva de puntos fixos nas condicions da bifurcacion de tipo sela nodo 2 clc , clear ; % limpar memoria da consola 3 % Para x positivos 4 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 5 x=0:0.01:20; % crea o vector x dende 0 ata 20 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 6 Y= = sqrt (x) ; % curva de puntos fixos 7 plot (x ,Y, 'k ' )% represento Y sobre x en cor negra 8 ylim([ = 5 5]) ; % limito ao eixe y 9 xlim([ = 5 25]) ; % limito ao eixe x 10 11 Z=sqrt (x) ; % curva de puntos fixos 12 plot (x ,Z, 'black ')% represento Z sobre x en cor negra 13 ylim([ = 5 5]) ; % limito ao eixe y 14 xlim([ = 5 25]) ; % limito ao eixe x 15 ax = gca ; % aspecto dos eixes localizados na orixe , debuxando as frechas e etiquetando os eixes 16 ax.XAxisLocation = ' origin ' ; 17 ax.YAxisLocation = ' origin ' ; 18 plot (max(x+4.7) , 0 , '>k ' ,'markerfacecolor ','k ' ) A.3. Representación de bifurcacións 53 19 plot (0 ,max(Z+0.4) , '^k ' ,'markerfacecolor ','k ' ) 20 xlabel ( 'r') ; ylabel ( 'x ' ) ; 21 22 % Para x negativos 23 clc , clear ; % limpar memoria da consola 24 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 25 x= = 20:0.01:0; % crea o vector x dende = 20 ata 0 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 26 Y= = sqrt( = x) ; % curva de puntos fixos 27 plot (x ,Y, 'k ' )% represento Y sobre x en cor negra 28 ylim([ = 5 5]) ; % limito ao eixe y 29 xlim([ = 25 5]) ; % limito ao eixe x 30 31 Z=sqrt( = x) ; % curva de puntos fixos 32 plot (x ,Z, 'black ')% represento Z sobre x en cor negra 33 ylim([ = 5 5]) ; % limito ao eixe y 34 xlim([ = 25 5]) ; % limito ao eixe x 35 ax = gca ; % aspecto dos eixes localizados na orixe , debuxando as frechas e etiquetando os eixes 36 ax.XAxisLocation =' origin ' ; 37 ax.YAxisLocation =' origin ' ; 38 plot (4.5 ,0 , '>k ' ,'markerfacecolor ','k ' ) 39 plot (0 ,4.8 , '^k ' ,'markerfacecolor ','k ' ) 40 xlabel ( 'r') ; ylabel ( 'x ' ) ; Programa A.8: Código da representación da bifurcación de tipo sela-nodo na Figura 2.7 (BifurSN.m) 1 % Representacion da bifurcacion de tipo sela nodo 2 clc , clear % limpar memoria da consola 3 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 4 x=0:0.01:20; % crea o vector x dende 0 ata 20 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 5 Y= = sqrt (x) ; % curva de puntos fixos 6 plot (x ,Y, ' == r ' )% represento Y sobre x en cor vermella con lina 54 A. Código de MATLAB discontinua ( inestable ) 7 ylim([ = 5 5]) ; % limito ao eixe y 8 xlim([ = 5 25]) ; % limito ao eixe x 9 10 Z=sqrt (x) ; % curva de puntos fixos 11 plot (x ,Z, 'b ' )% represento Z sobre x en cor azul ( estable ) 12 ylim([ = 5 5]) ; % limito ao eixe y 13 xlim([ = 5 25]) ; % limito ao eixe x 14 ax = gca ; % aspecto dos eixes localizados na orixe , debuxando as frechas e etiquetando os eixes 15 ax.XAxisLocation =' origin ' ; 16 ax.YAxisLocation =' origin ' ; 17 plot (max(x+4.7) ,0 , '>k ' ,'markerfacecolor ' ,'k ' ) 18 plot (0 , max(Z+0.4) , '^k ' ,'markerfacecolor ','k ' ) 19 xlabel ( 'r') ; ylabel ( 'x ' ) ; Programa A.9: Código da representación da bifurcación transcrítica na Figura 2.9 (BifurTC.m) 1 % Representacion da bifurcacion transcritica 2 clc , clear ; % limpar memoria da consola 3 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 4 x= = 10:0.01:10; % crea o vector x dende = 10 ata 10 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 5 U=x . * (x<0) ; % curva de puntos fixos definida por condicions 6 plot (x ,U, ' == r ' )% represento U sobre x en cor vermella con lina discontinua ( inestable ) 7 ylim([ = 10 10]) ; % limito ao eixe y 8 xlim([ = 10 10]) ; % limito ao eixe x 9 10 V=x . * (x>0) ; % curva de puntos fixos definida por condicions 11 plot (x ,V, 'b ' )% represento V sobre x en cor azul ( estable ) 12 ylim([ = 10 10]) ; % limito ao eixe y 13 xlim([ = 10 10]) ; % limito ao eixe x 14 ax = gca ; % aspecto do eixe y localizado na orixe , debuxando as frechas e etiquetando os eixes 15 ax.YAxisLocation =' origin ' ; 16 plot (max(x) = 0.3 ,0 , '>k ' ,'markerfacecolor ','k ' ) A.3. Representación de bifurcacións 55 17 plot (0 , max(x) = 0.3, '^k ' ,'markerfacecolor ','k ' ) 18 text (max(x) = 0.4 ,1 , ' r ' ,' color ' ,'k ' ) ; 19 ylabel ( 'x ' ) ; Programa A.10: Código da representación da bifurcación de tipo pitchfork na Figura 2.10 (BifurPF.m) 1 % Representacion da bifurcacion de tipo pitchfork 2 clc , clear ; % limpar memoria da consola 3 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 4 x=0:0.01:20; % crea o vector x dende 0 ata 20 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 5 y= = 20:0.01:0; % crea o vector y dende = 20 ata 0 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 6 U= = sqrt (x) ; % curva de puntos fixos 7 plot (x ,U, 'b ' )% represento U sobre x en cor azul ( estable ) 8 ylim([ = 10 10]) ; % limito ao eixe y 9 xlim([ = 10 15]) ; % limito ao eixe x 10 11 V=sqrt (x) ; % curva de puntos fixos 12 plot (x ,V, 'b ' )% represento V sobre x en cor azul ( estable ) 13 ylim([ = 10 10]) ; % limito ao eixe y 14 xlim([ = 10 15]) ; % limito ao eixe x 15 16 Z=0. * x ; % curva de puntos fixos 17 plot (x ,Z, ' == r ' )% represento Z sobre x en cor vermella con lina discontinua ( inestable ) 18 ylim([ = 10 10]) ; % limito ao eixe y 19 xlim([ = 10 15]) ; % limito ao eixe x 20 21 T=0. * y ; % curva de puntos fixos 22 plot (y ,T, 'b ' )% represento T sobre x en cor azul ( estable ) 23 ylim([ = 10 10]) ; % limito ao eixe y 24 xlim([ = 10 15]) ; % limito ao eixe x 25 ax = gca ; % aspecto do eixe y localizado na orixe , debuxando as frechas e etiquetando os eixes 26 ax.YAxisLocation =' origin ' ; 56 A. Código de MATLAB 27 plot (15 = 0.3 ,0 , '>k ' ,'markerfacecolor ','k ' ) 28 plot (0 ,10 = 0.3 , '^k ' ,'markerfacecolor ','k ' ) 29 text (15 = 0.5 , 1 , ' r ' ,' color ' ,'k ' ) ; 30 ylabel ( 'x ' ) ; Programa A.11: Código da representación da bifurcación de dobramento de período na Figura 2.11 (BifurDP.m) 1 % Representacion da bifurcacion de dobramento de periodo 2 clc , clear ; % limpar memoria da consola 3 hold on % conserva as representacions nos eixes actuais para que as novas representacions engadidas aos eixes non eliminen as existentes 4 x=0:0.01:20; % crea o vector x dende 0 ata 20 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 5 y= = 20:0.01:0; % crea o vector y dende = 20 ata 0 con paso 0.01 , e dicir , a diferenza entre un valor e o anterior e 0.01 6 U=sqrt (x) ; % curva de puntos fixos de f^2 7 plot (x ,U, 'b ' )% represento U sobre x en cor azul ( estable ) 8 ylim([ = 2 2]) ; % limito ao eixe y 9 xlim([ = 2 1]) ; % limito ao eixe x 10 11 V= = sqrt (x) ; % curva de puntos fixos de f^2 12 plot (x ,V, 'b ' )% represento V sobre x en cor azul ( estable ) 13 ylim([ = 2 2]) ; % limito ao eixe y 14 xlim([ = 2 1]) ; % limito ao eixe x 15 16 T=0 * y ; % curva de puntos fixos 17 plot (y ,T, 'b ' )% represento T sobre x en cor azul ( estable ) 18 ylim([ = 2 2]) ; % limito ao eixe y 19 xlim([ = 2 1]) ; % limito ao eixe x 20 21 R=0 * x ; % curva de puntos fixos 22 plot (x ,R, ' == r ' )% represento Z sobre x en cor vermella con lina discontinua ( inestable ) 23 ylim([ = 2 2]) ; % limito ao eixe y 24 xlim([ = 2 1]) ; % limito ao eixe x 25 ax = gca ; % aspecto do eixe y localizado na orixe , debuxando as frechas e etiquetando os eixes A.3. Representación de bifurcacións 57 26 ax.YAxisLocation =' origin ' ; 27 plot (1 = 0.05 ,0 , '>k ' ,'markerfacecolor ','k ' ) 28 plot (0 ,2 = 0.07 , '^k ' ,'markerfacecolor ' ,'k ' ) 29 text (1 = 0.1 ,0.18 , ' r ' ,' color ' ,'k ' ) ; 30 ylabel ( 'x ' ) ; Programa A.12: Código representación dos modelos 3 (chamadas.m) 1 % Diagrama Cobweb do modelo Beverton = Holt ( r=2) con 5 iteracions en [0 ,2] con 2 condicions i n i c i a i s x01=0.2 ( converxe pola esquerda ) e x02=1.8 ( converxe pola dereita ) 2 BevertonHolt2=@(x) 2. * x./(1.+x) ; 3 cobweb( BevertonHolt2 , 0.2 , 5 , 0 , 2) ; 4 cobweb( BevertonHolt2 , 1.8 , 5 , 0 , 2) ; 5 6 % Diagrama Cobweb do modelo Ricker ( r1=1.5 e r2=3) con 6 e 10 iteracions respectivamente (en [0 ,4] e en [0 ,2]) cos puntos i n i c i a i s x01=3 e x02=1.4 7 Ricker1=@(x) x. * exp(1.5 * (1. = x) ) ; 8 cobweb( Ricker1 ,3 ,6 ,0 ,4) ; 9 Ricker3=@(x) x. * exp(3. * (1. = x) ) ; 10 cobweb( Ricker3 ,1.4 ,10 ,0 ,2) ; A.3.1. Modelo de Beverton-Holt Programa A.13: Código representación do diagrama de bifurcación do modelo de Beverton-Holt na Figura 3.3 (BHbifur.m) 1 % Diagrama de bifurcacion do modelo de Beverton = Holt 2 Npre = 200; N = 100; % Npre e o numero de puntos da orbita que queremos descartar (numero de iteracions previas ) , e N, o numero de puntos da orbita que queremos representar para cada valor do parametro r (numero de iteracions ) 3 x = zeros (N,1) ; % crea unha matriz x de dimension (N x 1) con todos os valores 0 4 5 for r = 0.1:0.001:5.0 % crea un bucle dende r=0.1 ata 5 con paso 0.001 , todo o seguinte farase para cada valor de r 58 A. Código de MATLAB 6 x(1) = 0.5; % tomamos como valor i n i c i a l 0.5 , e dicir , o punto de inicio da orbita 7 8 for n = 1: Npre % calculamos as preiteracions pero non as almacenamos e sustituimos o valor i n i c i a l polo novo valor do punto da orbita do modelo de Beverton = Holt , enton x cambiaria a sua primeira columna da seguinte forma , x=(f ^(Npre) (x(1) ) ,0 ,... ,0) onde agora x (1)= f ^(Npre) (x(1) ) 9 x(1) = r * x(1) /(1 + x(1) ) ; 10 end 11 12 for n = 1:N = 1% calculamos as iteracions a partir do punto da orbita calculado no bucle anterior , e dicir , creamos o vector columna x =(x(1) , f (x(1) ) , . . . , f^N(x(1) ) ) 13 x(n+1) = r * x(n) /(1 + x(n) ) ; 14 end 15 16 plot ( r * ones (N,1) ,x , ' .b ' ,' markersize ' ,2) % creamos o vector columna de todos os valores r , para que tena a mesma lonxitude que o vector columna x , e representase os puntos de x sobre os puntos de r de cor azul , representamos mediante puntos e reducimos o tamano destes para que non se superponan 17 hold on ;% usamolo para reter todos os puntos representados mentres trazamos orbitas para os valores posteriores do parametro r 18 end % todo isto para cada valor de r ata que remate o bucle 19 xlabel ( 'r') ; ylabel ( 'x_n ' ) ; % anade etiquetas aos eixes 20 set (gca , ' xlim ' ,[0 5.0]) ; % facemos zoom nos parametros de r en [0 ,5] 21 plot ([1 5] ,[0 0] , ' : r ' ,'LineWidth ' ,1.5); % represento cunha recta de puntos vermellos (y=0 e x en [1 ,5]) a inestabilidade de Beverton Holt e aumento o seu tamano 22 hold off ;% deixamos de reter A.3.2. Modelo de Ricker Programa A.14: Código da función do diagrama de bifurcación do modelo de Ricker na Figura 3.5 (bifRickerdiag.m) 1 % Funcion diagrama de bifurcacion do modelo de Ricker con zoom A.3. Representación de bifurcacións 59 2 function bifRickerdiag (Npre , N, r0 , rf , rN) 3 x = zeros (N,1) ; % crea unha matriz x de dimension (N x 1) con todos os valores 0 4 dr = ( rf = r0 )/rN; 5 for r = r0 : dr : rf % crea un bucle dende r=r0 ata rf con paso dr , todo o seguinte farase para cada valor de r 6 x(1) = 0.5; % tomamos como valor i n i c i a l 0.5 , e dicir , o punto de inicio da orbita 7 8 for n = 1: Npre % calculamos as preiteracions pero non as almacenamos e sustituimos o valor i n i c i a l polo novo valor do punto da orbita do modelo de Ricker , enton x cambiaria a sua primeira columna da seguinte forma , x=(f ^(Npre) (x(1) ) ,0 ,... ,0) onde agora x(1)= f ^( Npre) (x(1) ) 9 x(1) = x(1) * exp( r * (1 = x(1) ) ) ; 10 end 11 12 for n = 1:N = 1% calculamos as iteracions a partir do punto da orbita calculado no bucle anterior , e dicir , creamos o vector columna x =(x(1) , f (x(1) ) , . . . , f^N(x(1) ) ) 13 x(n+1) = x(n) * exp( r * (1 = x(n) ) ) ; 14 end 15 16 plot ( r * ones (N,1) ,x , ' .b ' ,' markersize ' ,2) % creamos o vector columna de todos os valores r , para que tena a mesma lonxitude que o vector columna x , e representase os puntos de x sobre os puntos de r de cor azul , representamos mediante puntos e reducimos o tamano destes para que non se superponan 17 hold on ;% usamolo para reter todos os puntos representados mentres trazamos orbitas para os valores posteriores do parametro r 18 end % todo isto para cada valor de r ata que termine o bucle 19 plot ([2 4] ,[1 1] , ' : r ' ,'LineWidth ' ,1.5) % represento cunha recta de puntos vermellos (y=1 e x en [2 ,4]) un punto fixo inestable e aumento o seu tamano 20 plot ([0 4] ,[0 0] , ' : r ' ,'LineWidth ' ,1.5) % represento cunha recta de puntos vermellos (y=0 e x en [0 ,4]) un punto fixo inestable e aumento o seu tamano 21 xlabel ( 'r') ; ylabel ( 'x_n ' ) ; % anade etiquetas aos eixes