基于傅里葉函數(shù)擬合的滑坡變形監(jiān)測(cè)分析

董耀剛，頡 麗，姚正學(xué)，楊 軍，劉興榮，王翔宇

(1.甘肅省科學(xué)院 地質(zhì)自然災(zāi)害防治研究所，甘肅 蘭州 730000；2.中國(guó)石油天然氣股份有限公司 西南管道貴陽(yáng)輸油氣分公司，貴州 貴陽(yáng) 550000)

滑坡深部位移監(jiān)測(cè)是滑坡變形監(jiān)測(cè)的重要內(nèi)容,對(duì)準(zhǔn)確確定滑動(dòng)面的位置、評(píng)價(jià)滑坡活動(dòng)現(xiàn)狀及未來(lái)發(fā)展趨勢(shì)具有重要意義[1]。目前主要是根據(jù)位移-時(shí)間曲線確定滑坡的穩(wěn)定性，根據(jù)累積位移-深度曲線和鉆孔柱狀圖確定滑動(dòng)帶的位置[2]。具體應(yīng)用中，前期主要是根據(jù)滑坡深部位移變形曲線形態(tài)來(lái)定性分析滑動(dòng)面的位置和滑坡變形發(fā)展過(guò)程，但由于監(jiān)測(cè)結(jié)果是時(shí)間序列的隨機(jī)變量[3]，一些學(xué)者嘗試用非線性的方法，如Apriori算法[4]、Kalman濾波數(shù)據(jù)融合[5]、小波去噪后回歸[6]、信息維[7]、多重分形[8]、Verhulst-ARMA組合預(yù)測(cè)模型[9]、半?yún)?shù)模式的 GM(1,1)[10]等，進(jìn)行監(jiān)測(cè)數(shù)據(jù)挖掘，從而了解滑坡形變及其原因。這些方法的應(yīng)用取得了一定效果，但也存在諸多局限，如預(yù)測(cè)精度不高、適應(yīng)性不強(qiáng)、部分算法較為復(fù)雜等。研究計(jì)算簡(jiǎn)便、精度較高、運(yùn)算快捷的模型來(lái)分析滑坡監(jiān)測(cè)數(shù)據(jù)是研究者一直追求的目標(biāo)。本研究以2015年G316線稍子坡段(K2556+345～K2563+000)10#滑坡深部位移監(jiān)測(cè)數(shù)據(jù)為例，采用MATLAB軟件中的curve fitting tool工具對(duì)監(jiān)測(cè)數(shù)據(jù)進(jìn)行處理，具有計(jì)算簡(jiǎn)單方便、結(jié)果通俗易懂且精度較高的特點(diǎn)，可為更加準(zhǔn)確地分析滑坡的發(fā)展趨勢(shì)提供基礎(chǔ)數(shù)據(jù)。

1 滑坡特征

稍子坡段10#滑坡位于G316線K2560+526～K2560+777 km處，屬于老滑坡?；麻L(zhǎng)137 m、寬222 m、平均厚度8.5 m、面積25 400 m2、體積約29.4萬(wàn)m3，滑體由擾動(dòng)破碎的棕紅色含礫泥巖及少量砂巖組成，主滑方向?yàn)镹E43°。滑坡后緣坡陡，上方山體地形坡度較緩，滑體前緣滑坡舌推覆于河床之上，逼迫河道，前后緣相對(duì)高差180 m?；伦笥覂蓚?cè)的變形破壞有明顯差異：滑體左側(cè)的中后部地勢(shì)相對(duì)較高，前緣堆積在泥盆系褶皺地層上，并已被河流侵蝕呈陡坡；滑體右側(cè)在主滑之后至少又發(fā)生過(guò)3次小范圍的滑動(dòng)，現(xiàn)已完全滑落，后部地勢(shì)較低、坡面較緩，地下水埋深僅1.0 m左右，有濕地和3處泉水出露，有幾處小范圍淺層坡體因含水量高而明顯變形[11]，同時(shí)滑坡舌受河流侵蝕經(jīng)常發(fā)生坍塌。國(guó)道316線從滑坡中部通過(guò)，受滑坡左右兩側(cè)差異變形影響，通過(guò)滑體右側(cè)的公路損壞嚴(yán)重。2015年初相關(guān)單位邀請(qǐng)專家對(duì)該路段滑坡進(jìn)行考察會(huì)診時(shí)，考慮到該滑坡規(guī)模大、投資大和治理相對(duì)困難等，專家一致建議暫不治理，先進(jìn)行深部變形監(jiān)測(cè)，根據(jù)監(jiān)測(cè)結(jié)論再作決定。

