python ricker model 李亚普洛夫指数
怪人怪事 2023-05-04 14:31www.bnfh.cn怪人怪事
python ricker model 李亚普洛夫指数
我们的ricker model 定义为
其李亚普罗夫指数定义为
单个r的李亚普罗夫指数计算程序
import numpy as np
# 定义Ricker模型
def ricker(x, r):
return x np.exp(r (1 - x))
# 计算李雅普诺夫指数
def lyapunov(r, x0, niter):
x = x0
sum = 0
for i in range(niter):
x_next = ricker(x, r)
sum += np.log(abs(r (1 - 2 x)))
x = x_next
return sum / niter
# 运行计算李雅普诺夫指数的函数
r = 2.5
x0 = 0.1
niter = 1000
l = lyapunov(r, x0, niter)
print(f"Lyapunov exponent for Ricker model ith r={r:.1f} and x0={x0} is {l:.4f}")
参数r变化的李亚普罗程序
import numpy as np
import matplotlib.pyplot as plt
# 定义Ricker模型
def ricker(x, r):
return x np.exp(r (1 - x))
# 计算李雅普诺夫指数
def lyapunov(r, x0, niter):
x = x0
sum = 0
for i in range(niter):
x_next = ricker(x, r)
sum += np.log(abs(r (1 - 2 x)))
x = x_next
return sum / niter
# 定义绘制分支图的函数
def plot_bifurcation(x0, rmin, rmax, step, niter):
rs = np.arange(rmin, rmax, step)
xs = []
lyaps = []
for r in rs:
x = x0
for i in range(niter):
x = ricker(x, r)
lyap = lyapunov(r, x, niter)
lyaps.append(lyap)
for i in range(niter):
x = ricker(x, r)
xs.append([r, x])
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(8, 10))
ax1.plot(rs, lyaps, 'b-', l=2)
ax1.set_xlabel('r')
ax1.set_ylabel('Lyapunov exponent')
ax1.set_title('Lyapunov exponent for Ricker model')
ax2.scatter(np.array(xs)[:, 0], np.array(xs)[:, 1], s=0.1, c='b')
ax2.set_xlabel('r')
ax2.set_ylabel('x')
ax2.set_title('Bifurcation diagram for Ricker model')
plt.tight_layout()
plt.sho()
# 运行绘制分支图的函数
x0 = 0.1
rmin, rmax, step = 2.5, 4.0, 0.01
niter = 1000
plot_bifurcation(x0, rmin, rmax, step, niter)
分支图、李亚普罗夫指数图、cobeb图常用来分析一维差分方程的动力学行为。
logistic x(n+1)=ux(n)(1-x(n)) 可以同样分析。
其差分方程平衡点、稳定性分析请查阅相关文献。