0.   引言

随着科技的进步与发展，环境振动问题愈发受到关注，其对人体健康、结构安全和仪器正常工作等皆可造成不良影响，被列为世界七大环境公害之一^[1]。环境振动主要指由自然(气象、海洋、地下构造等)或人工(交通、施工等)活动引起的随机波动，工程实践中此类环境振动的振源往往不易确定，波场十分复杂，本文将此称为无源环境振动^[2](后简称无源振动)。无源振动的特点包括在整个工程考虑范围内没有明显振源的影响，场地振动表现出一种类似地脉动的性质，即振动都处于幅值较小且相近的量级，位移幅值一般在纳米到微米级别等。但地面的无源振动以远场瑞利波为主导，与地脉动仍有所差别。

现代科学技术迅猛发展，以同步辐射光源和自由电子激光为代表的大型高精密科学装置的建设也不断涌现。表 1给出了中国近20 a来一些大型高精密科学装置对基础振动幅值的控制要求，相较于规范^[3]中精密设备及仪器容许振动值，大型高精密科学装置振动分析更为关注振动位移均方根，且容许值普遍小于500 nm，较规范值更为严格，需开展准确定量研究。目前现场实测是研究环境振动对高精密科学装置影响的主要手段，如德国电子同步加速器研究中心(DESY)组织Amirikas等^[4]对全球近20个同步辐射光源实验室的地基振动进行了标准化的测试与分析，并形成了开放的数据库。其中，对于振动测试得到的离散位移数据{x_i}(i=1, 2, …, n)，振动位移均方根计算公式为

表  1  高精密科学装置对基础振动的控制标准^{[2, 8, 10]}

Table  1.  Control criteria of vibration of foundation of high-precision scientific equipments ^{[2, 8, 10]}

装置名称	1~100 Hz竖向振动位移均方根/nm
上海光源	< 150(安静时段)，< 300(嘈杂时段)
软X射线自由电子激光	< 500
硬X射线自由电子激光	< 150
北京高能同步辐射光源	< 25
注：安静时段系指0：00—4：00；嘈杂时段系指4：00—24：00。

下载: 导出CSV 

| 显示表格

现场实测证明了即便装置选址在无明显振源且环境振动较小的区域，亦存在振动超标的情况，上海光源(SSRF)便是典型案例之一。黄茂松等^[5-6]、周仁义^[7]在对上海光源环境微振动控制的研究中提出了无源振动的概念，将瑞利波随深度衰减规律用于无源振动分析，为较高环境振动水平的软土地基中成功建设上海光源做出了积极贡献，同时也为中国在精密装置基础微振动控制技术领域积累了宝贵经验。冯忠俊^[8]对上海自由电子激光项目(SXFEL)进行了详细的振动测试，并基于有限元方法对比了排桩和隔振沟两种减振方案的优劣性。岳建勇等^[9-10]对上海光源的微振动和微沉降控制进行了详细的介绍，提出了微振动监控方法和微变形控制技术，并将现场测试与计算结果进行了对比。

目前，新建高精密科学装置对环境微振动的控制要求越来越严苛，如正在建设的北京高能同步辐射光源工程(HEPS)对基础振动的控制要求为1~100 Hz的振动位移均方根应小于25 nm。动力有限元方法是实际重大工程振动分析的有力工具，本文依托中国“十三五”重大科技基础设施北京高能同步辐射光源工程筏板基础的环境微振动分析需求，在文献[5~7]首次提出的无源振动有限元分析方法的基础上，形成了现场实测、理论分析和有限元计算相结合的无源环境微振动分析方法，并对该方法进行了验证，研究成果为类似重大工程的环境微振动分析提供了重要参考。

1.   无源振动有限元分析方法

1.1   隔离体分析

常规的环境振动分析通常包括振源、地基(传播路径)、受振体(结构)3个部分。在无明显振源的条件下，振动分析需要解决振动输入问题。本文提出了一个无源振动分析模型，模型只包含由受振体和地基组成的隔离体，所有的振动都从模型侧边界输入，如图 1所示。该模型相当于在传播路径上将振源和受振体分开考虑，边界上的节点振动对模型而言即为振源，即将原有振源等效为其在模型边界上引起的振动场。显然，无源振动分析方法的关键是确定边界每个节点的振动输入。理论上可以通过大量的振动测试解决输入问题，然而实际工程中往往只对少量测点进行测试，因此需要通过合理的理论方法计算获得其他地面以及地下各点的振动时程。

