基于多元非線性回歸的極震區(qū)泥石流物源量估算模型

任志剛,裴向軍，顧文韜

(1.中國電建集團(tuán)成都勘測設(shè)計研究院有限公司，四川 成都　610072;2.成都理工大學(xué) 地質(zhì)災(zāi)害防治與地質(zhì)環(huán)境保護(hù)國家重點(diǎn)實(shí)驗(yàn)室,四川 成都　610059)

?

基于多元非線性回歸的極震區(qū)泥石流物源量估算模型

任志剛1,裴向軍2，顧文韜2

(1.中國電建集團(tuán)成都勘測設(shè)計研究院有限公司，四川 成都610072;2.成都理工大學(xué) 地質(zhì)災(zāi)害防治與地質(zhì)環(huán)境保護(hù)國家重點(diǎn)實(shí)驗(yàn)室,四川 成都610059)

汶川地震極震區(qū)高川鄉(xiāng)、清平鄉(xiāng)、映秀鎮(zhèn)的28組泥石流樣本，考慮溝域面積、相對高差、發(fā)震斷裂距離、溝域巖性及巖體結(jié)構(gòu)4個物源量影響因子，基于控制變量法找出單個因子與泥石流物源儲量的相關(guān)關(guān)系，通過MATLAB多元非線性回歸方法建立了綜合因素影響下的極震區(qū)泥石流物源量估算模型。最后，利用非擬合數(shù)據(jù)以外的極震區(qū)6條溝對模型進(jìn)行了驗(yàn)證，結(jié)果表明模型對于汶川地震極震區(qū)泥石流物源儲量的估算具有較強(qiáng)的準(zhǔn)確度及適宜性。

極震區(qū)；泥石流；物源；儲量

0　引　　言

“5.12”汶川地震具有震級大、震源淺、持時長、烈度高等特點(diǎn)。震后，在強(qiáng)大的地震力作用下，直接觸發(fā)了約15000處地質(zhì)災(zāi)害[1]，沿龍門山主斷裂帶展布的狹長極震區(qū)(地震烈度大于Ⅸ地區(qū))山體產(chǎn)生了大量崩滑災(zāi)害。大量崩滑堆積體及潛在震裂山體在泥石流溝域的劇增，使得降雨條件下極震區(qū)泥石流表現(xiàn)出明顯的“災(zāi)后效應(yīng)”，給災(zāi)區(qū)人民生命財產(chǎn)造成了嚴(yán)重的二次損害。

目前，國內(nèi)外關(guān)于泥石流物源特征的研究方法也不拘一格，主要是在定性描述物源分布規(guī)律的基礎(chǔ)上通過一定的方法對泥石流物源量進(jìn)行估算。對于泥石流物源量估算的主要方法有形態(tài)評價法[3]、比例統(tǒng)計法[4]、灰色理論法[5]、遙感解譯法[6]等，但這些方法極少與地震效應(yīng)進(jìn)行結(jié)合。因此，研究震后泥石流的物源量估算方法對于了解該類泥石流發(fā)展趨勢、評價危害程度及科學(xué)治理具有重要意義。

1　研究區(qū)泥石流概況

由于汶川地震極震區(qū)震后泥石流表現(xiàn)出規(guī)模大、危害強(qiáng)、反復(fù)發(fā)作等特征，本文主要選取汶川地震極震區(qū)震后泥石流物源量較多、啟動條件較低、危害性較大的地區(qū)作為研究區(qū)域進(jìn)行特大地震后泥石流物源特征研究。主要選取的研究區(qū)為地震烈度為Ⅺ度的汶川縣映秀鎮(zhèn)岷江流域、地震烈度為Ⅹ度的安縣高川鄉(xiāng)高川河流域、地震烈度為Ⅹ度的綿竹市清平鄉(xiāng)綿遠(yuǎn)河流域。其中映秀鎮(zhèn)岷江流域研究區(qū)選取11條泥石流溝作為研究樣本，高川鄉(xiāng)高川河流域研究區(qū)選取9條泥石流溝作為研究樣本，清平鄉(xiāng)綿遠(yuǎn)河流域研究區(qū)選取8條泥石流溝作為研究樣本。

2　極震區(qū)泥石流物源類型及發(fā)育特征

2.1物源類型

