连杆机构运动分析力分析Matlab语言m文件使用说明及算例(新)
基于matlab的四杆机构运动分析
基于matlab的四杆机构运动分析一、四杆机构基本概念四杆机构是一种通过变换连杆长度，改变机构运动形态的机械系统。
四杆机构通常由固定连杆、推动连杆、连接杆和工作连杆四个连杆组成，其中固定连杆和推动连杆固定不动，连接杆和工作连杆则沿固定轴线的方向做平动或旋转运动。
四杆机构的基本构造如下图所示：四杆机构的四个连杆的长度和构造参数，以及驱动机构的运动决定了机构的运动特性。
在进行四杆机构运动分析时，需要通过求解运动学关系式和动力学方程，得到连杆的运动规律和力学特性。
二、四杆机构运动学分析1.运动学基本方程四杆机构的运动学分析基本方程是连杆长度变化的定理，即:l₁²+l₂²-2l₁l₂cosθ₂=l₃²+l₄²-2l₃l₄cosθ₄其中，l₁，l₂分别为固定连杆和推动连杆长度；l₃，l₄分别为连接杆和工作连杆长度；θ₂，θ₄分别为推动连杆和工作连杆的夹角。
2.运动学求解方法根据四杆机构运动学基本方程，可以求解机构中任意连杆的角度和位置，从而分析机构运动规律。
在matlab程序中，运动分析可以采用分析法或图解法。
分析法通常采用向量法或坐标法，即将四杆机构中各连杆和运动副的运动量表示为向量或坐标，然后根据连杆长度变化的定理，求解四个未知角度θ₁、θ₂、θ₃、θ₄。
图解法则先通过画图确定机构的运动规律，在图上求解连杆的角度。
比如可以采用伯格（Bourgeois）图法或恰普利恩（Chaplygin）图法等。
四杆机构动力学分析基本方程包括平衡方程和力平衡方程。
平衡方程：当四杆机构处于平衡状态时，连杆的受力关系可以表示为：ΣF=0其中ΣF为各连杆受力的合力。
ΣF=m×a其中，m为每个连杆的质量，a为连杆的加速度。
四杆机构动力学求解方法以matlab为工具，可借助matlab的求解器完成求解。
具体可以利用matlab的优化工具箱、控制工具箱和系统动态学工具箱等，来实现机构模型的动态模拟、仿真和优化设计。
基于Matlab的机构运动分析教学实践
基于Matlab的机构运动分析教学实践1. 引言机构运动分析是机械设计中的重要环节之一，其研究涉及运动学和动力学等多专业知识。
对于机械工程专业的学生而言，学习机构运动分析并进行相关实践，可有效提高其综合素质和能力，为其未来工作打下坚实基础。
本文主要介绍基于Matlab的机构运动分析教学实践，探讨如何通过实践来巩固和深化相关知识，并为机械专业学生提供一定的参考。
2. 实践内容2.1 实验要求本次实践以单平面机构为例，要求学生使用Matlab软件对其结构和运动学特性进行分析，并绘制相关曲线图。
具体要求包括：1.根据机构结构，计算机构的连杆长度、摆动角度等运动参数。
2.使用Matlab编程，绘制机构的位置分析图和速度分析图。
3.对机构运动学分析结果进行总结，分析机构的运动特性和运动规律。
2.2 实验过程在实践过程中，学生首先需要了解单平面机构的构造和运动规律，掌握机构运动学计算的方法和Matlab编程技能。
接着，学生应根据实验要求设计程序，并输入对应参数，运行程序后绘制各种曲线图。
最后，学生应对实验结果进行总结和分析，深化对机构运动学的理解。
2.3 实验成果对于单平面机构的位置分析图和速度分析图，学生应绘制出各连杆的运动轨迹和速度曲线，并对其进行标注解释。
通过对实验结果的深入分析，学生能够更好地理解机构运动学的基本概念和运动规律，为今后从事与机械设计相关工作打下坚实基础。
3. 实践心得通过本次实践，笔者深刻认识到机构运动分析是机械设计中的重要环节之一，也是机械工程专业核心知识点之一。
学生应尽早掌握基本概念和计算方法，并通过实践来加深对机构运动规律的理解。
Matlab 软件在机构运动分析中具有广泛的应用价值，学生应熟练掌握其编程技能，以更有效地完成实验任务。
4. 结论通过本文的介绍，我们可以了解到基于Matlab的机构运动分析教学实践的相关内容和要求。
学生应通过实践来深化对机构运动学的理解，并为今后从事与机械设计相关工作打下坚实基础。
matlab平面连杆结构分析(机械原理课程设计)
优化参数：连杆 长度、角度、质 量等
优化结果：得到 最优的连杆结构 设计
感谢观看
汇报人：
平面连杆结构的应用范围
机械工程：用于设计、分析和优化机械设 备
生物医学：用于设计、分析和优化假肢、 康复设备等
航空航天：用于设计、分析和优化飞机、 火箭等航天器
机器人技术：用于设计、分析和优化机器 人关节、机械臂等
汽车工业：用于设计、分析和优化汽车底 盘、悬挂系统等
建筑工程：用于设计、分析和优化建筑结 构、桥梁等
03
平面连杆结构的运动学分析
平面连杆结构的运动学方程
平面连杆结构的运动学方程是描述连杆系统运动状态的数学模型 运动学方程包括位移方程、速度方程和加速度方程 运动学方程的建立需要知道连杆系统的几何参数和运动参数 运动学方程的求解可以通过数值积分方法或解析方法进行
平面连杆结构的运动学特性
运动学方程：描述连杆结构的运动状态 运动学参数：包括位移、速度、加速度等 运动学约束：限制连杆结构的运动范围 运动学仿真：通过计算机模拟连杆结构的运动过程
平面连杆结构的形状优化
优化目标：提 高连杆结构的 稳定性和刚度
优化方法：有 限元分析、拓
扑优化等
优化参数：连 杆的长度、宽
度、厚度等
优化效果：提 高连杆结构的 承载能力和使
用寿命
平面连杆结构的拓扑优化
拓扑优化：通过改变材料的分布和形状， 约束条件：结构的刚度、强度、稳定
以实现最优的结构性能
性等性能要求
目标函数：最小化重量或体积，同时 满足给定的性能要求
优化方法：遗传算法、粒子群算法、 模拟退火算法等
设计变量：材料的分布和形状
应用领域：汽车、航空航天、机械制 造等
基于matlab的小型快速压力机多连杆机构运动学分析
3
关于"3与"4的一阶导数和二阶导数运算，通过
MATLAB中的diff函数实现。
连杆4与连杆7的较接点,47的位移方程为：
"47="34+*4COS"4 %47$%34+*4Sin"4
/、 （ 10）
通过一阶求导、二阶求导的速度方程、加速度方
程为：
!"47=="34-*4Sin"4
（ 口）
!%47$=%34+*4COS"4 ."47=>"34-*4（ Sin"4*COS"4）
/(⑷、
!%78$=%47+*7COS"7
( 15)
.78$>%47+*7( COS"7-Sin"7)
( 16)
3多连杆杆系参数确定
该50kN小型快速伺服压LAB的软件环境，编制用于MATLAB 运动学仿真的M函数，搭建 连杆机构的参数化仿
11
第54卷
锻压装备与制造技术
%14$%34+*41Sin"4
将式（5）与式（2）联立，令 /二"14-rcos"、1二%14-
osin"，得到：
• 22+32+*32-*412
”2
2*3!（22+32）
3
/、
⑼
"4=arcSi•n--2--2-+--3- 2+*44112-*332---arc.tan2——
2*4&!（22+32）
文章编号：1672—0121（2019）05—0010—04
MATLAB程序设计M文件程序控制结构程序调试程序举例
y=input('Please input y=:');
[rho,the]=tran(x,y);
rho
the MATLAB中，函数能够嵌套调用，即一种函数能够调用别旳函数， 甚至调用它本身。一种函数调用它本身称为函数旳递归调用。
例6 分别建立命令文件和函数文件，将华氏温度f转换为 摄氏温度c。
21.1111
x=
21.1111
(5). 函数文件旳其他有关概念
①局部变量和全局变量：
局部变量只存在于单个函数工作空间，全局变量 经过global定义，可穿行于不同函数工作空间， 涉及基本工作空间workspace。
函数文件旳内部变量是局部旳，与其他函数 及MATLAB内存相互隔离； 而假如在若干函数中把某一变量定义为全局 变量，那么这些函数将公用这个变量，全局 变量旳作用域是整个MATLAB旳工作区，即 全程有效，全部函数都能够对其存取和修改。
if (nargin == 1) tol = max(size(x)) * max(s) * eps;
程序部分
end
r = sum(s > tol);
（2）命令m文件建立及其运营
建立 涉及下列环节： 进入m文件编辑器 输入程序 定义文件名，保存程序
命令M文件旳运营方式： 直接在命令窗口输入该文件旳文件名 在m文件编辑器中打开该文件后点击工具条中
② M文件模式
将matlab语句构成旳程序存储成以m为扩展名 旳文件，然后再执行该程序文件，这种工作模式 称为程序文件模式。
程序文件不能在命令窗口下建立，因为命令窗口 只允许一次执行一行上旳一种或几种语句。
(完整)基于matlab的四杆机构运动分析
1平面连杆机构的运动分析1。
1 机构运动分析的任务、目的和方法曲柄摇杆机构是平面连杆机构中最基本的由转动副组成的四杆机构，它可以用来实现转动和摆动之间运动形式的转换或传递动力。
对四杆机构进行运动分析的意义是:在机构尺寸参数已知的情况下，假定主动件（曲柄)做匀速转动，撇开力的作用,仅从运动几何关系上分析从动件(连杆、摇杆）的角位移、角速度、角加速度等运动参数的变化情况。
还可以根据机构闭环矢量方程计算从动件的位移偏差。
上述这些内容，无论是设计新的机械，还是为了了解现有机械的运动性能，都是十分必要的，而且它还是研究机械运动性能和动力性能提供必要的依据.机构运动分析的方法很多，主要有图解法和解析法。
当需要简捷直观地了解机构的某个或某几个位置的运动特性时,采用图解法比较方便，而且精度也能满足实际问题的要求。
而当需要精确地知道或要了解机构在整个运动循环过程中的运动特性时,采用解析法并借助计算机，不仅可获得很高的计算精度及一系列位置的分析结果，并能绘制机构相应的运动线图，同时还可以把机构分析和机构综合问题联系起来，以便于机构的优化设计.1。
2 机构的工作原理在平面四杆机构中，其具有曲柄的条件为:a.各杆的长度应满足杆长条件，即：最短杆长度+最长杆长度≤其余两杆长度之和。
b。
组成该周转副的两杆中必有一杆为最短杆,且其最短杆为连架杆或机架（当最短杆为连架杆时,四杆机构为曲柄摇杆机构;当最短杆为机架时，则为双曲柄机构）。
在如下图1所示的曲柄摇杆机构中，构件AB为曲柄,则B点应能通过曲柄与连杆两次共线的位置。
1.3 机构的数学模型的建立1。
3。
1建立机构的闭环矢量位置方程在用矢量法建立机构的位置方程时,需将构件用矢量来表示，并作出机构的封闭矢量多边形。
如图1所示,先建立一直角坐标系.设各构件的长度分别为L1 、L2 、L3 、L4 ，其方位角为、、、.以各杆矢量组成一个封闭矢量多边形，即ABCDA。
其个矢量之和必等于零。
机械原理4-23MATLAB平面连杆机构动力学分析
基于MATLAB/Solidworks COSMOSMotion的平面连杆机构动力学分析07208517王锡霖4-23在图示的正弦机构中,已知l AB =100 mm,h1=120 mm,h2 =80 mm,W1 =10 rad/s(常数)，滑块2和构件3的重量分别为G2 =40 N和G3 =100 N,质心S2 和S3 的位置如图所示，加于构件3上的生产阻力Fr=400 N，构件1的重力和惯性力略去不计。
试用解析法求机构在Φ1=60°、150°、220°位置时各运动副反力和需加于。
构件1上的平衡力偶Mb分别对三个构件进行受力分析如图：构件3受力图构件2受力图构件1受力图（1）滑块2：V S2 =L AB W1 ①a s2 = L AB W12②构件3：S=L AB sinΦ1 ③V3=L AB W1 COSΦ1 ④a3=-L AB W12 sinΦ1 ⑤(2)确定惯性力：F12=m2as2=(G2/g)LABW12 ⑥F13=m3a3=(G3/g)LABW12sinΦ1 ⑦(3)各构件的平衡方程：构件3：∑Fy=0,FR23 =Fr-F13∑Fx=0,FR4’=FR4∑MS3 =0,FR4=FR23LAcosΦ1/h2构件2：∑Fx=0,FR12x=F12cosΦ1∑Fy=0,FR12y=FR32-F12sinΦ1构件1：∑Fx=0,FR41x=FR12x∑Fy=0,FR41y=FR12y∑MA =0,Mb=FR32LABcosΦ1总共有八个方程，八个未知数。
归纳出一元八次方程矩阵：1 0 0 0 0 0 0 0 FR23 Fr-F130 1 -1 0 0 0 0 0 FR4’ 0-LAB COSΦ1/h20 1 0 0 0 0 0 FR40 0 0 1 0 0 0 0 FR12x = F12cosΦ1-1 0 0 0 1 0 0 0 FR12y -F12sinΦ10 0 0 -1 0 1 0 0 FR41x 00 0 0 0 -1 0 1 0 FR41y 0-LABCOSΦ1 0 0 0 0 0 0 1 Mb 0 AX=B进而可得：X=A\B。
二连杆正逆运动学求解matlab
二连杆正逆运动学求解matlab二连杆是机械系统中常见的一种连杆机构，由两个相互连接的连杆组成。
在机械工程中，正逆运动学是研究连杆运动的重要内容之一。
本文将使用Matlab来求解二连杆的正逆运动学问题。
正逆运动学是机械工程中的两个基本问题，正运动学问题是已知连杆的几何参数和运动关系，求解连杆末端的位置和速度；逆运动学问题则是已知连杆末端的位置和速度，求解连杆的几何参数和运动关系。
我们需要定义二连杆的几何参数，包括连杆的长度和初始位置。
假设连杆1的长度为L1，连杆2的长度为L2，连杆1的初始位置为(0,0)，连杆2的初始位置与连杆1的末端相连。
我们需要求解的是连杆2末端的位置和速度。
对于正运动学问题，我们可以根据连杆的几何关系和运动关系，通过几何推导来求解连杆末端的位置。
在Matlab中，我们可以使用向量和矩阵运算来实现这一过程。
我们可以根据连杆的几何关系，利用三角函数来求解连杆2末端的位置。
假设连杆1与x轴的夹角为θ1，连杆2与连杆1的夹角为θ2，则连杆2末端的位置可以表示为：x = L1*cos(θ1) + L2*cos(θ1+θ2)y = L1*sin(θ1) + L2*sin(θ1+θ2)其中，x和y分别表示连杆2末端的x坐标和y坐标。
接下来，我们可以根据连杆的运动关系，求解连杆的速度。
假设连杆1的角速度为ω1，连杆2的角速度为ω2，则连杆2末端的速度可以表示为：v_x = -L1*ω1*sin(θ1) - L2*(ω1+ω2)*sin(θ1+θ2)v_y = L1*ω1*cos(θ1) + L2*(ω1+ω2)*cos(θ1+θ2)其中，v_x和v_y分别表示连杆2末端的x方向和y方向的速度。
对于逆运动学问题，我们可以先根据连杆末端的位置，通过几何推导来求解连杆的几何参数。
然后，根据连杆末端的速度，通过运动学关系来求解连杆的角速度。
在Matlab中，我们可以使用符号计算工具箱来求解逆运动学问题。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
格式：格式：，格式：3x8. 作用有平衡力的构件力分析9. 平面连杆机构运动分析算例例1图示曲柄摇杆机构，已知l 1=150mm ，l 2=220mm ，l 3=250mm% % % % 绘制运动线图 theta1=theta1*180/pi; figure(1) subplot(3,1,1);plot(theta1,cx,'-',theta1,cy,':r'),grid on xlabel('曲柄转角( ° )');ylabel('位移(mm/s)'); legend('C 点x 方向位移','C 点y 方向位移');subplot(3,1,2);plot(theta1,vcx,'-',theta1,vcy,':r'),grid onxlabel('曲柄转角( ° )');ylabel('速度(mm/s)');legend('C点x方向速度','C点y方向速度');subplot(3,1,3);plot(theta1,acx,'-',theta1,acy,':r'),grid onxlabel('曲柄转角( ° )');ylabel('加速度(mm/s^2)'); legend('C点x方向加速度','C点y方向加速度'); figure(2)subplot(2,1,1);plot(theta1,theta3*180/pi),grid onxlabel('曲柄转角( ° )');ylabel('摇杆角位移( ° )'); subplot(2,1,2);plot(theta1,omiga3),grid onxlabel('曲柄转角( ° )');ylabel('摇杆角速度(rad/s)'); subplot(2,1,3);plot(theta1,alpha3),grid onxlabel('曲柄转角( ° )');ylabel('摇杆角加速度(/rad/s^2)');例2 图示曲柄滑块机构，已知l 1=150mm ，l 2=150mm ，e =25mm ，曲柄以n 1=955r/min 逆时针匀速转动，分析该机构的运动。
主程序% 曲柄滑块机构运动分析 clc,clearl1=50;% 曲柄长度 l2=150;% 连杆长度 e=25;% 偏距 n=955;% 曲柄转速m=1;% RRP II 级杆组装配模式系数omiga1=2*pi*n/60;alpha1=0;% 曲柄角速度、角加速度ax=0;ay=e;vax=0;vay=0;aax=0;aay=0;% A 点位置、速度、及速度 px=0;py=0;% 滑块导路上一定点（选为O 点）的坐标 vpx=0;vpy=0;apx=0;apy=0;% 滑块导路上一定点的速度、加速度 theta3=0;omiga3=0;alpha3=0;% 滑块导路角位置、角速度、角加速度 phi=0;% 曲柄的结构参数theta1=0:30:360;% 曲柄转角(每隔10°计算一次) theta1=theta1*pi/180;% 调用crank 函数，计算B 点运动参数 [bx,by]=p_crank(ax,ay,theta1,phi,l1);[vbx,vby]=v_crank(l1,vax,vay,omiga1,theta1,phi); [abx,aby]=a_crank(l1,aax,aay,alpha1,omiga1,theta1,phi); % 调用RRP 函数，计算BC 杆和滑块的运动参数 [cx,cy,sr,theta2]=p_RRP(bx,by,px,py,theta3,l2,m);[vcx,vcy,vr,omiga2]=v_RRP(bx,by,cx,cy,vbx,vby,vpx,vpy,theta2,theta3,l2,sr,omiga3); [acx,acy,ar,alpha2]=a_RRP(bx,by,cx,cy,px,py,abx,aby,apx,apy,theta3,vr,omiga2,omiga3,alpha3); % 绘制运动线图crankx=l1.*cos(theta1);cranky=e+l1.*sin(theta1); theta1=theta1*180/pi; figure(1) subplot(3,1,1); plot(theta1,cx,'r'),grid onxlabel('曲柄转角( ° )');ylabel('滑块位移(mm/s)');曲柄滑块机构subplot(3,1,2);plot(theta1,vcx,'r'),grid onxlabel('曲柄转角( °)');ylabel('滑块速度(mm/s)'); subplot(3,1,3);plot(theta1,acx,'r'),grid onxlabel('曲柄转角( °)');ylabel('滑块加速度(mm/s^2)');例3=840mmlCDe=0;%%%%plot(theta1,theta3*180/pi,'r'),grid onxlabel('曲柄转角( ° )');ylabel('摇杆角位移( ° )');subplot(3,1,2);plot(theta1,omiga3,'r'),grid onxlabel('曲柄转角( ° )');ylabel('摇杆角速度(rad/s)');subplot(3,1,3);plot(theta1,alpha3,'r'),grid onxlabel('曲柄转角( ° )');ylabel('摇杆角加速度(rad/s^2)');例4 图示六杆机构，已知l AB =80mm ，l BC =260mm ，l DE =400mmH 2=170mm 滑块5% e=0;% % % % [ex,ey]=p_crank(dx,dy,theta3,phi2,l31);[vex,vey]=v_crank(l31,vdx,vdy,omiga3,theta3,phi2); [aex,aey]=a_crank(l31,adx,ady,alpha3,omiga3,theta3,phi2); % 调用RRP 函数，计算F 点运动参数 [fx,fy,sr,theta4]=p_RRP(ex,ey,px,py,theta5,l4,m2);[vfx,vfy,vr,omiga4]=v_RRP(ex,ey,fx,fy,vex,vey,vpx,vpy,theta4,theta5,l4,sr,omiga5); [afx,afy,ar,alpha4]=a_RRP(ex,ey,fx,fy,px,py,aex,aey,apx,apy,theta4,vr,omiga4,omiga5,alpha5);% 绘制运动线图theta1=theta1*180/pi;subplot(3,1,1);plot(theta1,sr,'r'),grid onxlabel('曲柄转角( ° )');ylabel('滑块5位移(mm)'); subplot(3,1,2);plot(theta1,vr,'r'),grid onxlabel('曲柄转角( ° )');ylabel('滑块5速度(mm/s)'); subplot(3,1,3);plot(theta1,ar,'r'),grid onxlabel('曲柄转角( ° )');ylabel('滑块5加速度(mm/s^2)');10.平面连杆机构力分析算例H=0.65 m构件质量m4构件1%e=0;%omiga1=pi*n1/30;alpha1=0;% 曲柄角速度、角加速度% 曲柄AB的运动参数theta1=0:10:360; 曲柄转角(每隔10°计算一次)theta1=theta1*pi/180;% 调用crank函数，计算B点运动参数[bx,by]=p_crank(ax,ay,theta1,phi,lab);[vbx,vby]=v_crank(lab,vax,vay,omiga1,theta1,phi);[abx,aby]=a_crank(lab,aax,aay,alpha1,omiga1,theta1,phi);% 调用RPR函数，计算CD杆和D点运动参数[dx,dy,sr3,theta3]=p_RPR(cx,cy,bx,by,e,lcd,m1);[vdx,vdy,omiga3,vr3]=v_RPR(cx,cy,bx,by,dx,dy,vbx,vby,vcx,vcy,theta3);[adx,ady,alpha3,ar3]=a_RPR(cx,cy,bx,by,dx,dy,abx,aby,acx,acy,vr3,omiga3,theta3);% 调用RRP函数，计算E点运动参数[ex,ey,sr,theta4]=p_RRP(dx,dy,px,py,theta5,lde,m2);[vex,vey,vr,omiga4]=v_RRP(dx,dy,ex,ey,vdx,vdy,vpx,vpy,theta4,theta5,lde,sr,omiga5); [aex,aey,ar,alpha4]=a_RRP(dx,dy,ex,ey,px,py,adx,ady,apx,apy,theta5,vr,omiga4,omiga5,alpha5); % 调用crank函数，计算CD杆和DE杆质心点的运动参数[s3x,s3y]=p_crank(dx,dy,theta3,phi,lcd/2);[vs3x,vs3y]=v_crank(lcd/2,vcx,vcy,omiga3,theta3,phi);[as3x,as3y]=a_crank(lcd/2,acx,acy,alpha3,omiga3,theta3,phi);[s4x,s4y]=p_crank(dx,dy,theta4,phi,lcd/2);[vs4x,vs4y]=v_crank(lcd/2,vdx,vdy,omiga4,theta4,phi);[as4x,as4y]=a_crank(lcd/2,adx,ady,alpha4,omiga4,theta4,phi);% 力参数m1=0;m2=0;m3=20;m4=15;m5=62;% 各构件质量Js1=0;Js2=0;Js3=0.11;Js4=0.18;Js5=0;% 各构件绕其质心的转动惯量M1=0;M2=0;M3=0;M4=0;M5=0;% 作用于构件1、2、3、4、5的主矩F1=[0,0];F2=[0,0];F3=[0,0];F4=[0,0];% 作用于构件1、2、3、4的主矢fr=110;fr5=(omiga3>0)*fr;% 生产阻力for i=1:length(theta1)F5=[fr5(i),0];% 作用于构件5的主矢% 力分析% 调用RRP II级杆组力分析函数，计算运动副D、E中的约束反力[R34x(i),R34y(i),R45x(i),R45y(i),R56x(i),R56y(i),lcn]=F_RRP([dx(i),dy(i)],[ex(i),ey(i)],...[s4x(i),s4y(i)],[ex(i),ey(i)],m4,m5,Js4,Js5,M4,M5,F4,F5,theta5,...[as4x(i),as4y(i)],[aex(i),aey(i)],alpha4(i),alpha5);% 调用RPR II级杆组力分析函数，计算运动副B、C中的约束反力[R12x(i),R12y(i),R23x(i),R23y(i),R36x(i),R36y(i),lcn(i)]=F_RPR([cx,cy],[bx(i),by(i)],...[dx(i),dy(i)],[bx(i),by(i)],[s3x(i),s3y(i)],m2,m3,Js2,Js3,M2,M3,F2,F3,...[R34x(i),R34y(i)],theta3(i),[abx(i),aby(i)],[as3x(i),as3y(i)],alpha3(i));% 调用单杆力分析函数，计算运动副A中的约束反力及平衡力矩[R16x(i),R16y(i),Mb(i)]=F_Bar([ax,ay],[bx(i),by(i)],s1,m1,Js1,M1,F1,...[R12x(i),R12y(i)],as1,alpha1);endtheta1=theta1*180/pi;figure(1)subplot(3,1,1)plot(theta1,sr),grid onxlabel('曲柄转角（°）');ylabel('滑块位移（m）');subplot(3,1,2)plot(theta1,vr),grid onxlabel('曲柄转角（°）');ylabel('滑块速度（m/s）');subplot(3,1,3)plot(theta1,ar),grid onxlabel('曲柄转角（°）');ylabel('滑块加速度（m/s^2）');figure(2)plot(theta1,Mb),grid onxlabel('曲柄转角（°）');ylabel('平衡力矩（N.m）');。