基于相位调制的平稳非高斯桥面不平顺随机激励的模拟

2021-10-11陈水生李锦华夏钰桓

振动与冲击 2021年18期

陈水生，赵辉，，李锦华，夏钰桓

（1.华东交通大学土木建筑学院，南昌 330013；2.湖北恒大建设工程有限公司，湖北武穴 435400）

车辆在桥上行使时，高低不平的桥面状况不仅影响车辆行驶安全和乘坐舒适性，也是桥梁产生随机振动的主要激励源。在空间上处于车轮与桥梁间的桥面高低起伏状况具有很大的随机性，对车辆和桥梁的振动响应影响较大[1-2]。特别是最近十年，中国桥梁以每年3万多座的速度递增，目前我国公路桥梁数量达83.25万座，全国桥梁总数达100万座，已成为世界第一桥梁大国[3]，庞大的桥梁数目，桥面养护需求越来越大，但是，桥面养护的技术、材料、资金和人员不足，导致桥面养护工作相对滞后，并且，大量重载、超载车辆在桥梁上行驶，加速了桥面的老化和破损，很多桥梁的桥面出现了不同程度的病害，降低了桥面等级，桥面路况的恶化又加重了车辆对桥梁的冲击作用[4]。鉴于此，如何建立有效的、贴近实际的桥面不平顺随机激励输入模型，对研究车桥耦合系统的振动规律至关重要。

多年来，国内外学者对车桥耦合振动的研究从未间断，做了很多工作，由起初的单车荷载到目前的多车荷载，创造了很多有价值的成果，但有关桥面不平顺随机激励的数值模拟，学者们大多采用滤波白噪声法、谐波叠加法和逆傅里叶变换法等方法，把桥面不平顺看成是各态历经的平稳高斯随机过程[5-10]。然而，伴随着对各类随机过程的研究和认知的不断深入，人们认识到，很多自然现象和工程实际并非高斯随机过程，不仅建筑结构上的风压、土木工程与机械工程中的振动、水体的波浪、声波等具有非高斯特征，桥面不平顺也具有明显的非高斯特征。特别是车辆在桥上行驶时，其行驶速度不可能也无法保证严格匀速，行驶速度多变，其在桥上加速、减速或突然制动所引起的桥面不平顺激励已不再是平稳高斯过程，明显具有非高斯特征，但有关桥面不平顺非高斯随机激励的研究却很少。诚然，人们也在其他领域努力开发模拟非高斯过程的方法，非高斯过程的模拟是一个更具挑战性和实际重要性的问题。在过去的二十年里，平稳非高斯随机过程的生成算法也相继问世[11]，虽然生成的方式各种各样，但归纳起来有代表性的方法有二类:第一类，根据目标的特征统计参数（例如均值、方差、偏度与峰度）和目标功率谱密度（power spectral density，PSD）函数，先产生平稳高斯过程，利用非线性转换关系间接生成满足目标要求的平稳非高斯过程[12-15]，但该类方法在同时逼近目标PSD和目标统计特征参数时，两者之间很难协调；第二类，根据目标边际概率密度函数（probability density function，PDF）和目标PSD函数多次迭代产生平稳非高斯随机过程[16-18]，但该类方法目标PSD与目标PDF会出现不相容，需要反复迭代，计算量大，也可能出现不收敛的情况。在上述两类平稳非高斯生成方法的基础上，本文提出基于相位调制的平稳非高斯桥面不平顺随机激励生成方法，研究平稳非高斯桥面激励的特征，对比分析平稳非高斯桥面激励对车桥耦合系统的影响。

1 平稳非高斯随机过程生成理论

1.1 平稳高斯随机过程

采用IFFT的方法生成平稳高斯随机过程时，桥面不平顺激励的自功率谱必须转换成激励信号的幅值谱，并对幅值谱进行相位调制，然后才能进行IFFT（inverse fast Fourier transform）变换生成离散桥面激励信号。为研究方便，桥面不平度用q（n）表示，其幅值谱序列为Q（f），根据双边自谱的定义式

进而可以得到单边谱Gq（f）