(1)地震崩滑堆積型。強(qiáng)震時地震波在坡體內(nèi)產(chǎn)生的反復(fù)拉-壓-剪切作用，導(dǎo)致斜坡以不同的動力破壞形式失穩(wěn)，崩滑體進(jìn)入溝道后以直接參與或間接補(bǔ)給的方式參與泥石流運(yùn)動，是極震區(qū)震后泥石流物源量暴增的主要途徑。

(2)震裂山體失穩(wěn)型。地震作用下山體震裂損傷主要表現(xiàn)在兩個方面：一是在地震波反復(fù)作用下巖體劈裂、結(jié)構(gòu)面擴(kuò)張而導(dǎo)致巖體結(jié)構(gòu)劣化[7]，震裂山體在后期降雨作用及時效變形的作用下發(fā)生崩滑災(zāi)害；二是地震作用下坡表覆蓋層產(chǎn)生震動松弛現(xiàn)象[8]，使得密實(shí)度降低，后期暴雨條件下產(chǎn)生了大量土體滑坡。

(3)震裂坡表侵蝕型。震后，坡面的侵蝕現(xiàn)象主要表現(xiàn)在兩個方面：一是地震作用下坡表植被覆蓋率降低，土體密實(shí)度減小，暴雨條件下覆蓋層轉(zhuǎn)換成坡面侵蝕物源；二是震后松散的崩滑堆積體在強(qiáng)降雨作用下產(chǎn)生坡面侵蝕現(xiàn)象，成為坡面侵蝕物源。該類物源的產(chǎn)生與帶動多與泥石流水源的匯集同步發(fā)生。

2.2發(fā)育特征

極震區(qū)泥石流物源的發(fā)育特征是泥石流溝域斜坡自身發(fā)育特征以及地震作用、降雨作用等誘發(fā)因素的綜合疊加作用結(jié)果。

(1)坡體內(nèi)因。極震區(qū)泥石流溝道兩側(cè)的地形、侵蝕坡度及臨空條件是斜坡失穩(wěn)的重要因素，坡度的陡緩以及臨空程度的高低控制坡面的張力范圍以及坡腳的應(yīng)力集中帶[9]；斜坡巖土體的性質(zhì)及組合方式，是斜坡應(yīng)力集中、邊界發(fā)展貫通以及坡體失穩(wěn)的重要因素；泥石流溝道兩側(cè)坡體相對于河谷的高差大小將直接決定坡體地震加速度與振幅放大效應(yīng)的程度。

(2)地震作用。地震發(fā)生后，地震波在介質(zhì)中傳播，能量呈對數(shù)衰減趨勢衰減，根據(jù)震后地質(zhì)災(zāi)害調(diào)查，極震區(qū)地質(zhì)災(zāi)害沿主斷裂帶表現(xiàn)出明顯的“距離效應(yīng)”[10]；地震時逆沖發(fā)震斷裂的地震加速度峰值上盤衰減較慢，而下盤衰減較快，研究表明汶川地震所誘發(fā)的上盤的地質(zhì)災(zāi)害不僅在密度分布上，而且在規(guī)模上都大于下盤；發(fā)震斷裂的局部“鎖固段”在地震過程中被剪斷、破裂，將會釋放出較大的能量，產(chǎn)生強(qiáng)烈震動，形成“次級震源”，從而加強(qiáng)地表的地質(zhì)災(zāi)害的發(fā)生。

3　極震區(qū)泥石流物源量影響因子

在考慮地震效應(yīng)、巖土體結(jié)構(gòu)特征、溝域規(guī)模、地質(zhì)環(huán)境條件等情況下，將極震區(qū)泥石流物源量主控影響因子概括為：溝域面積、相對高差、離發(fā)震斷裂的距離、巖性及結(jié)構(gòu)特征。通過控制變量法找出各因子與物源儲量的相關(guān)關(guān)系后，進(jìn)行因子的多元非線性擬合。

3.1溝域面積

泥石流溝域面積對物源的影響主要表現(xiàn)在兩個方面：首先，溝域面積概化了溝的規(guī)模，同等條件下溝域面積與物源儲量呈線性增長(見圖1)；其次，溝域面積還影響靜儲量—動儲量的轉(zhuǎn)化，極震區(qū)泥石流溝多年儲量變化統(tǒng)計表明，溝域越大其匯水面積就越大，暴雨時坡表侵蝕量及溝道帶動量就越多，靜儲量—動儲量的轉(zhuǎn)化量就越高。

