第八次作业

作业

---

背景

一个单摆，在加入其它影响因素后，在一定的初始条件下展现出了混沌性质。我们来进一步观察这一现象。

摘要

我们可以用单摆来创造一个简单的混沌系统。用计算机编程，利用 Euler-Cromer method 来进行计算。绘出单摆的运动轨迹，进行定性分析。

正文

Exercise 3.20

Calculate the bifurcation diagram for the pendulum in the vicinity of FD = 1.35 to 1.5. Make a magnified plot of the diagram (as compared to Figure 3.11) and obtain an estimate of the Feigenbaum delta parameter.
计算在FD = 1.35 - 1.5附近的钟摆的分叉图。 绘制图的放大图（与图3.11相比），并获得Feigenbaum delta参数的估计

$$\omega _{i+1}=\omega _{i}-(\frac{g}{l}\sin \theta _{i}-q\omega _{i}+F_{D}\sin (\omega Dt_{i}))\Delta t$$
$$\theta _{i+1}=\theta _{i}+\omega _{i+1}\Delta t$$
$$t_{i+1}=t_{i}+\Delta t$$

代码：

python
P70 3.20
"""
import math
import pylab as pl
class Bifurcation_diagram:
    def __init__(self, F_D = 1.2, total_time = 60, init_theta = 0.2, q = 0.5, l = 9.8, g = 9.8, Omega_D = 2/3, time_step = 0.04):
        self.omega = [0]
        self.theta = [init_theta]
        self.t = [0]
        self.dt = time_step
        self.total_time = total_time
        self.q = q
        self.l = l
        self.g = g
        self.Omega_D = Omega_D
        self.F_D = F_D
    def calculate(self):
        n = 1
        self.tmp_theta = []
        self.tmp_F_D = []
        m = 100000
        self.total_time = 2 * 400 * math.pi / self.Omega_D
        self.dt = self.total_time / m
        for i in range(m + 1):
            self.omega.append(self.omega[i] + ((-self.g / self.l) * math.sin(self.theta[i]) - self.q * self.omega[i] + self.F_D * math.sin(self.Omega_D * self.t[i])) * self.dt)
            self.theta.append(self.theta[i] + self.omega[i + 1] * self.dt)
            self.t.append(self.dt * (i + 1))
            if self.theta[i + 1] < -math.pi:
                self.theta[i + 1] = self.theta[i + 1] + 2 * math.pi
            if self.theta[i + 1] > math.pi:
                self.theta[i + 1] = self.theta[i + 1] - 2 * math.pi
            #$\Omega_Dt=2n\pi$
            if  abs(self.t[i + 1] - 2 * n * math.pi / self.Omega_D) < self.dt / 2 :
                if i > 75000:
                    if round(self.theta[-1], 4) not in self.tmp_theta:
                        self.tmp_theta.append(round(self.theta[-1], 4))
                        self.tmp_F_D.append(self.F_D)
                n += 1
bif_F_D = []
bif_theta = []
len_bif_theta = 0
F_k = []
k = []
loop_i = True
step = 300
for i in range(step + 1):
    F_D = round(1.35 + 0.15 / step * i, 4)
    start = Bifurcation_diagram(F_D)
    start.calculate()
    bif_F_D += start.tmp_F_D
    bif_theta += start.tmp_theta
    print(start.tmp_theta, i, "/", step)
    temp_len_bif_theta = len(start.tmp_theta)
    if len_bif_theta != temp_len_bif_theta:
        F_k.append(F_D)
        k.append(temp_len_bif_theta)
        len_bif_theta = temp_len_bif_theta
pl.plot(bif_F_D, bif_theta, 'k.')
pl.xlabel('$F_D$', fontsize=20)
pl.ylabel('$\\theta$(radians)', fontsize=20)
pl.title('$\\theta$ versus $F_D$', fontsize=20)
pl.show()
for i in range(len(F_k)):
    print("F_k, k的值:", F_k[i], ",", k[i])

图像：

结论

对于不同的外界影响，产生混沌的图像有较大区别，且混沌的产生地方也显得无规律可循，不能单一的表示出来。
但是对于不同的外界影响来说，混沌都是在一定时间之后才产生的，但是不同的外界影响对混沌的范围有影响.

致谢

感谢张盛同学在程序方面的指导，感谢蔡浩老师的教导及《计算物理》。