通海臺垂直擺受建筑載荷影響的數(shù)值模擬研究

岳 沖 李利波 閆 偉 洪 敏

1 中國地震臺網(wǎng)中心，北京市三里河南橫街5號，100045 2 中國地震局地質(zhì)研究所地震動力學(xué)國家重點(diǎn)實(shí)驗(yàn)室，北京市華嚴(yán)里甲1號，100029 3 云南省地震局，昆明市北辰大道148號，650224

為精確測量地傾斜潮汐變化，目前常用的儀器有水管傾斜儀、水平擺傾斜儀和垂直擺傾斜儀，上述儀器的觀測精度可以達(dá)到10-10～10-9。以往研究結(jié)果顯示，地震之前往往存在趨勢轉(zhuǎn)折、破年變、短臨加速等異常變化[1-3]。柯昌安等[4]研究發(fā)現(xiàn)，2008年寧強(qiáng)5.7級地震前7 d，安康、寶雞、寧陜和西安臺垂直擺、水管儀和伸縮儀均出現(xiàn)短臨加速的異常現(xiàn)象；馬震等[5]研究發(fā)現(xiàn)，2016年門源6.4級地震前，周邊門源、湟源等形變臺站傾斜類觀測儀器出現(xiàn)短期趨勢轉(zhuǎn)折異常變化，各測項(xiàng)的異常時(shí)間較為同步，且各儀器出現(xiàn)異常的分量方位與震中方位較為吻合。地形變觀測儀器能夠直觀地反映區(qū)域構(gòu)造應(yīng)力場的調(diào)整變化，進(jìn)而捕捉地震前兆異常，但由于測量精度較高，觀測儀器易受到各種干擾因素的影響，如氣象因素、人為干擾、環(huán)境干擾、儀器更換或故障等[6-9]。山東泰安臺應(yīng)變及傾斜類觀測儀器受8 km外水庫抽水、補(bǔ)水等產(chǎn)生的載荷變化影響，出現(xiàn)2～3 d持續(xù)的短臨加速變化[10]；2014年魯?shù)?.5級地震前，距離震中不到50 km的昭通地震臺出現(xiàn)多個(gè)測項(xiàng)的短臨異常變化，研究發(fā)現(xiàn)該異常與降雨相關(guān)性較高[11]；Yue等[12]利用有限元模型分析發(fā)現(xiàn)，南北向長、東西向短的雙陽湖水位增加更易在湖周邊地表引起沿南北向形變，造成雙陽臺水管傾斜儀出現(xiàn)NS及EW方向差異性變形特征。因此，形變臺站觀測儀器出現(xiàn)變化后，除分析與地震相關(guān)的區(qū)域構(gòu)造應(yīng)力調(diào)整引起的變化外，也需要結(jié)合儀器不同測項(xiàng)的變化特征，分析其受到干擾因素的影響，進(jìn)而對儀器變化性質(zhì)做進(jìn)一步判定。

通海臺垂直擺NS和EW測項(xiàng)自2019-04起出現(xiàn)同步趨勢轉(zhuǎn)折變化，類似的趨勢轉(zhuǎn)折變化在2014年景谷6.6級地震前也出現(xiàn)過，但2019-04以來的趨勢轉(zhuǎn)折過程中存在建筑施工的影響。因此，本文針對通海臺垂直擺NS與EW測項(xiàng)趨勢轉(zhuǎn)折變化，構(gòu)建通海臺建筑載荷影響三維有限元模型，定量分析建筑載荷對通海臺垂直擺NS與EW測項(xiàng)影響的差異性特征，進(jìn)一步探討通海臺垂直擺趨勢轉(zhuǎn)折變化機(jī)理，為進(jìn)一步分析此次趨勢轉(zhuǎn)折變化對未來該區(qū)域應(yīng)力調(diào)整、地震危險(xiǎn)性評價(jià)影響等提供參考依據(jù)。

1 通海臺垂直擺

1.1 通海臺地質(zhì)概況