圖1　控制變量條件下溝域面積與總儲量關(guān)系

3.2相對高差

研究表明，地震過程中斜坡端部地震峰值加速度(PGA)總體隨高程增大[11]。為摒除泥石流溝底原始高差的影響，本文選擇溝頂與溝底的相對高差作為研究依據(jù)。高川鄉(xiāng)與映秀鎮(zhèn)研究區(qū)分段高差與震后崩滑堆積體方量的統(tǒng)計數(shù)據(jù)表明：映秀鎮(zhèn)和高川鄉(xiāng)相對高差分別在1 200 m和1 000 m后高程端部放大效應(yīng)開始突顯，崩滑體平均體積與總體積呈現(xiàn)陡增的趨勢，分段高差與崩滑體總體積總體呈指數(shù)關(guān)系增長(見圖2)。清平鄉(xiāng)研究區(qū)各崩滑災(zāi)害點(diǎn)的統(tǒng)計數(shù)據(jù)表明：崩滑體方量在相對高差800～900 m處出現(xiàn)陡增現(xiàn)象，崩滑體方量與相對高差總體呈指數(shù)關(guān)系增長(見圖3)。

圖2映秀鎮(zhèn)研究區(qū)崩滑體方量與高差關(guān)系圖3高川鄉(xiāng)研究區(qū)崩滑體方量與高差關(guān)系

3.3發(fā)震斷裂距離

對高川鄉(xiāng)、清平鄉(xiāng)各溝單位面積物源儲量與斷裂距離的統(tǒng)計結(jié)果表明：控制變量條件下，隨斷裂距離的增加，單位面積物源儲量總體呈對數(shù)關(guān)系衰減(見圖4)。

3.4巖性及巖體結(jié)構(gòu)

根據(jù)溝域巖體結(jié)構(gòu)特征采用各溝控制變量條件下的巖石堅硬程度、結(jié)構(gòu)面張開度、結(jié)構(gòu)面連通率、巖體裂隙度四個指標(biāo)評價溝域巖體質(zhì)量等級。

圖4　高川鄉(xiāng)、清平鄉(xiāng)研究區(qū)距離效應(yīng)擬合關(guān)系

(1)溝域巖體質(zhì)量等級評級系統(tǒng)。極震區(qū)研究區(qū)各溝域巖體質(zhì)量評價模型的建立主要運(yùn)用層次分析法(AHP)[14]。

本文巖體質(zhì)量評價模型所選擇的因子為：溝域巖體堅硬程度、溝域巖體結(jié)構(gòu)面張開度(e)、溝域巖體結(jié)構(gòu)面連通率(k)、溝域巖體節(jié)理平均間距(d)四個因子。

將歸一化后的矩陣進(jìn)行求行和后，再進(jìn)行列歸一化處理，得到Aij的特征向量w：

在評價因子和因子權(quán)重值的基礎(chǔ)上，結(jié)合地質(zhì)力學(xué)系統(tǒng)的理論方法(地質(zhì)力學(xué)系統(tǒng)，也稱巖體權(quán)重系統(tǒng)，即“巖體評分”)，建立巖體質(zhì)量等級評價系統(tǒng)(見表1)。

根據(jù)表1，溝域巖體的質(zhì)量等級指數(shù)N為：

aibi

式中a——為評價因子影響強(qiáng)度評分值；

b——為評價因子的權(quán)重值。

將各溝域巖體質(zhì)量等級評價因子的參數(shù)值代入巖體質(zhì)量等級評價系統(tǒng)，計算得出研究區(qū)各溝域巖體的質(zhì)量等級指數(shù)N值見表2。

表2　各溝域巖體質(zhì)量等級指數(shù)

注：N值越大則巖體質(zhì)量越差，越容易發(fā)生崩滑災(zāi)害。

(2)溝域巖體質(zhì)量與物源儲量的相關(guān)關(guān)系。本文對上述10條溝單位面積物源儲量進(jìn)行統(tǒng)計，結(jié)合各溝域巖體質(zhì)量等級評價指數(shù)，得出溝域巖體質(zhì)量與單位面積物源儲量整體呈線性關(guān)系(見圖5)。

圖5　溝域巖體質(zhì)量與單位面積儲量關(guān)系