图  1  无源振动分析模型示意图

Figure  1.  Schematic diagram of analysis model for vibration without apparent sources

下载: 全尺寸图片 幻灯片

1.2   地面振动确定方法

实际工程通常只在地面若干个点进行振动实测，因此地面上其它节点处振动需根据实测振动进行计算获得。该过程按实测测点的数量可分为两类，如图 2所示。若地面有多个测点，可将各个测点的信号进行傅里叶变换，得到对应的幅值谱和相位谱，然后将各点的幅值谱、相位谱按至需计算点的距离进行线性插值，得到需计算点的幅值谱和相位谱；将其相加得到傅里叶谱，再进行傅里叶逆变换，即可得到需计算点的地面振动时程。若地面仅一个测点，首先将地面振动时程信号进行傅里叶变换，得到对应的幅值谱和相位谱；其后保持幅值谱不变，将[-1, 1]范围内的随机数与相位谱相乘得到需计算点的相位谱；最后将原幅值谱与新相位谱相加后进行傅里叶逆变换，即可得到需计算点的地面振动时程。

图  2  地面振动确定方法的流程图

Figure  2.  Flow chart for determining ground surface vibration

下载: 全尺寸图片 幻灯片

利用上述地面振动确定方法可以得到地面任意位置的振动时程，其幅值与实测地面幅值大小接近，频域特性也基本一致，符合无源振动的特点。

1.3   振动随深度衰减模型

地下振动实测难度大且成本高，且通常测点数量较少。因此，无源振动分析采用理论方法建立振动随深度的衰减模型，再根据地面点的振动时程来计算地下各点的振动时程。

对于均质土层中瑞利波随深度衰减规律，可使用Richart等^[11]的弹性解析方法来获得。但实际地基通常为层状地基，且控制远场振动幅值的瑞利波在层状介质中具有弥散性，其在浅部软弱土层中可能会出现与均质土层不符的振动放大现象^[12]，传统的弹性解析方法不能反映层状地基振动特性，故本文采用Thomson-Haskell传递矩阵法^{[7, 13-14]}来求解瑞利波弥散曲线和各频率振动随深度衰减因子。

假定地基由均质、各向同性的弹性土层构成，层与层间的交界面平行于地基表面，根据弹性动力学理论，地面下的位移和应力可用标量势φ和矢量势ψ表达，将其代入拉梅方程并省略面波乘积因子，写成矩阵形式：

式中，u为x方向位移分量，w为z方向位移分量，${\sigma _z}$和τ_xz分别为z方向的正应力和切应力，i为虚数，k为波数且k=ω/c，ω为圆频率，c为相速度，上标+代表向上传播的波，-代表向下传播的波。

此处矩阵T仅是土体特性的函数，与频率无关，

式中，ρ为土层密度，r，r_α和r_β为相速度与拉梅常数有关的参数，具体可参考文献[6，12]。

由此，可将土层划分为n层，利用传递矩阵法可得到划分的第m层(m=1, 2, …, n)的解析表达式：

式中，z_m为第m层底部到地面的距离，G_m为层间传递矩阵，起到将矩阵S(z_m-1)向下传播的作用，E_m为相速度、拉梅常数、波数和土层厚度所构成的矩阵，具体可参考文献[6，12]。

根据土层在交界面上位移和应力的连续条件，由传递矩阵法可得：

将式(5)代入式(2)可得

由于土层内部没有振源，因此不存在向上传播的波。同时，地面为自由表面，其应力均为0。由此边界条件得

将式(7)的边界条件代入式(6)可得

式中，r_ij为矩阵R的元素。

式(8)可简化为

矩阵[LS]是瑞利波波速c和波长λ_R的函数。为保证地面存在振动，u和w不能同时为零，矩阵行列式必为零，即det[LS]=0，此即瑞利波的特征方程。由此，对于某一给定波长的瑞利波，可用式(9)求解出相应的瑞利波波速，通过对不同波长求解方程便可得到整个土层中瑞利波的弥散曲线。