通海臺位于通海盆地南緣，距離杞麓湖南側(cè)約3 km，通海盆地為斷陷盆地，呈長條形EW走向。通海臺始建于1971-07，位于寶秀-建水?dāng)嗔驯眰?cè)約1 km處[13]，地理坐標(biāo)為102.78°E、24.10°N，海拔高程為1 821 m，觀測山洞深38 m，全長54 m，觀測硐室覆蓋層厚23 m，硐室?guī)r性主要為白云巖、石灰?guī)r，山洞內(nèi)部巖體部分破碎。

通海臺地處川滇菱形塊體南部，該區(qū)域受NW向石屏-建水?dāng)嗔褞А⑶瓟嗔岩约敖麼S-NNE向小江斷裂帶南段等控制和影響(圖1)，上述區(qū)域歷史地震活動強(qiáng)烈，通海、峨山、石屏等地發(fā)生過10多次6.0級以上破壞性地震，震級最大的一次為1970-01-05通海M7.8 地震。自1970年以來(截至2021-01-01)，該地區(qū)發(fā)生多次5級左右地震，尤其是2018-08-13與2018-08-14兩次MS5.0地震均發(fā)生在通海盆地北緣，距離臺站僅約12 km。

1.2 通海臺垂直擺變化特征

目前通海臺硐室內(nèi)共布設(shè)有VS型和VP型垂直擺傾斜儀，儀器觀測墩類型為石灰?guī)r塊墩。其中，VS型垂直擺傾斜儀自2008-11-08開始觀測，截至2021-06，除NS測項(xiàng)于2010-01～02間存在部分缺數(shù)外，其他時(shí)段內(nèi)數(shù)據(jù)完整性較好，可為本文研究提供較長的背景數(shù)據(jù)資料。VP型垂直擺傾斜儀于2019-08才開始觀測，觀測時(shí)間較短，因此數(shù)據(jù)缺乏可靠性。本文重點(diǎn)結(jié)合VS型垂直擺2019年前數(shù)據(jù)變化的背景趨勢，對2019-04后NS和EW測項(xiàng)出現(xiàn)的同步變化特征進(jìn)行分析。

圖1 通海臺位置Fig.1 Location of Tonghai station

VS型垂直擺2009年以來的觀測數(shù)據(jù)資料如圖2(a)、(b)所示，圖中已對通海臺周邊歷史地震進(jìn)行標(biāo)注(標(biāo)注規(guī)則：震中距100 km內(nèi)MS5.0～5.5地震，200 km內(nèi)MS5.5～6.0地震，300 km內(nèi)MS6.0以上地震)。由原始觀測曲線可知，垂直擺NS測項(xiàng)2009～2014年一直持續(xù)N傾變化，2014年初N傾趨勢出現(xiàn)轉(zhuǎn)平后逐漸轉(zhuǎn)為S傾變化，2019-04后再次轉(zhuǎn)為N傾變化；垂直擺EW測項(xiàng)自觀測以來一直持續(xù)E傾變化，但E傾趨勢在2013年出現(xiàn)速率變緩，2019-04后E傾趨勢再次加速。

圖2 通海臺垂直擺NS和EW測項(xiàng)原始曲線、擬合曲線、去趨勢時(shí)序曲線Fig.2 Time series of original data, fitting data and de-trended data of vertical pendulum inclinometer NS and EW components of Tonghai station

單從垂直擺兩個(gè)測項(xiàng)時(shí)序曲線來看，兩個(gè)測線變化較為同步，均在2013年前后和2019-04出現(xiàn)趨勢轉(zhuǎn)折變化，而震例統(tǒng)計(jì)結(jié)果顯示，2013年出現(xiàn)趨勢轉(zhuǎn)折后，距離臺站259 km的景谷于2014-10-07發(fā)生MS6.6地震，2013年前后NS和EW測項(xiàng)的趨勢轉(zhuǎn)折變化是否與此次地震有關(guān)，還需進(jìn)一步深入研究。自2019-04起，距離硐室北東側(cè)約50 m處陸續(xù)開始建設(shè)醫(yī)院等建筑物(圖2中紅色陰影時(shí)段)，垂直擺NS和EW測項(xiàng)數(shù)據(jù)變化時(shí)間與建筑物載荷影響時(shí)間較為吻合，但該測項(xiàng)變化是否均由建筑物施工造成以及NS和EW測項(xiàng)變化量為何相差較大等問題仍需作進(jìn)一步分析研究。

