CN109459129A - 一种高效求解高速列车过桥振动随机特性的分析方法 - Google Patents

Abstract

本发明公开一种高效求解高速列车过桥振动随机特性的分析方法，将虚拟激励算法与自适应高斯积分技术相结合，自适应地识别出来随机激励对耦合系统具有关键性影响的频率区间段，并且智能化地细化关键性频率区间段以及粗画或摒弃非关键性频率区间，基于智能化处理后的频率区间段，寻找对应的频率点，代入车—桥系统耦合振动进行计算；基于各频点所得动态响应，运用高斯积分求得具有统计性及随机特性的计算结果，并且运用Monte Carlo方法对计算结果进行验证。本发明的方法有效的减少了计算成本，提高了计算效率。

Description

一种高效求解高速列车过桥振动随机特性的分析方法

技术领域

本发明涉及高速列车运行安全技术领域，具体为一种高效求解高速列车过桥振动随机特 性的分析方法。

背景技术

1964年日本建成了世界上第一条列车时速超过200公里的高速铁路，标志着现代高速铁 路时代的来临。随后，法国、德国和中国修建了大量的高速铁路，并取得了较好的经济效益 和社会效益。2008年京津城际铁路正式通车，这预示着我国高速铁路时代正式开启。随着高 速铁路趋于高速化和轻型化的发展，为保证列车高速运行的安全性、稳定性和乘客舒适性， 必须确保线路的高平顺性、稳定性及可靠性，因此高速铁路线路多采用桥梁结构代替常规的 路基结构。桥梁所占比例日益增大，长大跨度桥梁数量与日俱增。例如，京沪高速铁路中桥 梁总里程占据总里程数的80.5％。因此，在高铁工程建设中，桥梁结构被广泛的应用，车—桥 系统耦合振动的随机特性研究具有较广泛的应用前景。

由于轨道不平顺、机车车辆质量的不完全对中以及轮轨缺陷等因素，车辆与桥梁结构之 间具有自激特性，列车在桥梁结构上运行过程中，车辆与桥梁会发生相互作用，引起耦合振 动。轨道不平顺是车—桥耦合系统振动的主要激励源，轨道不平顺存在随机特性，车—桥系 统的耦合振动势必会体现出随机特性。以相同初始条件对车辆过桥的全过程进行多次独立的 数值模拟时，车辆和桥梁结构的动态响应结果具有较强的离散性。车辆上桥前的初始振动同 样具有很强的随机性，可作为随机变量，已有研究表明该部分振动对桥梁结构冲击系数的影 响达到40％左右。此外，桥梁结构和车辆结构的结构参数均具有一定的随机特性，也可看作 为随机变量。针对车—桥耦合系统随机响应，应采用随机振动理论，从概率和统计的角度， 对随机激励作用下的耦合系统的稳定性以及可靠性进行求解，并得到具有统计性的结果。

当以相同初始条件对车辆过桥的全过程进行多次独立的数值模拟时，车辆和桥梁结构的 动态响应结果呈现出较强的离散性，需要在随机振动理论框架下对车—桥系统耦合振动的随 机特性进行研究。传统的Monte Carlo方法是以概率统计理论为基础的一种重要的数值计算方 法，适用于车—桥耦合系统响应随机特性的模拟分析，但是需要大量的样本计算用于保证计 算结果的可靠性，计算成本过高，计算效率较低。虚拟激励算法能够有效地求解了耦合系统 随机振动的问题。但是，随着科学研究及工程技术的不断进步，车辆和桥梁结构趋于精细化 和复杂化，车—桥耦合系统的频率响应范围逐渐加宽，若仍采用常规积分形式与虚拟激励算 法相结合来分析系统响应的随机特性，会导致所需要计算的随机激励的离散频点数量过多， 计算工作量仍然较大。如果单纯为了提高计算效率而减少离散频点数，则很容易遗漏关键频 率点，导致计算结果无法满足精度要求。

到目前为止，诸多学者采用了虚拟激励算法对在不同类型随机激励作用下的车—桥耦合 系统的随机振动问题进行了研究，并取得了丰硕的成果。然而，车—桥耦合系统的参数矩阵 将会随着车辆在桥梁运行的位置不同而发生变化，只能够通过时域求解方法对动态响应进行 计算，进而在随机振动理论框架中分析系统耦合振动随机特性。随着桥梁结构复杂化和跨度 不断增大以及车辆模型精细化，耦合系统的激振频率区间变宽，计算频率点数增多，导致了 计算时间较长，计算效率不高。因此，在保证高计算精度要求下，急需发展能够进一步提高 车—桥系统耦合振动的随机特性计算效率的优化方法

