當量和爆高(埋深)反演的聲震分析方法

張亮永,張德志,肖衛國,梁旭斌,王同東,郭權勢,李 翱

(強脈沖輻射環境模擬與效應全國重點實驗室;西北核技術研究所：西安710024)

在國防和民用領域中,近地面爆源參數反演是非常重要的研究內容,快速預測爆炸當量和爆高(埋深)(height of burnst, HOB)等源參數,對爆炸事故和恐怖襲擊預警[1-4]、地震和火山噴發事件預報[5-7]、核查監測[8-10]及武器性能評估[11-12]等方面具有重要意義。采用聲學數據(近區為沖擊波)[13-21]或地震波數據(近區為地運動)[2-3, 22-24]預測爆源參數是常見的兩種近地面爆源參數反演方法,但采用這兩種方法反演的爆炸當量和爆高(埋深)分布在特定范圍內存在多個極值點(稱為折中關系[22-27]),給源參數估計造成了極大的不確定性。為提高源參數預測的精度,研究人員提出了聲震分析方法[22-27],該方法基于數據融合算法綜合地震波(地運動)數據和聲學(空氣超壓)數據反演爆源參數,通過同時利用聲震波動信息對源參數進行多重約束,減小折中關系的影響。聲震分析方法分為聲學模型、地震波模型和數據融合算法3個部分,本文將對這3部分目前的進展展開論述。

1 聲學模型

由于近地面爆炸能量會耦合到空氣中產生明顯的聲擾動,一般通過遠場聲學數據,利用超壓峰值、正向聲沖量及正向脈寬等波形特征量來預測地面爆炸當量[22-27],研究人員相繼提出適合快速計算,以超壓、聲沖量和持續時間為特征量的半經驗聲學模型及描述波形細節的全波形聲學反演方法。半經驗聲學模型包括ANSI模型[28, 29]、KG85標準模型[30]、BOOM模型[15, 31]和IPM聲學模型[19]等自由場聲學模型及地面爆炸聲學模型[2]和近地面爆炸聲學模型[22-26]。半經驗自由場聲學模型未考慮爆炸能量耦合到地下的影響,對地面爆炸事件通常假設兩倍等效爆炸當量進行爆源參數估計[20, 28-30, 32];而對于近地面爆炸事件,Zhang等[33-34]通過考慮地面反射、耦合作用及大氣環境因素,建立了基于自由場聲學模型的近地面聲學模型(簡稱為自由場近地面聲學模型),進行了爆炸當量和爆高估計。表1為不同的自由場聲學模型。地面爆炸聲學模型[2]是基于地面爆炸聲學數據建立的模型,考慮了地面爆炸能量耦合情況,用于地面爆炸事件的當量估計,但不適用于偏離地面位置的爆炸事件當量估計。而近地面爆炸聲學模型考慮了近地面聲震能量耦合特性,可用于近地面爆炸事件的當量和爆高(埋深)估計,但不適用比爆高為-2.4 ～-0.6 m·kg-1/3的淺埋情況。上述兩種模型需根據不同場地介質提前確定多個待定系數,適用的場地介質有限。Zhang等[25]在近地面爆炸聲來源機制[35-36]的基礎上,通過考慮氣體噴出地表和地沖擊耦合影響對近地面聲學模型進行改進,拓寬了比例爆高適用范圍,同時確定了花崗巖夾雜變質巖場地的模型待定系數。圖1為近地面爆炸聲來源機制[35-36]和改進前后近地面聲學模型[25]。

(a)Mechanism of near-surface explosion producing sound [35-36]

(b)Unimproved and improved near-surface acoustic models[25]

表1 自由場聲學模型

全波形聲學反演方法包括超壓波形模板匹配法[16]、經驗聲源模型法[20]和折合聲沖量法[14, 17, 24]等。超壓波形模板匹配法[16]以不同比距離下隨比時間變化的大量壓力波形數據為模板,將未知當量的全壓力波形和模板壓力波形進行匹配,對應匹配偏差最小的當量為估計當量,但該方法需大量不同當量不同比距離的壓力波形組成的樣本庫,但目前樣本數據庫還不夠充分;經驗聲源模型法[20]以G17HE為標準源模型,通過大量實測視壓力數據并假設標準高斯分布,獲取源模型的峰值壓力和持續時間平均值和方差(偏差分布),由此計算爆炸當量,但該模型聲測點數據需選取20～100 m·kg-1/3的比距離范圍,在未知當量時不好提前預設測點距離。圖2為超壓波形樣本庫[16]以及G17HE源模型和實測波形[20]。Kim等[14, 17]基于線性聲學理論反演爆源等效源時間函數,以KG85模型折合聲沖量作為參考函數,提出一種地面爆炸當量反演思路(本文稱為折合聲沖量法),并將其應用于剛性地面爆炸當量反演,具有較高的當量預測精度,但該方法只適用剛性地面爆炸當量估計,沒有考慮波傳播過程中地表結構的影響。圖3為等效源時間函數及預測和實測波形。張亮永等[37]對該方法進行了發展,通過考慮波傳播過程中地表結構的吸收衰減影響及近地面爆炸聲震能量分配過程,建立了多孔彈性場地近地面折合聲沖量模型,解決了多孔彈性場地的當量估計問題,給出當量相對偏差為7%,爆高相對偏差不超過19%。圖4為多孔彈性場地聲傳播示意圖和爆源參數反演結果。

