work_6th chapter 3 NO3.11

St: Wenya Wu 201530102016 WHU

3.11 For the three values of $F_D$ shown in Figeure 3.6, compute and plot the total energy as a function of the time and discuss the results.

综述

• 本题使用物理摆模型，在角加速度项除了回复力之外，附加了空气阻力和驱动力；
• 涉及方程为非线性方程;
• 计算过程为了保证能量值不发散采用euler cormer calculate;
• 偏离平衡位置的摆的角度限制在[$-\pi$,$\pi$]之间；
• 计算每一时刻的能量值并输出.

公式说明

euler cromer calculate method：

$$\omega_{i+1}=\omega_i-[(g/l)sin\theta_i-q\omega_i+F_Dsin(\Omega_Dt_i)]Δt$$ $$\theta_{i+1}=\theta_i+\omega_{i+1}Δt$$ $$\theta_i\in[-\pi,\pi]$$ $$H_i=\frac{1}{2}I\omega_i^2-lcos\theta_i$$
其中,$l$为摆线长度；$F_D$为策动力幅度，$\Omega_D$为策动力频率；$q$为阻力系数，$H_i$为当前能量值，$I$为转动惯量，$I=Ml^2$，$M$为摆质量.

计算过程及结果

• 设定各参数为$l=g=9.8;q=0.5;\Omega_D=0.67;dt=0.04$；
• 分别对比$F_D=0,0.5,1.2$时的 $\theta-t$图像；代码
  
• 画出相应的$\omega-t$图像（对$F_D>1$(大概就是1左右)，由于$\theta$此时出现了大于$\pi$,或者小于$-\pi$的情形，会出现导数不存在的点，这时绘制$\omega-t$图线应该用对$\theta$无限制的情形）；代码
  
  另:给出$F_D=1.2$时$\theta$角有限制和无限制的对比图
  (这里和课本上的有点不太一样，不知道为什么)代码
  
• 画出能量随时间的变化图(设置质量为$m=7kg$；零势能点为单摆悬点)代码
  (由于策动力项为$F_D=1.2$时，外界向单摆不断输入能量，所以能量值会大很多，画在同一张图上不能显示其他两种情况的细节，所以分开两幅图画)
  
  
  关于能量的讨论：
  
  1. 在无策动力的情形下动能很快耗散尽（由图中可见能量最终很快恒等于一个常数）
  2. 策动力较小时，能量表现出周期性，这是因为策动力与耗散一起作用的结果；
  3. 在策动力较大时，能量很大，此时策动力占据主导地位，并且处于数值极不稳定.
  4. 然后在$F_D=1.2$时，调整策动力频率为3，和调整阻力系数为0.999可以看到能量值变化规律变得稳定起来，这与经验相符合（增大策动力频率或者增大摩擦阻力，更利于运动的规则性）
     
     

总结

本次作业计算了物理摆在不同的策动力条件下能量随时间的变化情况，最终的可以看到能量的大小与策动力的幅值$F_D$关系密切，同时与策动力频率和摩擦阻力大小也有关.低频率的策动力更容易导致能量的不稳定.
另外，每次复制粘贴代码过于麻烦，准备学习一下类的使用，应该会更方便些.