在此设Gq（k）为桥面不平度激励数据的单边谱在频率点fk处的离散值，可以表示为

根据离散傅里叶变换公式

将式（4）代入式（1）可得

式（6）的离散表达式可以写成

结合式（2）、式（3）和式（6）推得

在给定桥面路谱的条件下，由式（8）可以计算桥面激励幅值谱离散值的模

由于信号的自功率谱不包含相位角信息，如果桥面激励幅值谱的模相同，相位角 ϕ（k）不同，经过IFFT变换，就可以得到平稳高斯桥面激励

1.2 平稳非高斯随机过程

高斯随机过程的高阶中心矩均为零，可以用功率谱密度来充分描述，而非高斯随机过程的高阶中心矩并非全都为零，通常采用三阶中心矩偏斜度S和四阶中心矩峭度K两个特征值来表示

如果随机过程的周期为T，则可以通过在T周期内对式（13）进行积分来近似前四阶中心矩，可得

进一步可求得

将式（15）、式（16）和式（17）代入式（11）和式（12），可得[19]

式（20）和式（21）中的相位角 ϕj，ϕk，ϕm，ϕn在[0，2π]内均匀分布，在保持所有幅值Aj，Ak，Am，An不变和部分相位角不变的情况下，只需要对少数相位角进行调制，就可以使两式中的某些项求和变大，从而使偏斜度和峭度增加或减小。

2 平稳非高斯桥面不平顺随机激励的生成算法

根据上文建立的峭度、偏度与相位角的关系，可以通过调制式（21）第二项的相位角ϕj，ϕk来改变峭度的大小，首先找出生成高斯激励的相位角序列中所有满足j=3k的多组相位角{ϕj1，ϕk1；ϕj2，ϕk2；…}，再对第一组相位角{ϕj1，ϕk1}中的ϕj1进行调制，用新的代替ϕj1，=3ϕk1，并采用逆傅里叶变换法生成桥面不平顺激励，最后计算生成信号的峭度值K1，若K1与目标峭度值的误差满足精度要求，则调制结束；若K1小于目标峭度值且误差不满足要求，则调制下一组相位角{ϕj2，ϕk2}中的 ϕj2，取新值=3ϕk2；若K1大于目标峭度值且误差不满足要求，则调制下一组相位角{ϕj2，ϕk2}中的 ϕj2，取新值=3ϕk2+π；再次计算峭度值K2进行判断，反复按照上述迭代方法进行调整，就可以生成满足目标峭度要求的平稳非高斯桥面不平顺随机激励。相同的方法，也可以对式（20）中的某些求和项进行不断调整，就可以生成满足目标偏斜度要求的平稳非高斯桥面不平顺随机激励。如果同时对峭度和偏斜度进行调整，就可以生成既满足目标偏斜度要求，又满足目标峭度要求的平稳非高斯桥面激励，其算法的流程图如图1所示。当然，为了降低同时调整峭度和偏斜度时相互干扰，峭度调整时只调制偶数位置的相位角，偏斜度调整时只调制奇数位置的相位角。

图1 平稳非高斯桥面不平顺生成算法流程图Fig.1 Flow chart of algorithm for generating stationary non-Gaussian deck irregularity

3 平稳非高斯桥面不平顺随机激励的数值模拟

通常情况下，假定车轮与桥面紧密接触，桥面不平度这一随机激励源将由两部分叠加组合而成，桥面本身因施工、材料、设计等因素引起的不平度和桥梁在车辆荷载作用下变形引起的附加不平度。为了验证相位调制的方法生成平稳非高斯桥面不平顺随机激励的有效性，此部分数值模拟的桥面不平度不包括桥梁变形的附加不平度，桥面本身存在的不平顺功率谱密度函数可以表示为

式中:当车辆以速度v匀速行驶时，桥面激励输入无论在频域还是在时域都是平稳随机过程；n0为空间参考频率；Gq（n0）为不同等级桥面的不平度系数；f为时间频率，当车辆行驶速度在36～108 km/h时，就可以保证时间频率的区间是0.33～28.3 Hz，这个时间频率范围包括了汽车振动的主要固有频率[20]。