(a)Pressure wave template[16] producing sound [35-36]

(b)Source model of G17HE with measured waveforms [20] acoustic models[25]

圖3 等效源時間函數(頂端黑色曲線)及預測(紅色曲線)和實測(黑色曲線)波形[14]

(a)Sound propagation diagram in three-layer medium model

(b)Yield and height of burst predicted by seismoacoustic analysis

2 地震波模型

近地面爆炸耦合的地震波(近距離稱為地沖擊或地運動)來源于爆炸直接耦合方式和空氣超壓(聲波)耦合方式,對近距離前者稱為直接地沖擊,后者稱為感生地沖擊[12, 38-39],而對遠距離,Albert等[40]將前者耦合的地震波稱為前驅地震波,后者耦合的地震波稱為聲耦合地震波。圖5為近地面爆炸耦合的地震波來源方式示意圖。

圖5 近地面爆炸耦合的地震波來源方式示意圖[40]

目前,近地面地震波模型以前驅地震波模型為主,隨著DM/DB[2],HRI-III[16, 23, 41],SAY[22, 42],FSE和SPE[24]等系列化爆試驗的開展,相繼建立了沖積土、沉積巖和花崗巖等場地介質的近地面前驅地震波爆炸當量估計模型[23-24, 26, 42]。其中,Templeton 等[42]基于HRI-III 和SAY 化爆試驗地震波數據對地震波模型進行了深入討論,對比和分析了沖積土和沉積巖兩種地質成分下的近地面爆炸地震波模型,建立了地震波數據分析和模型建立的基本方法,對近地面化爆地震波建模具有指導意義。圖6為不同場地介質的近地面地震波模型和實測數據。但上述化爆試驗場地介質類型主要是沖積土、沉積巖和花崗巖等,受試驗地質成分的限制,地震波模型局限在以上幾種地質成分,數據還不夠充分,模型還不夠完善[22, 24, 42-43],且未考慮爆源附近地質和其余位置地質不一樣的情況。

(a)Alluvium

(b)Sedimentary rock

Zhang等[25, 34, 44]基于試驗數據探討了花崗巖夾雜變質巖場地和土石混合場地的近地面地震波實測數據和模型,分別如圖7和圖8所示。

(a)Seismic data

(b)Seismic model

(a)Model one

(b)Model two

(c)Model three

(d)Model four

(e)Model five

為解決爆源附近地質和其余位置地質不一樣的混合場地地震波模型建模問題,給出了模型的待定系數,并通過引入放大系數。圖9為放大系數K不同取值對反演結果的影響。圖9中,綠色、紅色虛線分別為聲耦合地震波模型和地震波位移首峰值模型的折中曲線(最小偏差對應的當量和埋深),白色五角星為交叉和最小偏差位置。表2為土石混合場地不同模型的定義。

(a)K=2

(b)K=4

(c)K=6

(d)K=10

表2 土石混合場地模型的定義[44]

上述近地面前驅地震波模型,基于初至波首峰值建立當量和爆高(埋深)關系,適用中近區爆炸當量估計。對于遠區情況,Pasyanos等[10, 24]基于近地面爆炸地震波能量耦合關系,采用地震波尾波包絡法,建立近地面遠區地震波模型,進一步拓展了近地面地震波模型的適用范圍。圖10為尾波包絡法示意圖[46]和典型反演結果[24]。

(a)Schematic diagram of coda-wave method

(b)Typical inversion results

聲耦合地震波模型可解決單一地震計的近地面爆源當量估計問題,為爆源參數反演提供除前驅地震波外的新途徑,但相關研究較少。Zhang等[46-47]基于聲耦合系數建立聲耦合地震幅值和聲學幅值之間的關系,確定沙土場地的聲震耦合系數,并聯合近地面聲學模型建立聲耦合地震波模型。沙場地聲耦合地震波典型波形及耦合系數衰減特性分別如圖11和圖12所示[46]。表3為聲耦合地震波特征值的定義[46]。

(a)East

(b)North

(c)Verticle

(a)Z1

(b)ENZ1 1&2

(c)ENZ1 0&1

(d)Z2

(e)Zhalf

(f)ENZ2 1&2

(g)ENZhalf1 1&2

(h)ENZ2 0&1

(i)ENZhalf1 0&1

表3 聲耦合地震波特征值的定義[46]

3 數據融合算法

近地面爆炸能量會耦合到大氣和地介質中,產生大氣擾動(遠場為聲波)和地運動(遠場為地震波)。圖13為近地面爆炸產生的典型聲震信號[23]。

