地震災害風險普查中地震危險性計算效率

李得博 李井岡 趙艷南

1 中國地震局地震研究所,武漢市洪山側路40號,430071 2 中國地震局地震大地測量重點實驗室,武漢市洪山側路40號,430071

地震危險性分析是對工程場地未來可能遭遇的強地震動進行估計,是在工程場址周圍地震活動性和地震地質環境研究的基礎上估計工程場址的設計地震動參數。評定某個工程場址的地震危險性須考慮3個方面:1)確定該場址周圍未來一定時間內可能發生破壞性地震的地點、強度和性質;2)通過強震觀測數據和歷史地震資料,了解地震波傳播途徑的區域特性,確定地震動衰減關系;3)詳細勘察場地條件,分析局部場地條件對地震動的影響[1]。目前地震危險性分析方法大體上分為確定性方法和概率方法[2]兩大類,其中概率地震危險性分析方法由章在墉等[3]率先引入國內,經高孟潭[4]和時振梁等[5]改進。改進后的方法能夠反映我國地震活動的時空不均勻性,被廣泛應用于地震安全性評價和地震小區劃等工作中。

地震災害普查中要求按國家級30″×30″網格、省級6″×6″網格、市級3″×3″網格、縣級1″×1″網格劃分經緯度計算地震危險性,計算控制點數量達到千萬量級。若每個計算點耗時0.5 s,則省級計算點約1 000萬個,耗時約1 300 h,因此優化相關計算方法是一個重要的研究課題。大空間網格的計算結果顯示,各相鄰點之間的危險性計算結果變化較為平穩,若能夠通過大間隔網格控制點的危險性結果插值獲得可靠的小空間間隔控制點的地震危險性結果,將大大節省計算時間。基于此,本文先計算120″、60″、30″、15″、9″和6″經緯度劃分網格的地震危險性,然后插值得到較小網格的地震危險性結果,最后比較直接計算的小網格地震危險性結果和大空間網格插值得到的小網格地震危險性結果,評估由較大空間間隔插值所獲結果的可靠性、適用性及使用要求,以討論地震災害風險普查中兼顧可靠性和效率的計算地震危險性的最優插值網格。

1 概率地震危險性分析方法

概率地震危險性分析方法(PSHA)于20世紀80年代引入我國,之后逐漸形成了具有中國特色的考慮地震時空分布不均勻性的概率性地震危險性分析方法(CPSHA),GB 17741-2005《工程場地地震安全性評價》和GB 18036-2015《中國地震動參數區劃圖》(簡稱五代區劃圖)都采用CPSHA方法,其計算步驟如圖1所示,包含確定地震區帶及其地震活動性參數、選定地震動衰減關系和進行地震危險性計算等主要步驟。

圖1 概率地震危險性分析步驟Fig.1 Steps of probabilistic seismic hazard analysis

根據地震風險普查規范,區域的標準網格分為國家、省、市、縣4種,國家級網格為30″×30″(約1 km),省級網格為6″×6″(約200 m),市級網格為3″×3″(約100 m),縣級網格為1″×1″(約30 m)。其中30″標準網格文件可通過軟件平臺下載,6″標準網格嚴格按照標準制作,地震、地質、氣象等所有災種統一使用。標準網格基準點(西北角)為73°E、54°N,行列號向東及向南逐漸增加;標準網格制作中,行列號根據給定公式進行計算,作為對計算點的標簽,使用網格點中心點坐標(xcenter,ycenter)計算地震危險性;標準網格要求經緯度保留小數點后8位,空間位置與屬性信息字段的經緯度位置保持一致;分別計算50 a 63%、50 a 10%、50 a 2%和100 a 1%共4個概率水準的峰值加速度(PGA),不計算反應譜,給出基巖和場地調整的PGA,其中場地條件采用全國1∶100萬宏觀場地類別數據,場地調整系數參考《地震危險性圖編制規范FXPC/DZ P-01》,場地類別根據標準網格中心點進行判別,地震危險性分級方法按照《地震危險性圖編制規范FXPC/DZ P-01》以100 a 1%概率結果作為分級指標(不使用50 a 10%),分級排序采用降序(危險性最高為1級,最低為4級)。該策略的目的一方面是為了避免使用50 a 10%造成混淆誤解;另一方面,萬年一遇可反映計算點可能遭受的最大地震動。潛在震源區、地震活動性參數和衰減關系使用五代區劃圖資料,不考慮近幾年的潛源調整。

