随着人类对空间探索的不断深入和拓展，空间机器人被期望在太空在轨服务中执行高风险、复杂的操作任务，因此其相关动力学与控制技术受到各空间技术大国的广泛关注．考虑空间机器人的捕获操作过程是实现在轨服务最具挑战性的阶段，其研究已成为当前航天领域的热点[1-3]．目前，关于空间机械臂在轨捕获的研究已取得了一定的成果．文献[4]分析了空间机器人在轨插拔孔操作，并提出基于事件触发的滑模控制算法．文献[5] 采用伪谱法实现了空间柔性机器人的能耗最小控制．文献[6]考虑了近距离捕获目标的控制问题．文献[7]结合无差方法研究了柔性空间机器人的建模方法．以上研究关注的都是单臂捕获操作．针对双臂空间机器人系统，文献[8]研究了基于模糊神经网络的双臂捕获对接控制问题．文献[9]提出一种基于有效估计参数的双臂协同捕获方法．文献[10]针对捕获前双臂空间机器人协调轨迹规划问题提出基于粒子群的控制策略．文献[11]研究了一种空间人双臂捕获卫星的阻抗控制方案．由于当捕获具备高速、旋转特性的非合作卫星时，空间机器人会遭受巨大的碰撞冲击力矩，若冲击力矩过大，则很可能对最薄弱的关节处造成破坏，并最终导致整个任务的失败，因此针对双臂捕获操作的柔性控制有着重要的探索价值和意义．考虑到柔性机构——旋转型串联弹性执行器(RSEA)当地面机器人与外部环境发生碰撞时，对关节避免破坏具备良好的性能[12-13]，基于此，本研究针对空间机器人双臂捕获非合作卫星过程避免关节受冲击破坏的控制问题，将RSEA装置引入到空间机器系统．然而，RSEA装置的引入也会使系统关节具备柔性，使得配备了柔性机构的空间机器人在轨捕获的缓冲柔顺控制更为复杂．基于上述分析，本研究针对空间机器人在轨捕获柔性控制方法进行研究．采用奇异摄动法进行快慢变子系统分解，针对快变子系统，利用速度差值反馈控制器抑制关节柔性振动．考虑闭链混合体系统的动力学参数难以精确获得，因此针对慢变子系统，提出一种缓冲柔顺强化学习控制方法．所提强化学习控制策略通过与环境动态交互得到评价信号，并运用评价信号对系统代价函数进行估计，以实现系统控制的最大未来回报值[14-16]；同时，应用执行网络(ANN)单元对系统未知非线性函数进行逼近，并结合评价网络(CNN)对其进行优化，最终使控制器趋向最优状态，提高了空间机器人系统在轨操作的智能性．仿真结果表明了所提柔性控制的有效性．1 系统动力学及缓冲柔顺策略漂浮基双臂空间机器人系统由载体B0和刚性左右臂Bi(i=1,2,⋯,6)组成，如图1所示．选取O为惯性坐标系原点，建立惯性坐标系XOY，同时建立各分体的连体坐标系XiOiYi(i=0,1,⋯，10.13245/j.hust.240607.F001图1双臂空间机器人与目标卫星6)，图中：O0为载体质心；Oi为各杆关节铰中心；XtOtYt为卫星坐标系；Ot为卫星质心；ψ1和ψ2分别为O0O1和O0O4与X0的夹角；ri为各分体质心矢量；rc为系统总质心矢量；θL=[θ1,θ2,θ3]T为左臂各关节转角；θR=[θ4,θ5,θ6]T为右臂各关节转角；θm=[θ1m,θ2m,θ3m,θ4m,θ5m,θ6m]T为关节电机转角；θ0为载体姿态角．RSEA结构如图2所示，图中：R为扫臂半径；r为弹簧半径．当冲击力矩大于极限值时，设置关机阈值让所有电机关停．为避免电机频繁开关机影响性能，基于此，本研究设置电机开机阈值，当冲击力矩低于开机阈值时，电机再次开启．10.13245/j.hust.240607.F002图2RSEA结构图应用第二类拉格朗日建模方法，推导得捕获前空间机器人系统动力学方程为D(q)q¨+H(q,q˙)q˙=τc+JTF;Jmθ¨m+K(θm-θ)=τm;K(θm-θ)=τθ, (1)式中：D(q)∈R9×9为系统对称正定质量矩阵；H(q,q˙)q˙∈R9×1为系统包含科氏力和离心力项；q=[x0,y0,θ0,θLT,θRT]T，其中x0和y0为基座线位移；q˙与q¨为q对应的一阶及二阶导数，后续同类项亦为一阶导数及二阶导数；θ=[θLT,θRT]T；τc=[τaT,τ0,τLT,τRT]T，其中，τa=[0,0]T为载体位置控制力矩，τ0为载体姿态控制力矩，τL和τR∈R3×1分别为左右臂各关节输入力矩；τm=[τ1m,τ2m,τ3m，τ4m,τ5m,τ6m]T为各电机输出力矩；Jm=diag(J1m，J2m,J3m,J4m,J5m,J6m)为电机转子转动惯量；K=diag(k1,k2,k3,k4,k5,k6)为关节等效刚度，其中k1～k6为各关节等效刚度；J∈R6×9为机械臂末端对应的运动雅克比矩阵，由左右臂末端碰撞点运动雅克比矩阵JL和JR元素组成，其中JL和JR∈R3×6；F∈R6×1为机械臂末端所受作用力．采用牛顿-欧拉法获得捕获前卫星系统的动力学方程为Dtq¨t=JtTF'，(2)式中：Dt∈R3×3为被捕获卫星对称正定质量矩阵；qt=[xt,yt,θt]T为其广义坐标，其中，xt和yt为卫星线位移，θt为卫星姿态角；Jt∈R6×3为抓取点对应的速度雅克比矩阵；F'∈R6×1为末端执行器的反作用力，满足F'=-F．在捕获完成后，空间机器人左右臂满足J0Lθ˙L=J0Rθ˙R，(3)式中J0L和J0R∈R3×3分别为连体坐标系下左右臂末端速度雅克比矩阵．基于式(3)，定义qL=[x0,y0,θ0,θLT]T，其满足q˙=ATq˙L，(4)式中：A=[I6×6,A1T]T，其中A1=[03×3,J0R-1J0L]，I6×6为6阶单位矩阵，03×3为3阶零矩阵．对式(1)第一项进行积分，并整理得D(q)(q˙(t0+Δt)-q˙(t0))=JTP，(5)式中P=∫t0t0+ΔtFdt为捕获过程的碰撞冲量．考虑捕获碰撞时间Δt极小，则碰撞力可近似为F=P/Δt．(6)在捕获完成后，闭链混合体系统动力学方程[8]为Dbq¨L+Hbq˙L=τb-FbI，(7)式中：Db=Q+BJ1L-1JL，其中，Q=ADAT，B=AJT(JtT)+Dt；Hb=A(HAT+DA˙T)+BJ1L-1(J˙L-J˙1L∙J1L-1JL)；τb=Aτc；FbI=AJTFI．由矩阵运算可知FbI=04×1，则式(7)可改写为Dbq¨L+Hbq˙L=τb，(8)式中：τb=[τb1T,τb2T]T；τb1和τb2∈R3×1．为得到完全能控形式动力学方程，将式(8)改写为如下分块矩阵形式，即Db11Db12Db21Db22q¨aq¨θ+Hb11Hb12Hb21Hb22q˙aq˙θ=02×1τn, (9)式中：qa=[x0,y0]T；qθ=[θ0,θLT]T；τn=[τ0,τb2]T．由式(9)第一行可解得q¨a的表达式，有Dnq¨θ+Hnq˙θ=τn; Jmθ¨m+K(θm-θ)=τm;K(θm-θ)=τθ , (10)式中：Dn=Db22-Db21Db11-1Db12；Hn=Hb22-Db21Db11-1Hb12．2 控制器设计闭链混合体控制器包含快变子系统控制器及慢变子系统控制器，因此系统的总控制律为τm=τs+τf，(11)式中：τs∈R6×1为慢变子系统控制力矩；τf∈R6×1为快变子系统控制力矩．引入奇异摄动因子μ，并使其满足K=K1/μ2，(12)式中K1∈R6×6为对角、正定矩阵．将式(12)代入式(10)的下面两项，可得到快变子系统动力学方程为μ2τ¨θ =Jm-1K1(τm-Jmθ¨-τθ )．(13)为了抑制系统关节柔性振动，设计如下快变子系统控制律，即τf=-Kf(θ˙m-θ˙)，(14)式中Kf=K2/μ，其中K2∈R6×6为对角、正定矩阵．将式(11)和(14)代入式(13)，可得μ2Jmτ¨θ =K1(τs-Jmθ¨-τθ )-μK2τ˙θ ．(15)当μ→0时，K→∞，进而可得出慢变子系统的动力学方程为Dcq¨θ+Hcq˙θ=τcθ，(16)式中：Dc=Dn+Jmθ，其中Jmθ=[04×1,Jmθ1]T，Jmθ1=[03×1,(A2Jm)T]T；Hc为当θ˙=θ˙m时Hn对应的新矩阵；τcθ=[τ0,(A2τs)T]T．须要注意的是，混合体系统满足以下重要性质．性质1 Dc为对称、正定矩阵．性质2 D˙c-2Hc为交错矩阵，即满足xT(D˙c-2Hc)x=0,  ∀ x∈R4×1．为便于后续控制设计，选取x1=qθ，x2=q˙θ，则式(16)可表示为x˙1=x2;x˙2=-A3x2+B3, (17)式中：A3=Dc-1Hc；B3=Dc-1τcθ．定义轨迹跟踪误差为s1=x1-qθd，(18)式中qθd∈R4×1为系统期望轨迹．为评估系统动态交互过程控制策略的优化程度，选取系统代价函数为Q=∫t∞e(a-t)/λf(a)da，(19)式中：λ为正常数；f(t)为即时代价函数，满足f(t)=s1TΛ1s1+τcθTΛ2τcθ，(20)式中Λ1和Λ2∈R4×4为正定、对角矩阵．考虑当交互过程达到最小代价函数时，闭链混合体系统实现最优控制．鉴于代价函数无法直接获得，因此引入评价网络对代价函数进行逼近，即Q=WcTΦc(x)+ς，(21)式中：Wc∈Rn×1为RBF神经网络理想权值矩阵；Φc(x)=[Φ1,Φ2,⋯Φn]T，其中Φi(i=1,2,⋯，n)为基函数；ς为最优逼近误差．神经网络输入选取为x=s1．同时，代价函数估计值可有Q^=W^cTΦc(x)，(22)式中W^c为理想权值Wc的估计值，定义估计误差W˜c=W^c-Wc，且满足W˜˙c=W^˙c．结合式(19)，可得代价函数近似误差[14]为z=f-Q^/λ+Q^˙．(23)当常数值λ→∞时，则式(23)可表述为z=f+Q^˙=f+∇Q^x˙，(24)式中∇表示对x的梯度函数．选取CNN权值自适应律为W^˙c=-Kc∂Lc∂Wc，(25)式中：Kc为正常数；Lc=zTz/2．基于此，并结合式(24)，式(25)可改写为W^˙c=-Kcz(t)∂z∂W^c=-Kc(W^cTΓ+f)Γ，(26)式中Γ=∇Φcx˙-Φc/λ．选取虚拟控制变量η1，其满足η1=-λ1s1+q˙θd，(27)式中λ1为对角矩阵．同时，选取控制量s2，即s2=x2-η1．(28)对式(28)进行求导，可得s˙2=x˙2-η˙1=-Dc-1(τcθ-Hcx2)-η˙1．(29)选取如下形式李雅普诺夫函数，即V1=(s1Ts1+s2TDcs2)/2．(30)对式(30)进行求导，可得V˙1=s1Ts˙1+s2TDcs˙2+s2TD˙cs2/2=s1Ts2-s1Tλ1s1+s2T(τcθ-Dcη˙1). (31)基于以上设计，理想控制器模型可选取为τcθ=-s1-λ2s2+Dcη˙1，(32)式中λ2为对角矩阵．由于捕获过程混合体系统动力学模型难以精确获得，基于此，采用执行网络模块对其进行逼近，慢变子系统控制器选取为τcθ=-s1-λ2s2+W^aTΦa(xa)，(33)式中W^a为理想权值Wa的估计值．定义估计误差W˜a=W^a-Wa，且满足W˜˙a=W^˙a；同时，网络最优逼近误差为ςe．神经网络输入选取为xa=[qθT,q˙θT,qθdT,q˙θdT,q¨θdT]T．定义执行网络误差为E=ςa+λ3(Q^-Qd)，(34)式中：ςa=W˜aTΦa(xa)为拟合误差；λ3为对角矩阵；Qd=[0,0,…,0]T∈Rn×1为期望理想代价函数．定义La=ETE/2，选取ANN权值自适应律为W^˙a=-Ka∂La∂Wa，(35)式中Ka为正常数．结合La表达式，式(35)可写为W^˙a=-Ka(ςa+λ3Q^)Φa．(36)考虑ςa无法直接获得，因此式(36)改写为W^˙a=-Ka(W^aTΦa+λ3Q^)Φa．(37)定义如下李雅普诺夫函数V=(s1Ts1+s2TDcs2+W˜cTW˜c+W˜aTW˜a)/2．(38)当近似误差z→0时，有f=Q/λ-Q˙，可得f=-WcTΓ+ε1，(39)式中ε1=ς/λ+∇ςx˙，且满足ε1≤εc,max．对式(38)求导，结合式(39)，可得V˙=-s1Tλ1s1-s2Tλ2s2+s2T(W˜aTΦa-ςe)-KcW˜cT(W˜cTΓ+ε1)Γ-KaW˜aT(W^aTΦa+λ3Q^)Φa. (40)结合式(26)可知Q^TQ^≤2(W˜cTΦc)TW˜cTΦc+2(WcTΦc)TWcTΦc，代入式(40)，可得V˙≤-s1Tλ1s1-s2T(λ2-E)s2-(Ka-1)Φa2W˜a2/2-(KcΓ2-2Ka(λ3)2Φc2)W˜c2/2+KaΦa2⋅Wa2/2+Ka(λ3)2Φc2Wc2+ςe2/2+Kcεc,max2/2≤-cV+ρ, (41)式中：c=min{2λmin(λ1),2λmin(λ2-E),ΦA(Ka-1),ΦC(KcΓ2-2Ka(λ3)2Φc2)}，其中，Φa2≤ΦA，Φc2≤ΦC；ρ=KaΦa2Wa2/2+Ka(λ3)2Φc2Wc2+ςe2/2+Kcεc,max2/2．选取cρ/ϕ，则当V=ϕ时，V˙≤0为一个不变集，即当V(0)≤ϕ时，对t0，恒有V(t)≤ϕ．求解式(41)可得0≤V≤ρ/c+(V(0)-ρ/c)e-ct．结合李雅普诺夫稳定性定理，可知最终V(t)以ρ/c为界，系统轨迹跟踪误差最终一致有界．3 仿真结果与分析3.1　捕获过程系统缓冲性能模拟仿真所采用的空间机器人模型如图1所示，其动力学参数为：m0=150 kg，I0=40 kg∙m2，l0=1.062 m；mi=10 kg，Ii=5 kg∙m2，li=2 m，di=1 m (i=1,2,4,5)；mj=2.5 kg，Ij=2 kg∙m2，lj=0.5 m，dj=0.25 m (j=3,6)；mt=40 kg，It=18 kg∙m2，lt=0.5 m；ψ1=2.791 rad，ψ2=0.349 rad；Iim=0.05 kg∙m2，kim=1 kN/m(i=1,2,⋯,6)．关节等效刚度的计算公式[13]为K=2Km(3R2+r2)(2cos2φ-1)，式中：Km=diag(k1m,k2m,k3m,k4m,k5m,k6m)；R=0.1 m；r=0.01 m；φ为当空间机械臂末端施加载荷τ=[20 N∙m,20 N∙m,0 N,-20 N∙m，-20 N⋅m,0 N]T时扫臂的转角，选取φ=diag(3°，2°,1°,-3°,-2°,-1°)．采用配置/未配置RSEA装置的空间机器人对不同初速度(以空间机器人为参考)的卫星进行捕获模拟仿真，并计算得碰撞过程中关节所受冲击力矩，关节所受冲击力矩结果如表1所示．10.13245/j.hust.240607.T001表1缓冲性能模拟结果卫星初速度未配置RSEA所受最大冲击力矩/(N∙m)配置RSEA所受最大冲击力矩/(N∙m)冲击力矩最大降低百分比/%[0.5 m∙s-1，0.5 m∙s-1，0.0 rad∙s-1]T141.664.854.2[0.4 m∙s-1，0.0 m∙s-1，0.3 rad∙s-1]134.066.050.7[0.3 m∙s-1，0.3 m∙s-1，0.2 rad∙s-1]101.857.143.9由表1可看出：针对不同初速度卫星的捕获操作过程，配置RSEA装置后的空间机器人关节所受碰撞冲击力矩均得到有效降低，最大可减少54.2%的冲击力矩，最小也可减少43.9%，体现了其良好的缓冲性能．3.2　柔顺控制性能模拟为验证镇定运动过程柔顺控制策略的有效性，运用本研究所提控制方案进行数值仿真．系统的控制参数选取如下：K2=diag(5.5,5.5,5.5,5.5)，λ1=diag(1,1,1,1)，λ2=diag(400,400,400，400)，λ3=diag(0.1,0.1,0.1,0.1)；λ=0.3，ε=0.5，Kc=1，Ka=1，初始权值Wc=0100×4，Wa=060×4．初始q=[0.5 m,0.5 m,1  rad,2.27 rad，-1.40 rad,-1.22 rad,0.68 rad,1.87 rad,0.24 rad]T，当t0=0 s时对目标卫星进行捕获操作，此时目标卫星速度为vt=[0.3 m/s,0.3 m/s,-0.2 rad/s]T．捕获后期望位置选取为qθd=[π/18,2π/3,-π/3，-π/3,π/3,π/3,π/3]T rad．仿真时间选取为t=40 s．仿真结果如图3~8所示．10.13245/j.hust.240607.F003图3未采用缓冲柔顺控制的关节所受冲击力矩10.13245/j.hust.240607.F004图4采用缓冲柔顺控制的关节所受冲击力矩10.13245/j.hust.240607.F005图5评价网络范数10.13245/j.hust.240607.F006图6执行网络范数10.13245/j.hust.240607.F007图7载体姿态角轨迹10.13245/j.hust.240607.F008图8关节角轨迹图3给出了镇定过程未采用缓冲柔顺控制的关节所受冲击力矩τ．假设关节电机所能承受冲击力矩极限为50 N∙m，可发现此时的冲击力矩虽然因柔性机构的缓冲得到降低，但依然超出安全阈值，因此须结合主动开启、关闭电机策略进行控制．选取关机阈值为36 N∙m，开机阈值为3 N∙m，如图4中的虚线所示．图4为采用本研究所提缓冲柔顺控制时关节所受冲击力矩，对比图3与图4可知所提方法有效实现了对关节电机的保护．图5和图6分别为评价网络(CNN)与执行网络(ANN)的二范数NCNN与NANN，可发现其与环境动态交互过程，通过评价网络对控制信号不断评估，进而优化下一步控制策略，由于控制策略的优化，系统评价网络的值逐步降低，系统最终达到优化、稳定的状态．图7和图8为采用缓冲柔顺控制的镇定轨迹，分别采用所提强化学习(RL)算法、未采用评价网络的神经网络(NN)算法及比例-积分(PD)算法进行镇定控制，图中的绿色虚线为期望轨迹．PD算法选取比例参数为Kp=300，积分参数为Kd=300，神经网络参数选取与所提强化学习相同．PD算法具备结构简单、调整方便等特点，因此实时性好，且计算量小、存储数据量少，但也存在对外界环境变化的适用性及学习能力差的缺点；NN算法具备较好的自学习和自适应能力，但相较于PD控制，也存在实时性差、计算存储数据量大的特点；而RL算法能够在常规NN的基础上进一步改善控制性能，通过与外界的交互学习优化NN性能．对比RL算法与NN算法镇定轨迹可知：由于NN未采用评估网络对控制过程进行优化，因此其控制器相较于RL优化程度更低，进而收敛速度更慢．而对比RL算法与PD算法可知：由于RL算法具备较好自主学习能力，因此可更快实现对闭链混合体系统的镇定，有效提升了系统对环境的适应能力．4 结语为了避免冲击力矩对关节造成破坏，本研究在关节电机与机械臂之间设计了一种柔性机构，并提出了一种与之配合的主动开启、关闭电机的柔性控制方法．通过仿真试验可知：所提策略在捕获碰撞阶段最大可减少54.2%的冲击力矩，最小也可减少43.9%，体现了良好的缓冲性能．在镇定过程中，所提基于柔性机构的强化学习柔顺控制方法保证了关节所受冲击力矩限定在安全范围，从而实现了对关节的保护．与传统神经网络控制方法对比可知：所提强化学习控制方案可通过动态交互实现优化，从而减少空间机器人系统的控制成本，有效地延长在轨服务年限．此外，在仿真试验及计算机硬件条件达到一定条件后，可为空间机器人在轨捕获操作的柔顺控制实际操作和应用提供技术参考．