隧道掘进机(tunnel boring machine，TBM)是一种集机、电、液、传感、信息等技术于一体的大型现代化隧道施工装备．双护盾隧道掘进机因其岩层适应性强、掘进速度快、转弯半径小等优势而得到广泛应用[1-2]．在双护盾隧道掘进机的掘进过程中，为了计算出刀盘切口中心和盾尾中心相对设计轴线的偏差量，从而进行纠偏控制并确保隧道最终精确贯通[3]，须要准确测量前后各盾体的位姿．目前工程中广泛应用的盾体位姿测量方法主要包括棱镜法、激光标靶法等，均须配合全站仪的自动照准激光进行测量[4-5]．虽然双护盾隧道掘进机后盾体的位姿通过这些方法可以测量得到，但对于前盾体的位姿，由于前后盾铰接区域复杂的装配关系、空间狭小等原因，直接利用棱镜法或激光标靶法进行测量容易导致全站仪的照准激光被后盾体和管路等设备遮挡，从而导致测量失效[5-6]．因此，为了准确测量出双护盾隧道掘进机前盾体的位姿，如何在测量得到后盾体位姿的基础上精确测量出前后盾体间的相对位姿成为导向系统研究中的关键问题．关于双护盾隧道掘进机盾体间的相对位姿测量问题，已有研究主要包括基于单目相机的视觉测量方法和基于并联机构的位姿求解方法．例如，德国VMT公司采用在后盾上安装工业相机拍摄前盾上安装的二维码的方法进行盾体间相对位姿测量；Lin等[7]提出在后盾前端安装相机，拍摄前盾上多个发光二极管(LED)发光点的方法进行测量解算；黄喆等[8]在前盾体上设计并安装了一个异形主动发光标靶，利用后盾工业相机拍摄解算出相对位姿．这些方法均利用单目视觉技术结合相机投影算法求解出相机和拍摄的前盾参照物之间的位姿关系，结合初始标定参数，得到前后盾体之间的相对位姿．然而在复杂恶劣的现场工况环境下，隧道掘进机掘进引起的粉尘、雾气极易污染遮挡相机镜头及二维码、发光特征点等参照物，并干扰相机的通视空间，导致该类方法测量失效．此外，张强等[9]也提出类似并联机构位姿正解的方法，根据盾体间各个液压杆长约束条件建立非线性方程组，通过牛顿法迭代求解出相对位姿，然而数值法迭代求解依赖于初值的选取，容易不收敛，并且该方法只能用于铰接油缸交叉安装这一类构型的小型盾构机，对于油缸平行安装等其他构型的大中型隧道掘进机则不适用．基于对现有的盾体间相对位姿测量方法的调研和分析，本研究提出了一种基于多传感器信息的盾体间相对位姿测量方法，该方法可适用于任意构型的双护盾隧道掘进机，有助于解决恶劣工况环境下的盾体间相对位姿测量和导向问题．利用盾体间相对位姿反解多传感器信息相对容易的特点，本研究将复杂的盾体间相对位姿正向求解问题转化为一个多元函数的最优化问题．在此基础上，提出了一种动态自适应差分进化算法进行相对位姿解算，从而准确测量出双护盾隧道掘进机盾体间的相对位姿．经过仿真试验分析和试验台的实际测试验证，本文方法能够实现稳定可靠的盾体间相对位姿测量，满足工程应用中双护盾隧道掘进机的导向需求．1 基于多传感器信息的相对位姿正解1.1　机构分析及测量方案设计1.1.1　双护盾TBM机构分析如图1所示，以铁建重工公司某型号双护盾隧道掘进机为例，其组成包括前后两节圆柱形盾体，盾体之间通过多根铰接油缸相连且前后盾体的油缸铰接点在其铰接平面内的分布相同．零位状态下，前后盾体沿轴线对齐且没有相对滚动，各铰接油缸之间相互平行．因为前后盾体与油缸之间的连接为球铰副S，油缸柱塞和缸套之间的连接为移动副P，故双护盾隧道掘进机前后盾间的结构关系可等效为一个6自由度n-SPS并联机构．10.13245/j.hust.241206.F001图1双护盾隧道掘进机结构组成类比并联机构运动学问题，已知油缸长度求解相对位姿为运动学正解．然而双护盾隧道掘进机前后盾体上的铰接点分布相同，其对应的n-SPS并联机构存在着构型奇异问题，仅通过油缸杆长进行正向运动学计算容易因为运动分岔或多解问题导致相对位姿求解失败[10-11]．同时，在掘进过程中，双护盾隧道掘进机前后盾体组成的并联机构多处在奇异位置附近，此时其运动学雅各比矩阵的性质较差，采用数值法解算的相对位姿误差也很大[12-13]．因此，有必要采用安装附加传感器的方法[14-15]，通过前后盾体间多种传感器测量的信息条件进行充分约束，从而避免出现类似并联机构正解的奇异性问题，求解出盾体间相对位姿．1.1.2　多传感器测量方案设计综合考虑传感器在工业现场应用中的可靠性及抗干扰能力，选取在前后盾体铰接区安装四根油缸行程传感器、激光器及感光靶系统以及前后盾各一个倾角传感器的组合测量方案来对相对位姿进行充分约束．其中，油缸行程传感器采用磁致伸缩式，检测各个油缸的行程量Hi(i=1，2，3，4)；激光感光靶系统包括分别安装在前盾和后盾上的感光靶与激光发射器，激光光线在现场环境下能够穿透粉尘雾气照射在感光靶平面上，感光靶通过图像处理技术检测激光照射的光斑中心在其靶平面内的二维坐标（x'，y'）；前后盾上安装的倾角传感器用于测量不同的相对位姿状态下各自倾角之差∆γ．为了便于描述分析盾体间相对位姿和多传感器信息之间的定量对应关系，建立如图2所示的多传感器测量模型．其中，全站仪坐标系{B}在设站时建立，前盾坐标系{F}原点OF为前盾体轴线和铰接平面的交点，OF与前盾铰接点AF1连线为XF轴，盾体轴线为其YF轴，ZF轴由右手规则确定，后盾坐标系{S}根据同样的方式建立．前后盾体的轴线方向可由拟合圆法测出盾体前后端面的中心点确定[16]，铰接点坐标可由全站仪换站到盾体内部测量得到．10.13245/j.hust.241206.F002图2多传感器测量模型同时，在图2所示的测量模型中定义各传感器安装后的特征参数．油缸行程传感器的特征参数定义为对应的油缸铰接点AFi和ASi(i=1，2，3，4)在各自盾体坐标系下的坐标，由设计尺寸可得 FFi=ASSi=rcos θi,0,-rsin θiT，式中r和θi分别为铰接点所在分布圆的半径及分布角度．激光感光靶系统的特征参数分别定义为激光器的发射点E及其方向矢量L，以及感光靶二维平面的坐标原点M，X轴上点P，Y轴上点Q．同时，定义前、后盾倾角传感器的特征参数分别为其敏感轴矢量IF和IS．各传感器的特征参数值均可在盾体间零位状态下测量标定得到．由于盾体间相对位姿分量之间的相互耦合影响，根据多传感器的测量信息结合几何关系直接计算的相对位置和姿态角并不能准确反映前后盾坐标系间的相对位姿变换关系，因此本研究将各传感器信息作为相对位姿求解的约束条件，基于该测量模型中相对位姿和多传感器信息之间的正反解关系，设计合适的算法求解出解耦后的前后盾体间的相对位姿．1.2　相对位姿反解多传感器信息1.2.1　盾体位姿描述及坐标变换根据绕固定坐标轴X-Y-Z旋转的方式对前后盾体坐标系之间的位姿关系进行描述，定义坐标系{F}与坐标系{S}间的相对位姿向量Q=(x，y，z，α，β，γ)的含义为：在坐标系{F}与坐标系{S}初始重合的情况下，{F}先分别绕{S}的X，Y，Z轴转动α，β，γ角度，再分别沿{S}的X，Y，Z轴平移x，y，z长度．则该相对位姿向量代表的旋转变换矩阵FSR及平移向量FSP分别为：   RFS=cos γ-sin γ0sin γcos γ0001⋅ cos β0sin β010-sin β0cos β1000cos α-sin α0sin αcos α;FSP=[x,  y,  z]T.(1)因此，坐标系{F}中的点FP变换到坐标系{S}中即为 SP=RFS∙P F+PFS．(2)坐标系{F}中的方向矢量Fn变换到坐标系{S}中即为 Sn=RFS∙n F．(3)1.2.2　多传感器信息的反解计算假设已知前后盾相对位姿Q，按照式(2)和式(3)结合空间几何计算，可反解出各传感器的测量信息．a．油缸行程传感器测量的位移量Hi设安装有位移传感器的四根油缸的杆长为Li(i=1，2，3，4)，由铰接点坐标计算可得Li=FSR∙A FAFi+PFS-ASASi．(4)结合油缸零行程杆长L0，可得各油缸位移传感器测量的行程量Hi=Li-L0    (i=1,  2,  3,  4)．(5)b．感光靶测量的激光光斑中心坐标(x´，y´)根据激光发射点坐标 SE=[xE,yE,zE]T和激光方向矢量 SL=[xL,yL,zL]T，可得坐标系{S}中激光直线方程x-xExL=y-yEyL=z-zEzL．(6)将感光靶平面原点FM，X轴特征点FP及Y轴特征点FQ按照式(2)进行变换，可得其在坐标系{S}中的对应坐标FM，FP，FQ，由此三点可建立感光靶平面方程Ax+By+Cz+D=0．(7)联立式(6)和式(7)构建方程组，可得光斑点空间坐标SJ．同时易知感光靶二维平面水平方向单位向量 Su=(P S-M S)/  SP-M S和竖直方向单位向量 Sv=(Q S-M S)/  SQ-M S．根据空间投影关系，光斑点在感光靶平面内的二维坐标即为：x'=(J S-M S)u S;y'=(J S-M S)v S.(8)c．前后盾倾角传感器的倾角读数之差∆γ已知后盾姿态可得水平面法矢量SG，结合后盾倾角仪的敏感轴矢量SIS，可得其倾角读数γS=arcsin[(I SIS∙G S)/( SIS∙ SG)]．(9)同理，已知相对位姿可知前盾倾角仪敏感轴矢量 SF=RFS∙IFF，由此可计算前盾倾角仪的读数γF=arcsin[(I SIF∙G S)/( SIF∙ SG)]．(10)于是可反解前后盾倾角传感器的读数之差Δγ=γF-γS．(11)1.3　相对位姿正解优化问题建模若已知前后盾间的相对位姿Q=(x，y，z，α，β，γ)，基于多传感器测量模型，结合相对位姿和各传感器测量信息之间的反解关系式，能够计算出油缸行程Hi(i=1，2，3，4)、激光光斑点平面坐标(x'，y')以及前后盾倾角仪读数差∆γ等传感器信息，故可将对应的反解关系归纳为fi(x,y,z,α,β,γ)=Ii;                 Ii=Hi    (i=1,2,3,4),x'    (i=5),y'    (i=6),Δγ    (i=7), (12)式中：fi为根据相对位姿Q反解第i个传感器信息Ii的非线性映射函数；Ii为第i个传感器测量信息．因此，当正向求解盾体间相对位姿时，将各传感器测量信息Ii作为已知条件，同样可以联立式(12)所示的求解关系式，此时Ii为已知量，Q=(x，y，z，α，β，γ)为待求解量．然而直接求解该相对位姿正解关系复杂困难，基于上述建立的反解对应关系并利用反解过程相对容易的特点，可将该正解问题等价转化为下式所示的一个连续优化问题Q=(x,  y,  z,  α,  β,  γ);min F(Q)=∑i=17fi(Q)-Ii;s.t.  xmin≤x≤xmax ,  ymin≤y≤ymax,zmin≤z≤zmax ,  αmin≤α≤αmax,βmin≤β≤βmax ,  γmin≤γ≤γmax,(13)式中：fiQ-Ii为各传感器信息的反解值与测量值之间的绝对偏差；F(Q)为该相对位姿正解优化问题的目标函数．双护盾隧道掘进机的盾体间相对位姿正向求解问题转化为对于目标函数F(Q)的最小值优化问题，然而因为该类高维目标函数存在多峰多极值的情况，相比于传统数值方法依赖于计算初值的选取、计算效率低、求解易不收敛等问题，智能启发式算法具有启发式搜索、并行求解、全局收敛等优势[17-18]，故本研究采用智能优化类算法求解式(13)所示的高维连续函数最小值优化问题．2 基于差分进化算法测量相对位姿2.1　差分进化算法设计差分进化(differential evolution，DE)算法作为一种智能启发式算法，仿照自然界优胜劣汰和生物进化原则产生群体智能来指导优化搜索，采用浮点数编码和差分变异进化的方式，在多元变量的优化问题上已得到较多应用[19-20]，本研究将其用于解决式(13)所述盾体间相对位姿测量优化问题．为了保证盾体间相对位姿求解的精度和效率，须要综合算法的前期全局寻优能力和后期局部开发能力，提升寻优的收敛速度，为此本研究设计了一种自适应差分及动态更新策略对种群中的个体进行变异进化及迭代更新，算法具体内容如下．2.1.1　初始化种群设置规模为N、维度为D的相对位姿目标种群，在各相对位姿分量的可行域范围内，按照均匀概率对种群随机初始化，xij,0=xjmin+rand(0,1)∙(xjmax-xjmin)，(14)式中：xij，0为第i个相对位姿目标个体在第j维度的分量；xjmin和xjmax分别为xj下限值和上限值；rand(0，1)为(0，1)区间内随机选取的一个实数．2.1.2　自适应差分策略变异算子自适应差分变异操作为    vi,G+1=xr1,G+λ(G)(xbest,G-xr1,G)+F(G)(xr2,G-xr3,G) ,(15)式中：r1，r2，r3∈[1，2，⋯，N]，且r1≠r2≠r3≠i；G为种群代数；xbest，G为第G代种群中的最优相对位姿向量；λ(G)为组合因子函数；F(G)为变异因子函数；xbest，G- xr1,G为选择性差分矢量；xr2,G-xr3,G为随机性差分矢量．选择性差分矢量xbest，G- xr1,G的引入使种群中个体在变异过程中受到相对位姿最优个体xbest，G的引导，引导程度由组合因子λ(G)控制，λ(G)按下式随迭代轮数自适应递增，λ(G)=λmax-(λmax-λmin)exp(-kG)，(16)式中：λmin为组合因子的初值；λmax为组合因子上限值；参数k控制组合因子变化快慢．组合因子的变化趋势如图3所示，在迭代搜索前期λ取较小值，减弱xbest，G对种群中个体更新的指引，保证种群的多样性和全局探索能力，避免过早收敛到局部最优；在迭代中后期λ缓慢增大，使种群有更多的向xbest，G靠拢的趋势，加快算法向相对位姿最优解进行探索和收敛．10.13245/j.hust.241206.F003图3自适应组合因子λ变化曲线同时，随机性差分矢量xr2,G-xr3,G为个体在变异过程中添加的随机扰动量，由缩放因子F(G)控制，F(G)按下式随迭代轮数自适应递减，F(G)=Fmax-(Fmax-Fmin)(G/Gmax)2，(17)式中：Fmax为变异因子的初值；Fmin为变异因子终值；Gmax为算法设置的最大迭代次数．如图4所示，F(G)采用抛物线下降的形式，在相对位姿求解前期保持相对恒定值，保证种群的随机多样性；在搜索中后期F加快减小，起到减小步长的作用，并与逐步增大的λ相配合，加强种群在相对位姿最优解附近的局部开发，提高收敛速度和精度．10.13245/j.hust.241206.F004图4自适应变异因子F变化曲线2.1.3　交叉算子将变异种群中的相对位姿个体按一定概率与目标种群中的对应个体进行二项交叉，得到试验种群uij,G+1=vij,G+1    (rand(0,1)≤Cr或j=jrand);xij,G    (其他),(18)式中：Cr为交叉因子，控制相对位姿向量各维度发生交叉的概率；jrand为[1，2，⋯，D]范围内的随机整数，确保交叉后生成的试验个体至少有一维“基因”由变异个体提供．2.1.4　选择算子分别计算目标种群和试验种群中对应相对位姿个体的适应度值，根据优胜劣汰的原则选取适应度值较优的留下作为下一代目标种群中的个体，xi,G+1=ui,G+1   (fitness(ui,G+1)fitness(xi,G));xi,G   (其他),(19)式中fitness(  )为适应度函数．2.1.5　动态更新种群及最优个体策略传统的以群体为单位的进化更新方式可以保证种群内的多样性，但会导致最优解的进化更新存在一定的滞后性[21]，不利于提高相对位姿求解效率．因此，本研究采取在每轮迭代中对目标种群中的相对位姿个体进行连续的变异、交叉、选择操作的方式，使其进化为对应位置的下一代个体，随后立即更新相对位姿最优个体，以此方式遍历整个种群，更新方式对比如图5所示．这样在本轮迭代更新中，种群个体按式(15)进行的自适应差分操作将基于实时动态更新的选择性差分矢量xbest，G- xr1,G以及随机性差分矢量xr2,G-xr3,G，从而提高最优个体及整体种群更新的速度，加快得到盾体间相对位姿最优解．10.13245/j.hust.241206.F005图5种群及最优个体更新方式对比2.2　适应度函数设计适应度函数作为算法寻优的评价指标，用来衡量算法搜索的相对位姿和真实的相对位姿之间的偏离程度．如式(13)所示，理想情况下，对于一组绝对精确的传感器测量数据Ii，总存在一组相对位姿Q使得各项反解值fi(Q)与Ii完全相等，即目标函数的理论最小值为0，对应的相对位姿最优解Qbest即为相对位姿理论真值．然而实际情况中各传感器存在测量误差噪声并且其分辨力有限，对于一组带有误差噪声的传感器信息Ii，并不能总是找到一组相对位姿Q满足目标函数中各项反解偏差量|fi(Q)-Ii|均等于0，此时为了让最终求解的相对位姿误差尽可能小，须要对各传感器信息的反解偏差量|fi(Q)-Ii|在目标函数中所占的权重进行合理分配．本研究提出一种根据各传感器测量信息Ii的精度σIi，分别赋予目标函数中各项传感器信息偏差量|fi(Q)-Ii|一个与精度σIi大小成反比的影响因子ηi的方法．对于测量精度越高的传感器数据Ii，对应的影响权重因子越大，因此当适应度函数达到最小值时，其对应的|fi(Q)-Ii|相比于其他分量越小，即算法在寻优过程中更趋向于精度较高的传感器信息Ii进行搜索求解，从而提高相对位姿的求解精度．按照下式对各影响因子ηi进行归一化处理，得到算法的适应度函数fitness(Q)，ηi=(1/σIi)/∑i=17(1/σIi)  (i=1, 2, …, 7),∑i=17ηi=1;(20)fitness(Q)=∑i=17ηi∙fi(Q)-Ii，(21)式中σIi按照3σ原则根据对应传感器测量数据的精度进行取值．2.3　相对位姿测量算法流程在实际测量盾体间相对位姿的过程中，对于任意一组多传感器测量数据Ii(i=1，2，⋯，7)，结合测量模型的位姿反解关系fi及各影响因子ηi，均可按式(21)建立对应的适应度函数，从而利用本文差分进化算法快速迭代求解出盾体间相对位姿．算法测量盾体间相对位姿的步骤如下：步骤1 输入各传感器的测量数据Ii，建立适应度函数fitness()；步骤2 设置算法的初始参数及运行终止条件；步骤3 在前后盾相对位姿约束范围内随机初始化相对位姿向量目标种群；步骤4 计算初代目标种群的适应度值，找出初代种群中的相对位姿向量最优个体；步骤5 按本文自适应差分及动态进化策略更新目标种群中的单个体及最优个体；步骤6 根据种群适应度值的变化情况及种群的迭代次数判断是否达到寻优终止条件，是则继续下一步，否则转步骤5；步骤7 停止迭代寻优，输出种群中的相对位姿最优个体，即为该组多传感器测量数据所对应的盾体间相对位姿最优解．3 盾体间相对位姿测量仿真试验3.1　仿真试验数据生成以中国铁建重工公司生产的某型双护盾隧道掘进机为例，建立如图2所示的多传感器测量模型，经过标定取值得到各传感器的特征参数：油缸零行程杆长L0=1 500 mm，铰接点分布圆半径r=2 135 mm，铰接点分布角度θi (θ1=0°，θ2=85°，θ3=180°，θ4=265°)，激光发射点SE=[152.20 mm，933.47 mm，1 957.47 mm]T，激光方向矢量SL=[0.00，0.99，0.02]T，感光靶原点FM=[147.44 mm，-260.42 mm，1 970.99 mm]T，感光靶X轴点FP=[181.40 mm，-263.40 mm，1 939.03 mm]T，感光靶Y轴点FQ=[173.80 mm，-260.66 mm，1 995.70mm]T，前倾角仪敏感轴矢量FIF=[cos 41.81°，0，-sin 41.81°]T，后倾角仪敏感轴矢量SIS=[cos 39.32°，0，-sin 39.32°]T．因为在双护盾隧道掘进机的一个掘进周期内，后盾体利用撑靴支撑岩壁进行固定，前盾体通过推进油缸支撑后盾体向前掘进，所以仿真试验中将后盾体设为固定并已知水平面法矢量SG=[0.66，0.00，0.75]T．根据施工现场的某一典型工况数据，设定前盾体相对后盾体从相对位姿Q=(0 mm，1 750 mm，0 mm；0 mrad，0 mrad，0 mrad)的零位状态掘进到相对位姿Q=(100 mm，1 850 mm，100 mm；50 mrad，50 mrad，50 mrad)的状态，按均匀间隔选取该过程中的50组盾体间相对位置和姿态作为理论真实值，根据1.2节将反解得到的50组对应的多传感器信息作为相对位姿测量仿真试验的已知条件．因为实际测量中各传感器存在相应的精度和分辨力，所以根据各传感器的出厂参数对生成的各传感器信息Ii添加对应的高斯随机噪声，同时根据各传感器的分辨力将对应的测量数据保存相应的小数位数，其中油缸行程数据精度为1.0 mm，保留整数；光斑点二维坐标精度为0.5 mm，保留两位小数；前后盾倾角的精度为0.01°，保留三位小数，以此生成50组包含相对位姿理论真实值及对应的多传感器测量数据的仿真测试实例．3.2　仿真试验对比分析为了比较本文算法相比于其他智能启发式进化算法在解决盾体间相对位姿测量连续优化问题上的适用性和效果，引入文献[22-24]中多个变种改进的差分进化算法进行仿真试验对比测试．各算法的参数设定如下：a．文献[22]中的DE算法，变异率F=0.8，交叉率Cr=0.9，变异策略为DE/best/1，记为DE1；b．文献[23]中SDE算法，采用在[0.40，0.95]范围内呈现“S”型变化的变异率，交叉率Cr=0.9，变异策略为DE/rand/1，记为DE2；c．文献[24]中AD-DE算法，采用自适应变异率，双策略差分进化，交叉率Cr=0.9，记为DE3；d．本文差分进化算法，采用组合因子λmin=0.2，λmax=1，k=0.007，变异因子Fmax=0.5，Fmin=0.1，交叉因子Cr=0.9，记为DADE．算法的种群规模N均设置为50，个体维度D=6，终止条件设置为：若运行达到最大迭代次数Gmax=1 000，或者种群平均适应度值与最优个体适应度值之差小于ε=1×10-5(说明种群已几乎收敛)，则停止寻优并输出最优解．仿真程序在Python环境下编写和运行，计算机平台为Intel(R) Core(TM) i7-10875H CPU，主频为2.30 GHz，内存为16 GiB．图6为利用上述几种改进的差分进化算法及本文算法进行仿真试验时种群的最优适应度值及平均10.13245/j.hust.241206.F006图6算法适应度值变化对比适应度值的变化曲线对比．由图6可知：随着迭代次数的增加，DE1算法的最优及平均适应度值减小较为缓慢，适应度曲线的下降速率小于DE2、DE3及本文算法，并且有连续几十轮的进化代数中最优适应度值没有降低，说明DE1算法对于该优化问题不易跳出局部最优，容易出现搜索停滞的现象；DE2和DE3算法随着进化代数的增加，其最优适应度和平均适应度能够连续减小，最终达到收敛，但适应度曲线随迭代次数的下降速率明显低于本文算法；相比于另外三种算法，本文DADE算法的最优适应度及平均适应度下降速率更快，当种群进化代数达到65次左右时，本文算法已基本收敛，最优个体的适应度值小于1×10-5，说明本文算法的全局寻优能力较好，能够更快跳出局部最优并收敛．表1为上述算法和本文算法运行该50组仿真实例的平均迭代次数和最大迭代次数以及相应的计算时间对比，可知在同样的种群规模及运算条件下，DE1算法须平均迭代进化1 000次，所有算例均达到设置的进化代数上限值，平均计算时间长达15.2 s；DE2算法须平均迭代129次，平均计算时间为2.1 s；DE3算法须平均迭代153次，平均计算时间为2.6 s；而本文算法的平均迭代次数仅为65次，平均计算时间减小到1.1 s．相比于DE2算法和DE3算法，本文算法的平均求解效率分别提高约91%和136%．10.13245/j.hust.241206.T001表1算法迭代次数及计算时间对比算法平均迭代次数最大迭代次数平均时间/s最大时间/sDE11 0001 00015.216.1DE21291512.13.1DE31532002.63.8DADE65891.11.4表2对比了上述算法和本文算法运行该50组仿真实例的相对位姿测量误差的分布区间，以及最大绝对值和平均绝对值，可知：DE1算法求解该优化问题时存在较多陷入局部最优的情况，误差分布区间相比其他算法较大，相对位置和相对姿态误差的最大绝对值也分别达到了32.33 mm和10.96 mrad，超出了工程测量的规范要求[3]；而使用DE2算法、DE3算法及本文算法求解该50组仿真实例，则均能够较好地收敛于全局最优解附近，各相对位姿分量的平均绝对误差相近，且均能实现相对位置误差小于1 mm，相对姿态误差小于1 mrad，具有较高的相对位姿求解精度．10.13245/j.hust.241206.T002表2仿真试验各算法测量误差对比相对位姿误差误差指标DE1算法DE2算法DE3算法DADE算法Δx/mm分布区间[-32.33，13.05][-0.70，0.72][-0.71，0.70][-0.70，0.70]平均绝对值2.490.210.210.21最大绝对值32.330.720.710.70Δy/mm分布区间[-19.28，4.47][-0.51，0.55][-0.51，0.84][-0.51，0.51]平均绝对值1.370.210.200.20最大绝对值19.280.550.840.51Δz/mm分布区间[-7.34，7.89][-0.40，0.36][-0.40，0.40][-0.40，0.38]平均绝对值1.40.160.160.16最大绝对值7.890.40.40.4Δθx /mrad分布区间[-5.17，10.96][-0.42，0.35][-0.68，0.45][-0.30，0.37]平均绝对值1.130.130.160.12最大绝对值10.960.420.680.37Δθy /mrad分布区间[-6.50，4.28][-0.21，0.20][-0.21，0.20][-0.21，0.20]平均绝对值0.740.090.090.08最大绝对值6.50.210.210.21Δθz /mrad分布区间[-3.15，3.05][-0.34，0.34][-0.34，0.35][-0.32，0.34]平均绝对值0.640.130.120.12最大绝对值3.150.340.350.34综合分析仿真试验对比结果可知：当解决式(13)对应的相对位姿测量连续优化问题时，使用DE1算法容易陷入局部最优、迭代时间长、求解误差大，相比之下本文DADE算法具有更好的全局寻优性能，收敛精度高，所求相对位姿误差更小；同时，相比于DE2算法和DE3算法，在达到相同的收敛精度及求解误差近似的情况下，本文DADE算法具有更高的收敛速度和求解效率，能够大幅缩短相对位姿求解的时间周期．因此，本文DADE算法对于解决盾体间相对位姿测量高维空间连续优化问题具有更好的适用性和求解效果．3.3　算法精度分析及鲁棒性验证为了进一步测试本文算法在更多不同工况下对于不同传感器输入数据的相对位姿求解效果，从而分析验证其精度和鲁棒性，采用蒙特卡洛方法[25]在3.1节所述双护盾隧道掘进机的前后盾相对位姿范围内，随机抽样生成1 000组代表不同掘进工况的盾体间相对位姿，其中前盾相对后盾沿轴向的位置变化范围为1 500～2 000 mm，沿径向的位置变化范围为-200～200 mm，三个方向的姿态角变化范围为-90～90 mrad．按照1.2节反解生成各传感器测量数据，经添加噪声及保留分辨力对应的小数位数等预处理，得到对应的多传感器仿真测量数据，从而对本文算法进行求解精度及鲁棒性测试．算法的参数设置和终止条件与3.2节相同，最终得到盾体间相对位姿求解误差的统计结果如表4所示．根据表3可知：当用本文算法对不同工况下的盾体间相对位姿进行求解时，各分量误差几乎均符合均值为0的正态分布，未出现明显偏态的情况，根据大数定律及误差分布3σ准则[26]，可以认为本文算法的测量精度满足盾体间相对位置测量误差在±0.8 mm内，相对姿态测量误差在±0.5 mrad内．同时，当用本文DADE算法对约束范围内随机抽样选取的传感器数据进行求解时，均能够收敛于全局最优附近，没有出现相对位姿测量误差异常过大的情况．故本文算法对于不同工况下盾体间相对位姿测量优化问题的求解精度和鲁棒性满足要求．10.13245/j.hust.241206.T003表3算法精度及鲁棒性测试误差统计结果统计指标Δx/mmΔy/mmΔz/mmΔθx /mradΔθy /mradΔθz /mrad平均绝对值0.190.190.140.130.070.12最大绝对值0.870.760.590.630.290.77均值0.02-0.01-0.01-0.01-0.010.003倍标准差0.710.720.520.500.240.4610.13245/j.hust.241206.T004表4试验台应用测试相对位姿误差结果统计指标Δx/mmΔy/mmΔz/mmΔθx /mradΔθy /mradΔθz /mrad最大绝对值2.624.163.560.921.670.95平均绝对值1.680.952.140.341.010.54均值0.10-0.070.980.17-0.580.244 多盾体试验台应用测试试验为了进一步验证本文相对位姿测量方法在多盾体掘进机上的实际应用效果和可靠性，在铁建重工公司的多盾体掘进机试验台上进行了实物测试．试验台前后盾体的外径尺寸为6 m，盾体轴线长度约为3 m，零位状态下前后盾体间的轴距为1 500 mm，可实现最小50 m的转弯半径；试验中后盾体被焊接固定在地面支撑基座上，模拟实际掘进过程中后盾通过撑靴支撑岩壁固定的情况；通过液压系统调整控制前后盾铰接区域的油缸长度，可模拟掘进过程中前后盾间的相对运动．4.1　测试试验设计为了定量分析本文方法在多盾体试验台实际测试中的精度，须要知道准确的盾体间相对位姿作为基准值，因此设计多盾体试验台相对位姿测试试验流程如图7所示．其中三棱镜法[27-28]作为一种传统的刚体位姿测量方法，具有系统简单、成本低、可靠性高等优点．通过测量全站仪坐标系下三个棱镜中心点的坐标，结合刚体坐标系中该三点的坐标信息，即可准确计算出刚体的位姿．故在测试中每调整一次盾体间相对位姿并按照本文方法进行解算后，同时采用三棱镜法测量出对应状态下前后盾体间的相对位姿以作为基准值，进而对本文方法测量盾体间相对位姿的效果进行分析．10.13245/j.hust.241206.F007图7测试试验内容及流程4.2　测试结果及分析多盾体试验台测试过程中调整了19组不同的前后盾间相对位姿状态，按照上一节设计的试验流程进行测量计算和对比分析，汇总得到本文方法的盾体间相对位姿测量误差如表4所示．可知通过本文DADE算法结合试验台上各传感器读取的数据进行解算，得到盾体间相对位置测量误差的最大绝对值为4.16 mm，平均绝对值小于2.5 mm，均值小于1 mm；相对姿态测量误差的最大绝对值为1.67 mrad，平均绝对值小于1.5 mrad，均值小于1 mrad．故通过建立多传感器测量模型，盾体间相对位姿测量求解获得了充足的约束条件，避免了并联机构正解奇异性导致的求解误差过大等问题．虽然试验台测试时受到传感器特征参数标定误差及试验台铰接装配间隙等因素的影响，导致本文方法实际的测量误差比仿真分析的测量误差稍大，但测试结果表明，考虑这些因素后，本文方法最终能够实现盾体间相对位置测量误差小于5 mm，相对姿态测量误差小于2 mrad，各项误差均满足工程测量规范要求[3]．在试验台测试中，本文方法测量盾体间相对位姿平均所需计算时间为1.29 s，最大计算时间为2.36 s，也符合工程中的测量周期需求．可知本文基于差分进化算法的多传感器联合测量方法能够在实际应用中实现稳定可靠的盾体间相对位姿测量，满足工程应用中双护盾隧道掘进机的位姿测量和导向需要．5 结论本研究对双护盾隧道掘进机的盾体间相对位姿测量问题进行了分析，设计了测量方案并建立了多传感器测量模型，借鉴运动学反解容易的特点将复杂的相对位姿正解问题转化为一个多元连续优化问题，并提出了一种动态自适应差分进化算法进行求解．通过开展仿真试验及铁建重工公司的多盾体试验台测试试验，验证了本文盾体间相对位姿测量方法的精度和实际应用效果．得到主要结论如下：a．仿真试验表明，在相同条件下，本研究设计的差分进化算法在解决盾体间相对位姿测量连续优化问题时，相比其他同类型的改进或变种算法求解误差更小，效率更高，具有更好的适用性，同时也仿真验证了本文算法用于求解不同工况下盾体间相对位姿的精度和鲁棒性．b．铁建重工公司多盾体掘进机试验台的应用测试结果表明，本文方法可以测量出前后盾体间的相对位姿，实现相对位置测量误差小于5 mm，相对姿态测量误差小于2 mrad，并且本文算法测量求解所需的计算时间少于3 s，实际测量精度和计算周期均满足工程应用需求．c．本文测量方法为解决工程中油缸平行安装的大中型双护盾隧道掘进机的盾体间相对位姿测量问题提供了一种思路，在解决具有多体结构的不同连接构型的大型掘进装备的相对位姿测量问题上也具有较好的通用性．d．多盾体试验台应用测试中，受多传感器测量模型特征参数标定误差及铰接间隙等因素的影响，本文方法实际的测量误差比仿真试验中的误差稍大．后续研究可以对这些因素的影响大小进行具体分析，并考虑用误差补偿的方法来进一步提高盾体间相对位姿的实际测量精度．