工业机械臂建模

本文介绍如何在MATLAB中进行工业机械臂建模以及动力学仿真，翻译于Modeling an industrial arm
平台：MATLAB R2017a
机械臂模型由三部分组成：电机（Motor），减速机（Gear-box），臂结构（Arm structure），各部分都有各自的质量，分别绕着一个不受重力影响的轴转动。如图1所示。

图1：工业机械臂原理图
这个模型在两方面进行了简化：1.假设运动围绕一个不受重力影响的轴。一些机械臂（比如SCARA）由于结构的关系，建模时可以忽略重力项的影响，所以这种化简是合理的。2.令减速比r = 1。我们之后可以通过直接换算减速比来获得真实的物理参数，所以这种化简也是合理的。
具体的原理可以参考论文：E. Wernholt and S. Gunnarsson. Nonlinear Identification of a Physically Parameterized Robot Model. In preprints of the 14th IFAC Symposium on System Identification, pages 143-148, Newcastle, Australia, March 2006.
一． 建立状态空间方程
状态空间方程是整个现代控制理论的基础，让我们看看如何建立各个关节的状态空间方程。
状态变量：

对三个质量应用力矩平衡方程，得到状态空间方程：

其中J_m, J_g, J_a分别是电机，减速机，臂结构的惯量矩；d_g, d_a是“弹簧”阻尼系数，k_a是臂结构刚度。（弹簧代表刚体弹性）
机器人的输入为电机产生的转矩u(t)=tau(t)，电机的角速度y(t) = d/dt q_m(t)为测量输出（可以通过码盘读出）。质量被建模在齿轮箱之后和臂结构的末端，所以他们的角位置：q_g(t)和q_a(t)是不可测量的。齿轮箱的弹性由非线性弹簧建模，由弹簧扭矩tau_s(t)描述，它位于电机和第二个质量之间；而臂结构的灵活性则由线性弹簧在最后两个质量之间建模。系统的摩擦主要作用于第一个质量，这里用非线性摩擦力矩tau_f(t)来建模。
齿轮箱摩擦扭矩tau_f(t)，被建模为包括库仑摩擦、粘滞摩擦并考虑Stribeck效应（将静摩擦转换为动摩擦的过渡阶段视作自然指数增长）的模型：

其中Fv是粘滞系数，Fc是库仑摩擦系数，Fcs和是反映Stribeck效应的系数，beta是用来获得x3(t)速度从负向正平稳过渡的参数。
弹簧的扭矩tau_s(t)，用一个无平方项的三次多项式来描述：

其中k_g1和k_g3为齿轮箱弹簧刚度参数。
因为总体惯量J已知，故引入待定的比例因子a_m和a_g，进行参数重整化：

综上所述得到最终的状态空间结构，涉及13个不同的参数:Fv, Fc, Fcs, alpha, beta, J, a_m, a_g, k_g1, k_g3, d_g, k_a, d_a（其中前六个固定，后七个待求）

二． IDNLGREY机械臂模型对象
上面的模型结构被输入到一个名为robotarm_c的C（mex）文件中。(整个文件可以通过命令“type robotarm_c.c”查看)。在状态更新函数中，注意我们在这里使用了两个中间双变量，一方面是为了增强方程的可读性，另一方面是为了提高执行速度(tau在方程中出现了两次，但是只计算了一次)。

  /* State equations. */
  void compute_dx(double *dx, double *x, double *u, double **p)
  {
      /* Declaration of model pa