進(jìn)一步對垂直擺NS和EW測項(xiàng)作趨勢擬合處理，以消除儀器零漂等產(chǎn)生的影響，其中NS測項(xiàng)采用冪函數(shù)進(jìn)行擬合，相關(guān)系數(shù)為0.923；EW測項(xiàng)采用對數(shù)函數(shù)進(jìn)行擬合，相關(guān)系數(shù)為0.97，結(jié)果如圖2(c)、(d)所示。去趨勢后的NS和EW測項(xiàng)在2019-04后趨勢轉(zhuǎn)折更加突出，且NS測項(xiàng)趨勢變化量顯著大于EW測項(xiàng)，2019-04后NS測項(xiàng)N傾變化量達(dá)到約2 750 ms，而EW測項(xiàng)E傾變化量達(dá)到約580 ms，NS測項(xiàng)約為EW測項(xiàng)的4.7倍。

為排除觀測儀器更換、格值變化等對儀器觀測的影響[9]，計(jì)算2017～2020年NS和EW測項(xiàng)M2波潮汐因子(圖3)。由圖可知，2017～2020年NS測項(xiàng)M2波潮汐因子在0.7左右波動，其中2018-08～09出現(xiàn)較為突出的M2波潮汐因子增大現(xiàn)象，可能是因?yàn)槭艿酵ê纱蜯S5.0以及墨江MS5.9地震的影響；EW測項(xiàng)M2波潮汐因子4 a間較為穩(wěn)定，均在0.52左右波動。NS與EW測項(xiàng)在2019年的變化趨勢(圖中藍(lán)色曲線)與其他年份相差不大，表明2019-04至今儀器觀測穩(wěn)定，未發(fā)生儀器標(biāo)定等格值變化現(xiàn)象，而數(shù)據(jù)的轉(zhuǎn)折變化未能在潮汐頻段引起同步變化。

圖3 通海臺NS和EW測項(xiàng)M2波潮汐因子Fig.3 M2 wave tidal factor of vertical pendulum inclinometer NS and EW components of Tonghai station

2 模型設(shè)計(jì)

2.1 模型及邊界條件設(shè)置

2019年觀測硐室北東側(cè)約50 m的縣醫(yī)院住院部開始施工，該工程2019-09開始搭建主體建筑，整個(gè)樓層共設(shè)計(jì)12層，總面積為12 643.1 m2。整個(gè)模型尺寸設(shè)計(jì)為600 m×600 m，為避免模型邊界效應(yīng)等產(chǎn)生的影響，將通海臺設(shè)置于模型中央，距離模型左右邊界分別為300 m，模型深度設(shè)計(jì)為800 m。醫(yī)院位于通海臺東北部，其西南角距通海臺50 m，醫(yī)院模型尺寸設(shè)計(jì)為50 m×50 m。對醫(yī)院工程量載荷進(jìn)一步簡化，每層為0.5 m鋼混結(jié)構(gòu)(ρ=2 500 kg/m3)，共計(jì)12層，計(jì)算醫(yī)院對地面壓強(qiáng)為150 000 Pa(重力加速度取值10 m/s2)。地下介質(zhì)結(jié)構(gòu)參考距臺站30 km的玉溪市防震減災(zāi)局地震監(jiān)測井(DZKC1)，結(jié)合鉆井及通海臺山洞巖性特征，對模型沿深度方向進(jìn)行剖分，其中0～100 m為砂巖，100～200 m為粘土層，200～800 m為石灰?guī)r(圖4(a))。整個(gè)模型底部沿垂直方向進(jìn)行固定，模型表面自由，僅考慮建筑物產(chǎn)生的載荷壓強(qiáng)作用，不考慮重力作用，采用彈性體模型進(jìn)行模擬計(jì)算。整個(gè)垂直擺布設(shè)于石灰?guī)r塊墩上，折合擺長為10 cm，因此在實(shí)際模擬計(jì)算中將垂直擺簡化為模型中心點(diǎn)位。為保證計(jì)算結(jié)果的準(zhǔn)確性，對臺站及醫(yī)院區(qū)域網(wǎng)格進(jìn)行加密，其中模型邊界網(wǎng)格大小為50 m，臺站及醫(yī)院周邊網(wǎng)格間隔為1 m，網(wǎng)格類型為C3D8R，共計(jì)獲得466 760個(gè)網(wǎng)格單元(圖4(b))。