2 滑坡深部位移監(jiān)測(cè)方案

滑坡深部位移監(jiān)測(cè)是通過(guò)儀器測(cè)量地下巖土體相對(duì)于穩(wěn)定地層的位移量，據(jù)此確定移動(dòng)巖土體的滑移面和變形變化的規(guī)律[12]。其監(jiān)測(cè)原理示意和監(jiān)測(cè)點(diǎn)照片見(jiàn)圖1、2。

滑坡深部監(jiān)測(cè)孔的位置一般選在滑坡體主滑方向和具有代表性的主滑軸線上[12]。根據(jù)該滑坡的發(fā)育特征，在公路上方，沿其主滑斷面自上而下布設(shè)ZK10-9#、ZK10-8#、ZK10-5#三個(gè)監(jiān)測(cè)孔，以確定滑動(dòng)面深度。受人為活動(dòng)影響，ZK10-8#監(jiān)測(cè)孔孔口嚴(yán)重堵塞，因此在與ZK10-8#監(jiān)測(cè)孔斜距2.5 m處補(bǔ)設(shè)BZ8#監(jiān)測(cè)孔，以代替ZK10-8#監(jiān)測(cè)孔。此外，為便于繼續(xù)觀測(cè)該滑坡變形特征，在公路外側(cè)增設(shè)1個(gè)監(jiān)測(cè)孔BZ1#。各監(jiān)測(cè)孔基本參數(shù)如表1所示。監(jiān)測(cè)儀器采用成都世紀(jì)安巖儀器有限公司生產(chǎn)的AY-1B型固定式測(cè)斜儀，可直接監(jiān)測(cè)和讀取位移數(shù)據(jù)，并有微電腦時(shí)控開(kāi)關(guān)實(shí)時(shí)控制將數(shù)據(jù)定時(shí)存入電腦。

圖1 滑坡深部監(jiān)測(cè)原理示意

圖2 滑坡深部監(jiān)測(cè)孔照片

表1 滑坡深部監(jiān)測(cè)孔基本參數(shù)

參數(shù)各監(jiān)測(cè)孔基本參數(shù)ZK10-9#BZ8#ZK10-5#BZ1#孔口標(biāo)高/m1 711.061 700.901 690.231 672.23孔深/m25.823.525.023.0測(cè)深/m25.523.024.522.5總測(cè)深/m95.5

3 監(jiān)測(cè)數(shù)據(jù)成果

滑坡深部位移監(jiān)測(cè)可以提供大量的觀測(cè)數(shù)據(jù)，根據(jù)這些數(shù)據(jù)可以整理出各種曲線，從不同角度了解滑坡的性狀和發(fā)展趨勢(shì)。2015年稍子坡段10#滑坡各監(jiān)測(cè)孔監(jiān)測(cè)數(shù)據(jù)見(jiàn)表2。

表2 2015年各監(jiān)測(cè)孔累計(jì)位移量

注：其中“—”表示因施工覆蓋等原因?qū)е伦冃尾幻黠@，未測(cè)得數(shù)據(jù)。

4 監(jiān)測(cè)數(shù)據(jù)的處理及分析

4.1 累計(jì)位移-深度曲線

