基于SQP算法的概率積分法參數反演

代 陽

(1.安徽理工大學 測繪學院，安徽 淮南 232001；2.礦山采動災害空天地協同監測與預警安徽省教育廳重點實驗室，安徽 淮南 232001)

煤炭在我國能源供給結構中占有主體地位。在地下開采時，礦區巖體受力平衡遭到破壞，引發上覆巖層移動變形，導致地表沉降移動。且隨著開采范圍的擴大，周邊建筑與環境勢必受到影響，嚴重時可能會引發建筑物塌陷、地下水資源被破壞、地表塌陷等地質災害和次生災害，嚴重影響人們的生活安全[1-2]。因此，研究礦區開采引發的地表移動變形規律，在災害治理與防范中有著重要的意義。

概率積分法是預計礦區開采沉陷變化規律的積分模型，其模型的主要思想是基于LITWINISZYN提出的隨機介質理論，利用采空區實測數據擬合得到的預計參數，預測礦區地表移動變化。構建概率積分法模型的關鍵，在于保證輸入預計參數的準確性。由于該模型非線性程度高、參數之間相關性強、求解易受迭代初值影響，采用傳統的求解參數方法往往會帶來局部解甚至無解[3]。為解決上述問題，學者們引用了優化及智能優化算法反演模型參數。葛家新等以皖北礦區為例，采用模矢法進行反演參數，擬合預測的下沉效果較好，但易受到參數設計初值的影響[4]。查劍鋒等采用遺傳算法進行模型參數估計，并與最小二乘法和模矢法相比，證明遺傳算法的準確性與穩定性更優[5]。蘇軍明等采用模擬退火的方法，在搜索最優參數的過程中加入隨機因素，更好地得到全局最優解[6]。徐孟強等采用粒子群算法反演參數能夠得到準確的預計參數，且在一定范圍內受初值誤差影響較小[7]。智能優化算法通過模仿自然現象及生物進化規律，迭代逐步逼近參數最優值，然而在搜索的過程中，容易早熟且每次反演得到的參數具有一定的隨機性，為構建概率積分模型帶來影響[8]。因此，有必要引入新的算法反演概率積分法中的參數。

概率積分法參數反演可以看作是一種非線性多參數求解最優值的過程，提出采用二次序列規劃(Sequential Quadratic Programming，SQP)算法反演參數，以達到擬合參數精度最優的目的。SQP算法屬于最優化算法，通過二次規劃確定每次迭代的下降方向，并通過價值函數獲取步長，多次迭代逼近最優值[9]。在處理非線性問題中，SQP算法不依賴于初值，并具有良好的收斂與自我校正能力，結果的穩定性與可靠性強。目前，該算法已經被成功運用于各種非線性最優化問題中，例如在飛機數字化裝配[10]、陣列信號處理[11]、高密度電阻率等領域[12]。

根據概率積分模型構建目標函數，引用SQP算法進行參數反演，驗證SQP算法運用在開采沉陷預測參數的可能性。通過模擬監測數據受到不同程度的誤差及粗差的影響，進行大量仿真實驗，將SQP算法是反演結果與常用的優化及智能優化算法作比較，驗證SQP算法抗隨機誤差及粗差的性能。

1 概率積分法預計模型

概率積分法是一種預計由于礦區開采引起的地表沉陷與水平移動的模型[13](P116-120)。其中，預測地表任意點下沉值的函數模型為：

(1)

式中，(x,y)為地表任意點平面坐標；W(x,y)為對應坐標下的地表下沉值；W0(x)為傾向方向充分采動時，橫坐標值為x時在走向面的下沉值；W0(y)為走向方向充分采動時，縱坐標為y時在傾向面的下沉值；W0為地表最大下沉值，限于篇幅，W0、W0計算公式見參考文獻[13]，不再累述。

預測地表任意點的水平移動值的函數模型為：

(2)

式中，φ為由x軸正向逆時針旋轉到指定方向的方位角；U(x,y,φ)為對應坐標下的地表水平移動值；U0(x)為走向主斷面上的水平移動值；U0(y)為傾向主斷面上的水平移動值。限于篇幅，U0詳細計算見參考文獻[13]。

概率積分法模型反演參數屬于求解多參數最優化極值問題，目的是尋求一組變量，使預測的地表下沉值與水平移動值與實測數據的差值平方和最小。其目標函數的具體表現形式為：

(3)