同时，利用瑞利波的特征方程，也可以得到水平位移与垂直位移的比值。因此，可将地表垂向位移设为1.0，并结合式(7)，便可得到满足瑞利波特征方程的边界条件：

最后，联立式(4)，(9)和(10)，分别对垂向与水平向地面位移进行归一化，便可求得到任意深度的瑞利波位移衰减因子(位移衰减因子是指某一深度处的位移与地面处的位移之比)。

上述Thomson-Haskell传递矩阵法可快速且简单地获得瑞利波的弥散曲线和位移随深度的衰减因子。由此，便可利用地面实测振动来获得特定深度处的振动时程，如图 3所示。首先对地表瑞利波的位移时域曲线进行傅里叶变换，得到地表振动的傅里叶频谱。其次，根据不同频率的瑞利波位移随深度衰减因子，对地表位移傅里叶频谱按频率进行衰减，得到不同深度处的位移傅里叶频谱。最后对其进行傅里叶逆变换即可得到计算深度处的振动位移时程。

图  3  计算某一深度处位移时程的流程图

Figure  3.  Flow chart for calculating vibration at a certain depth

下载: 全尺寸图片 幻灯片

在求解出有限元模型所有边界节点振动时程并施加于相应的边界节点后，即可使边界位移场与实际振动位移场趋于一致，使有限元计算可较为准确地模拟工程场地实际振动响应。

2.   北京光源的无源环境微振动分析

2.1   工程背景

北京高能同步辐射光源(HEPS，简称北京光源)是我国重大科技基础设施建设“十三五”规划中提出优先建设、世界最高亮度的第四代同步辐射光源项目。北京光源位于北京市怀柔科学城内，南侧和西侧距离永乐大街与京加线约1 km，周边道路交通流量相对较少，属于典型的无源振动情况。北京光源效果图如图 4所示，主要实验设备位于装置区内，包括增强器和储存环等。其主体结构的基础形式采用筏板基础，并辅以素混凝土换填，部分设备机房采用柱下独立基础，其动力设备自身采取减振隔振措施^[15]。

图  4  北京光源效果图^[15]

Figure  4.  Schematic diagram of HEPS ^[15]

下载: 全尺寸图片 幻灯片

北京光源工程对地基基础微振动控制提出了极高的要求，包括各个方向频率在1~100 Hz的振动位移均方根(RMS)需小于25 nm，且地基基础需对地面振动位移幅值有较为明显的抑制作用。中国尚无同等水平的地基基础环境微振动的设计及评估经验。一旦施工后振动超出限值，将对装置的正常运行带来严重危害。因此，非常有必要对现有地基基础设计方案进行微振动分析，以评估其可行性及合理性。

2.2   有限元模型

图 5给出了北京光源工程有限元分析模型图。模型为圆柱体，其半径根据地面实测测点(即边界振动输入)到储存环中心位置的距离确定，为410 m。模型深度根据瑞利波影响深度确定，瑞利波影响深度在一倍波长左右^[16]，考虑到地表填土层会被挖除，故以第二层卵石土作为表面波主要传播介质。该层卵石土剪切波速为344 m/s，实测场地基频为3 Hz左右，一倍波长约为115 m，故模型深度确定为120 m。分析模型的结构部分由实验大厅、储存环及增强器组成，如图 5(b)所示。储存环隧道外环半径为237 m(含实验大厅)，内环半径为210 m，储存环隧道宽7.6 m，高5 m，隧道壁厚0.8 m，其基础筏板厚度为1 m，素混凝土换填厚度为3 m。增强器外环半径为72 m，内环半径为67 m，增强器隧道宽5 m，高4 m，隧道壁厚0.5 m，其下筏板厚度为0.7 m，素混凝土换填厚度为1.5 m。

图  5  北京光源工程有限元模型图

Figure  5.  Finite element model for HEPS

下载: 全尺寸图片 幻灯片

