主程序Ex714a.m global LO L1 L2 L3 th1 L0=20;L1=8;L2=25;L3=20;%输入杆长 theta1=input('当前角theta1=); theta3=input('对应于theta1的theta3近似值=)； %求当前输出角theta3 th1=theta1;theta3=fzero('ex714f,theta3); theta2 asin((L3*sin(theta3)-L1*sin(theta1))/L2); w1=input(w1=') w3 L1*w1*cos(pi/2-theta1+theta2)/... (L3*cos(theta3-pi/2-theta2))
主程序 Ex714a.m global L0 L1 L2 L3 th1 L0=20;L1=8;L2=25; L3=20; % 输入杆长 theta1=input('当前角theta1= '); theta3=input('对应于theta1的theta3近似值= '); %求当前输出角theta3 th1=theta1;theta3=fzero('ex714f',theta3); theta2 = asin(( L3*sin(theta3)- L1*sin(theta1))/L2); w1=input('w1= '); w3 = L1*w1*cos(pi/2-theta1+theta2)/ … (L3*cos(theta3-pi/2-theta2))
主程序ex714b.m global LO L1 L2 L3 th1 L0=20;L1=8;L2=25;L3=20;%输入基线及三根杆的长度L1,L2,L3 w1=input('杆1角速度w1=); theta1=linspace(0,2*pi,181); %杆1每圈分为180份，间隔2度。 theta3=input(('对应于theta1最小值的theta3(近似值)=)； dt=2*pi/180/w;%杆1转2度对应的时间增量 th1=theta1(1);:theta3(1)=fzero(ex714f,theta3g);%求初始theta:3 fori=2:181 th1=theta1(i); theta3()=fzero('ex714f,theta3(i-1);%逐次求theta3 end subplot(1,2,1),plot(theta1,theta3),ylabel('theta3),grid %画曲线 w3=diff(theta3)/dt; %求杆3的角速度 subplot(1,2,2),plot(theta1(2:length(theta1)),w3);grid
主程序 ex714b.m global L0 L1 L2 L3 th1 L0=20;L1=8;L2=25; L3=20;% 输入基线及三根杆的长度L1,L2,L3 w1=input('杆1角速度w1= '); theta1=linspace(0,2*pi,181); %杆1每圈分为180份，间隔2度。 theta3=input('对应于theta1最小值的theta3(近似值)= '); dt = 2*pi/180/w; % 杆1转2度对应的时间增量 th1=theta1(1);theta3(1)=fzero('ex714f',theta3g);%求初始theta3 for i=2:181 th1=theta1(i); theta3(i)=fzero('ex714f',theta3(i-1));% 逐次求theta3 end subplot(1,2,1),plot(theta1,theta3),ylabel('theta3'),grid %画曲线 w3 = diff(theta3)/dt; % 求杆3的角速度 subplot(1,2,2),plot(theta1(2:length(theta1)),w3);grid
子程序（函数程序）ex714f.m function y=ex714f(x) global LO L1 L2 L3 th1 y=L1.*cos(th1)+L2*sqrt(1-(L3*sin(x)- L1*sin(th1)).^2/L2/L2)-L3*cos(x)-LO; ·%在上述程序中注意h1上一个标量，而 theta1在ex714b中是一个数组，因为函数 ex714f中的theta1只能单点作为参数输入， 所以引入了h1作为其当前值
子程序(函数程序)ex714f.m function y=ex714f(x) global L0 L1 L2 L3 th1 • y=L1.*cos(th1)+L2*sqrt(1-(L3*sin(x)- L1*sin(th1)).^2/L2/L2)-L3*cos(x)-L0; • ％在上述程序中注意th1上一个标量，而 theta1在ex714b中是一个数组，因为函数 ex714f中的theta1只能单点作为参数输入, 所以引入了th1作为其当前值
程序运行结果 ·以L0=20,L1=8,L2=25,L3=20(及15)的条件下运行ex714b 的曲线见图7-1-5左。相应的角速度变化规律见图7-1-5右. 2.5 100 w3 L3=20 L3=15 50 2 0 1.5 L3=20 -50 3=15 -100 0.5 -150 0 5 theta1 10 5 theta1 10
程序运行结果 • 以L0=20,L1=8,L2=25,L3=20(及15)的条件下运行ex714b 的曲线见图7-1-5左。相应的角速度变化规律见图7-1-5右. 0 5 10 0.5 1 1.5 2 2.5 theta1 theta3 L3=20 L3=15 theta1 0 5 10 -150 -100 -50 0 50 100 L3=20 w3 L3=15 theta1
例7-1-5质点动力学例 考虑空气阻力时抛射体的质心飞行轨迹，计算质点飞行的轨迹 和距离，设空气阻力的方向与速度向量相反，大小与速度 的平方成正比，即Fz=-cw.如图7-1-6 解：◆建模： 考虑空气阻力后弹丸的运动方程如下，c为空气阻力系数。 dx dt 二Vx y二vy _dt dt 1 dv,=cv sino-mg --cw,-mg m dt
例7-1-5 质点动力学例 考虑空气阻力时抛射体的质心飞行轨迹,计算质点飞行的轨迹 和距离,设空气阻力的方向与速度向量相反,大小与速度 的平方成正比,即Fz=-cvv.如图7-1-6. 解：◆建模： 考虑空气阻力后弹丸的运动方程如下，c为空气阻力系数。 2 2 2 , cos sin x y x x x y y dx dy v v dt dt dv v m cv cv cvv dt v dv m cv mg cvv mg dt = = = − = − = − = − − = − −