CN104715138A - 钢轨扣件反力计算方法和地铁振动频域响应计算方法 - Google Patents

Abstract

本发明提供了钢轨扣件反力计算方法和地铁振动频域响应计算方法，具体涉及地铁隧道环境振动预测领域。本发明所提供的钢轨扣件反力计算方法通过在预先获取的第一模型的基础上，增加整体道床与隧道叠合梁模型，以确定第二模型，其中，第一模型是车辆-轨道耦合模型，第二模型是车辆-轨道-隧道耦合模型。并在第二模型中计算隧道指定部分的动参数随频域振动位移的变化量，其中，动参数包括动阻尼和动刚度，还进一步根据动参数随频域振动位移的变化量，计算在第二模型中的不同工况下的钢轨扣件反力，由于在建立的第二模型中考虑了相关技术中并未考虑到的轨道系统内高分子减震材料的动参数，因而使计算得出的扣件力更为准确。

Description

钢轨扣件反力计算方法和地铁振动频域响应计算方法

技术领域

本发明涉及地铁隧道环境振动预测领域，具体而言，涉及钢轨扣件反力计算方法和地铁振动频域响应计算方法。

背景技术

地铁列车在轨道上行驶时,由于车轮偏心,车轮与道岔、钢轨的碰撞以及线路不平顺等原因,引起车轮的振动,经钢轨→扣件→轨枕→道床→隧道结构→围护地层传至地面及建筑物。

若将整个隧道视为振源,可近似作为半空间内埋深面波源,振动以横波、纵波、瑞利表面波合成的复杂波动形式,通过大地介质向外传播,近场的振动波型主要以弯曲波形式传播,远场主要以表面波形式传播。对于地铁振动传播,上述主要影响参数中以列车质量、运行速度、轮轨、扣件、道床条件、隧道结构等因素对地铁振动源特性影响较大。通过地铁隧道环境振动预测能够便于对地铁的工作进行理论分析和指导地铁进行调整。

在地铁隧道环境振动预测方面，目前主要有4种方法，分别是现场测试、模型试验、经验公式与理论解析。现场测试虽然可以获取有限时段、有限断面、有限范围的振动数据，并掌握其在部分区域内沿线路纵向、横向与深度方向的衰减规律，但是地铁全线大范围的测试成本极高，不具推广价值，而且仅适用于既有线路的环境振动评估；模型试验在实验室内进行，操作方便且条件可控，但是实验条件过于理想，结果可靠性与准确性备受争议；经验公式是建立在大量测试数据基础上的，且仅适用于特定地点的地质环境，然而不同地区的地质类型及其分布特征差异极大，所以该方法不具有良好的普适性；理论解析可分为时域分析与频域分析两大类，时域预测的时间成本太高，不利于地铁全线环境振动的快速预测，然而尽管频域预测的效率较高，但是现有众多的频域分析模型中所有材料模型均基于弹性假设，也就导致了现有方法的预测精度并不高(特 别是频域内的预测精度方面)、适用范围也非常有限且均无法实现地铁隧道全线环境振动的可视化管理与评估。

综上，现有的地铁隧道环境振动预测方面，难以客观的获取进行地铁隧道环境振动预测所需要的数据，如钢轨扣件反力。

发明内容

有鉴于此，本发明实施例的目的在于提供钢轨扣件反力计算方法，以提高钢轨扣件反力计算方法计算的准确度。

第一方面，本发明实施例提供了钢轨扣件反力计算方法，包括：

在预先获取的第一模型的基础上，增加整体道床与隧道叠合梁模型，以确定第二模型，第一模型是车辆-轨道耦合模型，第二模型是车辆-轨道-隧道耦合模型；

在第二模型中计算隧道指定部分的动参数随频域振动位移的变化量，其中，动参数包括动阻尼和动刚度；

根据动参数随频域振动位移的变化量，计算在第二模型中的不同工况下的钢轨扣件反力。

结合第一方面，本发明实施例提供了第一方面的第一种可能的实施方式，其中，第二模型为：

$\begin{array}{c}{\stackrel{\·\·}{q}}_k^{\mathrm{bt}}\left(t\right)+\underset{i=1}{\overset{N^{\mathrm{bt}}}{\mathrm{\Σ}}}{C}_{\mathrm{fi}}{Y}_k^{\mathrm{bt}}\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right){\stackrel{\·}{q}}_k^{\mathrm{bt}}\left(t\right)+\frac{E^{\mathrm{bt}}{I}^{\mathrm{bt}}}{m^{\mathrm{bt}}}{\left(\frac{\mathrm{k\π}}{l^{\mathrm{bt}}}\right)}^4{q}_k^{\mathrm{bt}}\left(t\right)\\ +\underset{i=1}{\overset{N^{\mathrm{bt}}}{\mathrm{\Σ}}}{K}_{\mathrm{fi}}{Y}_k^{\mathrm{bt}}\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right)q_k^{\mathrm{bt}}\left(t\right)=\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{F}_{\mathrm{rsj}}\left(t\right)Y_k^{\mathrm{bt}}\left(x_j\right)\end{array},$