将虚拟激励算法(Pseudo Excitation Method，简称PEM)与自适应高斯积分技术(Self-adaptive Gauss integration，简称SGI)相结合(PEM-SGI)，并且有效地应用于求解车 —桥系统耦合振动随机特性。该方法能够自适应地自动化识别出来随机激励对耦合系统具有 关键性影响的频率区间段，并且智能化地细化关键性频率区间段并寻找对应的频率点，代入 车—桥系统耦合振动进行计算，基于各频点所得动态响应，积分求得具有统计性及随机特性 的计算结果，并且运用Monte Carlo方法对计算结果进行验证。此外，对所提出的优化算法与 传统计算方法进行比对，定量地表明优化算法所提高的计算效率。

发明内容

针对上述问题，本发明的目的在于提供一种能够大幅度提高车—桥系统耦合振动随机特 性的求解效率的优化自适应求解方法，能够自适应地识别出来随机激励对耦合系统具有关键 性影响的频率区间段，并且智能化地细化关键性频率区间段以及粗画或摒弃非关键性频率区 间，能够减少计算成本，提高计算效率。技术方案如下：

一种高效求解高速列车过桥振动随机特性的分析方法，包括以下步骤：

步骤1：建立车-桥耦合系统，有：

其中，m_c为车体质量，y^c为车体位移，y为梁体位移，r为轨道不平顺性，K和c分别 为结构刚度和阻尼；

对于车辆运行过程中的非平稳随机振动激励，用平稳过程乘以时变函数表示的非平稳过 程表示为：

x′(t)＝g(t)·x(t)

其中，x(t)为平稳过程，其功率谱密度函数为S_xx(ω)，g(t)为时变调制函数；

针对非平稳随机振动，构造虚拟激励：

其中，e^iωt表示所构造的简谐函数；

对所述虚拟激励进行时程积分得到非平稳虚拟响应y′^*(t)，则响应的功率谱密度函数为：

其中，表示为积分所得的非平稳虚拟响应；

步骤2：基于虚拟激励算法的基本理论，将车辆非平稳随机振动激励作用下，产生的响应 y(t)的方差表示为：

其中，小标y为响应结果；

假设公式的积分上下限截止频率分别为ω_l和ω_u，将区间[ω_l,ω_u]粗略等分为N个子区间 j＝1～N；

步骤3：基于相应结果，用高斯积分公式计算各子区间内的积分

步骤4：计算各子区间积分之和

步骤5：判断如果表示第j个子区间收敛，则将该子区间作为需要考察的 区间段，并将该子区间再细分为N个子区间段；否则，该子区间未收敛，不作细分；ε为预 给定的常数，即两次计算的误差值收敛值；

步骤6：针对未收敛的子区间令P＝N，i＝1～P，则将子区间划分为P个 子区间段[ω_i，ω_i+1]，则有计算；

步骤7：基于之前所给的误差收敛值ε，判断如果则需要对第i个子区 间段再细分，以此类推；

步骤8：如前后两次计算的结果已经收敛，则退出计算；反之，返回第步骤6继续计算；

步骤9：当所有的计算停止，所有的收敛结果相加，求得车辆非平稳随机振动激励作用下 相应的方差的最终结果

本发明的有益效果是：本发明通过自主地判断计算结果的收敛情况且功能性地分解频率 区间求解得到最终结果的过程，实现了自适应优化算法中自动化识别随机激励对系统随机特 性具有重要与不重要影响的频率区间并智能化细分重要频率区间的功能性；其中，需要进行 多次细化才能够得到最终收敛结果的频率区间段，即为随机激励对车—桥耦合系统随机振动 特性的重要区段；而无需细化的，即为不重要频率区间段，仅仅需要对其进行粗略分段或直 接摒弃，减少了计算成本，提高了计算效率；同时，研究人员可根据工程和科研项目的要求， 灵活地选择ε，调整误差目标值，从而找到符合实际需求的计算效率和精度的平衡点。

附图说明

图1为PEM-SGI求解车—桥耦合系统振动随机特性流程图。

图2为混凝土桥塔总体布置图；(a)为桥塔顺桥向构造图；(b)为桥塔横桥向构造图；(c) 为塔柱截面，即(b)中A-A处的横截面。

图3为车—桥耦合系统垂向加速度响应均值对比结果：(a)为车体；(b)为桥梁。

