import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import newton

def logic(x,r):
        return r*x*(1-x)

Caso r<3

a=[]
b=[]
N=20
x0=0.6
r=2.7
aux=x0
for i in range(N):
    a.append(aux)
    aux=logic(aux,r)
    b.append(aux)
    a.append(aux)
    b.append(aux)

plt.plot(a,b)
#plt.plot(a2,b2)
plt.plot(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),np.linspace(min(min(a),min(b)),max(max(a),max(b)),100))
plt.plot(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),logic(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),r))

[<matplotlib.lines.Line2D at 0x1cea921ce50>]

Caso r<3 comparativo

a=[]
b=[]
N=70
x0=0.6
r=2.6
aux=x0
for i in range(N):
    a.append(aux)
    aux=logic(aux,r)
    b.append(aux)
    a.append(aux)
    b.append(aux)

aux2=x0+0.001
a2=[]
b2=[]
for i in range(N):
    a2.append(aux2)
    aux2=logic(aux2,r)
    b2.append(aux2)
    a2.append(aux2)
    b2.append(aux2)
    

plt.plot(a,b)
plt.plot(a2,b2)
plt.plot(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),np.linspace(min(min(a),min(b)),max(max(a),max(b)),100))
plt.plot(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),logic(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),r))

[<matplotlib.lines.Line2D at 0x1cea930dfa0>]

Caso r>3 comparativo (caótico)

a=[]
b=[]
N=70
x0=0.6
r=3.6
aux=x0
for i in range(N):
    a.append(aux)
    aux=logic(aux,r)
    b.append(aux)
    a.append(aux)
    b.append(aux)

aux2=x0+0.001
a2=[]
b2=[]
for i in range(N):
    a2.append(aux2)
    aux2=logic(aux2,r)
    b2.append(aux2)
    a2.append(aux2)
    b2.append(aux2)
    

plt.plot(a,b)
plt.plot(a2,b2)
plt.plot(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),np.linspace(min(min(a),min(b)),max(max(a),max(b)),100))
plt.plot(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),logic(np.linspace(min(min(a),min(b)),max(max(a),max(b)),100),r))

[<matplotlib.lines.Line2D at 0x1cea9379880>]

Clase Matriz de Leslie

Conejos de Fibonacci

Matriz_c=[[1,1],[1,0]]
Matriz_c

[[1, 1], [1, 0]]

print(np.linalg.eig(Matriz_c)[0])

[ 1.61803399 -0.61803399]

(1+np.sqrt(5))/2

1.618033988749895

El máximo es el primer valor propio y su vector asociado sera:

print(np.linalg.eig(Matriz_c)[1].T[0])

[0.85065081 0.52573111]

$x_n= A_1 \lambda_1^n + A_2 \lambda_2^n$

$1=x_0=A_1 + A_2 $

$1=x_1=A_1 \lambda_1 + A_2 \lambda_2 $

Reemplazando $1=A_1 \lambda _1 + (1-A_1) \lambda_2=\lambda_2+A_1\sqrt{5}=\dfrac{1-\sqrt{5}}{2}+A_1\sqrt{5} $,

por lo que $A_1=\dfrac{1+\sqrt{5}}{2\sqrt{5}}$ y $A_2=\dfrac{\sqrt{5}-1}{2\sqrt{5}}$

Finalmente para $n$ grande $x_n\approx \dfrac{1+\sqrt{5}}{2\sqrt{5}}\left(\dfrac{1+\sqrt{5}}{2}\right)^n$

Comparacion con http://bandicoot.maths.adelaide.edu.au/Leslie_matrix/leslie.cgi?initial_pop[0]=0&initial_pop[1]=100&initial_pop[2]=100&initial_pop[3]=100&initial_pop[4]=0&initial_pop[5]=0&initial_pop[6]=0&initial_pop[7]=0&birth_rates[0]=0&birth_rates[1]=0&birth_rates[2]=0.8&birth_rates[3]=0.5&birth_rates[4]=0&birth_rates[5]=0&birth_rates[6]=0&birth_rates[7]=0&survival_rates[0]=0.95&survival_rates[1]=0.9&survival_rates[2]=0.739&survival_rates[3]=0&survival_rates[4]=0&survival_rates[5]=0&survival_rates[6]=0&survival_rates[7]=0&Submit+Leslie+Matrix=Submit+Leslie+Matrix