$\begin{array}{c}{F}_{\mathrm{rsj}}\left(t\right)=C_{\mathrm{pi}}\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_i\right)\stackrel{\·}{q}\left(t\right)+K_{\mathrm{pi}}\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_i\right)q_h\left(t\right)\\ -C_{\mathrm{pi}}\underset{h=1}{\overset{N{m}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right){\stackrel{\·}{q}}_h^{\mathrm{bt}}\left(t\right)-K_{\mathrm{pi}}\underset{h=1}{\overset{N{M}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right)q_h^{\mathrm{bt}}\left(t\right)\end{array};$

其中，mbt为道床-隧道叠合梁单位长度的质量(kg/m)；Ebt为整体道床-盾构隧道叠合梁的弹性模量(Pa)，Ibt为道床-隧道梁的截面惯量(m4)；lbt为道床-隧道简支梁长度(m)；NMbt为道床-隧道叠合梁所截取的模态阶数；Nbt为道床-隧道梁下离散支点数；Ykbt(xi)为道床-隧道叠合梁的振型；qkbt(t)为道床-隧道叠合梁的正则振型坐标；Kfi、Cfi为隧道底部地基沿纵向每米的支承刚度(N/m)与阻尼(N·s/m)，Frsj(t)为钢轨扣件反力(N)，Kpi、Cpi 分别为轨下胶垫刚度(N/m)与阻尼(N·s/m)，NM为钢轨所截取的模态阶数，Yh(xi)与qh(t)分别是钢轨振型和钢轨正则振型坐标。

结合第一方面，本发明实施例提供了第一方面的第二种可能的实施方式，其中，在第二模型中计算隧道指定部分的动参数随频域振动位移的变化量包括：

根据第二模型得到第二模型的矩阵表达式，  $\left[M\right]\left\{\stackrel{\·\·}{q}\right\}+\left[C\right]\left\{\stackrel{\·}{q}\right\}+\left[K\right]\left\{q\right\}=\left\{Q\right\}=\left[K_f\right]\left[Z_0\right],$ 其中，[M]为系统的质量矩阵，[C]为系统的阻尼矩阵，[K]为系统的刚度矩阵，[Kf]为系统的转换矩阵，{Q}为系统广义力向量，{q}为系统广义位移向量，{Z0}为不平顺位移向量；

根据第二模型的矩阵表达式分别确定隧道指定部分随频率、幅值与温度变化的阻尼矩阵的传递函数和随频率、幅值与温度变化的刚度矩阵传递函数；

计算第二模型系统的不同自由度响应的功率密度谱。

结合第一方面，本发明实施例提供了第一方面的第三种可能的实施方式，其中，根据第二模型的矩阵表达式分别确定隧道指定部分随频率、幅值与温度变化的阻尼矩阵的传递函数和随频率、幅值与温度变化的刚度矩阵传递函数包括：

令{Z0}＝[I]e-jωt，[q]＝[H(ω)]e-jωt，得出第二模型系统动参数随频率、幅值与温度变化情况下的传递函数，(-ω²[M]+iω[C_ω,A,T]+[K_ω,A,T])[H(ω)]＝[K_f][I]，其中，[C_ω,A,T]为随频率、幅值与温度变化的阻尼矩阵，[K_ω,A,T]为随频率、幅值与温度变化的刚度矩阵，ω是不同振动频率对应的角速度，A是不同频率对应的幅值，T是环境温度；

计算第二模型系统的不同自由度响应的功率密度谱包括，根据如下公式计算第二模型系统在不同自由度响应的功率谱密度S_qi(ω)＝|H(ω)|²S₁₁(ω)，其中，Sqi(ω)为第二模型系统第i个自由度下的位移响应功率谱密度，为系统第i个自由度下的加速度响应功率谱密度，S11(ω)为系统第1位轮对处轨道不平顺功率谱密度。

结合第一方面，本发明实施例提供了第一方面的第四种可能的实施方式，其中，在得出第二模型系统动参数随频率、幅值与温度变化情况下的传递函数后还包括：

在第二模型系统中，根据隧道指定部分材料在不同频率下，振幅小于预设数值的动参数，计算各频率振幅，以确定参考振幅；

根据参考振幅重新选择对应的动参数，进一步求解材料应变，并判断不同频率下动变形与动参数是否收敛至某一定值时停止；

若否，则再次执行步骤在第二模型系统中，根据隧道指定部分材料在不同频率下，振幅小于预设数值的动参数，计算各频率振幅，以确定参考振幅。

结合第一方面，本发明实施例提供了第一方面的第五种可能的实施方式，其中，根据参考振幅重新选择对应的动参数，进一步求解材料应变包括：采用1/3倍频程来分频带地统计应变，再利用不同频带应变大小与动参数的对应关系，重新调整系统刚度与阻尼矩阵。

第二方面，本发明实施例还提供了地铁振动频域响应计算方法，基于第一方面的钢轨扣件反力计算方法，还包括：

在预设的图层中，将不同土层分界线和隧道边界线分别离散成由多个两两间隔预定距离的坐标点所组成的坐标点集合，以确定土层分界点集和隧道边界点集；

沿里程方向，逐层比较纵向土层分界点集中，和/或隧道边界点集中，相邻点的深度坐标变化幅度；

若相邻点的深度坐标变化幅度超过预定变化幅度数值，则根据变化幅度超过预定数值的点的横坐标为参考坐标，沿竖直方向分别获取土层分界线和隧道边界线上与参考坐标相同的位置的点的参数信息作为建模信息；

根据获取到的建模信息建立隧道与地层二维有限元模型，并求解不同预测断面内地铁隧道环境振动频域响应。

结合第二方面，本发明实施例提供了第二方面的第一种可能的实施方式，其中，根据获取到的建模信息建立隧道与地层二维有限元模型，并求解不同预测断面内地铁隧道环境振动频域响应包括：

根据建模信息中，预测断面内隧道与土层的位置坐标，划分代表不同土层或隧道层的多个区域；

将各区域划分为4节点4单元的矩形网格，并使用线弹性材料模拟隧道、道床结构及远场土体，使用非线性粘弹性材料模拟近场土层；

根据稳态动力学分析方法进行频域求解，即可直接得到不同预测断面内地铁隧道环境振动频域响应。

结合第二方面，本发明实施例提供了第二方面的第二种可能的实施方式，其中，预定变化幅度数值为8％-15％。

结合第二方面，本发明实施例提供了第二方面的第三种可能的实施方式，其中，预定变化幅度数值为10％。

本发明所提供的钢轨扣件反力计算方法，与相关技术中的计算钢轨扣件反力的方法相比，由于考虑了动参数随频域振动位移的变化量，因而得出的钢轨扣件反力更为准确和贴合实际。

为使本发明的上述目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附附图，作详细说明如下。

附图说明

为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对范围的限定，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。

图1示出了本发明实施例所提供的钢轨扣件反力计算方法的基本流程图；

图2a示出了本发明实施例所提供的钢轨扣件反力计算方法的剪切模量随剪应变变化示意图；

图2b示出了本发明实施例所提供的钢轨扣件反力计算方法的阻尼随剪应变变化示意图；

图3a示出了本发明实施例所提供的钢轨扣件反力计算方法的弹性模量与阻尼随温度变化示意图；

图3b示出了本发明实施例所提供的钢轨扣件反力计算方法的弹性模量与阻尼随频率变化示意图；

图4示出了本发明实施例所提供的钢轨扣件反力计算方法的车辆/整体道床/隧道/地基垂向耦合非线性随机振动分析模型示意图；

图5示出了本发明实施例所提供的钢轨扣件反力计算方法的道床与盾构隧道叠合梁的受力分析模型示意图；

图6示出了本发明实施例所提供的地铁振动频域响应计算方法的隧道与周边地层的纵断面及其预测横断面的自动划分、建模信息的自动提取示意图；

图7示出了本发明实施例所提供的地铁振动频域响应计算方法的指定的某个预测横断面的二维有限元模型示意图；

图8示出了本发明实施例所提供的地铁振动频域响应计算方法的车辆-轨道-隧道-地基垂向耦合的非线性随机振动分析模型示意图，其中左图是整体道床与盾构隧道复合梁模型示意图，右图是整体道床与盾构隧道质量块模型示意图；

图9示出了本发明实施例所提供的地铁振动频域响应计算方法的整体道床与盾构隧道复合梁的受力分析模型示意图；

图10示出了本发明实施例所提供的地铁振动频域响应计算方法的整体道床梁与盾构隧道质量弹簧的受力分析模型示意图；

图11a示出了本发明实施例所提供的地铁振动频域响应计算方法的隧道地基与轨下胶垫动参数随频率、幅值非线性变化前后扣件力的功率谱密度中，动阻尼不变动刚度随频率、幅值非线性变化的示意图；

图11b示出了本发明实施例所提供的地铁振动频域响应计算方法的隧道地基与轨下胶垫动参数随频率、幅值非线性变化前后扣件力的功率谱密度中，动刚度不变动阻尼随频率、幅值非线性变化的示意图；

图12示出了本发明实施例所提供的地铁振动频域响应计算方法的地铁隧道上方地表振级的俯视色彩云图(不含建筑)；

图13示出了本发明实施例所提供的地铁振动频域响应计算方法的地铁隧道上方地表振级的俯视等势线图(不含建筑)；

图14示出了本发明实施例所提供的地铁振动频域响应计算方法的地铁隧道上方地表振级的俯视色彩云图(含建筑)；

图15示出了本发明实施例所提供的地铁振动频域响应计算方法的地铁隧道上方地表振级的空间曲面色彩云图(不含建筑)。

具体实施方式

下面将结合本发明实施例中附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件可以以各种不同的配置来布置和设计。因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明的实施例，本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明实施例提供了钢轨扣件反力计算方法，如图1所示，包括如下步骤：

S101，在预先获取的第一模型的基础上，增加整体道床与隧道叠合梁模型，以确定第二模型，第一模型是车辆-轨道耦合模型，第二模型是车辆-轨道-隧道耦合模型；

S102，在第二模型中计算隧道指定部分的动参数随频域振动位移的变化量，其中，动参数包括动阻尼和动刚度；

S103，根据动参数随频域振动位移的变化量，计算在第二模型中的不同工况下的钢轨扣件反力。

需要说明的是，相关技术中的车辆-轨道耦合是已经存在的技术。具体的，一辆地铁车可被简化为具有两系悬挂的、10个自由度的多刚体系统，分别包括车体的点头与沉浮运动、前后两个转向架的点头与沉浮运动以及四个车轮的沉浮运动；钢轨被近似视为离散点支承的经典Euler梁模型；轮轨非线性接触弹簧可近似处理为线性轮轨接触弹簧。需要说明的是，这些是车辆-轨道线性随机振动频域分析法的基本假定，属于已有的内容。

在车-轨-隧-地基耦合动力学频域分析方面，以往计算模型所有材料的参数均为常量，并未考虑关键材料参数的动力非线性特征，特别是隧道地基与轨道高分子减振材料的动力非线性特征。然而在实际情况下，它们的动力学参数(主要是动刚度与动阻尼)是随激励振幅的变化而非线性变化的(我国常见地基土的动模量、动阻尼与应变的非线性关系曲线如图2所示)，而且轨道系统内高分子减振材料的动参数还与激振频率、环境温度密切相关，如图3所示。

为了考虑这些因素，借鉴非线性动力问题的等效线性处理方法(即分段线性逼近法)，提出了车辆-轨道-隧道耦合模型(第二模型)。下面以目前地铁领域最为常用的长枕埋入式整体道床与盾构隧 道为例，详细阐述该模型的建立过程及其算法的实现过程，其理论模型如图4所示。

本申请核心技术主要是两个方面，分别是模型与算法的改进。模型的改进是指在既有车辆-轨道耦合模型(第一模型)下方增加了整体道床与隧道叠合梁模型，考虑了隧道与地基对轮轨激振荷载(在本专利中是指钢轨之下的扣件反力)的影响。算法的改进是指在创新建立的车辆-轨道-隧道-地基垂向耦合动力学频域分析模型中首次实现了动力非线性参数的等效线性算法，即非线性问题的分段线性逼近法。

模型的改进部分：道床与盾构隧道可近似被看作连续弹性离散点支承的有限长Euler叠合梁模型，如图5所示。其中，F_rsi为钢轨扣件反力(N)；Z_bt为整体道床与盾构隧道叠合梁的振动位移变量(m)。与钢轨振动微分方程类似，道床-隧道叠合梁振动微分方程可推导得出：

$\begin{array}{c}{\stackrel{\·\·}{q}}_k^{\mathrm{bt}}\left(t\right)+\underset{i=1}{\overset{N^{\mathrm{bt}}}{\mathrm{\Σ}}}{C}_{\mathrm{fi}}{Y}_k^{\mathrm{bt}}\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right){\stackrel{\·}{q}}_k^{\mathrm{bt}}\left(t\right)+\frac{E^{\mathrm{bt}}{I}^{\mathrm{bt}}}{M^{\mathrm{bt}}}{\left(\frac{\mathrm{k\π}}{l^{\mathrm{bt}}}\right)}^4{q}_k^{\mathrm{bt}}\left(t\right)\\ +\underset{i=1}{\overset{N^{\mathrm{bt}}}{\mathrm{\Σ}}}{K}_{\mathrm{fi}}{Y}_k^{\mathrm{bt}}\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right)q_k^{\mathrm{bt}}\left(t\right)=\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{F}_{\mathrm{rsj}}\left(t\right)Y_k^{\mathrm{bt}}\left(x_j\right)\end{array}---\left(1\right)$