利用累計(jì)位移-深度曲線可以進(jìn)行滑動(dòng)面判斷，滑動(dòng)面位置一般在累計(jì)位移-深度曲線相鄰一對(duì)正負(fù)曲率最大點(diǎn)之間，若為多滑動(dòng)面或滑動(dòng)帶則將出現(xiàn)幾對(duì)正負(fù)曲率極大值點(diǎn)。埋設(shè)在滑坡后部的監(jiān)測(cè)孔ZK10-9#無(wú)明顯的位移，說(shuō)明滑坡后緣處于相對(duì)穩(wěn)定狀態(tài)；BZ8#監(jiān)測(cè)孔在6.5～8.0 m處出現(xiàn)位移突變，說(shuō)明在這個(gè)深度范圍內(nèi)可能存在1個(gè)滑動(dòng)面(帶)；ZK10-5#監(jiān)測(cè)孔在7.0～11.5、16.5～17.5、21.0～23.5 m位置均出現(xiàn)位移突變，其中7.0～11.5 m處僅在1月28日出現(xiàn)位移，可認(rèn)為是儀器誤差等偶然事件導(dǎo)致，而在16.5～17.5、21.0～23.5 m處存在滑動(dòng)面的可能性較大。綜合3個(gè)監(jiān)測(cè)孔的監(jiān)測(cè)結(jié)果，可以在滑坡縱剖面上確定出滑動(dòng)面位置，如圖3所示。

4.2 位移-時(shí)間曲線

利用位移-時(shí)間曲線可以確定滑坡的穩(wěn)定性，但想通過(guò)簡(jiǎn)單的位移變化規(guī)律來(lái)判斷滑坡?tīng)顟B(tài)是困難的，這就需要采用一定的數(shù)學(xué)分析方法來(lái)預(yù)測(cè)滑坡深部位移的變形趨勢(shì)，其實(shí)質(zhì)就是一組數(shù)據(jù)的曲線擬合問(wèn)題。曲線擬合的方法很多，僅在MATLAB的curve fitting tool中就提供了12大類近40種擬合模型。早在19世紀(jì)早期，法國(guó)數(shù)學(xué)家傅里葉就提出“任何函數(shù)都能用無(wú)窮多個(gè)正弦和余弦函數(shù)的和表示”，說(shuō)明傅里葉函數(shù)具有廣泛的適應(yīng)性，特別適合對(duì)周期性函數(shù)進(jìn)行描述。考慮到滑坡位移數(shù)據(jù)具有周期性的特點(diǎn)，本研究選擇傅里葉函數(shù)進(jìn)行擬合，為簡(jiǎn)便起見(jiàn)，選擇一階傅里葉函數(shù)進(jìn)行擬合，其基本公式為

y(t)=a0+a1coswt+b1sinwt

(1)

式中：t為時(shí)間，d；y(t)為位移，mm；a0、a1、b1分別為待定系數(shù)；w為基頻，hz。