因本文主要關注計算效率和插值方法的有效性問題,故區域范圍統一采用計算點周圍150 km范圍截取潛在震源的方法確定,暫不考慮150 km范圍外遠大震的影響,PGA計算結果統一保留小數點后1位。

2 測試區域選取

測試區范圍為34.55°～35.88°N、114.80°～116.42°E(圖2),主要基于以下考慮:1)根據五代區劃圖,區域場地地震動分類較全,0.05～0.20g都有,便于比較不同等級分解處插值是否合理;2)地震安全性評價報告資料豐富,便于對比;3)區域大小合適,采用不同網格間隔進行劃分,網格點數跨越幾千到幾百萬量級,利于全面研究不同量級的計算效率和精度。

圖2 測試區位置及120″間隔網格點分布Fig.2 Location of test area and distribution of 120″ interval grid

計算區間及網格劃分參數為:全國網格劃分統一起點緯度為54°,起點經度為73°,網格間隔δ為3″ ,分別使用120″、60″、30″、15″、9″、6″及3″進行網格劃分;行號起點為10 001+int[((ycenter-54)×3 600-δ/2)/δ+0.5],列號起點為10 001+int[((xcenter-73)×3 600-δ/2)/δ+0.5],式中,int表示取整數。按照規范,不同級別的網格行列號的起點數不同,本文不作區分,只要能夠通過行列號找到對應的點即可。

單個控制點地震危險性分析計算輸出結果保存在以“點序號-列號-行號”命名的單獨文件夾中,其中包括5個周期(0 s,0.2 s,1.0 s,2.0 s,6.0 s)的危險性結果及相應潛源的貢獻、150 km范圍內五代區劃圖中潛源分布(圖3)、各潛源的空間分布函數和風險普查要求的4個概率水準PGA。前期不同周期的危險性結果和潛源貢獻分別存放在不同的文件中;后期為減少文件讀寫數目,盡可能提高計算效率,將危險性結果和潛源貢獻結果合并存儲在同一個文件中。實際計算表明,該措施對計算速度的影響不大,但可減少小文件的數目,提高計算機存儲空間的利用率。風險普查要求計算4個概率水準(50 a 63%、50 a 10%、50 a 2%和100 a 1%)下5個周期的譜加速度值,由于數據量巨大,且不同周期、不同場地條件下結果類似,本文僅展示基巖場地50 a超越概率10%的PGA比較結果。

圖3 測試區內某個計算點周邊150 km范圍內五代區劃圖潛源分布Fig.3 Distribution of potential sources within 150 km of a calculation point in the test area in five generations zoning map

3 提高計算效率的方法及相關結果

將測試區劃分出120″、60″、30″、15″、9″、6″和3″間隔的網格,計算所有網格點各種概率水準下的PGA。如表1所示,不同間隔劃分的總網格點數隨著劃分間隔的減小呈指數增長,相應的計算耗時也大幅增長。首先分別計算不同空間間隔的各網格點的地震危險性,然后根據較大空間間隔的計算結果,按照空間距離平方倒數加權的插值方式得到所有較小空間間隔點的PGA,最后將插值所得PGA與直接計算的網格點PGA(準確值)進行比較,計算誤差最大值、最小值及標準差(表2,單位Gal)。

表1 不同空間間隔的計算量和耗時Tab.1 The computation amount and time of different space intervals

表2 不同大網格插值結果的誤差最大值、最小值、標準差Tab.2 The maximum, minimum and standard deviation error of interpolation results with different larger grids

表1給出不同處理器的計算耗時,其中,CPU為i7-11700的臺式機計算了全部數據,其他幾種處理器只處理了15″以上間隔的數據,更小間隔的數據是根據單點速度乘以計算量估算出來的。

由于每個點地震危險性計算相互獨立,不存在數據交叉,可以完全并行計算而不沖突。程序采用多線程,總體來看,每個程序的線程數對速度影響并不大。單個程序運行時,CPU為i7-9850H的移動工作站計算速度最慢,計算效率約5 000點/h,CPU為i7-11700的臺式機計算速度最快,計算效率約10 000點/h。為進一步提高計算速度,測試同一臺計算機采用2、4、8個相同程序同時計算不同區域的方法。結果表明,計算速度總和有所提高,但無論采用多少個應用程序,單臺計算機總最高計算速度約為18 000點/h,這個速度上限可能是由計算機的CPU處理能力決定的。

