基于產(chǎn)狀Fisher分布的巖石節(jié)理最小樣本容量確定方法

楊梓楓，夏華宗，杜玉芳，王曉明，3

1.河北地質(zhì)大學 城市地質(zhì)與工程學院，河北 石家莊 050031；2.河北地礦建設工程集團有限責任公司，河北 石家莊 050000；3.河北省地下人工環(huán)境智慧開發(fā)與管控技術創(chuàng)新中心，河北 石家莊 050031

0 引言

節(jié)理是巖體中普遍發(fā)育的一種結構面，這些結構面將完整巖體切割成大小不等、形狀各異的巖石塊體，而且節(jié)理通常強度低、易變形，同時能夠為地下水提供存儲空間和運移通道，對巖體的完整性、力學特性、滲流特性及穩(wěn)定性起著控制作用。節(jié)理產(chǎn)狀是影響巖體結構及其工程性質(zhì)的關鍵參數(shù)之一，其平均值及離散性對邊坡、地下洞室等工程巖體的穩(wěn)定性及可靠性具有顯著影響[1，2]。因此，準確計算節(jié)理的平均產(chǎn)狀及其離散性對掌握巖體的工程特性進而開展工程巖體的評價、設計及施工具有重要意義。

節(jié)理平均產(chǎn)狀的準確性受樣本容量(即采樣數(shù)量)和離散性的影響，產(chǎn)狀離散性越大，則節(jié)理的采樣數(shù)量應越大，反之亦然。為保證統(tǒng)計結果的準確性，節(jié)理樣本應不少于某一數(shù)量，該數(shù)量即為最小樣本容量。在實際工程中，受采樣條件的限制，如節(jié)理露頭較少或測量困難而造成采樣數(shù)量較少，平均產(chǎn)狀的統(tǒng)計準確性勢必會受到影響。因此，確定節(jié)理的最小樣本容量對于節(jié)理統(tǒng)計分析具有現(xiàn)實意義。

巖體中的節(jié)理往往具有一定的方向性，通常會沿2～5個優(yōu)勢方向發(fā)育，產(chǎn)狀相近的節(jié)理應劃分為同一組。節(jié)理的產(chǎn)狀一般用傾向(或走向)和傾角表示，同一組節(jié)理的產(chǎn)狀通常服從一定的分布形式。常見的節(jié)理分布類型包括Fisher分布、Bingham分布和雙正態(tài)分布等[3-6]。其中，F(xiàn)isher分布由于僅包含1個參數(shù)且為可積分函數(shù)便于實現(xiàn)Monte Carlo模擬，得到了廣泛應用[7-12]。此外，一些研究表明，實測節(jié)理產(chǎn)狀能夠較好地擬合Fisher分布[13，14]。由于節(jié)理發(fā)育的隨機性和離散性，在工程實踐中通常將節(jié)理進行分組統(tǒng)計來確定每組節(jié)理的平均產(chǎn)狀[15-17]。

基于節(jié)理產(chǎn)狀服從Fisher分布，本文探討產(chǎn)狀離散性對節(jié)理采樣數(shù)量的影響。通過Monte Carlo隨機模擬技術產(chǎn)生具有不同離散性的產(chǎn)狀隨機數(shù)，在保證準確性的條件下確定節(jié)理的最小樣本容量，并建立最小樣本容量與產(chǎn)狀離散性的經(jīng)驗關系，最后通過算例驗證其適用性。

1 Fisher分布及產(chǎn)狀離散性

該分布于1953年由Fisher提出，是一種將節(jié)理產(chǎn)狀視為球面上關于平均產(chǎn)狀對稱的分布形式[3]。Fisher分布的概率密度函數(shù)為：

(1a)

(1b)

式中，θ′和φ′分別是當把坐標系z軸旋轉至節(jié)理的平均矢量方向(θm，φm)后在新坐標系中的方向角；參數(shù)κ反映了節(jié)理產(chǎn)狀分布的離散性或集中程度，κ越大表明節(jié)理的產(chǎn)狀圍繞平均產(chǎn)狀越集中，即離散性小。如果κ=0，則Fisher分布轉化為球面上的均勻分布。均質(zhì)巖體中的構造節(jié)理產(chǎn)狀通常較為集中，具有較大的κ值；而非均質(zhì)巖體中的節(jié)理產(chǎn)狀往往具有較大的離散性，其κ值通常較小。κ可按下式計算[3]：

(2)