图4为PEM-SGI和Monte Carlo法的车辆随机特性结果对比图：(a)为车体垂向位移均方差曲 线；(b)为车体垂向速度均方差曲线；(c)为车体垂向加速度均方差曲线；(d)为车体横向位移均 方差曲线；(e)为车体横向速度均方差曲线；(f)为车体横向加速度均方差曲线。

图5为PEM-SGI和Monte Carlo法的桥梁随机特性结果对比图：a)为桥梁垂向位移均方差曲 线；(b)为桥梁垂向速度均方差曲线；(c)为桥梁垂向加速度均方差曲线；(d)为桥梁横向位移均 方差曲线；(e)为桥梁横向速度均方差曲线；(f)为桥梁横向加速度均方差曲线。

具体实施方式

下面结合附图和具体实施例对本发明做进一步详细说明。本发明涉及到车辆动力学、桥 梁结构动力学、随机振动理论等多方面研究内容。将虚拟激励算法与自适应高斯积分技术相 结合(PEM-SGI)，并且有效地应用于求解车—桥系统耦合振动随机特性。该方法能够自适应 地识别出来随机激励对耦合系统具有关键性影响的频率区间段，并且智能化地细化关键性频 率区间段以及粗画或摒弃非关键性频率区间，基于智能化处理后的频率区间段，寻找对应的 频率点，代入车—桥系统耦合振动进行计算。基于各频点所得动态响应，运用高斯积分求得 具有统计性及随机特性的计算结果，并且运用Monte Carlo方法对计算结果进行验证。此外， 对所提出的优化算法与传统计算方法进行比对，定量地表明优化算法所提高的计算效率。该 算法不仅适用于求解车—桥耦合系统随机振动的相关问题，同样也适用于求解其他外荷载作 用下的结构或系统的随机振动问题。PEM-SGI方法的提出，采用由易到难的原则进行方法适 用性的验算，首先选取白噪声作用下的单质点线性振子对该方法进行初步的测评后；再应用 于求解随机外荷载作用下的结构随机振动分析(如风荷载作用下的桥塔或地震荷载作用下的 连续梁桥)当中，并与常规的虚拟激励算法结果进行比对，明确自适应优化算法的有效性和 高效性。最后，分析车—桥系统耦合振动随机特性，定量地分析计算效率的提高。

PEM-SGI求解车—桥耦合系统振动随机特性流程图如图1所示，高效求解高速列车过桥 振动随机特性的分析方法，包括以下步骤：

步骤1：建立车-桥耦合系统，有：

其中，m_c为车体质量，y^c为车体位移，y为梁体位移，r为轨道不平顺性，K和c分别 为结构刚度和阻尼；

对于车辆运行过程中的非平稳随机振动激励，用平稳过程乘以时变函数表示的非平稳过 程表示为：

x′(t)＝g(t)·x(t)

其中，x(t)为平稳过程，其功率谱密度函数为S_xx(ω)，g(t)为时变调制函数；

针对非平稳随机振动，构造虚拟激励：

其中，e^iωt表示所构造的简谐函数；

对所述虚拟激励进行时程积分得到非平稳虚拟响应y′^*(t)，则响应的功率谱密度函数为：

其中，表示为积分所得的非平稳虚拟响应；

步骤2：基于虚拟激励算法的基本理论，将车辆非平稳随机振动激励作用下，产生的响应 y(t)的方差表示为：

其中，小标y为响应结果；

假设公式的积分上下限截止频率分别为ω_l和ω_u，将区间[ω_l,ω_u]粗略等分为N个子区间 j＝1～N；

步骤3：基于相应结果，用高斯积分公式计算各子区间内的积分

步骤4：计算各子区间积分之和

步骤5：判断如果表示第j个子区间收敛，则将该子区间作为需要考察的 区间段，并将该子区间再细分为N个子区间段；否则，该子区间未收敛，不作细分；ε为预 给定的常数，即两次计算的误差值收敛值；

步骤6：针对未收敛的子区间[ω_j,ω_j+1]，令P＝N，i＝1～P，则将子区间[ω_j,ω_j+1]划分为P个 子区间段[ω_i，ω_i+1]，则有计算；

步骤7：基于之前所给的误差收敛值ε，判断如果则需要对第i个子区 间段再细分，以此类推；

步骤8：如前后两次计算的结果已经收敛，则退出计算；反之，返回第步骤6继续计算；

