2021iypt第七题珠子动力学（3）代码展示与文献链接

代码允许转载，但是放在附录时一定要记得把出处和制作人标清楚！！！！！！！！！！！！！！！！！！！！！！！！！！！

常见问题：

1.代码复制一定要复制到编辑器中，而不是命令行 

2.一定保证符号全英文

3.对代码修改后尤其是方程部分修改后，求解错误是正常现象，请一定要首先检查方程是否正确 

4.代码所用方法为欧拉法，精度并不高，个人思考难免有发现不了的错误，还请多多包涵。更好的方法欢迎评论区讨论

所用软件：MATLAB

代码：

%b站多一份问候

%iypt2021第七题，求解方法为欧拉法。

%方程求解代码，有错误请联系作者。

clear all

R=0.2; % m质量

w=80; % rad/s大环的角速度

T=1000;%计算时间

Q=zeros(T,1); % angle角度

Q(1)=0.01;%角度初始值，为零的话珠子不会抬起

u=0.4; % friction动摩擦系数

Q1=zeros(T,1); % angle derivative角速度

Q2=zeros(T,1); % angular acceleration角加速度

g=9.8;%重力加速度

t=0.001;%计算时间间隔

for i=1:T

Q2(i+1)=(R*sin(Q(i))*w^2-(R*cos(Q(i))*Q1(i)^2+g+u*2*Q1(i)*R*cos(Q(i))*w)*(sin(Q(i))+u*cos(Q(i)))/(cos(Q(i))-u*sin(Q(i)))-u*2*Q1(i)*R*cos(Q(i))*w+R*sin(Q(i))*Q1(i)^2)/(R*cos(Q(i))+R*sin(Q(i))*(sin(Q(i))+u*cos(Q(i)))/(cos(Q(i))-u*sin(Q(i)))) ;%方程

Q1(i+1)=Q1(i)+Q2(i+1)*t;%角速度

Q(i+1)=Q(i)+Q1(i+1)*t;%速度

end

r=R*sin(Q);%水平位移

h=R-R*cos(Q);%垂直偏离

hold on

plot(Q)%绘图

珠子动力学文献合集：

文献链接: https://pan.baidu.com/s/1_zlXjFJ1Ry4Lv273ltsuOg 提取码: grg9 复制这段内容后打开百度网盘手机App，操作更方便哦