圖4 模型設(shè)置及網(wǎng)格剖分Fig.4 Model setting and grid segmentation

2.2 模型參數(shù)

為對比模型參數(shù)對模擬結(jié)果的影響，構(gòu)建兩類模型參數(shù)進(jìn)行模擬計(jì)算：模型1中巖石楊氏模量參數(shù)取值為巖石參數(shù)下限，砂巖為2.5×107Pa，粘土層為1.0×108Pa，石灰?guī)r為5.0×1010Pa[14]；模型2中砂巖和粘土層楊氏模量取值為模型1的4～5倍，具體參數(shù)設(shè)置見表1。

2.3 模擬結(jié)果

通過施加建筑物載荷作用，模擬計(jì)算兩種模型地表位移結(jié)果(圖5)。由圖5(a)、(b)可知，沿NS向位移(N向位移為正，S向位移為負(fù))結(jié)果顯示，醫(yī)院(圖5中紅色虛線框)南部形成N向位移，醫(yī)院北部形成S向位移，通海臺位于N向位移影響范圍內(nèi)；沿EW向位移(E向位移為正，W向位移為負(fù))結(jié)果顯示，醫(yī)院西部形成E向位移，東部形成W向位移，通海臺位于E向位移影響范圍內(nèi)。為進(jìn)一步研究醫(yī)院建筑載荷對通海臺的影響，在模型表面作穿過通海臺沿NS向剖面1和沿EW向剖面2(剖面1由N向S，剖面2由W向E)，剖面位移場結(jié)果如圖5(c)、(d)所示。由圖可知，剖面1三維形變場以垂向沉降為主，最大沉降位移在剖面線225 m距離處，最大沉降量達(dá)40.45 mm，通海臺站處沉降量為-3.67 mm；剖面1沿EW向位移均為正值，表明該剖面均呈現(xiàn)E向運(yùn)動特征，在225 m距離處E向位移量達(dá)到最大10.09 mm，通海臺站處E向位移量為0.72 mm；剖面1沿NS向位移首先呈現(xiàn)S向位移運(yùn)動特征，NS向位移在225 m距離處變?yōu)?后逐漸轉(zhuǎn)為N向位移，通海臺站處N向位移量為2.21 mm。由于剖面2未穿過醫(yī)院建筑載荷區(qū)，該剖面三維形變場量級低于剖面1，三維形變場仍以垂向位移變化為主，在325 m距離處達(dá)到最大沉降量4.25 mm，EW向位移呈現(xiàn)E向運(yùn)動轉(zhuǎn)W向運(yùn)動特征，而NS向位移則均為N向運(yùn)動。

表1 模型參數(shù)設(shè)置及模擬結(jié)果

圖5 模擬位移結(jié)果及剖面1和剖面2三維形變場Fig.5 Simulated displacement results and 3D deformation field of profile 1 and profile 2

由此可知，位于通海臺東北部的醫(yī)院建筑物載荷在臺站附近產(chǎn)生N向以及E向位移加載，這與垂直擺NS和EW測項(xiàng)在2019-04出現(xiàn)的N傾和E傾趨勢轉(zhuǎn)折變化較為一致。模擬結(jié)果顯示，在不考慮地形等因素影響的條件下，模型1和模型2的巖石參數(shù)均顯示，建筑載荷對通海臺產(chǎn)生的N向位移量大于E向位移量，模擬得到的N向位移量約為E向位移量的3倍。

3 討 論