4　極震區(qū)泥石流物源量估算模型

泥石流物源儲量決定了泥石流的規(guī)模與危險等級，是研究震后泥石流的基礎(chǔ)。本文利用極震區(qū)三大研究區(qū)28個溝道隨機(jī)樣本(見表3)，在已分別建立的溝域面積、相對高差、斷裂距離、溝域巖體質(zhì)量與總儲量相關(guān)關(guān)系(見第3章)的前提下，利用多元非線性回歸方法建立四個因子與總儲量的關(guān)系，再建立總儲量與動儲量關(guān)系，得到極震區(qū)泥石流動儲量估算模型。

4.1基于MATLAB的多元非線性回歸模型

4.1.1回歸模型的建立

根據(jù)以上對控制變量條件下極震區(qū)物源儲量與溝域面積、相對高差、斷裂距離和溝域巖體質(zhì)量的關(guān)系分析可知：

表3　汶川地震極震區(qū)三大片區(qū)泥石流溝樣本物源參數(shù)總匯

(1)在控制其余變量條件下，溝域面積與總儲量呈線性正相關(guān)(見圖1)；(2)控制變量條件下，相對高差與總儲量呈指數(shù)正相關(guān)(見圖2、3)；(3)控制變量條件下，斷裂距離與總儲量呈對數(shù)負(fù)相關(guān)(見圖4)；(4)控制變量條件下，溝域巖體質(zhì)量與總儲量呈線性正相關(guān)(見圖5)。

根據(jù)以上分析及控制變量條件下的各因子擬合結(jié)果，極震區(qū)泥石流物源總儲量(V總)回歸模型可概化為：

V總=a1X1+a2ea3X2+a4ln(X3)+a5X4+a6

式中X1——為溝域面積,km2;

X2——為相對高差,m;

X3——為斷裂距離,km;

X4——為溝域巖體質(zhì)量指數(shù)(無量綱);

a1、a2、a3、a4、a5、a6——為特征系數(shù)。

4.1.2特征系數(shù)的確定

通過MATLAB多元非線性擬合運(yùn)算程序，將表5中Xi與Yi的28組數(shù)據(jù)輸入，將ai賦入初值，并進(jìn)行非線性回歸計算后，得到系數(shù)矩陣ai為：

ai={42.300 60.039 20.002 4

-0.026 90.458-1.986 2}T

將系數(shù)代入回歸模型，可得極震區(qū)泥石流總儲量計算模型：

V總=42.300 6X1+0.039 2e0.002 4X2-

0.026 9ln(X3)+0.458X4-1.986 2

上式計算模型回歸檢驗(yàn)結(jié)果如表4，r>r0.05(28)，F(xiàn)>F0.05(4，28)，檢驗(yàn)結(jié)果表明回歸模型相關(guān)性較好，擬合效果較好[16]。

4.1.3誤差分析

將表3中各溝域參數(shù)代入上式，可得回歸模型

表4　非線性回歸結(jié)果及檢驗(yàn)

注：r表示擬合相關(guān)性，F(xiàn)表示回歸效果顯著性檢驗(yàn)的統(tǒng)計量值，r0.05(28)、F0.05(4，28)檢驗(yàn)臨界值。

計算值，與實(shí)測值相比較，得出各點(diǎn)計算值與實(shí)測值的誤差關(guān)系如圖6。

圖6　實(shí)測值與回歸值關(guān)系

絕對誤差模型絕對誤差與相對誤差見圖7、8。

從以上誤差結(jié)果統(tǒng)計中可知，除樣本6、8、12外，其余樣本誤差均較小，符合相關(guān)性標(biāo)準(zhǔn)。其中，樣本6為高川鄉(xiāng)三叉溝，樣本8為高川鄉(xiāng)黃洞子溝，樣本12為清平鄉(xiāng)文家溝。根據(jù)許強(qiáng)等對汶川地震大型滑坡發(fā)育規(guī)律的分析[17]，清平—高川段為地震主斷裂的轉(zhuǎn)折和錯裂末端部位，具有“鎖固段效應(yīng)”，因此該區(qū)域?yàn)閿嗔褖簯?yīng)力高度集中區(qū)，地震時會釋放出更多能量。而研究區(qū)三叉溝、黃洞子溝和文家溝正好位于映秀—北川斷裂高川—清平段鎖固段內(nèi)，震時三溝內(nèi)分別誘發(fā)了火石溝滑坡、大光包滑坡和文家溝滑坡三個巨型滑坡，從而大大增加了此3個樣本的物源儲量。

