随着空间科学的不断发展^[1]，航天设备不断被送入太空，越来越多的废弃航天器停留在轨道上形成了空间垃圾，导致服役航天器与废弃航天器撞击的概率逐渐增加。2009年，美国的“铱-33”移动通信卫星与俄罗斯已废弃的“宇宙- 2251”军用通信卫星在西伯利亚北部上空相撞，引发了如何处理废弃卫星的问题^[2-4]。清除地球静止轨道空间碎片^[5]不但可以降低卫星相撞事故的概率，还可以节省宝贵的轨道位置资源，为今后的发展提供更多可能^[6-8]。

针对空间碎片移除问题，Starke等^[9]提出通过飞网和飞爪系统抓获废弃航天器碎片，并利用绳索拖拽离轨的ROGER系统，并对两种方案进行了分析。赵国伟等^[10]考虑了收集飞行器与废弃飞行器的姿态运动以及系统组合体的整体运动，建立了绳系拖拽离轨系统动力学模型，分析了摆动特性与平稳控制。朱仁璋等^[11]则讨论了稳定状态下绳系系统的运动中心及相关力学问题，确定了运动中心在绳系上的位置，并导出了绳系拉力的表达式。王东科等^[12]则针对空间绳系机器人对目标抓捕后的复合体姿态稳定控制问题进行了研究，分析了绳系连接到目标后，由于自旋和碰撞所带来的姿态不稳定问题。

考虑到太空实验的成本和难点，大量有关全物理仿真实验、半物理仿真实验的研究逐步开展。文献[13]利用气浮台来仿真卫星在外层空间的低摩擦运动；文献[14]基于滑台式气浮台, 采用回转中心与质心重合方法解决了气浮台力矩平衡难题；文献[15]以单轴气浮台作为卫星刚性主体，进行了大型挠性结构卫星全物理仿真试验的研究。全物理仿真可以将大型实物，或者较难实物试验的装置进行等效仿真试验。文献[16]利用物理建模的思想, 用SimMechanics建立了空间对接地面半物理仿真台的机械系统, 用Matlab/Simulink建立了控制系统模型, 搭建了虚拟空间对接地面半物理仿真台，验证了试验的可行性；文献[17]基于dSPACE实时仿真机、单轴气浮转台、星载计算机及陀螺和反作用飞轮的卫星姿态控制系统，设计了半物理仿真，并对仅用反作用飞轮的卫星姿态大角度机动控制模式进行了半物理仿真验证；文献[18]针对对接机构地面测试的空间对接半物理仿真系统，提出系统的动力学解算方案，构建动力学仿真大回路, 着重对大回路仿真系统的稳定性进行分析。在航天工程中，半物理实验常常结合全物理仿真与实物进行试验^[19]，进行更加可靠的可行性验证。

考虑到回收拖拽过程中系绳中张力测量误差以及等效误差，本文通过半物理实验的方法，对系绳拖拽过程进行仿真验证，进行了一种绳系卫星拖拽系统的地面半物理实验。首先，设计了由半物理仿真单元和全物理加载单元组成的机械结构，分析了等效误差；然后，基于绳系卫星拖拽动力学模型，考虑到实际试验过程中存在的张力测量误差扰动，设计了基于模型预测的控制策略；最后通过仿真，验证了基于所提控制策略的试验响应误差不大于5%，张力误差不大于5%，验证了半物理实验装置的可行性。

1 绳系卫星半物理机械系统设计 1.1 绳系拖拽过程动力学建模

废星绳系拖动控制动力学模型原理如图 1所示。

$ \dot{x}_{1}=v_{0}^{1} $

(1)

$ \dot{x}_{2}=v_{0}^{2} $

(2)

$ M_{1} \ddot{x}_{1}=F_{\mathrm{t}}+f_{\mathrm{s}} $

(3)

$ M_{2} \ddot{x}_{2}=f_{\mathrm{s}} $

(4)

$ l_{\mathrm{m}}=\frac{f_{\mathrm{s}}}{k} \times \frac{\left(x_{0}+x_{1}+x_{2}\right)}{0.1} $

(5)

式中: v₀¹为平台卫星初速度, v₀²为目标星初速度, M₁为平台卫星质量, M₂为目标星质量, F_t为飞行器推力, f_s为试验台张力传感器测试张力, l_m为张力控制单元位移指令, k为对应绳单元(10 cm)等效刚度，由实验测试出, x₀为飞行器间初始距离。