步骤9：当所有的计算停止，所有的收敛结果相加，求得车辆非平稳随机振动激励作用下 相应的方差的最终结果

2、涉及到的虚拟激励算法具体如下：

假设初始静止的n自由度时变系统在m维随机激励作用下的运动方程为

其中，[M],[C],和[K]分别为随时间变化的n自由度的结构质量、阻尼和刚度矩阵。{y}, 和分别表示为相应系统的位移、速度和加速度向量。[E]为表征作用力分布状况的n×m矩阵。{f(t)}为m维作用力向量，假设其功率谱密度函数为[S_ff(ω)]。

针对平稳随机振动过程，即激励为平稳随机过程，系统为线性系统，根据传统随机振动 理论，{y}的功率谱密度函数矩阵能够写成

[S_yy(ω)]＝[H(ω)]^*·[S_ff(ω)]·[H(ω)]^T (2.2)

其中，[H(ω)]频率响应函数的矩阵，[H(ω)]^*和[H(ω)]^T分别为其的共轭矩阵和转置矩阵。 通常来说，激励矩阵[S_ff(ω)]必为Hermitian矩阵，使用Cholesky分解对其进行处理，可以得 到

其中，[L]是下三角矩阵，{l}_k则是[L]的第k列矩阵。将式(2.3)代入式(2.2)得到

构造虚拟激励，如下式所示

针对每一个虚拟激励有它引起的位移响应向量必为

由式(2.4)到式(2.6)，可得

基于所得到的响应的方差值可以表示为

显然，的对角线元素就是结构位移的方差值。

假设结构阻尼矩阵满足正交性，通过振型分解的方法对上述过程进行计算，从而提高计算 效率。将式(2.5)的k阶虚拟激励作为输入激励，对多自由度运动方程进行分解，可得

其中，ω_j和ξ_j分别为结构第j阶模态的自振频率和阻尼比，u_j,k为在k阶虚拟激励作用下， 第j阶模态坐标响应。

其中，为结构第j阶的模态信息，u_j,k的稳态响应为

其中，第j阶频率响应函数为

基于振型叠加原理，在第k阶虚拟激励作用下的结构响应为

将式(2.13)代入式(2.7)，可得

令

因此，式(2.14)则可以写成

对于非平稳随机振动，即系统仍然为线性系统，激励为非平稳的随机过程。假设x(t)为零 均值高斯平稳随机过程的激励荷载，其功率谱密度函数为[S_ff(ω)]，x′(t)为非平稳随机过程的 激励荷载，g(t)为描述非平稳特征的时变调制函数，则非平稳随机过程为

x′(t)＝x(t)·g(t) (2.17)

由式(2.3)与式(2.17)，构造虚拟激励

对比式(2.5)和式(2.18)可知，非平稳随机过程的虚拟激励荷载多了一个时变调制函数，虽然 不能够直接采用频率响应函数直接得到响应结果，需要使用时程积分得到非平稳的虚拟响应 结果。

根据Duhamel积分，线性结构在非平稳随机激励x′(t)的作用下，结构响应y(t)为

其中，h(t-τ)代表了脉冲响应函数

因此，该结构响应y(t)的自相关函数可得

其中，R_xx(τ)＝R_xx(τ₂-τ₁)＝E[x(τ₁)·x(τ₂)]。

根据维纳-辛钦公式可以得到

其中，S_xx(ω)是x(t)的功率谱密度函数。将式(2.21)代入式(2.20)，得

其中

根据式(2.22)，响应y(t)的方差可以表示为

同时考虑

则，对比式(2.24)和(2.25)，可以得到y(t)的功率谱密度函数为

S_yy(ω,t)＝S_xx(ω)·I(ω,t)·I^*(ω,t) (2.26)

根据式(2.18)，激励的确定性响应为

则，y(t)的功率谱密度函数可以写为

上述公式，说明了虚拟激励算法的本质，即通过构造虚拟的简谐激励荷载，输入系统得到 虚拟的简谐响应；基于所得的虚拟简谐响应，利用积分方法得到响应的功率谱密度函数。从 而，采用数值积分技术计算得到响应结果的方差结果。通常，线性结构在零均值高斯激励的 作用下，其动力响应必然服从高斯分布且均值为零。因此，由式(2.25)求得响应结果的方差值 其结构响应y(t)服从分布。针对多自由度体系，该算法能够有效的提高计 算效率。

3、虚拟激励算法求解车—桥耦合系统振动的随机特性的具体步骤：

对于梁体而言，作用的动荷载为