為檢驗(yàn)曲線擬合效果，一方面利用均方根誤差、誤差平方和,以及擬合度等統(tǒng)計(jì)量來(lái)檢驗(yàn)，同時(shí)考慮到本次監(jiān)測(cè)數(shù)據(jù)有限，故僅留最后一組數(shù)據(jù)(ZK10-9#監(jiān)測(cè)孔)用作預(yù)測(cè)效果檢驗(yàn)，其余數(shù)據(jù)均參與建模，另一方面選擇三次多項(xiàng)式擬合和高斯擬合兩種模型作為效果對(duì)比分析模型。三次多項(xiàng)式擬合模型和高斯擬合模型分別為

y(t)=p1t3+p2t2+p3t+p4

(2)

(3)

上二式中：p1、p2、p3、p4、a1、a2、b1、b2、c1、c2均為待定系數(shù)；其余參數(shù)意義同上。

采用MATLAB軟件中的curve fitting tool擬合工具箱，對(duì)上述擬合函數(shù)進(jìn)行最小二乘法求解，3種擬合方法的擬合結(jié)果見(jiàn)表3。由表3可知，各擬合函數(shù)的擬合度均大于0.9，擬合效果均較好，其中：BZ8#監(jiān)測(cè)孔高斯擬合效果最好，擬合度最高，誤差平方和、均方根誤差最小，但預(yù)留數(shù)據(jù)組預(yù)測(cè)誤差率最小的是傅里葉擬合；ZK10-5#監(jiān)測(cè)孔擬合效果最佳的是三次多項(xiàng)式，擬合度最高，均方根誤差最小，但就誤差平方和、預(yù)留數(shù)據(jù)預(yù)測(cè)誤差率來(lái)說(shuō)，效果最佳的是高斯擬合；BZ1#監(jiān)測(cè)孔誤差平方和及均方根誤差最小、擬合度最高、預(yù)留數(shù)據(jù)預(yù)測(cè)誤差率最小的都是傅里葉擬合。總體來(lái)說(shuō)傅里葉擬合效果較好，因此本研究采用傅里葉擬合函數(shù)對(duì)滑坡深部位移變形進(jìn)行預(yù)測(cè)，預(yù)測(cè)結(jié)果見(jiàn)圖4—6。由圖4—6知，該滑坡深部位移仍處于變形擴(kuò)展階段，急需采取有效措施以遏制滑坡變形的進(jìn)一步擴(kuò)展。

5 結(jié)果與討論

以2015年G316線稍子坡段10#滑坡深部位移監(jiān)測(cè)數(shù)據(jù)為例，利用MATLAB軟件中的curve fitting tool曲線擬合工具箱中的傅里葉函數(shù)對(duì)3組實(shí)測(cè)數(shù)據(jù)進(jìn)行擬合分析，得出如下幾點(diǎn)認(rèn)識(shí)：

(1)滑坡深部位移監(jiān)測(cè)是研究滑坡深部位移特征行之有效的手段，利用累計(jì)位移-深度曲線上發(fā)生突變的位置，結(jié)合滑坡縱剖面圖，可以較準(zhǔn)確地確定出滑動(dòng)面的位置。

(2)利用曲線擬合的方法，如三次多項(xiàng)式、高斯函數(shù)、傅里葉函數(shù)等，可以對(duì)滑坡位移-時(shí)間曲線進(jìn)行深層次處理，進(jìn)而對(duì)滑坡變形發(fā)展趨勢(shì)進(jìn)行預(yù)測(cè)，方法簡(jiǎn)便且精度較高。預(yù)測(cè)結(jié)果顯示，該滑坡目前仍處于變形擴(kuò)展階段，急需采取有效措施來(lái)遏制變形的發(fā)展。

(3)滑坡位移-時(shí)間曲線擬合效果可以通過(guò)擬合度、誤差平方和、均方根誤差等統(tǒng)計(jì)量來(lái)檢測(cè)，擬合結(jié)果顯示：3種擬合函數(shù)均具有較佳的擬合效果，擬合度均大于0.9，但綜合來(lái)看以傅里葉函數(shù)擬合效果最佳，擬合度最高、誤差平方和、均方根誤差最小。

(4)鑒于傅里葉函數(shù)擬合的普適性及對(duì)周期性函數(shù)描述效果最佳的特點(diǎn)，可以利用傅里葉函數(shù)擬合的方法建立滑坡預(yù)測(cè)預(yù)報(bào)模型，但有關(guān)滑坡階段劃分、臨滑判據(jù)等需要進(jìn)一步研究。

表3 3種擬合模型擬合結(jié)果分析

注：*表示實(shí)測(cè)值。

圖4 BZ8#變形趨勢(shì)預(yù)測(cè)曲線 圖5 ZK10-5#變形趨勢(shì)預(yù)測(cè)曲線 圖6 BZ1#變形趨勢(shì)預(yù)測(cè)曲線