实验中全物理仿真单元通过飞轮配重等效惯量实现被动加载功能，等效卫星质量；半物理仿真单元模拟目标在绳系拖拽作用下的动力学状态，测量绳系张力作为输入，实时计算不同平台质量和绳系长度下收放过程载荷，整体半物理装置如图 2所示。

1.2 全物理仿真单元等效原理及机构设计

全物理仿真单元用于被动模式下对绳系收放试验产品进行加载测试，具备全物理负载模拟与绳系卷绕等功能。其结构如图 3所示。

大理石平台通过支撑调整组件安装于地面，通过气源气体节流后形成的承载气膜实现运动组件低摩擦的漂浮于大理石平台上。卷筒用于卷绕绳系，卷筒轴固定于卷筒上，与回转支撑轴通过调心球轴承实现低摩擦转动。惯量盘为全物理加载的模拟负载，其转动惯量模拟空间卫星等效质量。

为保证全物理仿真单元负载模拟的有效性和准确性，对全物理加载单元负载等效模拟进行计算。设全物理试验系统的传动效率为η，惯量盘的转动惯量和角速度分别为J_s和ω_s，卷筒转动惯量、直径和角速度分别为J_r、D_r和ω_r，卫星等效质量为M_s，则通过等效原理可得

$ \frac{1}{2} M_{\mathrm{s}} v_{t}^{2} \eta=\frac{1}{2} J_{\mathrm{s}} \omega_{\mathrm{s}}^{2} $

(6)

将上式整理可得

$ w_{\mathrm{s}}=\sqrt{\frac{M_{\mathrm{s}} v_{\mathrm{t}}^{2} \eta}{J_{\mathrm{s}}}} $

(7)

当星体等效质量M_s=500 kg、传动效率η=0.99、模拟转动惯量为40~260 kg·m²(其中气足安装板、气足及其他运动组件转动惯量为40 kg·m²)时，惯量盘的最大角速度为16.80 r/min。

同理，当星体等效质量为M_s=2 000 kg时，模拟转动惯量为220 kg·m²，其他不变时，惯量盘的角速度为13.18 r/min。则通过计算可得，在最大等效质量模拟工况下，惯量盘的最大直径为725 mm。

通过对等效质量模拟的极限情况分析可得，能够进行等效模拟的卫星的质量范围为：290~3 100 kg，相应的模拟惯量的旋转角速度范围分别为1.33~1.76 rad/s，经过换算得到角速度的表示方法为：76.4~100.8°/s和12.7~16.8 r/min。

1.3 半物理仿真单元设计及误差估计

半物理加载单元主要由支撑基座、绳系累计收放长度测量单元、绳系张力测量单元、张力控制单元、卷扬机构和直线移动单元组成，其结构如图 4所示。

图 4中，半物理加载单元通过基座整体固定于加载单元试验小车上；绳通过绳系收放测量单元测量累计收放长度，并且夹轮组将绳系由水平方向换向至铅锤方向；绳穿过绳系收放长度测量单元后与张力测量单元连接，测量即时绳内张力的大小；绳通过张力测量单元的定滑轮换向后与张力控制单元连接，张力控制单元的电动缸伸缩运动调节绳系中的张力大小，从而实现张力的精确闭环控制。

1.3.1 绳系张力测量单元

绳系张力测量单元用于对绳系张力进行高精度动态测量，其主要由高精度测力传感器和定滑轮组件组成，结构如图 5所示。其中，绳穿过定滑轮将运动模式和动力传递给张力控制单元。

1.3.2 张力测量误差分析

由于张力测量单元的测量精度受到张力测量传感器误差、滑轮摩擦力和惯量力等因素的影响，张力测量精度要从这三方面去考察。

设收放实验产品绳系的张力为F，通过张力测量单元夹轮组后的张力为F₁，通过张力测量单元定滑轮后绳系的张力为F₂，张力测量单元的测量误差及张力测量值分别为u_F和F′，则根据绳系传递过程可得到：

$ F^{\prime}=\frac{F_{1}+F_{2}}{2}\left(1 \pm u_{\mathrm{F}}\right) $

(8)

$ F_{1}-F=F_{I_{1}} $

(9)

$ F_{2}-F_{1}=F_{I_{2}} $

(10)

式中F_I1、F_I2分别为绳系长度测量单元夹轮组和张力测量单元定滑轮组引入的惯性力和摩擦力，绳系张力测量表达式如下：

$ u_{F}=F^{\prime}-F $

(11)

当绳系拖曳力的真值为10 N时，通过式(11)计算得到，张力测量值的误差为0.029 6 N；将大量程张力测量真值(10~50 N)带入式(11)中，可得到张力测量的最大误差为0.3%。