3.1 亚高斯桥面不平顺随机激励

以车辆行驶速度10 m/s为例，从式（22）可以分别得到不同等级桥面的PSD理论值，在此设定目标PSD为C级桥面功率谱密度、目标峭度Kp=2和目标偏度Sp=0，模拟精度控制在±0.1范围之内，车辆行驶距离60 m，采样频率1 000 Hz，采用相位角调制的方法生成亚高斯桥面不平顺随机过程，生成信号的峭度K=2.09、偏度S=0.07，与目标设定值非常接近，生成误差都小于0.1。从图2的功率谱密度对比可以看出，模拟得到的高斯和亚高斯桥面不平顺随机激励PSD与理论PSD吻合非常好；高斯和亚高斯桥面激励的均值都为零，方差分别为15.18 mm，14.86 mm，充分说明调制相位角不会改变功率谱密度；图3给出了模拟得到的高斯和亚高斯桥面激励，可以看出，亚高斯桥面激励的幅值小于高斯桥面激励，且亚高斯桥面激励完全在±3σ区域内，而高斯桥面激励部分时间点的幅值超过±3σ，σ为C级桥面不平顺的标准差，σ=15.23 mm；从图4的概率密度曲线也可以看出，高斯桥面激励有拖尾，而亚高斯桥面激励几乎无拖尾，幅值分布相对集中。

图2 高斯和亚高斯桥面不平顺功率谱密度Fig.2 Power spectral density of Gaussian and sub-Gaussian deck irregularity

图3 高斯和亚高斯桥面不平顺Fig.3 Gaussian and sub-Gaussian deck irregularity

图4 高斯与和亚高斯桥面不平顺概率密度Fig.4 Probability density of deck irregularity of Gaussian and sub-Gaussian

3.2 超高斯桥面不平顺随机激励

设定车辆行驶速度10 m/s，目标PSD为C级桥面的功率谱密度，目标峭度Kp=6和目标偏斜度Sp=0，模拟精度控制在±0.1范围之内，采用本文相位调制方法生成的超高斯桥面不平顺随机激励峭度K=5.91、偏度S=-0.02，满足设定精度要求。从图5可知，生成的超高斯桥面不平顺PSD与目标PSD吻合很好；虽然生成的高斯和超高斯桥面不平顺激励的PSD相同，方差值也很接近（分别为15.48 mm和15.6 mm），但从图6可知，超高斯桥面激励超出±3σ区域的幅值增多，超出部分的最大值达到77.7 mm，而高斯激励只有个别时间点的幅值超出±3σ区域；从图7的概率密度曲线也可以明显看出，超高斯桥面不平顺有长长的拖尾。

图5 高斯和超高斯桥面不平顺功率谱密度Fig.5 Power spectral density of Gaussian and super-Gaussian deck irregularity

图6 高斯和超高斯桥面不平顺Fig.6 Gaussian and super-Gaussian deck irregularity

图7 高斯与和亚高斯桥面不平顺概率密度Fig.7 Probability density of deck irregularity of Gaussian and super-Gaussian

相同的C级桥面PSD，相同行车速度10 m/s和相同偏斜度值S=0，图8给出了不同峭度K的超高斯桥面不平顺随机激励，可以看出，随着峭度值的增大，超高斯桥面不平顺随机激励幅值的最大值就越大。

图8 不同峭度的超高斯桥面不平顺Fig.8 Super-Gaussian deck irregularity with different kurtosis

图9为C级桥面PSD，车辆行驶速度10 m/s时，采用相位调制的方法生成的同一峭度值K=7，不同偏斜度值S的超高斯桥面不平顺随机激励，从图中可以看出，当偏斜度值为正值时，桥面激励的峰值分布偏向于桥面不平度的正值，即桥面凸起点多；当偏斜度值为负值时，桥面激励的峰值分布偏向于桥面不平度的负值，即桥面凹陷点多；且偏斜度值越大，桥面激励的峰值分布偏的就越多。