(a)Seismo-acoustic energy distribution

(b)Seismic data of HRI and HRII

(c)Acoustic data of HRI and HRII

采用聲學數據或地震波數據對源參數進行反演,當量和爆高(埋深)之間通常會存在明顯的折中關系使求解不適定。而聲震分析方法通過融合地震波(地運動)數據和聲學(空氣超壓)數據使得反演問題變得適定,從而得到當量和爆高(埋深)估計值。圖14為典型試驗 (FSE-4) 的折中關系和聲震分析結果[24]。圖14中,五角星為真實值,空心三角形為地表處的最優點,實心三角形為整個區域的最優點。

(a)Results by seismic data

(b)Results by acoustic data

(c)Results by seismic and acoustic data

聲震融合算法是聲震分析的關鍵組成,隨著聲震分析方法的不斷改進,聲震融合算法從相對偏差法[22]、對數偏差法[23]及折中曲線交叉法[23]到似然函數法[24]和貝葉斯MCMC方法[25-26, 48]再到機器學習方法[27],得到了不斷發展。表4為數據融合算法特點和計算公式。

表4 數據融合算法特點和計算公式

圖15為不同數據融合方法的典型反演結果。表4列出的前4種方法沒有考慮模型參數的先驗信息,源參數缺少進一步約束,且計算采用格點搜索法,效率較低。貝葉斯馬爾可夫鏈蒙特卡羅 (Markov chain Monte Carlo,MCMC)方法[25-26, 48]同時考慮了偏差信息和先驗信息,采用統計的方法對爆源參數進行估計可給出更合理的統計解釋,同時采用MCMC方法獲得較高的求解效率。Johannesson等[48]給出了聲震融合的貝葉斯MCMC反演框架,Zhang等[25]和Ford等[26]進一步發展了該方法,基于MCMC方法求解源參數,并應用于近地面爆源參數估計。圖16為貝葉斯MCMC方法求解流程。以上數據融合方法基于已有的聲震模型反演爆源參數和不同場地巖土類型具有不同的聲震模型。數據融合方法的反演精度和聲震模型密切相關。研究表明[42, 45],聲震模型和場地條件的不匹配將導致較大的反演偏差。為減小因聲震模型匹配程度帶來的反演偏差,提高反演模型在不同場地的適用性,Stroujkova等[27]采用人工神經網絡(artificial neural network,ANN)和貝葉斯正則化方法融合小樣本聲震數據來估計當量和爆高(埋深),嘗試將機器學習算法應用在爆源參數反演。圖17為基于深度學習的爆源參數反演流程。

(a)Method of relative deviation[22]

(b)Method of logarithmic deviation[23]

(c)Tradeoff curve crossing method[23]

(d)Likelihood function method[24]

(e)Bayesian MCMC method[25, 26, 48]

(f)Machine learning method[27]

圖16 貝葉斯MCMC方法求解流程[25]

圖17 基于深度學習的爆源參數反演流程[27]

4 總結與展望

經過十幾年的研究,發展了多種近地面聲學模型和地震波模型,適用沖積土、沙土、沉積巖、花崗巖和土石混合體等多種地質條件,聯合相對偏差法、折中曲線交叉法和貝葉斯MCMC等數據融合算法,或基于聲震數據采用機器學習方法,可有效解決多種場地近地面爆源的爆炸當量和爆高(埋深)反演問題。聲震分析方法雖得到了快速發展,但還需重點解決以下問題：

(1)近區聲/震模型修正完善。受到場地條件限制,不同地質、不同爆高和不同距離的聲震數據有限,還需更多的試驗數據對近地面聲/震模型進行補充、修正和完善,從而獲取適用不同地質成分的通用模型。目前,聲/震數據主要通過外場試驗獲取,需進一步發展基于微藥量化爆模擬及超重力離心機等平臺的實驗室模擬技術。

(2)遠區聲/震模型快速建模。區域較大時,受地質地形和大氣環境等因素的影響,近地面爆炸耦合的地震波和聲波表現出復雜的傳播特性,給遠區爆源參數快速準確反演帶來了巨大挑戰。目前,聲震分析技術以近區為主,需進一步發展大區域快速建模技術,實現遠區爆源參數的快速準確反演。

(3)基于機器學習的快速反演方法。機器學習可有效解決數據融合和快速計算問題,但該方法在聲震分析應用較少,需基于機器學習算法,進一步發展小樣本聲震數據反演技術及基于物理機制的大區域爆源參數反演技術。

(4)未知爆源位置的聲震分析方法。目前聲震分析研究是在爆源位置已知條件下開展的,未考慮爆源位置估計偏差的影響,需進一步開展基于聲震數據的爆源位置、當量和爆高(埋深)的綜合反演方法研究?？赏ㄟ^貝葉斯MCMC或機器學習方法建立爆源多參數的綜合反演框架。