式中，i為測站數據編號；N為測站的總數目；WS為測站實測下沉觀測值；US為測站實測水平移動值。

2 SQP算法的理論基礎

SQP算法作為一種應用于非線性函數參數反演的最優化算法，應用于概率積分法反演求參具體表現形式為[14]：

(4)

式中，x=[q,b,θ,tanβ,S1,S2,S3,S4]T為待求優化參數，其中q為下沉系數，b為水平移動系數，θ為開采影響傳播角，tanβ為主要影響角正切，S1,S2,S3,S4為對應的下山方向、上山方向、左方向和右方向的拐點偏移距。f(x)對應上文的公式(3)，Gj(x)對應的限制條件，根據礦區情況和需求進行改變。該算法采用拉格朗日求極值法將參數反演轉化為求解二次近似函數最優值的問題：

(5)

通過線性化轉化后得到QP子問題，其目標函數變換為：

(6)

xk+1=xk+akdk

(7)

式(6)中，d為全變量搜索方向，為梯度，k為迭代次數，Bk為拉格朗日函數的海森矩陣，通過半正定牛頓近似法(Dense quasi-newton approximation，BFGS)進行計算得到[15]。式(7)中，ak為第k次迭代時的步長。

SQP算法反演概率積分模型預計參數的主要步驟為：

1)給定SQP算法反演參數x的初值、為算法反演提供初始解，同時令迭代次數k=0；

2)將反演參數值帶入，計算目標函數f(x)及約束函數Gj(x)，差分計算梯度向量值，求解二次規劃子問題，計算搜索方向dk；

3)如果dk=0或者xk已經符合精度檢驗要求，則迭代終止，得到的xk即為該非線性函數的最優解；否則，從點xk沿著dk方向搜索下一步迭代值，并確定步長ak，令xk+1=xk+akdk；

4)通過BFGS算法確定海森矩陣Bk，并轉化為正定矩陣。令k=k+1，轉至步驟(2)。

3 仿真實驗與結果統計分析

3.1 設計工作面概況

依據工作面觀測點的設計原則，參考于文獻[8]中的工作面，制定模擬礦區的觀測線。具體設計為，煤層開采厚度m為4 m、煤層傾角α為3°、工作面的走向長L和傾向長D分別為670 m和300 m、開采深度為260 m，管理頂板為全部垮落法進行管理[8]。概率積分參數的真值分別設計為：下沉系數q為0.64，水平移動系數b為0.3，開采影響角正切值tanβ為2.37，開采影響傳播角θ為88.5°，上山方向、下山方向、左方向與右方向的拐點偏移距S1、S2、S3、S4分別設計為52 m。并在礦區上方沿著走向線布設25個監測點(S1～S25)，沿著傾向線布設44個監測點(E1～E44)，各個監測點設計的點位布設，如圖1所示。

圖1 工作面及對應點位模擬圖

3.2 SQP算法參數反演準確性分析

以上述S1～S25、E1～E44共69個監測點的下沉與水平移動值作為反演數據的原始數據，結合礦區的基本信息，本文設計了五組不同的參數初始值，將反演結果與真值對比，對算法的穩定性進行評估。同時，通過比較不同初始參數設置下反演結果的差異，驗證SQP算法預計參數的結果是否受到初始設計值的影響。

表1 SQP算法反演概率積分法準確性分析

通過分析表1中的反演結果可以得到：1)對比五組不同參數設置下的SQP算法反演結果，隨著設置的初始值不斷增大，q、b、tanβ、θ與四個方向的拐點偏移距的預計參數并沒有發生明顯的變化。2)將反演結果與真值對比，參數q、b、tanβ和θ與四個方向的拐點偏移距的預計結果均與真值基本保持一致，相對誤差不到1%。根據上述實驗結果，驗證了SQP算法受到初始值設置的影響很小，且可靠性與穩定性高。

3.3 SQP算法對隨機誤差、粗差的抗干擾能力

在礦區觀測中，觀測數據必然會受到隨機噪聲和粗差的影響，因此測試算法抗隨機誤差與粗差的干擾能力，是評價算法性能的重要依據。本文模擬多組在不同程度的隨機誤差及粗差影響下的監測數據，將SQP算法、模矢法與遺傳算法進行對比實驗，具體步驟如下：