2 绳系卫星半物理仿真控制设计 2.1 控制电机建模 2.1.1 驱动电机模型

驱动电机是扭振系统及加载系统的运动输出执行环节，驱动电机一般由高速小惯量的直流电动机和减速环节等组成。驱动电机在经过减速器后，通过扭矩传感器分别与扭振系统及加载电机连接，驱动电机的模型包括电机开环模型(如图 6所示)，PWM驱动器模型及电机常用的速度、电流双闭环模型。

2.1.2 开环模型

从系统输入U到输出θ_r的开环传递函数：

$ \frac{\theta_{\mathrm{r}}(s)}{U(s)}=\frac{K_{\mathrm{PWM}} K_{\mathrm{T} 1}}{J_{\mathrm{r}} L_{\mathrm{r}} s^{3}+\left(R_{\mathrm{r}} J_{\mathrm{r}}+B_{\mathrm{r}} L_{\mathrm{r}}\right) s^{2}+\left(B_{\mathrm{r}} R_{\mathrm{r}}+K_{\mathrm{T1}} K_{\mathrm{e1}}\right) s} $

(12)

式中：U为直流伺服电机电枢电压, R_r为直流伺服电机电枢回路总电阻, L_r为直流伺服电机回路总电感, θ_r为直流伺服电机角位移, K_e1为直流伺服电机反电动势系数, K_T1为直流伺服电机转矩系数, J_r为等效转动惯量, B_r为直流电机阻尼系数。

2.1.3 PWM驱动环节

PWM脉宽调节装置由脉宽调制器和PWM变换器组成，由于PWM的变换电压要到下一个周期才能改变，因此PWM脉宽调节装置相当于延时环节。

设PWM脉宽调节装置的开关频率为f_PWM，对应的周期T_PWM即为延时时间，其传递函数为：

$ G_{\mathrm{PWM}}=K_{\mathrm{PWM}} \mathrm{e}^{-T_{\mathrm{PWM}}s} $

(13)

$ K_{\mathrm{PWM}}=\frac{u_{\mathrm{d}}}{u_{\mathrm{c}}} \dot{x}_{1}=v_{0}^{1} $

(14)

式中：K_PWM为放大系数, u_d为PWM变换器的输出电压, u_c为脉宽调制器的控制电压。

由于上式中包含指数函数，使得系统为非最小相位系统。为便于分析，将其按泰勒级数展开，则

$ \begin{aligned} &G_{\mathrm{PWM}}=K_{\mathrm{PWM}} \mathrm{e}^{-T_{\mathrm{PVM}}{s}}=\\ &\frac{K_{\mathrm{PWM}}}{1+T_{\mathrm{PWM}} s+\frac{\left(T_{\mathrm{PWM}} s\right)^{2}}{2 !}+\frac{\left(T_{\mathrm{PWM}} s\right)^{3}}{3 !}+\cdots} \approx \frac{K_{\mathrm{PWM}}}{T_{\mathrm{PWM}} s+1} \end{aligned} $

(15)

由于PWM周期很小，常用的开关频率为10 kHz，因此忽略式(15)中的高阶项，PWM脉宽调节装置近似为一阶惯性环节，其传递函数为

$ G_{\mathrm{PWM}}=K_{\mathrm{PWM}} \mathrm{e}^{-T_{\mathrm{PWM}} {s}} \approx \frac{K_{\mathrm{PWM}}}{T_{\mathrm{PWM}} s+1} $

(16)

2.1.4 转速、电流双闭环电机模型

电机采用转速、电流反馈双闭环控制驱动，其闭环数学模型如图 7所示。其中K_PWM为PWM驱动环节放大系数，T_PWM为PWM驱动环节延时周期，R_r为直流伺服电机电枢回路总电阻，L_r为直流伺服电机回路总电感，θ_r为直流伺服电机角位移，K_e1为直流伺服电机反电动势系数，K_T1为直流伺服电机转矩系数，J_r为等效转动惯量，B_r为直流电机阻尼系数。

2.2 位置模式控制

整个系统的控制可以等效为测试加载机构的张力控制和绳系拖拽机构的位置控制实现。张力控制通过力传感器引入张力值与设定值的差值控制绳子中的张力，为位置控制提供输入。位置内环控制器由实际位置、张力测得位置变化值及给定位置三者同时确定输入值。为提高张力控制模式响应速度，系统设计一种基于位置内环的高精度张力控制单元，如图 8所示。

