基于四元数的桥梁短线节段误差分析方法

李光泉

(中铁十四局集团有限公司 山东济南 250101)

1 引言

短线法预制拼装技术，目前在中等跨径预应力混凝土箱梁施工中得到了广泛应用，其关键技术在于线形控制和姿态调整相关的几何计算方法。

短线法采用分段预制拼装的方法，利用已建造的节段与相邻将要浇筑节段匹配。在施工过程中，不仅会有测量仪器产生的测量误差、匹配段定位误差等，还会有施工中振捣、混凝土收缩等原因产生的预制误差，包括:梁长误差、转角误差、扭转误差、错台误差等。这些误差会导致已浇筑的节段偏离设计线形，并且由于已浇筑节段无法更改，从而使得误差不断积累。此时需要识别出已有的误差，在已浇筑节段作为匹配段时，通过调整其位置和姿态参数，使后续节段修正该误差。识别和修正预制误差的几何计算，存在多种不同的方法。

王侃、李国平[1]以阶段接缝和横坡为基准建立坐标系，采用6个控制点形成水平控制线和高程控制线，来分别识别和修正误差。侍刚[2]等建立三维整体坐标系和节段局部坐标系，通过二者间的坐标变换研究其几何关系，并引入非线性最小二乘法，根据测量点坐标结果来建立方程组并求解。刘海东[3]等研究了在三维坐标系中考虑所有误差的误差修正方法。周凌宇[4]等对三维空间整体和局部坐标系间的空间位置不重合进行了分析。時学军[5]研究了适合于任意旋转角度的坐标变换方法，并通过非线性最小二乘法来识别误差。

从以上研究成果可见，误差分析和修正的方法是通过三维坐标变换，把整体坐标系和节段局部坐标系进行匹配，对节段测点的浇筑阶段和匹配阶段测量结果进行分析比对，通过非线性最小二乘法来建立方程组计算出已产生的预制误差[6]，然后通过整体坐标系的拼装坐标来确定误差调整方法[7]。但目前的方法仍存在一些弱点，几种误差仍是分别识别和修正，没有统一起来。三维变换方法使用的是旋转矩阵方法，使用欧拉角来作为误差调整目标，这就使调整的方式与旋转角的顺序有关，在非线性最小二乘法进行误差分析时用旋转角作为计算目标，计算过程较为复杂[8]。目前误差修正采用的大多为直接调整法[9]，如用分步调整法，目前的方法实现起来难度较大。

采用四元数进行三维几何计算，同样可以完成三维坐标变换[10]，并且可以把位置和姿态定义、空间的状态和变换、误差识别与调整统一起来。

2 四元数和坐标旋转

2.1 四元数的定义

定义四元数:

式中:a、b、c、d均为实数。 同时有:

四元数也可以表达为:

式中:v为向量，即任意向量可表示为d等于0的四元数[11]。

四元数的模:

q对应的单元四元数为。

q的共轭四元数q*，对于单位四元数，其共轭四元数也是其逆四元数q－1。

四元数乘法公式，可依据规则由式(1)推导:

2.2 坐标旋转

单位四元数可用于定义任意空间旋转，对于如下单位四元数:

q1＝pqp－1，表示对向量q绕轴v旋转2t角度后得到q1。对于坐标系的旋转，可以通过对其x、y、z轴单位向量的旋转通过计算得到。

坐标系的旋转变换可以通过旋转矩阵来定义。单位四元数可以与旋转矩阵互相转换，转换公式较为复杂不再列出。

定义旋转的另外一个方法是欧拉角，通过定义顺序施加的三个旋转角来定义空间旋转。欧拉角可以有多种不同定义方法，视实际工程应用需要来使用。本文所用旋转角与标准欧拉角正好相反，即对于坐标系D1，其坐标轴为三个单位向量(x1，y1，z1)，本文定义旋转规则:先绕z1轴旋转ψ，再绕旋转后y1轴旋转θ，最后绕x1轴旋转φ，最终与坐标系D2(x2，y2，z2)重合。

旋转角转换为四元数:

四元数转换为旋转角:

3 短线法几何计算流程