图9 不同偏斜度的超高斯桥面不平顺Fig.9 Super-Gaussian bridge deck irregularity with different skewness

4 平稳非高斯桥面激励的重构

为了更好的认识桥面不平顺的特征，课题组采用LXBP-5型公路连续式八轮平整度仪对江西上饶鄱余公路的乔木湾乐安河大桥进行了桥面不平度测试。在桥梁的两个行车方向上采集多组数据，因为车辆行驶速度凭借观察测速仪显示屏来控制，车速控制较难，故选用其中外场环境最理想和车速控制最好的一组数据，车辆以2 m/s速度匀速行驶，采样频率50 Hz，采样点数为3 000，如图10所示。通过计算该桥面不平度随机激励的三阶和四阶中心矩可知，其左右轮迹桥面不平度的峭度分别为3.72，3.67，偏斜度分别为-0.37，-0.29，实测结果表明，在车辆匀速行驶的理想外场环境下，桥面不平顺也表现出了较强的非高斯特征，可见，将桥面不平顺视为平稳高斯随机过程就难以吻合实际桥面路况。

图10 实测桥面不平顺Fig.10 Actual measurement of bridge deck irregularity

根据实测桥面不平度数据的PSD，采用逆傅里叶变换法生成左轮迹的高斯激励，并通过不断调制相位角重构桥面不平顺非高斯激励，生成信号如图11所示，重构非高斯信号的峭度值为3.719 4，偏斜度值为-0.388 3，与实测桥面不平度很接近。图12为重构高斯和非高斯桥面不平度的功率谱密度和现场实测不平度的功率谱密度对比，重构信号的PSD与实测数据的PSD吻合的非常好，保证了方差的一致性。图13为重构桥面不平度的概率密度和现场实测数据的概率密度对比，由图可以看出，重构非高斯桥面不平度与实测值的PDF非常接近，拖尾长度一致，比高斯信号表现出更明显的非高斯特征。

图11 重构的桥面不平度数据Fig.11 Reconstructed bridge deck roughness data

图12 功率谱密度曲线Fig.12 Power spectral density curve

图13 概率密度曲线Fig.13 Probability density curve

同时，课题组也对多座在役桥梁的桥面平整度进行了测试，根据实测桥梁桥面路况的恶化程度，采用本文相位调制的方法分别进行了重构，如图14所示，从图可以看出，桥面路况恶化越严重，其峭度值就越大，且偏斜度值偏向于负值一侧，与桥面凹陷较多的实际情况相符；基于相位调制生成的桥面不平度与实测不平度吻合较好，桥面不平度的最大幅值也接近。

图14 桥面不平度模拟值与实测值对比Fig.14 Comparison of simulated and measured values of deck unevenness

可见，基于相位调制的平稳非高斯桥面不平顺生成方法有效可靠，与现在大量使用的平稳高斯桥面不平顺生成方法相比，采用本文方法生成的平稳非高斯桥面不平顺更切合实际桥面路况，实用性更强，计算方便且精度满足要求。生成的平稳非高斯路况比传统平稳高斯路况可以在保持功率谱密度不变的同时，提供更强烈的桥面激励，更好的再现了桥面路况等级降低和路况恶化这一实际情况，为研究行车舒适度和桥梁振动提供了更好的桥面激励生成方法，实用性更广泛。

5 平稳非高斯桥面激励对车桥耦合系统的影响

为研究平稳非高斯桥面激励对车桥耦合系统的影响，以图15所示的车桥耦合模型为研究对象，桥梁简化为一30 m跨径的简支梁桥，移动车辆简化为两自由度的弹簧质量阻尼系统，y（x，t），z1，z2分别为梁的动挠度、车体的动位移和悬架动位移。其中，车体质量为m1，悬架系统质量为m2，悬架系统的弹簧刚度为k1、阻尼为c1，轮胎刚度为k2、阻尼为c2，各参数如表1所示。

图15 车桥耦合模型Fig.15 Coupling model of vehicle and bridge

表1 桥梁结构参数及车辆参数Tab.1 Bridge structure parameters and vehicle parameters