m=[0,0,0.8,0.6]
s=[1,0.95,0.9,0.75]
l=[]
f=[]
aux=1
for i in range(len(s)):
    aux=aux*s[i]
    l.append(aux)
    f.append(l[i]*m[i])
def F(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))
    return aux2
def Fp(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))*(-(i+1))
    return aux2
def F2(r):
    return F(r)-1

r=np.log(1.02)

N=len(s)
M=np.zeros((N,N))
M[0]=m
for i in range(N-1):
    M[i+1][i]=s[i+1]
print("Original Numpy")
print(np.log(max(abs(np.linalg.eig(M)[0]))))
print("\n")
print("Original Pagina")
print(r)
print("\n")
print( "Newton")
print(newton(F2, 0))
print("\n")
print( "Aproximado")
print((1-F(0))/Fp(0))
print("\n")
print( "Explota")

Original Numpy
0.019802032589678185


Original Pagina
0.01980262729617973


Newton
0.0198020325896784


Aproximado
0.019145084934558642


Explota

m=[0,0,0.2,0.3,0.5,0.1,0]
s=[1,1,0.95,0.95,0.7,0.4,0.2]
l=[]
f=[]
aux=1
for i in range(len(s)):
    aux=aux*s[i]
    l.append(aux)
    f.append(l[i]*m[i])
def F(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))
    return aux2
def Fp(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))*(-(i+1))
    return aux2
def F2(r):
    return F(r)-1

r=np.log(0.949)

N=len(s)
M=np.zeros((N,N))
M[0]=m
for i in range(N-1):
    M[i+1][i]=s[i+1]
print("Original Numpy")
print(np.log(max(abs(np.linalg.eig(M)[0]))))
print("\n")
print("Original Pagina")
print(r)
print("\n")
print( "Newton")
print(newton(F2, 0))
print("\n")
print( "Aproximado")
print((1-F(0))/Fp(0))
print("\n")
print( "Se extingue")

Original Numpy
-0.05208964748385038


Original Pagina
-0.052346480372209236


Newton
-0.05208964748385024


Aproximado
-0.05854175316452894


Se extingue

m=[0,0,0.8,0.5]
s=[1,0.95,0.9,0.739]
l=[]
f=[]
aux=1
for i in range(len(s)):
    aux=aux*s[i]
    l.append(aux)
    f.append(l[i]*m[i])
def F(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))
    return aux2
def Fp(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))*(-(i+1))
    return aux2

def F2(r):
    return F(r)-1

r=np.log(0.999976)

N=len(s)
M=np.zeros((N,N))
M[0]=m
for i in range(N-1):
    M[i+1][i]=s[i+1]
print("Original Numpy")
print(np.log(max(abs(np.linalg.eig(M)[0]))))
print("\n")
print("Original Pagina")
print(r)
print("\n")
print( "Newton")
print(newton(F2, 0))
print("\n")
print( "Aproximado")
print((1-F(0))/Fp(0))
print("\n")
print( "Converge")

Original Numpy
-2.337279352769967e-05


Original Pagina
-2.4000288004632086e-05


Newton
-2.3372792509843604e-05


Aproximado
-2.3373717084514114e-05


Converge

Ejercicio 1.25

l=[1.0000,0.8335, 0.7313, 0.5881 ,0.4334 ,0.2928 ,0.1813 ,0.1029 ,0.0535,0.0255]
m=[0.0000,0.6504,2.3939,2.9727,2.4662,1.7043,1.0815,0.6683,0.4286,0.300]
f=[]
for i in range(len(l)):
    f.append(l[i]*m[i])
def F(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))
    return aux2
def Fp(r):
    aux2=0
    for i in range(len(f)):
        aux2=aux2+f[i]*np.exp(r*(-(i+1)))*(-(i+1))
    return aux2
print("El R_0 es:",F(0))
print("El r_1 aproximado es:",(1-F(0))/(Fp(0)))
print( "Newton:",newton(F2, 0))

El R_0 es: 5.90430658
El r_1 aproximado es: 0.20543241493649214
Newton: 0.4836376675564584

s=[l[0]]
for i in range(1,len(l)):
    s.append(l[i]/l[i-1])
N=len(s)
M=np.zeros((N,N))
M[0]=m
for i in range(N-1):
    M[i+1][i]=s[i+1]
print("Original Numpy")
print(np.log(max(abs(np.linalg.eig(M)[0]))))

Original Numpy
0.48363766755645826