式中，N為節(jié)理的數(shù)量；rN為節(jié)理的和向量，其方向即為節(jié)理的平均產(chǎn)狀，其計算過程將在2.2節(jié)中詳述。Watson的研究認為當κ>3時，公式(2)的計算結果具有足夠精度[18]。眾多研究表明[10-12]，同一組節(jié)理的κ一般大于3，因此公式(2)具有較好的適用性。

2 研究方法

2.1 Fisher分布產(chǎn)狀的隨機生成

節(jié)理的方向是以其單位法向量表示的，一般取指向上方的法向量。在直角坐標系中(圖1)，假設OP為某節(jié)理的單位法向量，點N為P點在xOy平面上的垂直投影，θ為OP與z軸正方向的夾角，φ為x軸正方向沿逆時針方向與ON的夾角。

Fisher分布是當把z軸旋轉至節(jié)理的平均矢量方向(θm，φm)之后的形式，此時它的兩個變量θ′和φ′是獨立的，其概率密度函數(shù)見公式(1)。其中φ′服從(0，2π)之間的均勻分布，可直接利用公式(3)產(chǎn)生，u為(0，1)范圍內(nèi)的均勻分布隨機數(shù)。

φ′=2πu

(3)

θ′的密度函數(shù)是一個可積分函數(shù)，其累計分布函數(shù)為：

(4)

這樣可利用逆變換法產(chǎn)生Fisher分布的隨機數(shù)，即：

(5)

得到：

(6)

如果稱旋轉前(z軸鉛直向上)的坐標系為舊坐標系，它的變量為(φ，θ)；稱旋轉后(z軸與節(jié)理的平均矢量一致)的坐標系為新坐標系，它的變量為(φ′，θ′)，則兩個坐標系的變化關系為[7，9]：

式中：(l，m，n)和(l′，m′，n′)分別為對應(φ，θ)和(φ′，θ′)的方向余弦。l，m，n分別為：

(8)

地質(zhì)學中，節(jié)理的方向是用傾向和傾角(α，β)表示的，而非(φ，θ)。如果采用圖3中的坐標系，(α，β)和(φ，θ)的關系可以用下式計算：

(9)

綜合上述分析，產(chǎn)生一組N個Fisher分布的隨機節(jié)理產(chǎn)狀的過程如圖2所示。

圖2 Fisher分布隨機節(jié)理產(chǎn)狀的產(chǎn)生流程Fig.2 Generation flowchart of random joint orientation obeying Fisher distribution

2.2 節(jié)理平均產(chǎn)狀計算

節(jié)理平均產(chǎn)狀是指節(jié)理組中每一條節(jié)理面的單位法向量之和所在的方向。假設一組節(jié)理包含N條節(jié)理，每一條節(jié)理的傾向和傾角分別為φi和θi，通過公式(8)可以得到每條節(jié)理的單位法向量。由于節(jié)理產(chǎn)狀的離散性，節(jié)理產(chǎn)狀的分布范圍有時較大，當采用上半球赤平極射投影時可能出現(xiàn)“超越上半球投影的產(chǎn)狀”(OEBEUHP)[9，13]。圖3為200個隨機產(chǎn)生的服從Fisher分布的節(jié)理產(chǎn)狀，其平均產(chǎn)狀為45°∠60°，κ=10，其中11個節(jié)理極點為OEBEUHP。

圖3 平均產(chǎn)狀為45°∠60°的節(jié)理產(chǎn)狀極點圖(κ=10)Fig.3 Pole diagram of joint orientations with mean orientation of 45°∠60°(κ=10)

由于存在OEBEUHP，同一組節(jié)理向上的法向量之間的夾角可能大于90°，所以在計算節(jié)理平均產(chǎn)狀時，不能將節(jié)理的法向量ni簡單相加，應考慮采用ni還是-ni。鄭俊等提出了一種識別和調(diào)整OEBEUHP的策略，通過計算ni與合向量rN的夾角來識別OEBEUHP，并用-ni重新計算rN[9，13]。實際上，初始rN的計算已經(jīng)存在由OEBEUHP引起的問題，因此這種策略需要適當調(diào)整。本文在此基礎上，提出了一種計算rN的修正方法，具體過程見圖4。根據(jù)得到的rN，可利用公式(8)和(9)將其轉換為節(jié)理的平均產(chǎn)狀。節(jié)理產(chǎn)狀的隨機生成和平均產(chǎn)狀的計算均通過Python程序實現(xiàn)。

圖4 節(jié)理和向量rN的計算過程Fig.4 Calculation process of joint and vector rN

3 最小節(jié)理采樣數(shù)量的確定