为保证计算精度，模型的单元尺寸需小于瑞利波最小波长的1/10^[17-18]。考虑土体剪切波速取344 m(由于填土层已被换填)，计算频率取场地基频3 Hz，则结构区域外土体的最大单元尺寸取10 m，储存环结构区域内土体单元尺寸为2~5 m，实验大厅单元尺寸为8 m，储存环隧道单元尺寸为0.8~2.5 m，增强器隧道单元尺寸为0.5~2.5 m。有限元计算采用ABAQUS软件，单元采用线性实体单元C3D8。模型底部固定，共包含291941个单元和314738个节点。

北京光源工程地基以卵石层为主要持力层，从上而下依次为填土层、卵石层、全风化基岩、强风化、中风化、微风化基岩等。无源环境微振动一般处于纳米到微米级别，可将土体视为线弹性材料。根据工程勘察报告^[15]，土体参数取值见表 2。结构的本构模型均采用线弹性模型。其中，隧道及筏板的弹性模量设为30 GPa，换填的素混凝土弹性模量设为22 GPa，两者泊松比均为0.167。模型未施加材料阻尼。

表  2  北京光源工程地基土层分布及参数

Table  2.  Soil profile and parameters at HEPS site

土层	厚度/m	密度 /(kg·m-3)	压缩波速/(m·s-1)	剪切波速/(m·s-1)	泊松比	弹性 模量 /MPa
①素填土	1	1910	636	173	0.460	166.9
②1卵石	11	2200	952	344	0.425	742.0
②2卵石	6	2200	1451	492	0.435	1528.4
④全风化基岩	8	2250	1962	777	0.407	3822.5
⑤强风化基岩	4	2350	2253	1011	0.374	6600.7
⑥中风化基岩	21	2590	3616	1547	0.388	17206.8
⑦微风化基岩	49	2680	4302	1936	0.373	27583.3

下载: 导出CSV 

| 显示表格

2.3   振源输入

北京光源场地振动实测由中国科学院高能物理研究所完成，其地面测点安静时段和嘈杂时段振动位移RMS值分别约为17，50 nm。本次分析选取了嘈杂时段中一个地面测点的实测振动位移时程曲线作为地面振动输入，如图 6所示。其3个方向的位移峰值均在160 nm左右，相应的RMS值在45 nm左右。

图  6  北京光源工程场地实测振动位移时程

Figure  6.  Time histories of measured displacement at HEPS site

下载: 全尺寸图片 幻灯片

图 7，8分别给出了该工程场地瑞利波的弥散曲线和各频率对应的位移随深度衰减因子。由图 7可知，在高频短波情况下，分层地基的瑞利波速接近于上层土的瑞利波波速，当波长增大(频率减小)，波速随之不断增大，分层地基的瑞利波波速趋近于下层土的瑞利波波速。这表明当高频短波时，瑞利波受上层土影响大；当低频长波时，瑞利波受下层土影响大。同时，由图 8可知，位移随深度衰减因子的分布规律和频率有着密切的关系，频率越高，波长越短，衰减越快，传播范围越浅，反之亦然。

图  7  北京光源分层地基瑞利波弥散曲线

Figure  7.  Dispersion curve of Rayleigh waves at HEPS site

下载: 全尺寸图片 幻灯片

图  8  北京光源地基各频率位移随深度衰减因子

Figure  8.  Attenuation factors of Rayleigh waves at HEPS site

下载: 全尺寸图片 幻灯片

利用北京光源地基各频率位移随深度衰减因子，结合各节点地面振动时程，即可计算得到对应节点各深度处的振动位移时程，如图 9所示。随着深度增加，振动幅值迅速减小，且竖向振动相较于水平振动衰减较慢。对比图 6和图 9(a)中地下12 m处振动幅值可知，竖向振动在地下12 m处较地面处出现了略微的放大，这是因为地面实测测点振动主频率约为10 Hz，而地基中10 Hz竖向振动位移随深度衰减因子在地下12 m处大于1.0，即该频率的竖向振动在此位置被放大。而其他深度的位移衰减因子均小于1.0，且其随深度的增加而迅速减小，因此其他深度的振动位移均随深度增加而迅速衰减。

图  9  各深度处的计算位移时程

