航天器质量辨识方法研究

航天器质量辨识方法研究

0关于质量特性的轨辨识方法为了进行高精度的轨道和姿态控制，必须正确识别控制对象的质量特征参数（包括质量、强度和质量中心位置）。因此，在卫星发射前，地面通过各种方法计算（基于ct技术）测量质量特征。然而,由于所建模型不可能与实际完全一致(如某些装配件常被省略),基于CAD技术的计算结果可能有较大误差;而在地面进行测量又无法完全克服重力、测试机构本身带来的影响,导致测试结果不会很准确。另外,航天器在轨运行中,由于燃料的长期消耗、航天器构型的变化(如太阳帆板或某些大型天线的展开、宇航员的移动等)、载荷的在轨捕获与释放(如从宇宙飞船或航天飞机上释放卫星;空间机器人捕获目标、清除轨道垃圾等)、与其他航天器的对接等,都会导致系统的质量特性发生变化,很多变化是剧烈的(如捕获、释放卫星,与目标的对接等)、无法确知的(如对非合作目标的操作、轨道垃圾的清理等)。若不解决质量特性的在轨辨识问题,控制系统很难保证在其发生大的变化后还能进行精确的轨道、姿态控制。目前为止,已有一些学者针对航天器质量特性在轨辨识问题开展了研究。Bergmann等提出了一种使用高斯二阶滤波的辨识方法,该方法很像扩展卡尔曼滤波器,但算法更复杂、计算量也更大;而且辨识模型忽略了ω×Iω项(ω×为ω的叉乘操作数),因此要求航天器运动速度足够慢,否则辨识精度满足不了要求。Tanygin和Williams提出了一种针对自旋卫星的最小二乘辨识算法。Wilson和Rock提出了基于指数加权递归最小二乘算法,利用速率陀螺的测量数据,可辨识系统惯量和质心位置,以及每个推力器的特性参数,然而无法辨识航天器的质量。后对该方法进行了改进,并增加了加速度计以提供航天器线加速度的测量,由此辨识航天器的质量、惯量、质心位置。在该辨识方法中,分别对质量、惯量、质心位置列近似线性方程,人为切割了参数内部的耦合关系。虽然采用“递归”的方法可逐步减少去耦合后导致的误差,但当航天器运动速度较快时,该误差仍然不小。Kim等利用递归最小二乘法辨识三轴航天器仿真器的惯量、质心位置,并实现仿真器重力平衡的自适应控制。其存在的问题与Wilson等的方法相同。国内在该方面的研究则很少,文献讨论了单轴气浮台转动惯量的测试,唐江河等研究了惯导平台参数的辨识方法,王书廷等则利用类似于文献的递归最小二乘法辨识航天器的惯量和质心位置。耿立辉等采用多变量两阶段辨识方法对卫星模拟装置的姿态动力学模型进行辨识。这些方法都没有很好解决航天器参数辨识中的两个主要问题:(1)能同时辨识质量、惯量、质心位置;(2)辨识中需要将惯量、质心位置等的耦合考虑进去,以提高辨识精度。另外,刘宇等提出了基于角动量守恒的动力学参数辨识方法,但是假定各刚体质心位置已知。针对如上问题,本文提出了两种航天器质量特性在轨辨识方法,能辨识航天器所有质量特性参数;而且所采用的方法,不需要假设航天器运动足够慢,惯量与质量、质心辨识的耦合关系直接考虑包含其中,克服了以往方法需要将其耦合关系强行分离的缺点。1坐标系定义和动力学建模1.1坐标系的关系不失一般性,假设所研究的航天器系统以4个飞轮(三正交+一斜装)和12个推力器作为执行机构,并装配了三轴速率陀螺和加速度计,分别提供航天器角速度和线加速度的测量(其它的常规测量设备如星敏感器、太阳、地球敏感器等不在此专门指出)。航天器采用整星零动量的控制方式。为方便讨论,定义如下坐标系:(1)轨道坐标系OoXoYoZo:原点位于系统质心,OoZo轴由系统质心指向地心,OoXo轴位于轨道平面且垂直于OoZo轴,沿飞行方向,OoYo由右手定则确定;(2)本体坐标系ObXbYbZb:原点位于星箭对接环中心,对地模式下ObXb轴指向飞行方向,ObZb轴指向地心,ObYb由右手定则确定。本体坐标系作为卫星上的几何参考基准,其它星体固连坐标系(如部件的安装坐标系、测量坐标系等)均以此为参考;(3)质心坐标系OcmXcmYcmZcm:原点位于系统质心(在本体系中的位置用矢量rcm表示),各轴指向与本体坐标系相同;(4)推力器i(i=1,…,12)坐标系OJiXJiYJiZJi:原点位于推力器i的作用点(在本体系中的位置用矢量ri表示),OJiXJi轴为推力器推力矢量方向;(5)加速度计测量坐标系OsXsYsZs:为加速度计测量的参考坐标系(原点在本体系中的位置用矢量rs表示)。各坐标系之间的关系如图1所示。卫星在发射前,地面将通过各种手段确定各个坐标系的原点和指向。除质心坐标系外,其它坐标系都很容易确定,且在卫星发射入轨工作后也不会发生变化。而质心坐标系的原点与质量分布有关,入轨后由于星体构型的变化、燃料的长期消耗等,质心位置会发生变化。1.2卫星姿态动力学方程卫星相对于其质心的角动量为Ηb=Jbω(1)其中,ω为卫星旋转角速度,Jb为卫星(含飞轮)相对于其质心的惯量矩阵,具有如下形式:Jb=[J11J12J13J12J22J23J13J23J33](2)则卫星姿态动力学方程表示为:Jb˙ω+ω×Jbω=τ(3)其中,ω×为ω的叉乘操作数,即ω×=[0-ωzωyωz0-ωx-ωyωx0](4)τ为系统质心的主矩。当反作用飞轮运动时,需考虑飞轮角动量,此时卫星系统的角动量为Ηsum=Ηb+Ηw(5)上式中,Hw为飞轮相对于本体的角动量,对于3正交+1斜装(与各轴等倾角)的配置,飞轮组的角动量为Ηw=[JwxΩx+1/√3JwsΩsJwyΩy+1/√3JwsΩsJwzΩz+1/√3JwsΩs](6)其中Jwx,Jwy,Jwz,Jws分别为x、y、z及斜装飞轮的转动惯量,Ωx,Ωy,Ωz,Ωs分别为x、y、z及斜装飞轮的旋转角速度。2识别仪器质量特征的方法2.1参数解算的最小二乘法识别方法2.1.1不同转速时转动惯量的辨识由于航天器采用整星零动量的方式,因此,使用飞轮进行姿态机动的过程中,系统的角动量近似为0。根据式(5)有Ηsum=Ηb+Ηw=0(7)将(1)代入上式,并化简得Jbω=-Ηw(8)其中Jbω=[J11J12J13J12J22J23J13J23J33][ωxωyωz]=[J11ωx+J12ωy+J13ωzJ12ωx+J22ωy+J23ωzJ13ωx+J23ωy+J33ωz]=[ωx00ωyωz00ωy0ωx0ωz00ωz0ωxωy][J11J22J33J12J13J23](9)令A=[ωx00ωyωz00ωy0ωx0ωz00ωz0ωxωy](10)x=[J11,J22,J33,J12,J13,J23]Τ(11)则根据(8)、(9),有Ax=-Ηw(12)式(12)为关于x的线性方程组,有三个方程六个未知数。对于给定采样时刻ti,星体旋转角速度ω和飞轮角动量均可测出。通过多个时刻的采样,可得到多组数据,利用最小二乘法即可辨识转动惯量。对于ti时刻,方程(12)表示为Aix=-Ηwi(13)假设有N1组数据,1≤i≤N1,则可组成如下的方程组:[A1A2⋮AΝ1]x=[-Ηw1-Ηw2⋮-ΗwΝ1](14)将式(14)表示成如下形式Κx=b(15)其中,Κ=[A1A2⋮AΝ1]∈R3Ν1×6,b=[-Ηw1-Ηw2⋮-ΗwΝ1]∈R3Ν1(16)则方程(15)的最小二乘解为:x=(ΚΤΚ)-1ΚΤb(17)由此可见,利用飞轮实现姿态机动,获取机动过程中基座姿态角速度(在本体坐标系中的分量)、飞轮角速度数据后,即可采用最小二乘法辨识出系统的转动惯量。上述辨识过程不受航天器质量、质心位置辨识结果的影响(虽然是基于系统对质心的角动量守恒方程辨识系统相对于质心的惯性矩阵,但辨识过程中不依赖于具体的质心位置),故是“解耦”的。2.1.2加速度计量值上述方法利用角动量守恒方程已经辨识了航天器的转动惯量,而要辨识出质量、质心位置,则必须通过喷气对航天器进行控制,根据作用效果,辨识模型参数。在喷气对航天器进行控制过程中,飞轮不参与控制,其转速为0。如上所述,航天器上安装有加速度计,其安装在rs的位置。加速度计可测出rs处的加速度˙vs。由刚体动力学可知,rs处的线速度与质心处的线速度有如下关系:vs=vcm+ω×(rs-rcm)(18)其中vs为加速度计测量坐标系原点的线速度,vcm为航天器质心的线速度,rcm为质心位置矢量(在本体中的表示)。对式(18)两边进行微分,有˙vs=˙vcm+˙ω×(rs-rcm)+ω×(ω×(rs-rcm))(19)而˙vcm=n∑i=1fi/m(20)式中,fi为第i个喷气的作用力(矢量),m为航天器质量,n为参与控制的推力器的个数(1≤n≤12)。从(20)可知,若系统合力为0,即n∑i=1fi=0,质量m将不出现在方程(19)中,则无法辨识航天器的质量。因此,在所设计的机动中,需要采用非平衡作用力,该结论与文献给出的一致。对式(19)进行进一步的处理,有˙vs=n∑i=1fi⋅m-1+˙ω×(rs-rcm)+ω×ω×(rs-rcm)=[n∑i=1fi,(˙ω×+ω×ω×)][m-1rs-rcm](21)将式(21)写成如下形式:By=z(22)上式中,B=[n∑i=1fi,(˙ω×+ω×ω×)]∈R3×4(23)y=[m-1rs-rcm]∈R4(24)z=˙vs∈R3(25)方程(22)的系数矩阵B可根据控制器的输出、角速度的测量结果来确定,z直接通过加速度计测出。因此,类似于2.1.1节的方法,可根据机动过程中的历史数据,采用最小二乘法求解y,从而辨识出航天器的质量m、质心位置rcm。2.2在pso下，po探测器的质量特征识别方法2.2.1从总体性能指标上辨识上节中,分别利用飞轮、喷气的控制作用,获取航天器运动数据,由此辨识出航天器的惯量,质量、质心位置。上述方法的好处是转动惯量与质量、质心位置的辨识真正解耦,对于只关心转动惯量而不关心质量、质心位置的情况,如航天器只需要姿态控制,而不需轨道控制的情况,仅采用2.1.1的方法即可满足要求。若需辨识所有质量特性参数,则要分别利用飞轮、喷气对航天器进行控制,分两步实现辨识(转动惯量,质量、质心位置)。本节将研究仅利用推力器的控制作用实现对所有质量特性参数(包括转动惯量、质量、质心位置)的整体辨识。与以往方法不同的是,没有人为切割参数辨识中的相互耦合关系,也无需附加任何假设(如要求航天器运动足够慢),而是将其整体考虑,将参数辨识问题转换为非线性系统优化问题,从总体性能指标上辨识所有参数。根据星体的姿态动力学方程(3),当采用推力器对星体控制时,有如下关系:Jb˙ω+ω×Jbω=n∑i=1(ri-rcm)×fi(26)其中,ri为推力器i作用点的位置矢量。另一方面,根据式(19),有:m(˙vs-(˙ω×+ω×ω×)(rs-rcm))=n∑i=1fi(27)式(26)和(27)可写为:Ι˙ω+ω×Ιω-n∑i=1(ri-rcm)×fi=0(28)m(˙vs-(˙ω×+ω×ω×)(rs-rcm))-n∑i=1fi=0(29)由上面两式知,若合力矩n∑i=1(ri-rcm)×fi=0,则利用喷气实现的机动无法辨识航天器的惯量矩阵I;若仅有平衡力作用,即n∑i=1fi=0,则根据式(29),无法辨识质量。所以在设计航天器的机动时,既要有姿态机动,也要有轨道机动。假设辨识出的参数分别记为ˆm、ˆΙ和ˆrcm,其中ˆΙ=[ˆΙ11ˆΙ12ˆΙ13ˆΙ12ˆΙ22ˆΙ23ˆΙ13ˆΙ23ˆΙ33],ˆrcm=[ˆrxˆryˆrz](30)根据(28)和(29),要求ˆΙ˙ω+ω׈Ιω-n∑i=1(ri-ˆrcm)×fi→Ο3×1(31)ˆm(˙vs-(˙ω×+ω×ω×)(rs-ˆrcm))-n∑i=1fi→Ο3×1(32)其中,O3×1为3×1的零矢量。对于第k个采样时刻,式(31)、(32)左边分别记为:gk(Ι^,r^cm)=Ι^ω˙k+ωk×Ι^ωk-∑i=1n(ri-r^cm)×fik=[gxgygz]Τ(33)hk(m^,r^cm)=m^(v˙sk-(ω˙k×+ωk×ωk×)(rs-r^cm))-∑i=1nfik=[hxhyhz]Τ(34)令x=[m^,Ι^11,Ι^22,Ι^33,Ι^12,Ι^13,Ι^23,r^x,r^y,r^z]Τ(35)则式(31)、(32)等价于Jk(x)=gx2+gy2+gz2+hx2+hy2+hz2→0(36)若有N2组数据,1≤k≤N2,则:J(x)=∑k=1Ν2Jk(x)→0(37)经过上面的推导,质量特性的辨识问题可转化为非线性系统(31)、(32)的优化问题,待优化的指标函数为式(37)。2.2.2pso算法描述微粒群算法(ParticleSwarmOptimization,PSO)是由Kennedy和Eberhart等于1995年开发的一种进化计算技术,来源于对简化社会模型的模拟。PSO算法同遗传算法类似,是一种基于群体搜索的方法,群体成员信息共享,依据“确定性”和“可能性”准则加快搜索速度。两者都随机初始化种群,而且都使用适应度值来评价系统。但PSO算法不像遗传算法那样对个体使用演化算子,如交叉(crossover)、变异(mutation)等,而是在解空间追随最优的微粒进行搜索,所以PSO算法结构简单,计算效率更高。PSO算法相比GA算法在非线性连续函数优化的问题上有很高的计算效率,因此本文采用PSO算法来优化式(37)。PSO算法将每个个体看作D维搜索空间中的一个没有体积的微粒,在搜索空间中以一定的速度飞行。这个速度根据它本身的飞行经验以及同伴的飞行经验进行动态调整。第i个微粒表示为Xi=(xi1,xi2,…,xiD),它经历过的最好位置(有最好的适应值)记为Pi=(pi1,pi2,…,piD),也称为个体极值pbest。在群体所有微粒经历过的索引号用符号g表示,即Pg,也称为群体极值gbest。微粒i的速度用Vi=(vi1,vi2,…,viD)表示。在找到这两个极值时,其第d维(1≤d≤D)根据如下方程来更新自己的速度和新的位置:vid=wvid+c1rand()(pid-xid)+c2rand()(pgd-xid)(38)xid=xid+vid(39)其中,w为惯性权重(inertiaweight),使微粒保持运动惯性,使其有扩展搜索空间的趋势,有能力搜索新的区域。c1和c2为加速度常数,分别代表将每个微粒pbest和gbest位置的统计加速项的权重。rand()和Rand()为两个在[0,1]范围内变化的随机函数。此外,微粒的速度Vi被一个最大速度Vmax所限制。如果某一维更新后的速度超过设定的Vmax,那么这一维速度被限定为Vmax。PSO算法的一个很有吸引力的地方是只需要调节少量的参数。函数(37)的最优化问题通过PSO算法来解决。首先,第i个微粒定义为:Xi=xi=(m^i,Ι^11i,Ι^22i,Ι^33i,Ι^12i,Ι^13i,Ι^23i,r^xi,r^yi,r^zi)(40)粒子数为np,且最大迭代次数为Nmax,迭代次数用k表示。那么,算法的实现流程如下:(1)令k=1,并初始化一群微粒,即{x1(k),x2(k),⋯,xnp(k)},1≤k≤Νmax(41)(2)通过式(37)计算每个微粒的适应度函数Fi,即Fi(k)=J(xi(k)),i=1,2,⋯,np(42)(3)计算每个微粒的最优适应度和最优位置Fi-best={Fi(k),ifk=1Fi(k),elseifFi(k)＜Fi-bestFi-bestelse(43)Ρi={ai(k),ifk=1ai(k),elseifFi(k)＜Fi-bestΡielse(44)(4)计算全局最优的适应度Fg-best和全局最优位置Pg,即Fg-best={min1≤i≤np(Fi(k)),ifk=1min1≤i≤np(Fi(k)),elseifmin1≤i≤np(Fi(k))＜Fg-bestFg-best,else(45)Ρg={ai-min(k),ifk=1ai-min(k),elseifmin1≤i≤np(Fi(k))＜Fg-bestΡg,else(46)(5)根据式(38)和(39)分别改变速度和位置;(6)如达到结束条件(如指标函数小于给定阈值,或迭代次数超过Nmax),则算法结束;否则,返回步骤(2)。算法收敛到最优值时的变量x*即是辨识的结果。3模拟研究3.1航天器的结构布局为对所提出的方法进行仿真验证,在Matlab/Simulink中建立了航天器的动力学模型,理论质量、惯量、质心位置如表2、表3所示;飞轮最大角动量15Nms;最大转速5000rpm;惯量2.86×10-2kg·m2。航天器装配了12个推力器,可实现6DOF(自由度)的控制,其布局如图2所示,图中的箭头方向为喷气矢量方向,产生的推力方向与此相反。各推力器的幅值均为2N,安装位置及推力矢量方向如表1所示。3.2参数解算器最小平方法的识别方法3.2.1转动惯量模型的辨识利用飞轮实现航天器从当前姿态机动到指定姿态,按一定采样周期记录机动过程中的基座姿态角速度、飞轮角速度,然后根据2.1.1节的方法辨识航天器的转动惯量。航天器的初始姿态和期望的终止姿态分别为:Ψ0=Τ(47)Ψd=[6°-7°8°]Τ(48)姿态机动过程如图3所示,机动过程中航天器本体角速度(在本体系中的表示)和飞轮角速度(斜装飞轮未使用)如图4所示。基于动量守恒辨识的转动惯量及其辨识误差如表2所示,从中可见,辨识结果是非常理想的。需要指出的是,由于所给出的结果仅保留了4位小数点,因此表中的“理论值”与“辨识值”看起来似乎一样(下同)。3.2.2推力器使用的顺位如上所述,要辨识航天器的质量,施加的作用力合力不能为0。因此,本节设计了每次仅开启一个推力器的机动方式(飞轮组关闭),推力器使用的次序依次为f1→f7→f2→f8→f3→f11→f4→f12→f5→f9→f7→f10,每个推力器的作用时间均为t=5s,总的机动时间60s。按一定采样周期记录机动过程中的基座姿态角速度、飞轮角速度,然后根据2.1.2节的方法辨识航天器的质量、质心位置,辨识的结果如表3所示,