5.1 车桥耦合振动方程

依据车辆系统振动原理，建立两自由度车辆振动方程

式中:Mv，Cv，Kv分别为车辆的质量矩阵、阻尼矩阵和刚度矩阵；为车辆振动过程中作用于车辆各自由度的惯性荷载列向量；z=[z1z2]T为车辆自由度列向量。

为减少车桥耦合计算矩阵维数，采用广义坐标离散的方法建立桥梁的振动方程

桥梁模型引入模态综合法，把相互耦联的节点运动方程解耦，转化为互相独立的模态方程，式（24）可以变为

考虑平稳非高斯桥面不平顺随机激励的影响，车轮作用在桥梁上的荷载可以写成

式中:N为车轮作用处位移场的插值函数；Nx为N对车辆行驶方向x坐标的导数；q为车轮的平稳非高斯桥面不平度随机激励值。

将式（26）代入式（25），然后联立式（23）可得车桥耦合振动方程

式中:Mbv，Cbv，Kbv分别为车桥相互作用模型的质量、阻尼和刚度矩阵；Fg为车辆自重引起的车桥耦合竖向荷载向量；Fw为桥面不平顺随机激励引起的车桥耦合竖向荷载向量；u为广义坐标向量，u=[y1y2…ynz1z2]T。将前文模拟得到的平稳非高斯随机激励样本q输入到式（29），通过编制基于Newmark-β算法的MATLAB程序求解车桥系统的振动响应。

5.2 车桥振动响应分析

5.2.1 确定性桥面不平顺激励

车辆以10 m/s的速度行驶在C级桥面上，分别生成亚高斯、高斯和超高斯桥面激励样本，样本偏斜度S=0，峭度值分别为K=2，K=3，K=4。图16给出了各激励样本作用下的桥梁跨中竖向振动和车体竖向振动时程曲线，由图可得，桥梁和车辆在超高斯桥面激励下的振动响应大于平稳高斯激励和亚高斯激励所产生的振动响应，且曲线波动大于平稳高斯激励和亚高斯激励。图17为车桥耦合系统振动响应最大值随桥面激励峭度的变化关系，由图可得，车桥振动响应随桥面激励峭度的增大波动较大，两者并非成线性比例关系；不同的峭度值，桥梁振动响应最大差值为0.82 mm，车辆振动响应最大差值为5.79 mm，说明峭度的增大对车辆振动响应的影响较大，对桥梁振动响应的影响并不显著。

图16 桥梁和车辆的振动响应Fig.16 Vibration response of bridge and vehicle

图17 峭度与车桥振动响应最大值的关系Fig.17 Relationship between kurtosis and maximum value of vehicle bridge vibration response

设定峭度值K=5，不同桥面激励偏斜度值与车桥振动响应最大值的关系，如图18所示。从图可以看出，桥梁振动响应最大差值为0.75 mm，车辆振动响应最大差值为1.25mm，偏斜度的改变对车辆振动响应的影响大于对桥梁振动响应的影响。

图18 偏斜度与车桥振动响应最大值的关系Fig.18 Relationship between skewness and maximum value of vehicle bridge vibration response

5.2.2 随机性桥面不平顺激励

为进一步应用随机振动理论对桥面不平顺随机激励下的车桥振动响应进行统计分析，采用传统的蒙特卡罗法对2 000个桥面随机激励样本进行统计，计算车桥耦合振动响应的均值和标准差，分析平稳非高斯桥面不平顺的随机性对车桥振动响应的影响。图19为车辆以10 m/s速度通过桥梁时，桥梁跨中竖向振动位移和车体竖向振动位移的均值，从图中可以看出:高斯激励与非高斯激励的桥梁竖向振动位移均值基本重合，充分说明桥梁的振动响应主要还是由车辆重力这一确定性荷载激励引起，平稳非高斯桥面激励对桥梁振动响应的均值影响很小。但是，车辆振动响应均值对桥面随机激励较敏感，平稳高斯激励的车体振动响应均值大于亚高斯激励的车体振动响应均值，小于超高斯激励的车体振动响应均值。