$\begin{array}{c}{F}_{\mathrm{rsj}}\left(t\right)=C_{\mathrm{pi}}\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_i\right)\stackrel{\·}{q}\left(t\right)+K_{\mathrm{pi}}\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_i\right)q_h\left(t\right)\\ -C_{\mathrm{pi}}\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right){\stackrel{\·}{q}}_h^{\mathrm{bt}}\left(t\right)-K_{\mathrm{pi}}\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right)q_h^{\mathrm{bt}}\left(t\right)\end{array}---\left(2\right)$

在公式1中，m^bt为道床-隧道叠合梁单位长度的质量(kg/m)；Ebt为整体道床-盾构隧道叠合梁的弹性模量(Pa)，由轨枕、道床和盾构管片强度计算；I^bt为道床-隧道梁的截面惯量(m4)；l^bt为道 床-隧道简支梁长度(m)；NM^bt为道床-隧道叠合梁所截取的模态阶数；Nbt为道床-隧道梁下离散支点数；Y_k ^bt(xi)为道床-隧道叠合梁的振型；q_k ^bt(t)为道床-隧道叠合梁的正则振型坐标；K_fi、C_fi为隧道底部地基沿纵向每米的支承刚度(N/m)与阻尼(N·s/m)，其大小随不同频域振动位移非线性变化，即因荷载频域响应不同，隧道底部不同频域对应的地基支承刚度与阻尼也不相同。F_rsj(t)为扣件力(N)，由公式2计算。在公式2中，K_pi、C_pi分别为轨下胶垫刚度(N/m)与阻尼(N·s/m)，这些参数在本模型中可考虑其随频率、幅值与温度变化。NM为钢轨所截取的模态阶数，Yh(xi)与qh(t)分别是钢轨振型和钢轨正则振型坐标。

算法的改进部分：基于上一段所描述的模型，借鉴非线性问题的等效线性思想(即分段线性逼近原则)，可方便实现非线性动参数随频率、幅值、温度的非线性变化特征。计算过程如下：

