Undifferenced Precise GPS Processing Software PLAOD and Its Performance Analysis
-
摘要:
PLAOD (precise location and orbit determination)是西南交通大学自主研发的GPS非差精密定位定轨软件. 该软件目前能够处理GPS观测数据,具有高精度定位、大气反演和低轨卫星轨道确定功能. 本文在介绍PLAOD相关理论的基础之上对其性能进行了分析. 利用多组GPS观测数据进行测试,通过与IGS (international GNSS service)等相关机构精密产品对比,PLAOD具备如下性能:静态定位精度可达毫米级;动态定位精度可达厘米级;静态对流层天顶延迟的估计精度通常优于1 cm,动态对流层天顶延迟的估计精度在2 cm以内;电离层天顶延迟估计结果与IGS的IONEX (ionosphere exchange)产品具有较高的吻合度;与JPL (jet propulsion laboratory)提供的简化动力学低轨卫星轨道相比,几何法确定轨道的准确度可达厘米级.
Abstract:PLAOD (precise location and orbit determination) is an undifferenced precise GPS data processing software for positioning and orbit determination fully developed by Southwest Jiaotong University. Currently, the software can process GPS observation data, and has the functions of precise positioning, atmospheric inversion and orbit determination for low earth obiter. Following the introduction to basic theory of PLAOD, its performance is focused. Multiple sets of GPS observation data were tested. Compared with the precise product provided by IGS (international GNSS service) and other institutes, PLAOD has shown following features. The static positioning accuracy can reach millimeter level, the kinematic positioning accuracy can reach centimeter level, the estimation precision of the tropospheric zenith delay in static mode is usually better than 1 cm, and the one in kinematic mode is better than 2 cm. There is a high consistence between the ionospheric zenith delay results from PLAOD and the IONEX (ionosphere exchange) product of IGS. In comparison with the reduced dynamic orbit of low earth orbiter provided by JPL (jet propulsion laboratory), the accuracy of orbit determined by geometry method is usually at the centimeter level.
-
钢轨波磨是铁路行业尚未解决的难题,在地铁线路上尤为严重. 任何一条地铁线路上,只要曲线轨道半径低于350 m,则该轨道内轨上出现波磨的几率接近100%;然而在大半径曲线轨道(R≥650~800 m)和直线轨道上,很少出现钢轨波磨:这是地铁线路上钢轨波磨的普遍规律[1-2]. 当列车通过波磨轨道时,轮轨界面上的剧烈振动加剧了车辆和轨道部件的疲劳断裂,甚至危害行车安全[3-4].
目前钢轨波磨的形成机理可分为3类:1) 轮轨表面初始不平顺激励下,轮轨界面的瞬时动力作用引起摩擦功波动导致钢轨波磨[5-10];2) 轮轨系统的粘-滑振动引起了钢轨波磨[11-13];3) 饱和轮轨蠕滑力激发了轮轨摩擦耦合系统的自激振动导致钢轨波磨[14-19]. 第一类理论模型为它激振动模型,该理论可以较好地解释地铁线路上的部分波磨现象,但不足以解释地铁线路上波磨主要出现在小半径曲线内轨,在大半径曲线和直线轨道上波磨发生率较低这一普遍现象;第二类理论模型预测得出的波磨频率约为20~80 Hz[12-13],远小于地铁线路上绝大多数波磨的通过频率;第三类波磨理论由陈光雄教授提出,该理论认为当车辆通过小半径曲线轨道时,导向轮对上的轮轨蠕滑力饱和,引发了轮轨系统的摩擦自激振动,致使钢轨表面出现波磨,该理论可以合理地解释地铁线路上规律性的波磨现象,正逐渐被学者们认可[20-23].
梯形轨枕因其具有良好的减振降噪效果被铺设在减振要求高的地段[24]. 然而,在小半径曲线梯形轨枕轨道上出现了严重的钢轨波磨[25],目前对于梯形轨枕轨道波磨成因的研究还不充分. 本文基于轮轨系统摩擦自激振动导致钢轨波磨的观点对小半径曲线上的梯形轨枕轨道波磨展开研究,建立了导向轮对-梯形轨枕轨道系统有限元模型. 该模型采用实体单元对扣件系统进行建模,考虑了弹条的扣压作用以及轨距挡块的横向限位作用,建模方式更加符合实际工况. 应用复特征值分析和瞬时动态分析分别求解了该摩擦耦合系统的运动稳定性和时域动态响应,揭示了小半径曲线上梯形轨枕轨道波磨的成因.
1. 数值模型和理论方法
1.1 导向轮对-梯形轨枕轨道系统有限元模型
当车辆通过小半径曲线梯形轨枕轨道时,导向轮对与钢轨的接触位置如图1所示. 外轮的轮缘与高轨的轨距面接触,内轮的踏面与低轨的轨顶面接触. 图1中:δL、δR分别为外轮、内轮与钢轨间的接触角,FL和NL分别为外轮与高轨接触界面上的横向蠕滑力和法向接触力,FR和NR分别为内轮与低轨接触界面上的横向蠕滑力和法向接触力,FSVL和FSVR为作用在导向轮对左、右轴端的垂向悬挂力,FSLL和FSLR为作用在导向轮对左、右轴端的横向悬挂力,KC为弹条刚度,KDL和CDL分别为横向缓冲垫的刚度和阻尼,KDV和CDV分别为垂向减振垫的刚度和阻尼.
扣件系统由轨下垫板、弹条、轨距挡块、铁垫板等组成,轨距挡块约束着钢轨的横向运动,弹条产生的扣压力将钢轨紧紧地压在轨道上,钢轨的纵向运动通过钢轨与轨下垫板间的摩擦力来约束,螺旋道钉将铁垫板固定在梯形轨枕上. 为了简化模型,本文仅对扣件系统中的轨下垫板和铁垫板进行实体单元建模.
依托于图1所示的轮轨系统接触模型,使用ABAQUS有限元分析软件建立了导向轮对-梯形轨枕轨道系统有限元模型,如图2所示. 车轮廓形为LM磨耗型踏面,轮对滚动圆半径为0.42 m;钢轨廓形为60 kg/m标准钢轨截面,模型中单根钢轨的长度为36 m;轨下垫板离散地分布在钢轨底部,相邻轨下垫板/铁垫板间的距离为0.625 m;梯形轨枕由横向钢管和预应力混凝土纵梁组合而成,相邻横向钢管间的距离为2.5 m;导向轮对、钢轨、轨下垫板、铁垫板和梯形轨枕通过实体单元建模,整个有限元模型采用C3D8I实体单元划分,模型中共包含
422968 个单元,606942 个节点. 轮轨间的接触以及钢轨和轨下垫板间的接触采用面-面接触算法和有限滑移准则进行建模. 轮轨间的摩擦系数为0.40,钢轨与轨下垫板间的摩擦系数为0.75[26]. 通过垂向压缩弹簧来模拟弹条的扣压作用,通过“Equation”约束模拟轨距挡块的横向限位作用,通过“Tie”约束模拟轨下垫板和铁垫板以及铁垫板和梯形轨枕的连接关系. 模型中采用接地弹簧和阻尼单元模拟侧向缓冲垫和垂向减振垫对梯形轨枕的缓冲减振作用,相邻缓冲减振垫间的距离为1.25 m. 缓冲减振垫的刚度为2 × 107 N/m,阻尼为7 × 104 N·s/m[27]. 各部件的材料参数见表1.表 1 材料参数Table 1. Material parameters部件 密度/
(×103 kg·m−3)弹性模量/MPa 泊松比 钢轨/轮对 7.80 210000 0.30 轨下垫板 1.19 12 0.45 铁垫板 7.80 173000 0.30 横向钢管 7.80 210000 0.30 混凝土纵梁 2.50 36000 0.20 1.2 轮轨摩擦耦合系统自激振动的有限元分析方法
采用复特征值分析来求解导向轮对-梯形轨枕轨道系统的运动稳定性. 在小半径曲线轨道段,导向轮对上的轮轨蠕滑力达到饱和,约等于滑动摩擦力. 在此情形下,轮轨摩擦耦合系统的运动方程为
M¨u+C˙u+Ku=0, (1) 式中:M、C、K分别为系统的质量、阻尼和刚度矩阵,u为节点位移向量.
式(1)的通解为[28]
u(t)=z∑i=1φieλit=z∑i=1φie(αi+jωi)t, (2) 式中:z为拟提取的特征值阶数;${\lambda _i} = {\alpha _i} + {\mathrm{j}}{\omega _i}$,为第i阶特征值;${{\boldsymbol{\varphi}} _i} $为特征向量;t为时间.
由式(2)可以看出,当特征值实部为正时,节点位移随时间指数式增长,这意味着轮轨系统发生了不稳定振动. 通常采用等效阻尼比$ \xi $来判断系统发生自激振动的趋势,如式(3).
ξ=−2αi|ωi|. (3) 当$ \xi <$0时,轮轨摩擦耦合系统发生了自激振动. 而且,$ \xi $越小,轮轨系统发生自激振动的趋势越强,钢轨波磨更易发生.
复特征值分析是一种频域分析方法,而瞬时动态分析是一种时域非线性分析方法,可以获取摩擦自激振动发生时系统的动态响应. 本文使用ABAQUS/Standard隐式求解器计算了轮轨系统的动态响应. 该求解器采用Newmark方法计算每个增量步结束时刻系统的平衡方程,进而获取节点的位移和速度[29],如式(4)、(5).
u|t+Δt=u|t+Δt˙u|t+Δt2((12−β)¨u|t+β¨u|t+Δt), (4) ˙u|t+Δt=˙u|t+Δt((1−γ)¨u|t+γ¨u|t+Δt), (5) 式中:Δt为积分步长,β = (1–α2)/4,γ = 0.5–α,−0.5≤α≤0.
2. 计算结果分析
2.1 小半径曲线段导向轮对-梯形轨枕轨道系统摩擦自激振动分析
复特征值分析的步骤为:
步骤1 在导向轮对轴端施加垂向悬挂力和横向悬挂力;
步骤2 引入轮轨间的摩擦耦合作用;
步骤3 提取轮轨系统的固有频率;
步骤4 进行复特征值分析,提取轮轨系统的摩擦自激振动.
通过定义轮对的横向滑动速度引入导向轮对和钢轨间的摩擦耦合作用. 应用复特征值分析提取了50~
1200 Hz内导向轮对-梯形轨枕轨道系统的摩擦自激振动,如图3所示. 可以看到,在该频率范围内系统共存在7个自激振动,对应频率分别为150.35、163.07、422.16、447.81、471.05、493.97 Hz和519.94 Hz.图4为7个自激振动的振动模态,自激振动主要发生在低轨和内轮上,与钢轨波磨发生在小半径曲线梯形轨枕轨道内轨上相吻合. 根据等效阻尼比越小,自激振动发生趋势越强的原则,可知频率为150.35 Hz和493.97 Hz的自激振动更易发生.
图5为曲线半径为350 m的梯形轨枕轨道上的内轨波磨,其波长$\lambda $为60~80 mm. 列车通过该轨道的平均速度$v $为37 km/h,由$\lambda = v/f$可得波磨通过频率为130~180 Hz[27]. 综上可得:当车辆通过小半径曲线段梯形轨枕轨道时,发生在导向轮对内侧轮轨界面上频率为150.35 Hz的自激振动致使内轨上出现了波长为60~80 mm的波磨.
2.2 轮轨法向接触力变化
应用瞬时动态分析求解轮轨系统的动态响应分2个步骤完成:第一步,在导向轮对的轴端施加垂向悬挂力和横向悬挂力;第二步,设置轮对运行速度为40 km/h,运行时间为0.15 s,采样间隔为5 × 10–5 s. 由于未施加任何外部激励且钢轨表面无初始不平顺,因此,计算所得的轮轨系统的响应为自激振动响应. 图6所示为轮轨界面上的法向接触力的时间历程曲线. 从图中可以看出,内侧轮轨界面上的法向接触力波动幅值要大于外侧轮轨界面上的,其中外轮与高轨间法向接触力的波动幅值约为20 kN,而内轮与低轨间法向接触力的波动幅值约为30 kN.
图7为轮轨法向接触力的功率谱密度分析结果. 可以看到,内外侧轮轨界面上法向接触力的振动主频为160 Hz. 由$\lambda = v/f$可得波长约为69 mm,该预测结果与线路上梯形轨枕轨道波磨的波长(60~80 mm)接近.
2.3 侧向缓冲垫阻尼对轮轨系统摩擦自激振动的影响
不同侧向缓冲垫阻尼CDL,轮轨系统自激振动频率的分布如图8所示. 从图中可以看出:轮轨系统自激振动的数量不随侧向缓冲垫阻尼变化,轮轨系统产生了7个自激振动;侧向缓冲垫阻尼对400~600 Hz内的自激振动的影响很小,而对150~200 Hz内的2个自激振动有较大影响,2个自激振动的等效阻尼比随侧向缓冲垫阻尼值的增大而增大,这说明系统发生自激振动的趋势减弱. 因此,增大侧向缓冲垫阻尼可以在一定程度上抑制梯形轨枕轨道上的内轨波磨.
2.4 垂向减振垫刚度对轮轨系统摩擦自激振动的影响
图9所示为不同垂向减振垫刚度KDV,轮轨系统自激振动频率的分布. 从图9可以看出,当垂向KDV=10 kN/mm时,轮轨系统有6个自激振动,最小等效阻尼比为–
0.00789 ;当KDV=20 kN/mm时,轮轨系统有7个自激振动,最小等效阻尼比为−0.00506 ;当KDV=30 kN/mm时,轮轨系统有6个自激振动,最小等效阻尼比为−0.00679 ;当KDV=40 kN/mm时,轮轨系统有6个自激振动,最小等效阻尼比为−0.00798 . 由以上分析可知,当KDV=20 kN/mm时,最小等效阻尼比高于其余3种工况,这意味着系统发生自激振动的趋势减弱. 过小或过大的垂向减振垫刚度均不利于抑制梯形轨枕轨道上的内轨波磨.2.5 梯形轨枕结构对轮轨系统摩擦自激振动的影响
图10所示为4种不同结构的梯形轨枕. 每块轨枕上横向钢管的数量N=2根时,对应结构的梯形轨枕相邻横向钢管间的距离为5.000 m;N=3根时,对应结构的梯形轨枕相邻横向钢管间的距离为2.500 m;N=5根时,对应结构的梯形轨枕相邻横向钢管间的距离为1.250 m;N=9根时,对应结构的梯形轨枕相邻横向钢管间的距离为0.625 m. 需要说明的是,横向钢管间距的设计值为2.500 m.
图11所示为4个不同结构的导向轮对-梯形轨枕轨道系统自激振动频率的分布. 从图中可以看出:梯形轨枕结构对轮轨系统的摩擦自激振动有显著影响,并非每块梯形轨枕上的横向钢管数量越多,轮轨系统就越稳定. 如蓝色椭圆圈内所示,当每块梯形轨枕上安装9根横向钢管时,轮轨系统产生自激振动的等效阻尼比反而更小,轮轨系统更不稳定,钢轨波磨更易出现. 红色椭圆圈内所示为低频段的2个自激振动,当每块梯形轨枕上安装2根横向钢管时,等效阻尼比更低,轮轨系统更易趋于不稳定. 当每块梯形轨枕上安装5根横向钢管时,等效阻尼比最大,轮轨系统发生自激振动的趋势减弱. 因此,在地铁小半径曲线地段,铺设横向钢管间距为1.250 m的梯形轨枕可以一定程度抑制梯形轨枕轨道波磨.
3. 结 论
本文研究了地铁小半径曲线梯形轨枕轨道上的内轨波磨,建立了导向轮对-梯形轨枕轨道系统有限元模型,模型中采用实体单元对扣件系统进行建模,应用复特征值分析和瞬时动态分析分别求解了轮轨系统的运动稳定性和时域动态响应. 研究了缓冲减振垫参数及梯形轨枕结构对轮轨系统摩擦自激振动的影响,得到以下结论:
1) 饱和轮轨蠕滑力引起的频率为150 Hz的轮轨系统共振振动是地铁线路小半径曲线段梯形轨枕轨道上波长为60~80 mm内轨波磨的成因.
2) 增大侧向缓冲垫阻尼可以一定程度上抑制梯形轨枕轨道上的内轨波磨,过小或过大的垂向减振垫刚度均不利于抑制钢轨波磨.
3) 改变单块梯形轨枕上横向钢管的数量实际上是改变了导向轮对–梯形轨枕轨道系统的刚度矩阵,进而影响了轮轨系统发生自激振动的趋势以及钢轨波磨的发展趋势. 铺设横向钢管间距为1.250 m结构的梯形轨枕可以抑制梯形轨枕轨道上的钢轨波磨.
-
表 1 误差改正
Table 1. Error correction
表 2 PLAOD参数处理方案
Table 2. Parameter processing strategy of PLAOD
参数 处理方案 测站坐标 静态时过程噪声为 0,动态时先验值源自伪距单点定位 模糊度 过程噪声为 0 对流层 随机游走,过程噪声默认为 4 cm2/h,动态时顾及历元间高差引起的影响 对流层梯度 随机游走 电离层 随机游走 电离层梯度 随机游走 接收机钟差 先验值源自伪距单点定位 卫星钟差 20 min弧段二阶多项式拟合 卫星坐标 10 阶滑动拉格朗日插值 表 3 测试站地理点位置
Table 3. Geographic location of testing stations
站点 经度/(°) 纬度/(°) CHAN(长春) 125.443 43.790 BJFS(北京房山) 115.892 39.609 BJNM(北京) 116.224 40.245 HKSL(香港小冷水) 113.928 22.372 HKWS(香港黄石) 114.335 22.434 JFNG(九峰) 114.491 30.516 LHAZ(拉萨) 91.104 29.657 TWTF(桃园) 121.164 24.954 表 4 站点解算坐标偏差
Table 4. Coordinate bias of solved stations
站点 年积日 dN
/mmdE
/mmdU
/mm收敛时间
/minCHAN 001 1.9 −9.2 −11.3 32.0 002 2.4 −5.7 −7.5 19.0 003 2.0 −7.1 −9.2 43.0 BJFS 001 2.0 1.2 −1.8 30.0 002 2.5 0.7 −0.5 21.0 003 2.5 −0.2 2.8 23.0 BJNM 001 2.6 −3.6 −0.8 29.0 002 4.5 0.4 1.4 50.0 003 4.0 −4.5 4.3 46.5 HKSL 001 2.0 0.5 −3.9 30.5 002 3.8 0.1 1.7 25.0 003 4.4 −2.5 2.0 20.5 HKWS 001 1.0 1.5 2.2 30.0 002 3.3 1.9 −0.1 52.5 003 3.5 −0.3 1.9 104.5 JFNG 001 2.0 3.2 0.6 20.0 002 2.8 4.4 9.0 13.0 003 3.1 3.4 8.1 19.0 LHAZ 001 3.6 −1.5 −3.8 42.0 002 4.1 1.5 −0.3 57.0 003 3.5 −1.5 −4.4 44.0 TWTF 001 0.8 3.9 −9.3 26.5 002 4.4 1.0 −5.4 21.5 003 4.0 0.1 −6.3 19.0 -
[1] ZUMBERGE J F, HEFLIN M B, JEFFERSON D C, et al. Precise point positioning for the efficient and robust analysis of GPS data from large networks[J]. Journal of Geophysical Research, 1997, 102(B3): 5005-5017. doi: 10.1029/96JB03860 [2] KOUBA J, HÉROUX P. Precise point positioning using IGS orbit and clock products[J]. GPS Solutions, 2001, 5(2): 12-28. doi: 10.1007/PL00012883 [3] GAO, Y, SHEN, X. A new method for carrier phase based precise point positioning[J]. Navigation, 2002, 49(2): 109-116. doi: 10.1002/j.2161-4296.2002.tb00260.x [4] 张小红,李星星,李盼. GNSS精密单点定位技术及应用进展[J]. 测绘学报,2017,46(10): 1399-1407. doi: 10.11947/j.AGCS.2017.20170327ZHANG Xiaohong, LI Xingxing, LI Pan. Review of GNSS PPP and its application[J]. Acta Geodaetica et Cartographica Sinica, 2017, 46(10): 1399-1407. doi: 10.11947/j.AGCS.2017.20170327 [5] GENG Jianghui. Rapid integer ambiguity resolution in GPS precise point positioning[D]. Nottingham: University of Nottingham, 2010. [6] SIMON B, RICHARD B. L, MARCELO C. S. The precise point positioning software centre: an insight into online PPP services[R]. Buenos Aires: IAG, 2009. [7] IBANEZ D, ROVIRA G A, SANZ J, et al. The GNSS laboratory tool suite (gLAB) updates: SBAS, DGNSS and global monitoring system[C]//9th ESA Workshop on Satellite Navigation Technologies and European Workshop on GNSS Signals and Signal Processing. Noordwijk: IEEE, 2018: 1-11. [8] TAKASU T. RTKLIB: open source program Package for RTK-GPS[R]. Tokyo: Tokyo University of Marine Science and Technology, 2009. [9] GENG J H, CHEN X Y, PAN Y X, et al. PRIDE PPP-AR: an open-source software for GPS PPP ambiguity resolution[J]. GPS Solutions, 2019, 23(4): 1-10. [10] BANVILLE S, GENG J H, LOYER S, et al. On the interoperability of IGS products for precise point positioning with ambiguity resolution[J]. Journal of Geodesy, 2020, 94(1): 1-15. doi: 10.1007/s00190-019-01332-z [11] ASHBY N. Relativity in the global positioning system[J]. Living Reviews in Relativity, 2003, 6(1): 1-42. doi: 10.12942/lrr-2003-1 [12] PETIT G, LUZUM B. IERS Conventions (2010) [EB/OL]. [2021-7-26]. https://www.iers.org/IERS/EN/Publications/TechnicalNotes/tn36.html. [13] WU J T, WU S C, HAJJ G A, et al. Effects of antenna orientation on GPS carrier phase[J]. Manuscripta Geodaetica, 1993, 18(2): 91-98. [14] WELCH G, BISHOP G. An introduction to the Kalman filter[EB/OL]. [2021-07-26]. https://www.cs.unc.edu/~welch/media/pdf/kalman_intro.pdf [15] GUO J Y, YUAN Y D, KONG Q L, et al. Deformation caused by the 2011 eastern Japan great earthquake monitored using the GPS single-epoch precise point positioning technique[J]. Applied Geophysics, 2012, 9(4): 483-493. doi: 10.1007/s11770-012-0360-2 [16] KOUBA J. Implementation and testing of the gridded Vienna mapping function 1 (VMF1)[J]. Journal of Geodesy, 2008, 82(4/5): 193-205. [17] BOEHM J, NIELL A, TREGONING P, et al. Global mapping function: a new empirical mapping function based on numerical weather model data[J]. Geophysical Research Letters, 2006, 33(7): 304-317. [18] BOEHM J, HEINKELMANN R, SCHUH H. Short note: a global model of pressure and temperature for geodetic applications[J]. Journal of Geodesy, 2007, 81(10): 679-683. doi: 10.1007/s00190-007-0135-3 -