其中，m_c为质量块质量，y^c为质量块位移。显然，与轨道不平顺性无显式相关。

对于质量块而言，有

其中，y为梁体位移，r为轨道不平顺性，K和c分别为弹性系数和阻尼。式(2.30)为

右边激励第二项本质为随机激励。如已有完善的车桥振动程序，可以方便地用虚拟激励 解决。对于一个线性时不变系统，输入信号的功率谱密度为S_xx(ω)，则输出信号的功率谱密 度为

S_yy(ω)＝|H(ω)|²S_xx(ω) (2.32)

其中，H(ω)为线性时不变系统的频率响应函数。其物理意义为输入信号为简谐波时，输 出信号也为同频率的简谐波，但其振幅放大倍数为H(ω)

由此，构造确定性虚拟激励

由于式(2.33)为一简谐波，则由频率响应函数H(ω)的物理意义可知，其对应的输出信号 也为同频率的简谐波，只是振幅放大H(ω)，即

则，比较式(2.32)和(2.34)，不难看出

虚拟激励算法的实质就是，通过构造虚拟的简谐激励，输入系统得到虚拟的简谐响应， 通过求虚拟简谐响应的模得到响应的功率谱。这一算法过程对于单质点振子无实质性的意义， 但是，对于多自由度体系，有明显的效率提高。

以上过程是基于平稳随机过程的，对于非平稳随机振动，也即系统仍然是线性的，但是 激励是非平稳的，如果激励能够写出一个平稳过程乘以一个时变函数，即

x′(t)＝g(t)·x(t) (2.36)

其中，x(t)为平稳过程，其功率谱密度函数为S_xx(ω)，g(t)为时变调制函数，x′(t)为调制后 的非平稳过程。构造虚拟激励

由于式(2.37)所示的激励多了一个时变调制函数，所以不能采用式(2.34)的频率响应函数 直接得到响应，需要通过时程积分得到非平稳虚拟响应y′^*(t)。类似地，响应的功率谱密度函 数为

根据虚拟激励算法的基本思想。对于梁体而言

其中，m_c为质量块质量，y^c为质量块位移。对于质量块而言，有

上述两个方程是耦合的。从原理上分析，如果采用迭代求解，对于式(2.39)，由于首先不 知道假设其为0，则式(2.40)为移动荷载过桥问题，得到y(t)和将其代入式(2.40)， 则车辆振动的激励分为两部分，为一确定性激励，其产生y^c响应的均值注 意此时均值是随时间变化的；为随机激励部分，由于式(2.40)是一线性系统，则由这 部分随机激励产生的y^c响应与轨道不平顺性r同频率。则

其中，α为一系数。将上式带入式(2.40)

对于上式，右边第一项仍然是对应了确定性的均值部分；第二项 (m_c·α·r)·δ(x-vt)，从数学上分析，实际上是平稳过程r乘以非平稳时间调制函数(m_c·α)·δ(x-vt)，所以其为一类似式(2.36)的非平稳过程。那么，其求解可以采用式(2.36)～(2.38)的思路。求得新的梁体响应后，注意此时由于激励非平稳，梁体响应也为非平稳， 反馈给式(2.40)，迭代循环，直至收敛，这一收敛计算过程只是前面分析的两步的重复。

上述分析可以看出，(1)车桥系统是非平稳的；(2)这种非平稳是由(m_c·α)·δ(x-vt)这 个调制函数引起，这个调制函数是隐式的，而非显式的。可以利用式(2.36)～(2.38)的思路，直 接进行时程积分计算，反映上述非显式的调制函数。

轨道激励有三种不同的类型：方向轨道不平顺Y_s(x)、高低轨道不平顺Z_s(x)、水平轨道 不平顺θ_s(x)。该系统可假设为均匀调制演变随机过程，即

Y_s(x)＝g_y(x)r_y(x) (2.43)

Z_s(x)＝g_z(x)r_z(x) (2.44)

θ_s(x)＝g_θ(x)r_θ(x) (2.45)

其中，r_y(x)、r_z(x)、r_θ(x)分别为各激励的零均值平稳随机过程。其功率谱密度函数分 别为S_y、S_z、S_θ；g_y、g_z、g_θ；分别为各激励确定性慢变调制函数。通过关系式x＝vt可 以将轨道谱由空间域转为时间域，即

考虑列车各车轮所受轨道激励的相位差，方向、高低、水平三类轨道激励为三维异相位 非平稳随机激励。假定三类轨道激励之间没有相关性，则其功率谱矩阵可以写为