第1步，根据车辆模型、轨道模型、新增的道床与隧道叠合梁模型以及轮轨线性耦合关系，即可得到车辆-轨道耦合系统的振动微分方程的矩阵表达式：

$\left[M\right]\left\{\stackrel{\·\·}{q}\right\}+\left[C\right]\left\{\stackrel{\·}{q}\right\}+\left[K\right]\left\{q\right\}=\left\{Q\right\}=\left[K_f\right]\left[Z_0\right]---\left(3\right)$

在公式3中，[M]为系统的质量矩阵，[C]为系统的阻尼矩阵，[K]为系统的刚度矩阵，[K_f]为系统的转换矩阵，{Q}为系统广义力向量，{q}为系统广义位移向量，{Z₀}为不平顺位移向量。

令{Z₀}＝[I]e^-jωt，则得出[q]＝[H(ω)]e^-jωt，再将其代入系统二阶线性微分方程(3)可求解系统振动传递函数矩阵[H(ω)]：

(-ω²[M]+iω[C]+[K])[H(ω)]＝[K_f][I]       (4) 

不难看出，公式4可计算轮轨系统动参数随频率、幅值与温度变化情况下的传递函数，见公式5。

(-ω²[M]+iω[C_ω,A,T]+[K_ω,A,T])[H(ω)]＝[K_f][I]      (5) 

在公式5中，[C_ω,A,T]为随频率、幅值与温度变化的阻尼矩阵，[K_ω,A,T]为随频率、幅值与温度变化的刚度矩阵，ω是不同振动频率对应的角速度，A是不同频率对应的幅值，T是环境温度。由此不难看出，非线性动参数随温度变化比较好实现，但是其随频率与幅值变化则需进行循环迭代运算，以保证其能够收敛。

又由线性随机振动理论可知，利用上面计算的传递函数，即可获得系统不同自由度响应的功率谱密度：

S_qi(ω)＝|H(ω)|²S₁₁(ω)          (6) 

${\stackrel{\·\·}{S}}_{\mathrm{qi}}\left(\mathrm{\ω}\right)={\mathrm{\ω}}^4{S}_{\mathrm{qi}}\left(\mathrm{\ω}\right)---\left(7\right)$

在公式6和7中，S_qi(ω)为系统第i个自由度下的位移响应功率谱密度，为系统第i个自由度下的加速度响应功率谱密度，S₁₁(ω)为系统第1位轮对处轨道不平顺功率谱密度(即为输入激励)。

其中，为实现上面公式5中非线性动参数随频率与幅值的循环迭代运算，首先依据材料不同频率下较小振幅(地基应变下限为10^-6，高分子减振材料应变下限为10^-3)对应的动参数计算各频率振幅，然后根据新计算的振幅大小重新选择对应的动参数，再一次求解材料应变。如此周而复始，直到不同频率下动变形与动参数收敛至某一定值时停止。

需要说明的是，输入的是轨道不平顺功率密度谱，输出的则是系统振动随机响应功率密度谱，但是频域内幅值大小与采样频率、 采样点数有关，因此其幅值大小没有实际意义，不过不同频带内的功率(也就是能量)却是不变的，它等于功率密度谱与横坐标响应频率围成的面积。基于以上分析，一方面考虑到振动系统中不同频带内的能量守恒，另一方面方便结合我国现有规范实现分频振动级的评价。本申请所提供的钢轨扣件反力计算方法中，采用1/3倍频程来分频带地统计应变能，再利用不同频带应变能大小与动参数的对应关系，重新调整系统刚度与阻尼矩阵，继而在随机振动响应算法中近似反映材料动参数随频率与幅值的非线性变化特征。

本申请所提供的车辆-轨道-隧道耦合模型(第二模型)不仅可以考虑材料参数的非线性问题，使计算结果更符合实际情况，而且由于该算法属于理论解析法，因此计算效率也比数值积分的时域求解法快，5秒左右即可完成一个计算工况，便于进行地铁隧道全线各类轨道、各类隧道与各类地基条件下地铁扣件系统的频域荷载谱(参见图11a和图11b，其中图11a中的实线为动阻尼为常量，动刚度为常量的结果，虚线为动阻尼为常量，动刚度随频率与复制变化的结果；图11b中，实线为动刚度为常量，动阻尼为常量的结果，虚线为动刚度为常量，动阻尼随频率与幅值变化的结果)的快速预测。从图中可清晰地看出，隧道地基与扣件胶垫动参数随频率、幅值变化前后地铁扣件系统频域荷载谱的差异，可充分说明在车辆-轨道-隧道-地基耦合动力系统中考虑这些因素的必要性。

具体的，1、车辆-轨道-隧道-地基耦合动力学的频域计算模型中，

基于车辆-轨道耦合动力学原理(第一模型的原理)，通过推导了车辆-轨道-道床与隧道叠合梁-地基垂向耦合非线性随机振动分析模型(简称“第一类模型”)。然后又在这个基础上，发展出了车 辆-轨道-道床与隧道复合梁-地基垂向耦合非线性随机振动分析模型(简称“第二类模型”)以及车辆-轨道-整体道床梁-管片环质量块-地基垂向耦合非线性随机振动分析模型(简称“第三类模型”)。这些模型的算法相同，只是模型结构逐步复杂且更接近实际情况。图4展示了车-轨-隧-地基垂向耦合非线性随机振动预测模型的“第一类模型”，图8(左部)为“第二类模型”，图8(右部)是“第三类模型”。

其中，1，第二类模型的振动微分方程

在第二类模型中，整体道床与盾构隧道可被近似看作是两个连续弹性离散点支承的有限长Euler梁，如图9所示。其中Zb为整体道床梁的振动位移变量(m)，Z_t为盾构隧道叠合梁的振动位移变量(m)。与第一类模型相比，第二类模型将整体道床与盾构隧道分开建模，可用于分析整体道床与盾构隧道的振动能量分担比及其它们之间的相互作用关系，比如盾构隧道由于下沉导致道床脱空的几何非线性问题。

由第一类模型的振动微分方程，不难得出第二类整体道床与盾构隧道复合梁模型的振动微分方程，如公式8和公式9所示。为了便于理解公式中参数的含义，需要说明的是，这些参数的物理意义相同，所不同的是其所指的对象不同，即参数的右上角标，参数右上角标“b”代表道床，参数右上角标“t”代表隧道。