基于位置内环的张力控制策略采用用力传感器构建力外环，把力传感器的输出作为内环位置控制器输入的修正值。其控制原理示意图如图 8所示，由力传感器测量实际张力F_t，并F_t将与期望力F_d进行比较，经过力外环控制器得到张力系统位置修正量ΔX_d，其表达式为

$ \Delta X_{\mathrm{d}}=K_{\mathrm{fp}}\left(F_{\mathrm{d}}-F_{\mathrm{t}}\right)+K_{\mathrm{fl}}\left(\dot{F}_{\mathrm{d}}-\dot{F}_{\mathrm{t}}\right) $

(17)

式中: K_fp、K_fI分别为力外环控制的比例和积分系数。从而修正后的位置控制量为

$ X_{\mathrm{d}}=X_{\mathrm{p}}+\Delta X_{\mathrm{d}}-X_{\mathrm{t}} $

(18)

式中，X_t为位置反馈量，修正后的位置指令为X_d。

由于位置内环控制较为简单，采用PID控制器就可以简单的实现效果较好的稳定。控制律如下：

$ u_{X}=K_{\mathrm{p}}\left[X_{\mathrm{d}}+\frac{1}{T_{i}} \int\limits_{0}^{t} X_{\mathrm{d}} \mathrm{d} t+T_{\mathrm{d}} \frac{\mathrm{d} X_{\mathrm{d}}}{\mathrm{d} t}\right] $

(19)

2.3 张力模式控制

基于张力控制的精度要求，在静态过程中，简单的PID控制器可以较有效的控制并输出力误差，用于修正位置。然而，在实际拖动过程中，由于张力测量单元的误差，若要使绳内张力始终保持在较小的动态变化中，从而使位置变化量在较小的精度范围内稳定控制，需要用到更复杂的模型预测控制。

模型预测控制采用阶跃响应模型，节省了控制器的设计时间，也同时为动态控制提供了条件。考虑到张力控制系统为SISO的定常稳态系统，通过单位阶跃输入下的输出响应{0, s₁, s₂, …, s_N, s_N+1, …}，假设系统N步后达到稳态，则系统输出可以表示为

$ y(k)=\sum\limits_{l=1}^{N-1} s_{l} \Delta u(k-l)+s_{N} u(k-N) $

(20)

式中，Δu(k-l)=u(k-l)-u(k-l-1)。

这样在时刻k就可以预测到未来任意时刻的输出值：

$ \bar{y}(k+1 \mid k)=\bar{y}_{0}(k+1 \mid k)+s_{1} \Delta u(k) $

(21)

在控制时域M，预测时域N时，P时刻预测值为

$ \begin{aligned} &\bar{y}(k+P \mid k)=\bar{y}_{0}(k+P \mid k)+s_{P} \Delta u(k)+ \\ &s_{P-1} \Delta u(k+1 \mid k)+\cdots+s_{P-M+1} \Delta u(k-M+1 \mid k) \end{aligned} $

(22)

式中，M≤P≤N。

在张力控制中，设定性能指标为

$ J(k)=\sum\limits_{i=1}^{P} \omega_{i} e^{2}(k+i \mid k)+\sum\limits_{j=1}^{M} r_{j} \Delta u^{2}(k+j-1 \mid k) $

(23)

式中：ω_i、r_j都是非负的标量；e(k+i|k)=y_s(k+i)- y(k+i|k)为张力误差量；y_s(k+i)为未来输出的设定值。式中第1项为最小化张力误差量，第2项为抑制过于剧烈的控制增量，防止张力控制出现震荡造成控制不当。通过最小化性能指标函数J(k)，求解每一步的控制增量Δu，滚动优化到下一时刻中去。这样，由模型预测控制算法所得设计的控制框图如图 9所示。

2.4 稳定性分析

考虑系统的输入输出，现将系统等效为二阶振荡环节做稳定性分析。

$ G(s)=\frac{k_{0}^{2}}{s^{2}+2 a k_{0} s+k_{0}^{2}} $

(24)

系统经采样和零阶保持器后的Z传递函数为：

$ G(z)=\frac{m(z)}{p(z)}=\frac{m_{1} z^{-1}+m_{2} z^{-2}}{1+p_{1} z^{-1}+p_{2} z^{-2}} $

(25)

式中，p₁=-2ecos 2πT，p₂=e²，e=exp(-2πTs₂)。根据张力模式控制中给出的控制率可知，y₀(k+j)为k时刻对于Δu(k)=0的原始输出，那么令a_i=1-eⁱ(s₂sin 2πTi+cos 2πTi)为系统阶跃响应采样值，在不考虑模型失配时可得系统闭环传递函数：