由于節(jié)理產(chǎn)狀的離散性，節(jié)理平均產(chǎn)狀的計算結果受樣本大小的影響。產(chǎn)狀離散性越大，即κ值越小，節(jié)理的樣本容量應越大才能保證平均產(chǎn)狀計算的準確性。本文以平均產(chǎn)狀為45°∠60°的節(jié)理為例，通過隨機生成不同數(shù)量、不同κ值的Fisher分布產(chǎn)狀來研究樣本容量對平均產(chǎn)狀計算準確性的影響，κ值分別取5，10，15，20，25，30，35，40，45，50，55，60，70，80，90和100。對于每一個κ值，共進行1 000次抽樣，樣本容量為ns。圖5繪制了一次隨機實現(xiàn)時不同κ值的節(jié)理極點分布圖，樣本容量ns為250。

若每次抽樣的平均產(chǎn)狀與45°∠60°的夾角小于5°，則認為滿足準確性要求。將滿足以上要求的抽樣次數(shù)記為np，那么準確率p可以定義為：

(10)

當κ≤40時，每次抽樣的樣本容量ns為5的倍數(shù)；當κ>40時，ns為自然數(shù)。

圖6反映了不同κ值時樣本大小對平均產(chǎn)狀準確率的影響。樣本容量對平均產(chǎn)狀的統(tǒng)計準確性影響顯著，樣本容量越大，平均產(chǎn)狀的統(tǒng)計準確率越高；產(chǎn)狀離散性(即κ值)對統(tǒng)計準確率也有明顯影響，離散性越大即κ值越小，達到一定準確率所需的樣本數(shù)越多。一般情況下，準確率p為0.95可滿足工程要求，因此將p=0.95對應的樣本容量定義為最小樣本容量。當樣本大小超過最小樣本容量時，可以保證產(chǎn)狀的統(tǒng)計準確率不低于0.95。圖7為具有不同κ值的節(jié)理產(chǎn)狀所應滿足的最小樣本容量ns，min，可以看出，ns，min隨著κ值的增大迅速減小，當κ值達到40以后，ns，min的變化很小且趨于水平。為方便在工程實踐中確定節(jié)理的最小樣本容量，采用反比例函數(shù)(公式(11))對ns，min和κ值的關系進行擬合，得到了圖7中的擬合曲線。公式(11)的確定系數(shù)R2為0.994，表明該式能準確反映ns，min與κ值的反比例關系，可用于確定節(jié)理的最小樣本容量。

圖6 不同κ值時樣本容量對統(tǒng)計準確率的影響Fig.6 Effect of sample size on statistical accuracy with different κ values

圖7 不同κ值對應的最小樣本容量Fig.7 Minimum sample size corresponding to different κ values

(11)

實際上，在對節(jié)理產(chǎn)狀進行采樣之前無法預知節(jié)理組的κ值。但隨著采樣過程的進行，會得到越來越多的產(chǎn)狀數(shù)據(jù)，可以利用這些數(shù)據(jù)對κ值進行估算，從而判斷最終的采樣數(shù)量ns是否滿足公式(11)的要求。如果ns不小于ns，min，則采樣數(shù)量滿足要求，平均產(chǎn)狀的準確率能夠得到保障。本文的研究基于產(chǎn)狀服從Fisher分布，實際巖體中大多數(shù)構造節(jié)理均成組發(fā)育且服從Fisher分布，因此公式(11)對構造節(jié)理具有較好的適用性。對于孤立的隨機節(jié)理，其發(fā)育沒有規(guī)律性，不適合采用公式(11)進行計算。

4 驗證算例

為驗證公式(11)的適用性，對三峽工程地下廠房的上游邊墻和烏東德水電站平洞PD5-2節(jié)理產(chǎn)狀進行了抽樣統(tǒng)計分析。三峽水電站是世界上最大的水利樞紐之一，位于中國湖北省宜昌市夷陵區(qū)。地下電站位于長江右岸白巖尖山體內(nèi)，主廠房洞室為圓拱直墻型，全長 311.3 m，軸線走向 223.5°，最大跨度 32.6 m，最大高度 87.3 m。廠房圍巖為堅硬的花崗巖，巖體完整～較完整，整體穩(wěn)定性較好。地下電站廠房洞室的開挖揭露了708條節(jié)理，分組后的節(jié)理極點分布情況見圖8，每組節(jié)理產(chǎn)狀的統(tǒng)計參數(shù)見表1。烏東德水電站平洞PD5-2位于金沙江烏東德壩址左岸850 m高程，為順江支洞。平洞洞形為城門洞形，斷面尺寸約2 m×2 m。本平洞共揭露節(jié)理354條，分組后的節(jié)理極點分布情況見圖9，每組節(jié)理產(chǎn)狀的統(tǒng)計參數(shù)見表2。