$\begin{array}{c}{\stackrel{\·}{q}}_k^b\left(t\right)+\underset{i=1}{\overset{N^b}{\mathrm{\Σ}}}{C}_{\mathrm{si}}{Y}_k^b\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_i\right){\stackrel{\·}{q}}_k^b\left(t\right)+\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{C}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_j\right){\stackrel{\·}{q}}_h^b\left(t\right)+\frac{E^b{I}^b}{m_r^b}{\left(\frac{\mathrm{k\π}}{l^b}\right)}^4{q}_k^b\left(t\right)\\ +\underset{i=1}{\overset{N^b}{\mathrm{\Σ}}}{K}_{\mathrm{si}}{Y}_k^b\left(x_i\right)\underset{h=1}{\overset{N{M}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_i\right)q_k^b\left(t\right)\stackrel{}{+\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{K}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{N{M}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_j\right)q_h^b\left(t\right)}\\ -\underset{i=1}{\overset{N_b}{\mathrm{\Σ}}}{C}_{\mathrm{si}}{Y}_k^b\left(x_i\right)\underset{h=1}{\overset{N{M}_t}{\mathrm{\Σ}}}{Y}_h^t\left(x_i\right){\stackrel{\·}{q}}_h^t\left(t\right)-\underset{i=1}{\overset{N_b}{\mathrm{\Σ}}}{K}_{\mathrm{si}}{Y}_k^b\left(x_i\right)\underset{h=1}{\overset{N{M}^t}{\mathrm{\Σ}}}{Y}_h^t\left(x_i\right)q_h^t\left(t\right)\\ =\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{C}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_j\right)\stackrel{\·}{q}\left(t\right)+\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{K}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_j\right)q_h\left(t\right)\end{array}---\left(8\right)$

盾构隧道梁振动微分方程：

$\begin{array}{c}{\stackrel{\·}{q}}_k^t\left(t\right)+\underset{i=1}{\overset{N^t}{\mathrm{\Σ}}}{C}_{\mathrm{fi}}{Y}_k^t\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^t}{\mathrm{\Σ}}}{Y}_h^t\left(x_i\right){\stackrel{\·}{q}}_k^t\left(t\right)+\underset{j=1}{\overset{N^b}{\mathrm{\Σ}}}{C}_{\mathrm{sj}}{Y}_k^t\left(x_j\right)\underset{h=1}{\overset{{\mathrm{NM}}^t}{\mathrm{\Σ}}}{Y}_h^t\left(x_j\right){\stackrel{\·}{q}}_h^t\left(t\right)+\frac{E^t{I}^t}{m_r^t}{\left(\frac{\mathrm{k\π}}{l^t}\right)}^4{q}_k^t\left(t\right)\\ +\underset{i=1}{\overset{N^t}{\mathrm{\Σ}}}{K}_{\mathrm{fi}}{Y}_k^t\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^t}{\mathrm{\Σ}}}{Y}_h^t\left(x_i\right)q_k^t\left(t\right)+\underset{j=1}{\overset{N^b}{\mathrm{\Σ}}}{K}_{\mathrm{sj}}{Y}_k^t\left(x_j\right)\underset{h=1}{\overset{N{M}^t}{\mathrm{\Σ}}}{Y}_h^t\left(x_j\right)q_h^t\left(t\right)\\ =\underset{j=1}{\overset{N^b}{\mathrm{\Σ}}}{C}_{\mathrm{sj}}{Y}_k^t\left(x_j\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_j\right){\stackrel{\·}{q}}_h^b\left(t\right)+\underset{j=1}{\overset{N^b}{\mathrm{\Σ}}}{K}_{\mathrm{sj}}{Y}_k^t\left(x_j\right)\underset{h=1}{\overset{N{M}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_j\right)q_h^b\left(t\right)\end{array}---\left(9\right)$

在公式8和公式9中，K_si、C_si为整体道床沿纵向每米的支承刚度(N/m)与阻尼(N·s/m)m^b(m^t)为道床梁(隧道梁)的单位长度质量(kg/m)；l^b(l^t)为道床(隧道)简支梁长度(m)；NM^b(NM^t)为道床梁(隧道梁)所截取的模态阶数；N^b(N^t)为道床梁(隧道梁)支点数；Z^b(x,t)(Z^t(x,t))为道床梁(隧道梁)的振动位移；E^b(E^t)为道床梁(隧道梁)的弹性模量(Pa)；I^b(I^t)为道床梁(隧道梁)的截面惯性矩(m4)。

2，第三类模型的振动微分方程

第三类模型与第二类模型相比，主要区别是在盾构隧道结构的模拟方面。第三类模型将盾构隧道看作是由1.2m宽的管片质量环组成，并可有效考虑管片环间的剪切作用，其受力模型如图10所示。因此，第三类模型不仅可考虑整体道床与盾构隧道的动力相互作用关系，还可用于分析管片环间的受力情况。

由第一类与第二类模型的振动微分方程，不难得出第三类模型中整体道床梁的振动微分方程(见公式10)以及盾构隧道的振动微分方程(见公式11)。

$\begin{array}{c}{\stackrel{\·}{q}}_k^b\left(t\right)+\underset{i=1}{\overset{N^b}{\mathrm{\Σ}}}{C}_{\mathrm{si}}{Y}_k^b\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_i\right){\stackrel{\·}{q}}_k^b\left(t\right)+\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{C}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_j\right){\stackrel{\·}{q}}_h^b\left(t\right)+\frac{E^b{I}^b}{m_r^b}{\left(\frac{\mathrm{k\π}}{l^b}\right)}^4{q}_k^b\left(t\right)\\ +\underset{i=1}{\overset{N^b}{\mathrm{\Σ}}}{K}_{\mathrm{si}}{Y}_k^b\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_i\right)q_k^b\left(t\right)+\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{K}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_h^b\left(x_j\right)q_h^b\left(t\right)\\ -\underset{j=1}{\overset{N^t}{\mathrm{\Σ}}}\underset{i=(j-1)n-1}{\overset{\mathrm{jn}}{\mathrm{\Σ}}}{C}_{\mathrm{si}}{Y}_k^b\left(x_i\right){\stackrel{\·}{Z}}_{\mathrm{tj}}-\underset{j=1}{\overset{N^t}{\mathrm{\Σ}}}\underset{i=(j-1)n+1}{\overset{\mathrm{jn}}{\mathrm{\Σ}}}{K}_{\mathrm{si}}{Y}_k^b\left(x_i\right)Z_{\mathrm{tj}}\\ =\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{C}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_j\right)\stackrel{\·}{q}\left(t\right)+\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{K}_{\mathrm{pj}}{Y}_k^b\left(x_j\right)\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_j\right)q_h\left(t\right)\end{array}---\left(10\right)$

在公式10中：Z_ti为盾构隧道管片环的振动位移；为盾构隧道管片环的振动速度；n为单个管片环上道床底部支点数，N^tn＝N^b。