圖7樣本絕對誤差關(guān)系圖8樣本相對誤差關(guān)系

4.2極震區(qū)泥石流物源動儲量估算模型

調(diào)查及實(shí)際研究表明，極震區(qū)泥石流物源總儲量與動儲量存在一定的相關(guān)性。將三大研究區(qū)28個泥石流樣本中實(shí)測物源總儲量與動儲量進(jìn)行相關(guān)性擬合[18]，得到相關(guān)關(guān)系如圖9。

圖9　極震區(qū)泥石流物源總儲量與動儲量擬合關(guān)系

得到極震區(qū)三大片區(qū)泥石流物源總儲量V總與動儲量V0擬合關(guān)系如下：

V0=0.295V總-5.192 (R2=0.855)

將上式代入極震區(qū)泥石流物源總儲量計算模型，可得極震區(qū)泥石流物源動儲量估算模型：

V0=12.479A+0.011 56e0.002 4H-

0.007 936lnL+0.135 11N-5.777 93

式中V0——為極震區(qū)泥石流溝域物源動儲量;

A——為泥石流溝域面積;

H——為泥石流溝域相對河流高差;

L——為溝域離發(fā)震斷裂的垂直距離;

N——為溝域巖體質(zhì)量等級指數(shù)(可根據(jù)溝域巖體質(zhì)量等級評級系統(tǒng)求得)。

4.3模型準(zhǔn)確性及適宜性檢驗(yàn)

為檢驗(yàn)?zāi)Ｐ偷倪m宜性，隨機(jī)選擇三大研究區(qū)以外的極震區(qū)6條泥石流溝作為適宜性檢驗(yàn)樣本，根據(jù)所得的極震區(qū)泥石流動儲量預(yù)測模型，代入相應(yīng)參數(shù)(見表5)進(jìn)行計算及誤差分析。

表5　極震區(qū)泥石流物源預(yù)測模型適宜性檢驗(yàn)樣本參數(shù)匯總

將適宜性檢驗(yàn)樣本參數(shù)代入評價模型，得出泥石流物源動儲量計算值，與實(shí)測值比較，其誤差關(guān)系見圖10。

圖10　極震區(qū)泥石流物源預(yù)測模型適宜性檢驗(yàn)誤差關(guān)系

由圖10可知，除北川縣城以外，其余樣本相對誤差及絕對誤差均較小。北川縣西山坡溝溝內(nèi)發(fā)育有一個110萬m3的大型滑坡[19]，具有難以預(yù)測的襲擾因素。摒除襲擾因素，本模型適宜于評價、估算極震區(qū)未發(fā)育有巨型滑坡泥石流溝域的物源動儲量。

5　結(jié)論與建議

本文主要大量現(xiàn)場調(diào)查資料及統(tǒng)計數(shù)據(jù)對汶川地震后極震區(qū)泥石流物源的類型、分布及發(fā)育特征、影響因素及估算方法等方面展開論述研究。結(jié)合前人研究，初步概括了汶川地震極震區(qū)泥石流物源的分布和發(fā)育特征，初步研究了泥石流物源量的影響因子，并建立了極震區(qū)物源估算模型對極震區(qū)泥石流物源量進(jìn)行預(yù)測研究。主要概括為以下幾點(diǎn)：

(1)根據(jù)前人研究成果分析了極震區(qū)物源的形成機(jī)制，通過現(xiàn)場調(diào)查總結(jié)了極震區(qū)震后泥石流物源類型，并從坡體內(nèi)因以及地震斷裂效應(yīng)兩方面分析了震后泥石流物源的形成機(jī)制及發(fā)育特征；

(2)在物源發(fā)育特征分析的基礎(chǔ)上總結(jié)了影響物源發(fā)育的主控因子，其中包括溝域面積、相對高差、發(fā)震斷裂距離、溝域巖性及巖體結(jié)構(gòu)等。通過大量統(tǒng)計數(shù)據(jù)分析了在控制變量作用下各因子與物源量的函數(shù)關(guān)系，在此過程中建立了溝域巖體質(zhì)量等級評價模型；