從表2和部分大網格插值得到的小網格結果的誤差分布(圖4～7)可以看出:1)利用120″網格插值結果的誤差最大值都達到2.0 Gal,這可能是由大網格相鄰2點危險性計算值的差值決定的;60″網格插值結果的最大誤差為1.9 Gal;30″網格插值結果的最大誤差為0.7 Gal;15″ 網格插值結果的最大誤差為0.4 Gal;9″ 網格插值結果的最大誤差為0.3 Gal;6″ 網格插值結果的最大誤差為0.2 Gal。2)所有大網格插值結果的誤差的標準差都較小,說明絕大多數點的誤差較小,插值結果誤差分布也說明了這一點。3)對于國家級30″ 網格,采用60″ 網格的計算值進行插值,誤差小于1 Gal;對于省級6″ 網格,采用30″、15″ 和9″ 網格的計算值進行插值,誤差小于1 Gal;對于市級3″ 網格,采用30″、15″、9″ 和6″ 網格的計算值進行插值,誤差小于1 Gal;同理,對于縣級1″網格,有理由推測采用30″、15″、9″、6″和3″ 網格的計算值進行插值,誤差也小于1 Gal。

圖4 不同大網格插值得到15″網格PGA的誤差分布Fig.4 The error distribution of 15″ grid PGA interpolated from different larger grids

圖5 不同大網格插值得到9″網格PGA的誤差分布Fig.5 The error distribution of 9″ grid PGA interpolated from different larger grids

圖6 不同大網格插值得到6″網格PGA的誤差分布Fig.6 The error distribution of 6″ grid PGA interpolated from different larger grids

圖7 不同大網格插值得到3″網格PGA的誤差分布Fig.7 The error distribution of 3″ grid PGA interpolated from different larger grids

為進一步研究插值誤差產生的原因,將誤差較大的點投影到地圖上。圖8給出由60″網格地震危險性計算值插值得到3″網格地震危險性結果時誤差大于0.5 Gal和大于1.0 Gal點的分布,其他大網格計算值插值得到較小網格結果的情況與之類似。圖8(a)顯示,誤差較大的點主要位于分區邊界處、較復雜處及計算區域邊界處;圖8(b)顯示,邊界處誤差較大,這主要是由于待插值的點位于計算值網格的外側,而插值約束點非常少,此時采用空間距離反距加權方法不能反映變化的趨勢性。可以通過擴大計算區范圍,將所有需要插值的點都作為內點進行插值,或對于計算區域外邊界點采用趨勢外推插值的方法進行插值。由此可見,目前地震風險普查實踐中一般要求計算點取行政邊界外圍擴大10 km的標準格網,以減少邊界效應對插值結果影響的做法是合適的。

圖8 由60″網格插值得到3″網格的PGA誤差分布Fig.8 The PGA error distribution of 3″ grid interpolated from 60″ grid

4 結 語

本文基于五代區劃圖的潛在震源區劃分方案,采用全國地震災害風險普查統一的危險性計算方法,分別獲得了120″、60″、30″、15″、9″ 和6″ 經緯度網格插值得到較小網格的地震危險性結果,并將其與直接計算的小網格地震危險性結果進行比較。結果表明,30″ 網格計算值插值得到各種小網格的地震危險性結果誤差都小于1 Gal,是兼顧精度和速度比較合理的網格間隔。

由于誤差較大的點主要位于地震動區劃圖分區邊界處、較復雜處及計算區域邊界處,通過擴大計算區范圍,將所有需要插值的點都作為內點進行插值,可以解決計算區域邊界處誤差較大的問題。當然,在不需要等待很長時間的情況下,還是選擇盡可能小的網格進行計算,特別是潛在震源區比較復雜的區域。

在實際地震災害風險普查的危險性計算中,根據精度要求和時間因素進行一些前期計算實驗,綜合評估精度和效率,再確定計算方案是很有必要的。值得注意的是,本文主要是研究地震災害危險性普查中地震危險性計算的效率問題,對地震危險性計算中的地震潛在震源、地震動衰減關系分析等都進行了簡化,因此本文結果不能直接應用于測試區的地震安全性評價。