其中：n代表随机变量的维数；λ_p和分别为S(ω)的第p阶特征值和特征向量；上标* 为矩阵的复共轭。

λ₁＝S_y，λ₂＝S_θ，λ₃＝S_z，

根据线性时变系统非平稳随机响应的多维虚拟激励法的基本概念，即

其中：G_r(t)＝diga([g_r1(t-t₁)g_r2(t-t₁)…g_rm(t-t₁)])为确定性时间包络函数矩阵， m为相干点数(对于车桥问题即为车轮数)。

通过离散频率域的方式来确定计算频点个数n；在第j个频点构造系统虚拟激励如上述 F_ri(ω_j,t)(j＝1,2,...,n)；根据公式生成虚拟激励的时程样本，将其当做轨道不平顺激励带入运 动方程当中，进行求解；通过杜哈梅积分求解得到对应的虚拟响应Y_i(ω_j,t)(j＝1,2,...,n)；判 断计算结果是否满足精度要求，如满足则证明计算收敛，退出程序；不满足，继续进行迭代 计算。根据所得虚拟响应，将所感兴趣部分的响应结果提取出来，计算得到响应的功率谱和 方差。

其中，△ω为频率增量。

4、关于高斯积分

在计算空间单元的荷载列车及刚度矩阵时，积分函数式较为复杂，且大多具有隐性特征， 即使能够写出积分显式，其积分也极其复杂，往往采用数值积分来代替函数的定积分。数值 积分是指，在积分区域内按照一定规则选出一些点，称为积分点。基于算选出的积分点，求 解得到积分函数对应的结果，乘以相应的加权系数并求和，作为近似的积分值。常用的数值 积分方法有：矩形法、梯形法和辛普森法，该三种求解方法，都是限定采用所选区间内的等 分点作为积分点，用不同的k次多项式近似求解被积函数 y(x)＝a_mx^m+a_m-1x^m-1+a_m-2x^m-2+…+a₀。矩形法为常数近似，梯形公式具有一阶代数精 度且线性近似，辛普森公式则具有三阶代数精度且抛物线近似。高斯积分并不需要采用等分 点作为积分点数，具有相同的积分点数达到较高精度或用较少的积分点数达到同样精度的特 征。高斯积分的基本公式为

式中，A_j(j＝1～n)为不依赖f(x)的求积系数，x_j(j＝1～n)为积分点。适当地选取A_j和x_j可 使式(2.53)具有2n-1次代数精度。

在已有的文献中给出了在[-1,1]积分区间下高斯积分公式的积分点x_j和权重因子A_j值。 对于n个积分点，只要选取恰当的权重因子及积分点位置，能够使式(2.53)在被积分函数为不 超过(2n-1)次多项式时精确成立。式(2.53)所示的任意区间下的积分，只需将其变换到[-1,1] 区间即可求解。显然，对于式(2.25)所示的积分进行全区间高斯数值积分很难获得满意的结果。 更为合理的做法是，将积分区间划分为若干子区间，对每个子区间单独进行高斯数值积分。 因而，如何划分高斯积分子区间，对于提高计算准确性和效率至关重要。

(1)单质点线性振子的随机振动分析

为了验证自适应高斯积分结合虚拟激励算法的有效性和高效性，以受平稳高斯白噪声激 励作用的单质点线性振子振动为例，其运动方程为

式中，ω₀和ξ为振子的圆频率和阻尼比，分析中取二者的值分别为1.0rad/s和0.03。f(t) 是高斯白噪声激励，强度S₀取为0.001。虚拟激励分析时，按单边谱进行分析，ω₀积分上限 取20rad/s，积分下限取0rad/s。

根据经典随机振动理论，单质点振子在平稳高斯白噪声激励作用下响应方差的精确理论 解为对于本节内容，均方差σ_y的解析解计算为0.2288。分别采用梯形积分和自 适应高斯积分求解响应均方差，自适应高斯积分计算中，ε取0.0001，N取10，计算结果如 表1所示。由表1可得，当采用梯形积分进行数值求解时，要达到小数点后4位的计算精度， 积分区间需要划分为1000等份子区间，也即是要对1000个频率离散点计算功率谱密度S_yy(ω)， 再对其进行梯形积分求和才可得到响应方差。当采用自适应高斯积分技术时，达到小数点后 4位的计算精度只需划分成100个积分子区间，积分区间数是梯形积分区间数的10％。而整 个计算过程共进行了136次区间划分，由于每个高斯积分区间有两个频率离散点，PEM-SGI 共需对272个频率离散点计算功率谱密度S_yy(ω)，计算效率提高了3.62倍。