1)設計SQP算法、模矢法與遺傳算法的實驗參數。模矢法的參數設計為：參數終止步長為0.0001；遺傳算法的參數設計為：種群個體數為100，遺傳代數為200代，交叉概率為0.9，變異概率為0.001；SQP算法的參數設計為：最大迭代次數200，終止誤差為0.0001。每個算法的初始參數的取值范圍均為：q∈[0.61,0.67]；b∈[0.25,0.35]；tanβ∈[1.97,2.77]；θ∈[80°,90°]；S1、S2、S3、S4的取值范圍均為S∈[32,72]。

2)將設計的采礦條件與沉陷預計參數真值帶入概率積分模型中，得到設計工作面觀測點的下沉與水平移動值的真值，在真值的基礎上分別加入中誤差為5 mm、10 mm與20 mm的隨機誤差，并在不同的隨機誤差情況下，各仿真20組實驗數據，并統計由SQP算法、模矢法與遺傳算法反演的各組參數的均值。

3)在實驗步驟(1)中的誤差為10 mm的仿真數據的基礎上，隨機選取0～5個觀測站的下沉與水平觀測值，在此基礎上添加5σ-20σ的粗差，并統計由SQP算法、模矢法與遺傳算法反演參數的均值。

SQP算法、模矢法與遺傳算法反演的參數均值統計值結果，如表2所示。

表2 SQP算法、模矢法與遺傳算法反演參數計算結果

表2給出了不同組別的仿真實驗，分析反演參數的結果可知：1)在不同誤差情況下，SQP算法對于q、b、tanβ和θ的反演結果并未有明顯變化，且與真值相比反演結果基本一致，最大誤差不超過1.2%；而對于S1、S2、S3、S4反演結果具有一定的隨機性，但與真值相比最大誤差仍不超過8%，仍具有一定的可靠性。2)在隨機誤差10 mm的基礎上加入了一定數目的粗差時，SQP算法反演結果并未發生改變，與之前相比變化不超過2.8%。3)隨著誤差的不斷增大，SQP算法、遺傳算法與模矢法反演參數的整體精度隨著誤差的增長有所降低，但結果與真值相差不大，仍具有一定的可靠性。4)對于q、b、tanβ和θ參數的反演結果，SQP算法受誤差因素影響較小，結果與模矢法和遺傳算法相比，均更加靠近真值。且在5 mm誤差的情況下，SQP算法優勢更加明顯，反演結果幾乎可以認為是真值。5)對于S1、S2、S3、S4參數的反演結果，SQP算法的精度隨著誤差的增大而有所降低。在5 mm與10 mm隨機誤差的影響下，SQP算法對于S1、S2、S3、S4的反演結果均優于遺傳算法與模矢法；在20 mm隨機誤差的影響下，SQP算法對于S1、S2的反演結果均優于模矢法與遺傳算法，但對于S3、S4反演的結果精度稍差于遺傳算法與模矢法；在10 mm粗差的影響下，SQP算法對于S1、S2、S3、S4的反演結果均優于模矢法，但在S1、S2、S3方向上預計結果與遺傳算法相比稍差，但考慮到算法的穩健性，其仍具有一定的可靠性。6)在監測數據受到不同程度的隨機誤差與粗差的影響下，SQP算法反演參數的整體精度最高，遺傳算法次之，而模矢法較差。因此，表2證明了SQP算法對觀測數據中含有隨機誤差和粗差具有良好的抗干擾能力。

4 結論

通過分析常用反演參數算法存在的缺陷，提出使用SQP算法反演概率積分法模型中的參數，進行大量的仿真實驗，得到如下結論：

1)SQP算法不受初始值設置影響，且反演參數的結果與真值基本一致，準確度高。

2)SQP算法對隨機誤差和粗差具有一定的抗干擾能力，其反演參數結果在隨機誤差與粗差的影響下，并未產生明顯的偏差，依然可靠。

3)在不同的誤差環境的影響下，SQP算法反演參數的精度基本高于模矢法與遺傳算法，但在少數情況下，個別方向上的拐點偏移距的預測精度可能稍低于遺傳算法，考慮到SQP算法的穩健性，其結果依然可靠。

概率積分法預測結果受參數誤差與模型誤差影響，本文針對如何降低參數誤差影響進行了改進。但在現實礦區中，實測的下沉值和水平移動值還會受到松散層厚度、斷層等其他因素影響，從而影響到模型自身精度。如何精化概率積分法的數學模型，使得函數預計曲線更符合實際的下沉曲線和水平移動曲線，值得進一步研究和探索。