$\begin{array}{c}{M}_{\mathrm{ti}}{\stackrel{\·\·}{Z}}_{\mathrm{ti}}\left(t\right)+(\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{\mathrm{sj}}+C_{\mathrm{fi}}+2C_{\mathrm{wi}}){\stackrel{\·}{Z}}_{\mathrm{ti}}\left(t\right)+(\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{K}_{\mathrm{sj}}+K_{\mathrm{fi}}+2K_{\mathrm{wi}})Z_{\mathrm{ti}}\left(t\right)\\ -\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{\mathrm{sj}}{\stackrel{\·}{Z}}_{\mathrm{bj}}\left(t\right)-\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{K}_{\mathrm{sj}}{Z}_{\mathrm{bj}}\left(t\right)-C_{\mathrm{wi}}{\stackrel{\·}{Z}}_{t(i+1)}\left(t\right)-K_{\mathrm{wi}}{Z}_{t(i+1)}\left(t\right)\\ -C_{\mathrm{wi}}{\stackrel{\·}{Z}}_{t(i+1)}\left(t\right)-K_{\mathrm{wi}}{Z}_{t(i-1)}\left(t\right)=0\end{array}---\left(11\right)$

在公式11中，K_wi、C_wi为管片环纵向剪切刚度(N/m)与阻尼(N·s/m)；M_t为盾构单个管片环质量(kg)；n为一个管片环上道床底部支点数。其中，

$Z_{\mathrm{bj}}\left(t\right)=\underset{k=1}{\overset{{\mathrm{NM}}^b}{\mathrm{\Σ}}}{Y}_k^b\left(x_j\right)q_k^b\left(t\right)---\left(12\right)$

进一步，在确定了钢轨扣件反力后，便可以进一步得出地铁扣件系统频域荷载谱。并可以进一步进行多断面的隧道-地层二维有限元模型的建立。也就是地铁振动频域响应计算方法包括如下步骤：

在预设的图层中，将不同土层分界线和隧道边界线分别离散成由多个两两间隔预定距离的坐标点所组成的坐标点集合，以确定土层分界点集和隧道边界点集；

沿里程方向，逐层比较纵向土层分界点集中，和/或隧道边界点集中，相邻点的深度坐标变化幅度；

若相邻点的深度坐标变化幅度超过预定变化幅度数值，则根据变化幅度超过预定数值的点的横坐标为参考坐标，沿竖直方向分别获取土层分界线和隧道边界线上与参考坐标相同的位置的点的参数信息作为建模信息；

根据获取到的建模信息建立隧道与地层二维有限元模型，并求解不同预测断面内地铁隧道环境振动频域响应。

具体的，为了实现地铁全线众多断面的隧道-地层二维有限元模型的循环建立与自动预测，需要利用地铁沿线地层与隧道纵、横断面的勘察设计CAD图及其基于VBA语言的二次开发技术，来自动划分预测断面与自动提取各断面的建模信息。所需的建模信息包括隧道埋深(埋置深度)、土层类别、土质参数及其对应深度。该部分的专利技术(地铁振动频域响应计算方法)的目的是依据地铁沿线土层与隧道建模信息的变化幅度(综合考虑预测效率与预测精度发现建模信息的变化幅度被设定在10％比较合理)，自动进行隧道-地层预测横断面的科学划分，同时还可为各断面自动读取土层的埋深、层厚、类型、参数以及隧道净空、隧道埋深的建模信息，便于实现地铁隧道全线不同里程所有预测横断面的批量化建模，并在建模的基础上进行频域响应的计算。

具体实现的过程如图6所示。首先，在地层与隧道纵断面CAD图内，将不同土层分界线(图6中的水平方向细线)与隧道边界(图 6中的水平方向粗线)设置为不同的图层，同时将这些分界线离散成1m间隔的坐标点，然后沿里程方向逐层比较纵向土层(或隧道顶层、底层)分界线上相邻离散点(图6中的圆圈)深度坐标的变化幅度，如果任一层分界线上沿纵向的离散点深度坐标的变化幅度超过10％，则自动插入竖直线(图6中的竖直方向的直线线)，每段红色竖直线即代表不同的预测横断面。在插入竖直线(即预测横断面)之后，利用CAD中图层可逐个显示的功能，循环显示各类土层与隧道层，每显示一图层就能自动读取不同竖直线(即预测横断面)与土层(或隧道层)上下界面的交点坐标(图6中的圆圈)，同时还可方便地将各类土层的物理力学指标与这些交点坐标形成对应，最终生成可供调用的隧道与地层建模信息文件。因此，这个建模信息文件既包括不同里程所有预测横断面内各类土层、隧道层分界线的位置坐标，还包括土层物理力学参数(这些参数包括土层的非线性动模量、非线性动阻尼以及剪切波速、密度、泊松比，并可根据预测断面内各土层分界线的位置坐标形成一一对应)。

进一步，在此基础上，也就是在地铁沿线预测断面自动划分与其建模信息自动提取的技术基础上，借助Python脚本语言，在ABAQUS有限元软件中自动读入隧道与地层的建模信息文件(在地铁振动频域响应计算方法中已经获得)，然后逐个断面地自动完成几何布置、网格划分、材料赋值与边界处理的重复性工作，这些工作的具体含义如图7所示。

从图7可以看出，在每个预测横断面的二维有限元模型的建立过程中，首先需要利用建模信息文件中预测断面内隧道与土层的位置坐标，划分不同的区域(其中不同颜色的各区域代表不同类型的土层或隧道层，如图7所示)。然后，将各区域划分为4节点4单元 的矩形网格(为了避免人为的低通滤波而造成的预测误差，需要将近场土层的网格边长取为地层剪切波长的1/12，地层剪切波长又等于剪切波速除以所关心的最大振动频率，其中剪切波速可从地质勘查资料中获取，最大关心频率一般是人体的最大敏感频率80Hz)，同时使用线弹性材料模拟隧道、道床结构及远场土体，近场土层则需使用非线性粘弹性材料来模拟。再之后，完成模型边界的设置，即在模型左端边界设置水平约束的对称边界，模型下部边界设为固定约束，模型上部边界与隧道的远端边界(即模型右端边界)为自由无约束边界。最后，运用直接稳态动力学分析方法(之所以选择直接稳态动力学分析法，是因为该求解方法的计算效率与计算精度均较高，一般在5分钟左右就可完成1个断面的计算，非常适合对地铁全线多断面的环境振动响应进行快速预测)进行频域求解，即可直接得到不同预测断面内地铁隧道环境振动频域响应。显然，这些工作在各个预测断面内具有相似性与重复性，因此可通过有限元软件的二次开发语言将上述过程编成程序，从而通过循环调用该程序来实现所有二维有限元预测断面的循环建立与自动求解，继而提高地铁隧道全线环境振动的预测效率。