表1单质点振子随机振动分析结果

(2)风荷载作用下的超高桥塔振动随机特性分析

以某斜拉桥的超高桥塔为例，进行顺桥向脉动风作用下的结构风致随机特性分析。桥塔 构造布置如附图2所示。该桥塔总高度为398.16m，桥塔塔柱沿高度方向为线性变截面形式， A-A断面以下的塔柱截面壁厚为2m，即t1＝t2＝2m；A-A断面以上的塔柱截面壁厚为1.5m， 即t1＝t2＝1.5m。桥塔采用C50混凝土，弹性模量取值为3.45×10⁴MPa，泊松比采用0.2，阻 尼比为0.02，密度为2549kg/m³。

计算采用Davenport提出的风荷载功率谱，其表达式如下

式中，为10m高度处平均风速，本次计算取值为20m/s；为高度z处平均风速；k为地面粗糙度系数，本文取值为0.03；ω为激励的圆频率，进行虚拟激励分析时，采用单边功率谱，ω积分上限取30rad/s，积分下限取0rad/s；

为作用在结构上z高度处的平均风荷载，表达式如下

式中，ρ是空气密度，取值为1.25kg/m³；A(z)为高度Z处迎风面积；C_d＝2.2311为阻力 系数，本文所采用的数据使用风洞试验数据。

平均风速沿高度变化采用指数方程，和的关系表达式为

式中，α是底面粗糙度指数，本文取值0.12。

z₁与z₂高度两点的风荷载的互谱密度定义为

式中，ρ(ω,z₁,z₂)为相干函数，可以写成

式中，C₁取值7.0。

分别采用梯形积分和自适应高斯积分技术对结构随机振动求解桥塔顶部的横向位移响应 方差，计算结果如表2所示。由表可得，如需达到小数点后两位有效数字的计算精度，梯形 积分至少需将计算区间划分为1500段积分子区间，即要对1500个频率离散点计算功率谱密 度S_yy(ω)。运用自适应高斯积分技术，仅仅需要将计算区间划分为55段积分子区间即可满足 小数点后两位的计算精度，是梯形积分区间数的3.7％。整个自适应过程的区间划分，需要计 算110个频率离散点，是梯形积分计算次数的7.3％。如需达到小数点后4位的计算精度，两 种方法得到的塔顶位移均方差均收敛于0.5334m。梯形积分需要的积分子区间数为2400段， 也即功率谱计算离散频点数为2500次。自适应高斯积分总共需要127次区间划分，对应的功 率谱计算离散频点数为254次，计算次数仅为梯形积分的10.6％。

表2桥塔风致随机振动分析结果

将自适应高斯积分方法引入随机振动的虚拟激励算法。通过自适应迭代将积分区间细分 为若干个合理的积分子区间，对每个积分子区间单独进行高斯数值积分，使得积分的计算效 率显著提高。通过数值算例可以表明，在同样精度要求下，自适应高斯积分技术能显著减少 计算次数，进一步提高了虚拟激励算法的计算效率。对于大型实际工程结构的随机振动分析， 该方法具有显著的工程实践意义。

(3)车—桥耦合系统振动的随机特性分析

为了系统地研究车-桥耦合系统的在随机激励作用下的动态特性，对于列车过桥的全过程 进行模拟，车-桥耦合系统的随机振动响应分为确定性响应和随机性响应。采用PEM-SGI研 究列车通过3跨简支梁桥时，车-桥耦合系统的随机振动响应。三种轨道不平顺的平稳随机过 程的功率谱采用德国高速线路不平顺轨道谱密度(低干谱)。模拟列车过桥全过程，分别采用 PEM-SGI和Monte Carlo统计方法对车-桥耦合随机振动问题进行求解。车辆与桥梁的垂向加 速度均值结果，如附图3所示。将PEM-SGI的均值结果与500次MonteCarlo算法的统计结 果进行比较，结果表明，车体和桥梁结构垂向加速度吻合度较高，所得结果能够满足计算要 求。