“十五”數(shù)字化改造后，由前兆數(shù)據(jù)觀測原理和數(shù)據(jù)產(chǎn)出過程可知，絕大多數(shù)地殼傾斜、應(yīng)變的初始產(chǎn)出物理量為電壓，通過與儀器“格值”相乘后可得到相應(yīng)的物理觀測量，如應(yīng)變量、傾斜量等[9,15]。由此可知，儀器觀測結(jié)果并非地殼“真實(shí)”形變量，這會進(jìn)一步增大本文模擬結(jié)果與實(shí)際觀測結(jié)果直接對比的不確定性。通海臺垂直擺擺體折合擺長為10 cm，垂直擺NS向自2019-04開始趨勢轉(zhuǎn)折，N傾轉(zhuǎn)折變化量達(dá)到約2 750 ms。假設(shè)2019-04趨勢轉(zhuǎn)折前垂直擺儀器處于絕對垂直狀態(tài)，N傾轉(zhuǎn)折變化量轉(zhuǎn)換為N向位移僅為0.001 3 mm，而該變形量遠(yuǎn)遠(yuǎn)小于模型1與模型2中NS向位移量。由于儀器硐室?guī)r性主要為白云巖、石灰?guī)r，而上述模型結(jié)果為建筑載荷對地表較軟砂巖和粘土層的影響。為探究上述差異是否由模型參數(shù)設(shè)置所引起，本文對模型參數(shù)進(jìn)一步調(diào)整，將模型中砂巖、粘土層參數(shù)均設(shè)置為石灰?guī)r參數(shù)，模擬建筑物直接作用在石灰?guī)r地層所產(chǎn)生的影響，得到通海臺N向位移量為0.031 mm，E向位移量為0.004 4 mm，變化量仍大于儀器觀測的結(jié)果。

建筑載荷作用到沉積層后在傳遞至儀器所在巖體的過程仍需更加深入的研究，但受到該地區(qū)精細(xì)的地下波速結(jié)果、鉆孔信息等資料限制，本文對地下介質(zhì)連續(xù)性、各向異性、物性參數(shù)、地表起伏變化等已進(jìn)行必要的簡化，且未考慮地下流體等因素的影響。通過上述簡化模型模擬不同巖石參數(shù)下建筑載荷對通海臺NS和EW測項(xiàng)的影響，結(jié)果均顯示建筑載荷對通海臺產(chǎn)生的形變場存在不均勻性，建筑載荷在NS向位移變化顯著大于EW向位移變化，這與實(shí)際儀器中NS測項(xiàng)變化量大于EW測項(xiàng)變化量的結(jié)論一致。

4 結(jié) 語

針對通海臺垂直擺NS和EW測項(xiàng)自2019-04出現(xiàn)同步趨勢轉(zhuǎn)折變化，利用三維有限元模型模擬建筑載荷對通海臺NS向和EW向位移影響的差異特征，主要結(jié)論如下：

1)對垂直擺NS、EW測項(xiàng)作趨勢擬合處理，以消除儀器零漂等產(chǎn)生的影響，結(jié)果表明，2019-04后NS測項(xiàng)N傾變化量達(dá)到約2 750 ms，而EW測項(xiàng)E傾變化量達(dá)到約580 ms，NS測項(xiàng)趨勢變化量顯著大于EW測項(xiàng)，NS測項(xiàng)約為EW測項(xiàng)的4.7倍。

2)計(jì)算2017～2020年NS和EW測項(xiàng)M2波潮汐因子，NS測項(xiàng)M2波潮汐因子在0.7左右波動，EW測項(xiàng)M2波潮汐因子4 a間較為穩(wěn)定，均在0.52左右波動；2019-04至今儀器觀測穩(wěn)定，未發(fā)生儀器標(biāo)定等格值變化現(xiàn)象，而數(shù)據(jù)的轉(zhuǎn)折變化未能在潮汐頻段引起同步變化。

3)位于通海臺東北部的醫(yī)院建筑載荷在臺站附近產(chǎn)生N向及E向位移加載，與垂直擺NS和EW測項(xiàng)在2019-04后出現(xiàn)的N傾和E傾趨勢轉(zhuǎn)折變化較為一致。在不考慮地形等因素的影響下，模型1和模型2結(jié)果均顯示，建筑載荷對通海臺產(chǎn)生的N向位移量均大于E向位移量，N向位移量約為E向位移量的3倍。

致謝：感謝中國電建集團(tuán)北京勘測設(shè)計(jì)研究院有限公司李善飛高級工程師在軟件使用過程中提供的幫助。