与以往地铁环境振动评估技术的最大不同之处在于，本专利致力于可视化呈现地铁长区间、大范围的环境振动响应，并首次创造出地铁隧道环境振动的地图式管理理念。为了使振动地图更有针对性，利用Matlab语言编制了我国各类环境振动控制标准中的不同评价方法(如1/3倍频程计权振动加速度级等)，并将理论预测结果处理成相应的评价指标，然后再将其绘制成图。(目前主要编制了《城市区域环境振动标准》(GB10070-1988)、《城市轨道交通引起建筑物振动与二次辐射噪声限值及其测量方法标准》(JGJ/T170-2009)、 《城市轨道交通(地下段)列车运行引起的住宅建筑室内结构振动与结构噪声限值及测量方法》(DB31/T470-2009)、《住宅建筑室内振动限值及其测量方法》(GB50355-2005))。另外，在振动地图可视化呈现方式研究方面，通过比选确立了3种适宜的呈现方式，分别是地铁环境振动级沿线路纵向与横向分布的俯视色彩云图(见图12)、振动级的等势线图(见图13)、振动级的空间曲面色彩云图(见图15)。同时，为了有效反映地铁线路两侧建(构)筑物的分布与其周边地铁环境振动水平的对应情况，还在振动地图中精确定位出地铁沿线的敏感建(构)筑物(见图14)，并需要相应存储它们的基础类型、结构型式、材料力学参数的数字化信息(这些信息需单独处理与逐步积累，是本专利中最费时的工作)，有助于进一步用以预测与评估地铁沿线建(构)筑物的二次结构噪声。具体的，图14中，第一部分代表华侨大厦(混9)，第二部分代表江东商业大厦(混8)，第三部分代表浦发银行(混9)，第四部分代表中医院江东门诊部(混4)，第五部分代表商业银行江东支行(混6)，第六部分代表江东华侨大楼(混6)。

下面以一个具体的实施案例来说明本发明的价值：以某城市地铁隧道为例，应用前面发明的地铁隧道环境振动地图技术，可得到该城市地铁隧道周边的环境振动水平。具体计算条件是，地铁车辆为A型车，车速为80km/h；轨道结构包括60kg/m钢轨、弹性分开式DT VI2型扣件(扣件胶垫材质是聚氨酯)与短枕埋入式整体道床；隧道是外径6.2m的盾构隧道，隧道地基以粉质黏土、细砂以及砂质黏土为主。另外，由于我国尚无地铁轨道不平顺标准谱，因此这里选用美国5级轨道不平顺谱(轨道不平顺的计算波长范围为0.1m至100m)作为输入激励进行预测计算。

根据上述计算条件，可得到隧道地基与轨下胶垫动参数随频率、幅值变化前后扣件力的功率谱密度，如图11a和图11b所示。从图中可以看出，如果不考虑隧道地基与轨下胶垫动参数随频率、幅值变化的真实力学行为，将会导致地铁荷载与实际情况严重不符。通过应用图11a和图11b的荷载谱，并与自动循环建立的隧道-地层二维有限元模型相结合，即可快速得到该城市地铁隧道上方地表纵向350m与横向45m范围内的振级色彩云图(见图12与图14)、振级等势线图(见图13)与振级空间曲线色彩云图(见图15)。需要说明的是，在图12至14中，横坐标轴为隧道纵向中心线在地面上的投影(单位：m)，纵坐标为距离隧道中心线的水平距离(单位：m)，地图最右边的色谱等级是我国环境振动常用评价指标——VL_z竖向振动加速度级(单位：dB)。在图15中，横坐标轴同样为隧道纵向中心线在地面上的投影(单位：m)，纵坐标亦为距离隧道中心线的水平距离(单位：m)，竖向坐标则为VL_z竖向振动加速度级(单位：dB)。

与相关技术相比，本专利具有以下至少一种优点：

1、与现有预测方法相比，地铁隧道环境振动地图技术的经济成本较低，预测效率与预测精度均较高，可适用各种车辆、各类轨道、各式隧道与多种地质环境，并可快速实现数十甚至数百公里范围地铁隧道环境振动的可视化评价。

2、所提出的车-轨-隧-地基大系统耦合的非线性动力学频域预测模型，可有效考虑轨道系统内高分子材料与隧道地基动参数(主要是刚度与阻尼)随激振频率、激振幅值、环境温度非线性变化的实际情况。

3、综合考虑预测效率与预测精度，设定了预测断面建模信息沿线路纵向的变化幅度限值为10％，一旦超过该限值，可自动插入隧道-地层环境振动预测横断面，同时还可为地铁全线所有预测横断面自动提取相对应的地层与隧道建模信息(包括隧道埋深、土层类别、土质参数及其对应深度)，为地铁沿线大量的隧道-地层二维有限元模型的建立提供了便利，避免了机械性的人工劳动，节省了大量时间，从而极大地提高了预测效率。

4、在调用地铁全线预测横断面的建模信息文件基础上，将有限元模型建立过程中几何布置、网格划分、材料赋值、边界处理及模型求解的重复性工作编制成程序，然后通过循环调用该程序来实现地铁沿线多组隧道-地层环境振动二维有限元模型的自动建立与求解，从而可为相关从业人员减少机械性的重复劳动量，节省了大量的时间成本。

5、首次提出了地铁隧道全线周边环境振动总体水平及其分频振动水平分布状况的可视化展示技术，并比选出了3种适宜的图像表达方式，可用以直接识别不同车辆、不同轨道、不同隧道以及地质环境下地铁周边环境振动的敏感区与安全区。

6、地铁隧道环境振动地图集成了我国各类环境振动评价规范，它们的评价对象、评价方法、评价指标及其限值均存在一定的差异，因此可方便且快捷地找出不同评价体系下地铁隧道沿线环境振动超限的具体方位及其严重程度。

7、利用地铁环境振动地图以及地铁周边建(构)筑物信息(包括基础类型、结构型式、材料力学参数及其与隧道的位置关系)，可进一步用以分析地铁隧道周边建(构)筑物本身的振动水平及其引发的二次结构噪声。

以上，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。

Claims (10)

1.钢轨扣件反力计算方法，其特征在于，包括：

在预先获取的第一模型的基础上，增加整体道床与隧道叠合梁模型，以确定第二模型，所述第一模型是车辆-轨道耦合模型，所述第二模型是车辆-轨道-隧道耦合模型；

在第二模型中计算隧道指定部分的动参数随频域振动位移的变化量，其中，动参数包括动阻尼和动刚度；

根据动参数随频域振动位移的变化量，计算在第二模型中的不同工况下的钢轨扣件反力。

2.根据权利要求1所述的钢轨扣件反力计算方法，其特征在于，所述第二模型为：

$\begin{array}{c}{\stackrel{..}{q}}_k^{\mathrm{bt}}\left(t\right)+\underset{i=1}{\overset{N^{\mathrm{bt}}}{\mathrm{\Σ}}}{C}_{\mathrm{fi}}{Y}_k^{\mathrm{bt}}\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{tb}}\left(x_i\right){\stackrel{.}{q}}_k^{\mathrm{bt}}\left(t\right)+\frac{E^{\mathrm{bt}}{I}^{\mathrm{bt}}}{m^{\mathrm{bt}}}{\left(\frac{\mathrm{k\π}}{l^{\mathrm{bt}}}\right)}^4{q}_k^{\mathrm{bt}}\left(t\right)\\ +\underset{i=1}{\overset{N^{\mathrm{bt}}}{\mathrm{\Σ}}}{K}_{\mathrm{fi}}{Y}_k^{\mathrm{bt}}\left(x_i\right)\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right)q_k^{\mathrm{bt}}\left(t\right)=\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{F}_{\mathrm{rsj}}\left(t\right)Y_k^{\mathrm{bt}}\left(x_j\right)\end{array},$