(3)在各因子顯著性分析的基礎(chǔ)上，將以上建立的各因子的函數(shù)關(guān)系概化為多元非線性方程，通過MATLAB程序語言將極震區(qū)三大片區(qū)泥石流物源因子及物源量的28組數(shù)據(jù)代入方程進(jìn)行回歸分析，從而確定回歸方程的系數(shù)，最終通過相關(guān)性驗(yàn)證及誤差分析；

(4)通過極震區(qū)擬合數(shù)據(jù)以外的另外6條溝對模型適宜性進(jìn)行檢驗(yàn)，檢驗(yàn)結(jié)果顯示模型適宜性較強(qiáng)。

[1]殷躍平.汶川八級地震地質(zhì)災(zāi)害研究[J].工程地質(zhì)學(xué)報,2008,16(4): 433-444.

[2]馮文凱,黃潤秋,許強(qiáng).地震波效應(yīng)與山體斜坡震裂機(jī)理深入分析[J]. 西北地震學(xué)報.2011,33(1):20-25.

[3]周必凡.泥石流防治指南[M].北京:科學(xué)出版社,1991.

[4]蔣忠信.泥石流固體物質(zhì)儲量變化的定量預(yù)測[J].山地研究,1994,12(3):155-162.

[5]喬建平，黃棟，楊宗佶，孟華君. 汶川地震極震區(qū)泥石流物源動儲量統(tǒng)計方法討論[J]. 中國地質(zhì)災(zāi)害與防治學(xué)報，2012，6(2):1-6.

[6]Chuan Tang, Zhu Jing. Approach for urban debris flow risk assessment [J]. Advances in Water Scient ,2006,17(3):67-71.

[7]唐川，朱靜.城市泥石流災(zāi)害風(fēng)險評價探討[J].水科學(xué)進(jìn)展，2006,17(3):67-71.

[8]Newmark N M. Effects of earthquakes on dams and embankments [J].Geotechnique 1965，15(2)：139-160.

[9]張倬元,王士天,王蘭生,黃潤秋,許強(qiáng),陶連金.工程地質(zhì)分析原理(第三版)[M].北京:地質(zhì)出版社,2008.

[10]黃潤秋,李為樂.汶川大地震觸發(fā)地質(zhì)災(zāi)害的斷層效應(yīng)分析[J].工程地質(zhì)學(xué)報,2009,17(1):19-27.

[11]鄒威.強(qiáng)震作用下均質(zhì)斜坡動力響應(yīng)的大型振動臺試驗(yàn)研究[D].成都理工大學(xué),2011.

[12]于海英,王棟,楊永強(qiáng),盧大偉,解全才,張明宇,周寶峰,江紋鄉(xiāng),程翔,楊劍.汶川8.0級地震強(qiáng)震動特征初步分析[J].震災(zāi)防御技術(shù).2008,3(4):321-336.

[13]Huang Runqiu ,Li Weile.Fault effect analysis of geohazards triggered by Wenchuan Earthquake ,China[J]. Bull. Of Engineering Geology,2009, 17(1):19-88.

[14]許樹柏.實(shí)用決策方法-層次分析法原理[M].天津:天津大學(xué)出版社.1988.

[15]袁進(jìn)科.斜坡震裂巖體結(jié)構(gòu)特征與震后崩塌識別體系研究[D].成都理工大學(xué),2012.

[16]崔圣華,裴向軍,黃潤秋.直線型斜坡滾石運(yùn)動速度特征研究[J].工程地質(zhì)學(xué)報.2013,21(6):912-919.

[17]許強(qiáng),李為樂.汶川地震誘發(fā)大型滑坡分布規(guī)律研究[J].工程地質(zhì)學(xué)報.2008,3(4):818-826.

[18]顧文韜,裴向軍,裴鉆,張雄,李天濤.極震區(qū)泥石流物源特征研究[J].自然災(zāi)害學(xué)報.2015,24(2):107-114.

[19]唐川,梁京濤.汶川震區(qū)北川9.24暴雨泥石流特征研究[J].工程地質(zhì)學(xué)報,2008,12(6):751-758.

2016-03-07

任志剛(1972-)，男，四川仁壽人，高級工程師，從事水力水電工程勘測設(shè)計工作。

P642.23

B

1003-9805(2016)03-0056-06