表1 三峽地下廠房上游邊墻節(jié)理產(chǎn)狀統(tǒng)計結果Table 1 Statistical results of the orientation of joints in the upstream side wall of the Three Gorges underground powerhouse

表2 烏東德水電站平洞PD5-2節(jié)理產(chǎn)狀統(tǒng)計結果Table 2 Statistical results of the orientation of joints in the adit PD5-2 of the Wudongde hydropower station

圖8 三峽地下廠房上游邊墻節(jié)理極點分布圖Fig.8 Distribution diagram of joint poles in the upstream side wall of the Three Gorges underground powerhouse

圖9 烏東德水電站平洞PD5-2節(jié)理極點分布圖Fig.9 Distribution diagram of joint poles in the adit PD5-2 of the Wudongde hydropower station

根據(jù)公式(11)，計算出三峽地下廠房上游邊墻每組節(jié)理的最小樣本容量應達到97、97和105。接下來，對每組節(jié)理產(chǎn)狀進行模擬抽樣，抽樣次數(shù)均為1 000次，樣本大小為最小樣本容量。最后利用公式(10)計算抽樣結果的準確率。經(jīng)程序計算，每組節(jié)理平均產(chǎn)狀的準確率分別為0.988、0.992和0.994，均滿足準確率要求。該研究區(qū)每組節(jié)理的采樣數(shù)量均大于相應的最小樣本容量，滿足采樣要求。因此，三峽工程地下廠房的上游邊墻的節(jié)理產(chǎn)狀具有統(tǒng)計學意義。

根據(jù)公式(11)，計算得到烏東德水電站平洞PD5-2每組節(jié)理的最小樣本容量應達到87、100和104。假設每組節(jié)理采樣數(shù)分別為80、90和95，未達到相應組別的最小樣本容量。利用公式(10)，計算出每組抽樣結果的準確率分別為0.955、0.948和0.946。結果表明，除第一組采樣能達到準確率以外，其余兩組未能滿足準確率要求。抽樣過程具有一定的隨機性和偶然性，在抽樣數(shù)量不足時，仍有較小可能達到準確率，但是為保證抽樣結果具有較高的可信度，節(jié)理抽樣應達到最小樣本容量的要求。本節(jié)通過正反兩個算例驗證了公式(11)的可靠性，表明節(jié)理最小樣本容量ns，min與κ值的經(jīng)驗關系適用于確定節(jié)理產(chǎn)狀的最小樣本容量。

5 結論

本文基于節(jié)理產(chǎn)狀服從Fisher分布，實現(xiàn)了節(jié)理產(chǎn)狀的隨機生成。為解決“超越上半球投影的產(chǎn)狀”(OEBEUHP)引起平均產(chǎn)狀計算不準確的問題，提出了一種計算平均產(chǎn)狀的修正方法。在此基礎上，隨機生成不同數(shù)量、不同離散性(κ值)的Fisher分布產(chǎn)狀來研究節(jié)理樣本容量對平均產(chǎn)狀準確性的影響，并建立了最小樣本容量與κ值的經(jīng)驗關系。研究結果表明：節(jié)理樣本容量對平均產(chǎn)狀的統(tǒng)計準確性影響顯著，樣本容量越大，平均產(chǎn)狀的統(tǒng)計準確率越高；產(chǎn)狀離散性對統(tǒng)計準確率也有明顯影響，離散性越大即κ值越小，達到一定準確率所需的樣本數(shù)越多。為保證平均產(chǎn)狀的統(tǒng)計準確率不低于0.95，節(jié)理的最小樣本容量隨著Fisher分布常數(shù)κ值的增大迅速減小；二者之間具有明顯的反比例關系，反比例系數(shù)為1 007.23，經(jīng)驗公式的確定系數(shù)為0.994。對三峽地下廠房上游邊墻和烏東德水電站平洞PD5-2的節(jié)理統(tǒng)計分析從正反兩方面驗證了本文方法的適用性，節(jié)理采樣具有一定的隨機性，在滿足最小樣本容量的條件下，產(chǎn)狀統(tǒng)計分析結果具有較高的可靠性。本文提出的方法適用于確定節(jié)理產(chǎn)狀的最小采樣數(shù)量，判斷節(jié)理的采樣數(shù)量是否具有統(tǒng)計意義，能有效指導節(jié)理的采樣工作。