$\begin{array}{c}{F}_{\mathrm{rsj}}\left(t\right)=C_{\mathrm{pi}}\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_i\right)\stackrel{.}{q}\left(t\right)+K_{\mathrm{pi}}\underset{h=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{Y}_h\left(x_i\right)q_h\left(t\right)\\ -C_{\mathrm{pi}}\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right){\stackrel{.}{q}}_h^{\mathrm{bt}}\left(t\right)-K_{\mathrm{pi}}\underset{h=1}{\overset{{\mathrm{NM}}^{\mathrm{bt}}}{\mathrm{\Σ}}}{Y}_h^{\mathrm{bt}}\left(x_i\right)q_h^{\mathrm{bt}}\left(t\right)\end{array};$

其中，m^bt为道床-隧道叠合梁单位长度的质量(kg/m)；E^bt为整体道床-盾构隧道叠合梁的弹性模量(Pa)，I^bt为道床-隧道梁的截面惯量(m⁴)；l^bt为道床-隧道简支梁长度(m)；NM^bt为道床-隧道叠合梁所截取的模态阶数；N^bt为道床-隧道梁下离散支点数；Y_k ^bt(x_i)为道床-隧道叠合梁的振型；q_k ^bt(t)为道床-隧道叠合梁的正则振型坐标；K_fi、C_fi为隧道底部地基沿纵向每米的支承刚度(N/m)与阻尼(N·s/m)，F_rsj(t)为钢轨扣件反力(N)，K_pi、C_pi分别为轨下胶垫刚度(N/m)与阻尼(N·s/m)，NM为钢轨所截取的模态阶数，Y_h(x_i)与q_h(t)分别是钢轨振型和钢轨正则振型坐标。

3.根据权利要求2所述的钢轨扣件反力计算方法，其特征在于，所述在第二模型中计算隧道指定部分的动参数随频域振动位移的变化量包括：

根据第二模型得到第二模型的矩阵表达式， $\left[M\right]\left\{\stackrel{..}{q}\right\}+\left[C\right]\left\{\stackrel{.}{q}\right\}+\left[K\right]\left\{q\right\}=\left\{Q\right\}=\left[K_f\right]\left[Z_0\right],$ 其中，[M]为系统的质量矩阵，[C]为系统的阻尼矩阵，[K]为系统的刚度矩阵，[K_f]为系统的转换矩阵，{Q}为系统广义力向量，{q}为系统广义位移向量，{Z₀}为不平顺位移向量；

根据第二模型的矩阵表达式分别确定隧道指定部分随频率、幅值与温度变化的阻尼矩阵的传递函数和随频率、幅值与温度变化的刚度矩阵传递函数；

计算第二模型系统的不同自由度响应的功率密度谱。

4.根据权利要求3所述的钢轨扣件反力计算方法，其特征在于，所述根据第二模型的矩阵表达式分别确定隧道指定部分随频率、幅值与温度变化的阻尼矩阵的传递函数和随频率、幅值与温度变化的刚度矩阵传递函数包括：

令{Z₀}＝[I]e^-jωt，[q]＝[H(ω)]e^-jωt，得出第二模型系统动参数随频率、幅值与温度变化情况下的传递函数，(-ω²[M]+iω[C_ω,A,T]+[K_ω,A,T])[H(ω)]＝[K_f][I]，其中，[C_ω,A,T]为随频率、幅值与温度变化的阻尼矩阵，[K_ω,A,T]为随频率、幅值与温度变化的刚度矩阵，ω是不同振动频率对应的角速度，A是不同频率对应的幅值，T是环境温度；

所述计算第二模型系统的不同自由度响应的功率密度谱包括，根据如下公式计算第二模型系统在不同自由度响应的功率谱密度S_qi(ω)＝|H(ω)|²S₁₁(ω)，其中，S_qi(ω)为第二模型系统第i个自由度下的位移响应功率谱密度，为系统第i个自由度下的加速度响应功率谱密度，S₁₁(ω)为系统第1位轮对处轨道不平顺功率谱密度。

5.根据权利要求3所述的钢轨扣件反力计算方法，其特征在于，在所述得出第二模型系统动参数随频率、幅值与温度变化情况下的传递函数后还包括：

在第二模型系统中，根据隧道指定部分材料在不同频率下，振幅小于预设数值的动参数，计算各频率振幅，以确定参考振幅；

根据所述参考振幅重新选择对应的动参数，进一步求解材料应变，并判断不同频率下动变形与动参数是否收敛至某一定值时停止；

若否，则再次执行步骤在第二模型系统中，根据隧道指定部分材料在不同频率下，振幅小于预设数值的动参数，计算各频率振幅，以确定参考振幅。

6.根据权利要求5所述的钢轨扣件反力计算方法，其特征在于，所述根据所述参考振幅重新选择对应的动参数，进一步求解材料应变包括：采用1/3倍频程来分频带地统计应变，再利用不同频带应变大小与动参数的对应关系，重新调整系统刚度与阻尼矩阵。

7.地铁振动频域响应计算方法，基于如权利要求1-6所示的钢轨扣件反力计算方法，其特征在于，还包括：

在预设的图层中，将不同土层分界线和隧道边界线分别离散成由多个两两间隔预定距离的坐标点所组成的坐标点集合，以确定土层分界点集和隧道边界点集；

沿里程方向，逐层比较纵向土层分界点集中，和/或隧道边界点集中，相邻点的深度坐标变化幅度；

若相邻点的深度坐标变化幅度超过预定变化幅度数值，则根据变化幅度超过预定数值的点的横坐标为参考坐标，沿竖直方向分别获取土层分界线和隧道边界线上与参考坐标相同的位置的点的参数信息作为建模信息；

根据获取到的建模信息建立隧道与地层二维有限元模型，并求解不同预测断面内地铁隧道环境振动频域响应。

8.根据权利要求7所述的地铁振动频域响应计算方法，其特征在于，所述根据获取到的建模信息建立隧道与地层二维有限元模型，并求解不同预测断面内地铁隧道环境振动频域响应包括：

根据建模信息中，预测断面内隧道与土层的位置坐标，划分代表不同土层或隧道层的多个区域；

将各区域划分为4节点4单元的矩形网格，并使用线弹性材料模拟隧道、道床结构及远场土体，使用非线性粘弹性材料模拟近场土层；

根据稳态动力学分析方法进行频域求解，即可直接得到不同预测断面内地铁隧道环境振动频域响应。

9.根据权利要求7所述的地铁振动频域响应计算方法，其特征在于，所述预定变化幅度数值为8％-15％。

10.根据权利要求9所述的地铁振动频域响应计算方法，其特征在于，所述预定变化幅度数值为10％。