对于任何一个坐标系，可以用原点位置坐标和姿态四元数来定义其坐标系参数。为此可以做如下约定:以大地坐标系作为整体坐标系，记为G0，以固定端模顶面为基准的坐标系为U0，对于每个节段，其局部坐标系为G1，见图1。以其前端接缝线为y轴，中点为原点，x轴在顶面内且与y轴垂直，z轴与xy平面垂直，注意这种定义方法是为了适应曲线梁和有横坡的情况，此时节段前后端中心线不与前后端中心连线重合。当阶段n作为匹配段时，节段局部坐标系记为U1，此时其父坐标系为U0，各测点的测量坐标中包含了节段误差，因此阶段n＋1会偏离设计坐标。计算的目标就是识别节段n的U1坐标系参数即其位置坐标和姿态四元数，并通过计算节段n＋1测点的整体坐标系，再变换回U0坐标系，计算作为匹配段的U1坐标系参数来修正之前的误差。

图1 节段局部坐标系

首先根据设计线形，可计算出节段各关键点的大地坐标，即全局坐标系G0参数和节段局部坐标系G1参数。

首节段浇筑完成，测量各测点U0下的坐标，通过坐标变换可得到其全局坐标。

对任意已经浇筑完成的梁段n，假定其已完成浇筑且得到各测点U0坐标并计算出实际全局坐标。依照梁段n＋1的G1坐标系，可根据全局到局部坐标的变换计算梁段n匹配段的坐标系参数，开始浇筑梁段n＋1。梁段n＋1浇筑完成后，测量梁段n＋1的各测点坐标和梁段n各测点坐标，显然后者包含了预制误差，根据本文算法识别出U1坐标参数和误差四元数，可根据局部坐标系到全局坐标系的变换计算出梁段n＋1的全局坐标。此坐标已经偏离了设计全局坐标，需要根据梁段n＋2的G1坐标系参数来计算需要修正的梁段n＋1匹配段参数，这样可以继续下一段的施工和计算。

4 节段预制误差四元数的识别

梁段n作为匹配段在n＋1段浇筑后，其前端局部坐标系U1的位置和姿态可以用表示位置的四元数p(x，y，z，0)和表示姿态的单位四元数q(a，b，c，d)来表示，共计7个参数，应采用其6个测点的U0坐标系下的浇筑后坐标v1(x1，x2，x3，0)和作为匹配段的测量坐标v2(x2，y2，z2，0)来进行识别。 由此我们可以得到如下方程:

利用四元数乘法可以得到:

同时引入约束条件:

利用式(3)可以定义:

这样利用式(2)就可以建立18个方程，加上式(4)，根据非线性最小二乘法，可以得到最优化目标函数:

由于目标函数是4次多项式函数，可以使用梯度下降法求解[12]。求解需要对7个参数求导，得到梯度函数:

而f函数对a、b、c、d的偏导如下:

5 实例验证

某3×40 m曲线连续梁，采用短线法拼装7个节段，每个节段长2.9 m。为了能有效验证本文方法的正确性，以桥梁实际数据为基础，发生的预制误差预先随机产生，为匹配段偏离原位置x、y、z，姿态误差旋转角转为四元数，同时引入正负0.3 mm的随机测量误差，然后用每一个施工步骤的测量数据来做输入，产生的计算结果与预设误差对比。对比结构见表1、表2。表1为位置误差对比，单位mm；表2为计算结果的四元数转换为旋转角后的对比结果，单位10－3弧度。

表1 位置误差对比结果

表2 姿态误差对比结果 10－3rad

从对比结果来看，位置误差小于0.2 mm，姿态误差小于10－4弧度。

6 结论与展望

(1)使用四元数，不仅可以完成三维坐标变换，结合非线性最小二乘法，可以正确、可靠地分析出带有测量误差的预制误差。

(2)基于四元数的非线性最小二乘法，由于目标函数为多项式函数，可以直接求导得出梯度函数向量，迭代算法更容易收敛。

(3)采用四元数的误差分析方法，对测点位置、数量，旋转角度都没有严格的限制，可更加方便推广到其他装配式桥梁施工方法中。

(4)三维仿真已经成为短线法软件的必要部分，能够很好地避免使用过程中的人为错误，而四元数结果可以直接用于三维建模和动画。

(5)四元数可以同时用于表示坐标、对象姿态、对象姿态的变换，可以很方便地进行插值，因此可用于分步调整法，这是后续应展开研究的方面。

(6)本文方法未考虑混凝土收缩的梁长误差，这需要多引入一个缩放因子作为参数，并重新推导目标函数，这是后续需研究的工作。