$ F_{\mathrm{s}}(z)=\frac{d_{\mathrm{s}} m(z)}{p^{*}(z)} $

(26)

式中，${d_i} = \frac{{{a_i}}}{{\sum\limits_{j = 1}^p {a_j^2} }}, {d_s} = {d_1} + {d_2} + \cdots + {d_p}$。

当不考虑控制权时，有下式：

$ p^{*}(z)=1+p_{1}^{*} z^{-1}+p_{2}^{*} z^{-2} $

(27)

其中：$p_1^* = {b_{20}} - 1 - 2e\cos 2\pi T, p_2^* = {e^2}\left( {{b_{20}} - } \right.\left. {{d_{p0}}{a_{p + 1}}} \right), {b_{20}} = \frac{{\sum\limits_{j = 1}^p {{a_j}} {a_{j + 1}}}}{{\sum\limits_{j = 1}^p {a_j^2} }}, {d_{p0}} = \frac{{{a_p}}}{{\sum\limits_{j = 1}^p {a_j^2} }}$。

根据式(26)、(27)，系统稳定性取决于分母多项式p^*(z)的系数，其稳定的充要条件为

$ \left\{\begin{array}{l} 1+p_{1}^{*}+p_{2}^{*}>0 \\ 1-p_{1}^{*}+p_{2}^{*}>0 \\ 1>p_{2}^{*} \end{array}\right. $

(28)

根据p₁^*，p₂^*定义可知，上式成立。则在式(23)的性能指标设定下，不考虑控制权值时，不论优化时域p取何值，系统总是稳定的。证毕。

3 仿真校验 3.1 位置模式仿真分析

通过编码器测量输出端位置，利用Simulink搭建仿真模型，调节PID控制器各参数，分析系统频率响应特性及Bode图如图 10~12所示。

通过正弦跟踪可以看出，通过PID校正，系统的位置模式控制可以达到良好的跟踪效果。系统的相角裕度P_m=22.8°，可以看出系统的稳定的。同时，系统在10 Hz时相角位移为1.84°，系统在10 Hz以内有较好的频率响应特性。

3.2 张力模式仿真分析

外环张力控制时，由于系统仿真救援是一个弹簧系统，拖拽过程中会产生一定频率的误差；同时由于张力测量误差的存在，可将这类干扰视为有界干扰噪声。仿真过程中，在系统输入信号中加入有界随机噪声作为抖动干扰，测试模型预测控制对干扰的抑制作用。可以看到系统频率响应如图 13、14所示。

3.3 模型预测控制控制器仿真

模型预测控制器设计中，采样间隔取0.4 s，预测时域长度为30 s，控制时域长度取10 s，输入权重为0.3，输出权重为1.0，暂不考虑绳系抖动因素，不加入约束条件。则控制器设定时相关参数曲线如图 15所示。

可以看出，通过设计模型预测控制器，即使加入了噪声信号，如图 16所示，频率响应曲线仍可以较好的跟踪输入，为力误差反馈提供了较好的条件。

3.4 误差分析

通过实验测量，得到如图 17、18所示张力误差曲线及联合电机阶跃响应曲线。可以看出，在实际工况中，力误差降低到了需求的5%以内，响应曲线也可以较好的满足试验要求。

4 结论

1) 本文针对空间废弃卫星绳系系统的回收过程，采用半物理仿真实验技术设计了一种地面仿真实验装置。通过对绳系拖拽过程进行动力学建模，设计了实验装置。实验装置分为全物理加载单元和半物理仿真单元，分别用来模拟卫星等效质量和卫星在绳系拖拽作用下的动力学状态。通过对两部分实验装置的参数分析和误差估计，实现了绳系拖拽系统的机械设计。

2) 为解决绳系拖拽半物理仿真实验的控制问题，建立了基于伺服系统的双电机模型。针对绳系拖拽系统的验证问题，将整个系统的控制等效为了测试加载机构的张力控制和绳系拖拽机构的位置控制，通过双闭环电机的仿真实验，验证了所设计控制系统的良好频率响应特性。

3) 考虑到实际绳系拖拽过程中会出现的张力测量误差扰动，通过将张力控制性能指标设为优化函数，并进行稳定性分析，设计了可抑制扰动的模型预测控制控制器。通过等效张力误差为有界干扰噪声的方法，进行了一组仿真验证。仿真结果表明，力误差降低到了需求的5%以内，实现了较好的控制效果。