Figure  9.  Time histories of calculated displacement at different depths

下载: 全尺寸图片 幻灯片

最后，由上述方法计算得到北京光源模型边界共3616个节点的振动位移时程曲线，将其以位移边界条件的形式输入到模型对应的节点，进行模型的动力有限元分析。

2.4   分析方法验证

为验证本文提出的有限元分析方法的合理性，将图 5中有限元模型的结构部分去除，按表 2中土层情况建立纯地基模型，并保持网格的划分不变，以便和场地实测、有结构模型进行对比。

首先，通过模态分析计算了上述有限元地基模型的特征频率，结果表明该地基模型的特征频率为3.20 Hz，与场地实测基频3.11 Hz非常接近^[15]，初步验证了几何参数与材料模型参数选取的合理性。此外，为进一步验证无源振动分析方法的合理性，计算了模型在每个边界节点输入按上述方法计算的振动时程后的振动响应。图 10给出了纯地基模型沿直径方向土层表面各点振动响应计算结果。由图可见，不同点各方向的位移峰值及RMS值均在一个较小的范围波动，X，Y，Z方向位移峰值的平均值分别为169，161，156 nm，对应的RMS值的平均值分别为43，45，43 nm，与实测振动响应水平非常吻合，结果符合无源振动条件下整个场地地面振动大小基本相同的理论预期，进一步表明本文提出的无源振动有限元分析方法的合理性。

图  10  纯地基模型地表各节点位移振幅

Figure  10.  Amplitudes of surface vibration displacement at different positions for model without structures

下载: 全尺寸图片 幻灯片

2.5   含结构模型的计算分析

图 11给出了含结构模型整体及不同结构区域处(图 5(b))的地表各点振动位移的RMS值。总体而言，结构区域的位移RMS值较场地振动水平明显降低，表明该地基基础设计方案中筏板的设置及素混凝土换填对环境微振动的控制效果显著。同时，筏板外的土层振动位移RMS值明显增大，这说明波在遇到换填的素混凝土垫层与钢筋混凝土筏板后，有一大部分能量被反射，从而导致结构范围外的土层振动增大。同时，由于部分能量被垫层反射后，传至结构的能量将减少，且结构的刚度也较大，从而使结构的振动幅值减小。此外，增强器隧道的振动大于储存环隧道及实验大厅，主要原因可能有两个：①增强器隧道的底板厚度及垫层厚度均小于储存环隧道及实验大厅；②增强器隧道范围较小，刚度变化十分剧烈，储存环隧道则有较大范围的筏板基础(包括实验大厅)进行刚度过渡，其刚度变化或波动反射相对平缓。

图  11  北京光源工程场地和基础振动位移RMS值

Figure  11.  RMS values of vibration of soil and foundation

下载: 全尺寸图片 幻灯片

由图 11可知，采用嘈杂时段实测振动输入条件下，结构区域X和Y方向的振动位移RMS值小于25 nm的控制要求，而Z方向绝大部分区域的振动位移RMS值满足控制要求，但在变刚度区域存在轻微的超标现象，如实验大厅最外侧、储存环隧道和增强器隧道相邻处和增强器隧道内侧等。因此，可考虑增加一定的刚度过渡范围来减小振动的反射。

3.   结论

本文系统介绍了一种无明显振源条件下环境微振动的有限元分析方法，并将其应用于大科学装置北京高能同步辐射光源工程地基基础的环境微振动分析，得出3点结论。

(1) 结合地表振动实测，通过地面振动确定方法和由传递矩阵法得到的各频率振动位移随深度衰减因子，能计算场地内地表和地下任意一点的振动时程曲线，为无源振动分析中边界输入提供了理论支撑。

(2) 本文的无源振动有限元分析方法可使有限元模型地表任意处的振动位移与场地实测结果相吻合，成功模拟了北京高能同步辐射光源工程场地的实际振动情况。

(3) 计算结果表明北京光源地基基础结构区域的位移峰值较场地振动水平明显降低，表明现有地基基础设计方案对环境微振动有较好的抑制作用，但在变刚度区域振动会出现放大的现象，建议增加一定的刚度过渡范围。