针对于车—桥耦合系统动态响应的均方差结果，车辆模型则选取了车体的横向及竖向振 动形态作为目标值进行统计，结果如附图4所示，竖坐标为车体各个响应的计算结果，横坐 标则表示了列车运行位置，0m处为桥梁端部，车辆由距离桥梁20m处开始运行(-20m)采 用500个样本计算所得的均方差曲线比50个样本计算所得的曲线更为平缓，与PEM-SGI所 求的曲线更为吻合。随着Monte Carlo方法的样本增加，其计算所得的统计均方差更为接近虚 拟激励算法所得的计算结果。这就反映了Monte Carlo法的重要特性，即计算结果的精确性与 样本数具有密切关系。为了求得精确的结果，需要选取足够多的样本数量。以车体垂向加速 度为例，Monte Carlo算法取50个样本时，与PEM-SGI所得结果的最大偏差达到22.9％；当 样本数量为500个时，最大误差只有5.1％。针对车体的横向运动状态，以车体的横向加速度 结果为例，Monte Carlo算法的样本个数为50时，与PEM-SGI计算结果的最大误差为23.23％， 随着样本数量增加到500，其最大误差随之减小到7.8％。

桥梁模型主要选取了第二跨跨中的横向及竖向的位移、速度和加速度作为考察目标进行 归纳，均方差曲线如附图5所示，竖坐标为车体各个响应的计算结果，横坐标则表示了列车 运行位置。在车辆驶入第二跨之前，第二跨桥梁结构无任何振动响应；在车辆行驶至40m处， 桥梁结构出现动态响应。因此，桥梁横坐标为40m～140m。由图可见，采用500个样本的Monte Carlo方法计算所得结果与PEM-SGI计算所得结果吻合度较高。以桥梁垂向加速度为例，当 样本数量为50时，与PEM-SGI计算结果的最大误差为30.2％；当Monte Carlo样本达到500 个时，则相对最大误差为11.3％。此外，对于桥梁横向加速度而言，随着样本数量由50增加 到500，Monte Carlo算法与PEM-SGI所得计算结果的误差由34.2％降低到9.65％。随着Monte Carlo样本数的增加，Monte Carlo法所得系统响应结果曲线逐渐趋于PEM-SGI所得计算结果， 吻合度较好。采用PEM-SGI所需的样本数为73个，而采用传统PEM则需要246个样本才能 够达到相同精度。这也说明了，PEM-SGI所花费的计算时间仅为梯形积分的29.7％，效率提 高了将近4倍，节约了426.7分钟。PEM-SGI能够在保证更高精度的情况下，较之于Monte Carlo 方法，计算效率提高了近7倍，节约了1053.3分钟。

Claims (1)

1.一种高效求解高速列车过桥振动随机特性的分析方法，其特征在于，包括以下步骤：

步骤1：建立车-桥耦合系统，有：

其中，m_c为车体质量，y^c为车体位移，y为梁体位移，r为轨道不平顺性，K和c分别为结构刚度和阻尼；

对于车辆运行过程中的非平稳随机振动激励，用平稳过程乘以时变函数表示的非平稳过程表示为：

x′(t)＝g(t)·x(t)

其中，x(t)为平稳过程，其功率谱密度函数为S_xx(ω)，g(t)为时变调制函数；

针对非平稳随机振动，构造虚拟激励：

其中，e^iωt表示所构造的简谐函数；

对所述虚拟激励进行时程积分得到非平稳虚拟响应y′^*(t)，则响应的功率谱密度函数为：

其中，表示为积分所得的非平稳虚拟响应；

步骤2：基于虚拟激励算法的基本理论，将车辆非平稳随机振动激励作用下，产生的响应y(t)的方差表示为：

其中，小标y为响应结果；

假设公式的积分上下限截止频率分别为ω_l和ω_u，将区间[ω_l,ω_u]粗略等分为N个子区间j＝1～N；

步骤3：基于相应结果，用高斯积分公式计算各子区间内的积分

步骤4：计算各子区间积分之和

步骤5：判断如果表示第j个子区间收敛，则将该子区间作为需要考察的区间段，并将该子区间再细分为N个子区间段；否则，该子区间未收敛，不作细分；ε为预给定的常数，即两次计算的误差值收敛值；

步骤6：针对未收敛的子区间令P＝N，i＝1～P，则将子区间划分为P个子区间段[ω_i,ω_i+1]，则有计算；

步骤7：基于之前所给的误差收敛值ε，判断如果则需要对第i个子区间段再细分，以此类推；

步骤8：如前后两次计算的结果已经收敛，则退出计算；反之，返回第步骤6继续计算；

步骤9：当所有的计算停止，所有的收敛结果相加，求得车辆非平稳随机振动激励作用下相应的方差的最终结果