剖分算法建模分析

1、剖分算法建模分析     引言地球重力場信息在地球構(gòu)造研究、國土資源調(diào)查、礦產(chǎn)勘探以及國民經(jīng)濟(jì)等領(lǐng)域中發(fā)揮著重大作用重力的正問題與反問題是解釋重力資料的核心研究內(nèi)容，其中前者主要研究不同形狀、產(chǎn)狀和場源密度等場源體或地質(zhì)體所引起的重力異常特征及分布等，其計算過程為正演計算；后者主要根據(jù)重力異常的分布來計算地質(zhì)體或場源體的密度分布，其計算過程為反演計算目標(biāo)地質(zhì)體的重力建模建立了地質(zhì)體的剩余密度值與地面重力異常值的關(guān)系，既可利用重力正演公式計算地質(zhì)體所引起的重力異常，也可以應(yīng)用于求解重力反問題因此，重力建模是重力正反演問題的關(guān)鍵環(huán)節(jié)，建模方法的優(yōu)劣則直接影響

2、到重力正反演的計算精度隨著對重磁資料解釋的深入研究、重力測量精度的提高及計算機(jī)的快速發(fā)展，重力正反演研究由過去的二維逐步發(fā)展為三維對于三維連續(xù)密度變化的目標(biāo)體，目前最常用的建模方法是塊體方法，即將目標(biāo)地質(zhì)體剖分為規(guī)則的構(gòu)造單元如長方體，計算所有單元產(chǎn)生的異常并疊加便得到目標(biāo)地質(zhì)體所產(chǎn)生的異常值樓海等給出了基于矩形網(wǎng)格模型的三維連續(xù)密度分布目標(biāo)地質(zhì)體的重力解析公式，并在實際應(yīng)用中得到較好的結(jié)果，但當(dāng)?shù)刭|(zhì)形狀比較復(fù)雜時，其計算結(jié)果并不理想，如果加密構(gòu)造單元，則又會增加計算量及存儲空間盛國平采用水平截面法對非規(guī)則幾何形狀的目標(biāo)體進(jìn)行分割，并假設(shè)每一個質(zhì)面層的密度值為常數(shù)，給出了較嚴(yán)密的數(shù)學(xué)解析表達(dá)式

3、張嶺等針對截面為任意形狀的二度體問題，利用二維剖分方法，將截面分割成若干三角形，二度體被分解成若干三棱柱的組合來計算二度體的重力異常值，這種方法表達(dá)了二維的連續(xù)密度變化，但僅限制于二度體的重力計算問題對于任意形狀和變密度的三度體，現(xiàn)有重力建模方法均存在不足，為此，本文提出基于剖分算法的重力建模方法采用剖分算法將三維目標(biāo)地質(zhì)體分解為若干變密度四面體體元，推導(dǎo)基于四面體體元的重力正演公式；比較分析常規(guī)塊體方法和本文方法應(yīng)用于重力正演的計算效果，并采用共軛梯度法加密度約束條件對非規(guī)則形狀變密度的地質(zhì)體進(jìn)行反演計算，驗證本文重力建模方法的正確性和有效性剖分的重力建模地質(zhì)體的建模地質(zhì)體的建模步驟為：首先

4、根據(jù)地質(zhì)體的三維坐標(biāo)信息及物性信息等相關(guān)資料將目標(biāo)地質(zhì)體離散化為空間離散點；其次，采用合理的建模算法組織空間離散點，比如傳統(tǒng)重力建模時采用的塊體方法，地學(xué)領(lǐng)域、及剖分方法等；最后將所有體元結(jié)構(gòu)（例如四面體、三棱柱體、六面體和長方體等）組合起來，根據(jù)體元頂點屬性特征采用顏色或紋理對模型進(jìn)行可視化對于規(guī)則形狀的地質(zhì)體建模，采用傳統(tǒng)塊體方法即可對地質(zhì)體進(jìn)行快速建模，并真實逼近模型結(jié)構(gòu)及地質(zhì)體的物性特征，但是對于非規(guī)則形狀變密度地質(zhì)體的建模，如褶皺、斷層等，傳統(tǒng)建模方法則會改變地質(zhì)體的結(jié)構(gòu)、密度分布和其他物性特征，所以本文引進(jìn)剖分方法來逼近真實的地質(zhì)體剖分算法由二維三角剖分算法演化而來，年為了限定二維

5、平面離散點的有效作用范圍，首先定義了二維平面上的圖，年由將圖演化出了更易于分析應(yīng)用的三角網(wǎng)圖和三角網(wǎng)為目前普遍接受和廣泛采用的分析研究區(qū)域離散數(shù)據(jù)的有力工具，已應(yīng)用到石油勘探、地質(zhì)、礦業(yè)、城市規(guī)劃和環(huán)境監(jiān)測等領(lǐng)域?qū)崿F(xiàn)了網(wǎng)生成算法，對邊界恢復(fù)算法進(jìn)行改進(jìn)，陳曉勇、李清泉等研究了四面體格網(wǎng)結(jié)構(gòu)（）模型的生成算法剖分算法的剖分結(jié)果為四面體體元，四個頂點為空間離散點，包含坐標(biāo)信息和物性特征四面體格網(wǎng)數(shù)據(jù)模型是二維三角形網(wǎng)（，）數(shù)據(jù)結(jié)構(gòu)在三維空間上的擴(kuò)展，與之相應(yīng)的數(shù)據(jù)結(jié)構(gòu)為四面體格網(wǎng)結(jié)構(gòu)（）模型以三角剖分為基礎(chǔ)，將目標(biāo)空間用緊密排列但不重疊的非規(guī)則四面體的組合來表示四面體應(yīng)具有如下特點：生成的四面體之

6、間無重疊部分；所形成的四面體組合可以覆蓋整個目標(biāo)體的三維空間；四面體的外接球不包含空間離散點集中除該四面體四個頂點外的任一點這些性質(zhì)保證了組成四面體的三角形近似等邊或等角，及四面體體元的組合更逼近真實的目標(biāo)地質(zhì)體為了比較常規(guī)塊體方法及剖分算法對重力建模的影響，本文設(shè)計了兩個目標(biāo)地質(zhì)體目標(biāo)體為連續(xù)密度變化的矩形體，上頂面埋深為，范圍為××（長×寬×高），八個角點的剩余密度值依次為、和，如圖所示；目標(biāo)體為傾斜臺階組合模型體，八個角點的剩余密度值與目標(biāo)體相同，具體參數(shù)參見圖對于目標(biāo)體的建模，首先將目標(biāo)體離散化為若干空間點位，再分別采用常規(guī)塊體剖分方法與四面體剖

7、分方法對這些離散點位建模在、和方向取分塊數(shù)為××，空間點數(shù)為××，對目標(biāo)體的剖分結(jié)果見圖比較圖和圖可知：兩種方法的建模結(jié)果均能準(zhǔn)確表達(dá)地質(zhì)體，理論上得到的重力正演計算結(jié)果應(yīng)該相差不大采用傳統(tǒng)塊體方法對目標(biāo)體的建模結(jié)果見圖當(dāng)采用剖分算法對目標(biāo)體建模時，如果直接利用空間離散點進(jìn)行四面體構(gòu)網(wǎng)，則會出現(xiàn)“跨越”現(xiàn)象，即左右傾斜臺階會出現(xiàn)錯位結(jié)構(gòu)，所以本文首先采用算法建立左右傾斜臺階接觸面的數(shù)據(jù)結(jié)構(gòu)，作為組合模型的特征面，然后進(jìn)行整體四面體構(gòu)網(wǎng)，其建模結(jié)果如圖所示從圖可以看出，具有約束的算法可以避免“跨越”現(xiàn)象，從而得到合理的三維模型比較圖和圖可知：塊體方法的建模

8、結(jié)果與真實地質(zhì)體有一定的差異，而剖分算法的建模結(jié)果能更好地逼近真實地質(zhì)體，因此理論上兩種方法的重力正演計算結(jié)果應(yīng)該存在較大差異；塊體方法的建模結(jié)果與分塊數(shù)密切相關(guān)，而剖分算法受分塊數(shù)的影響較小，隨著圖分塊數(shù)的增加，塊體方法的建模結(jié)果將逐漸逼近目標(biāo)體的形狀及密度分布，但計算量也將大幅度增加重力計算模型四面體體元為剖分算法的基本結(jié)構(gòu)單元，非規(guī)則形狀變密度地質(zhì)體對地面點產(chǎn)生的重力異常為所有四面體體元產(chǎn)生的重力異常之和任意一個四面體體元（見圖）對地面點產(chǎn)生的重力異常為式中，表示在四面體內(nèi)任意點（，）的剩余密度，由四個頂點的剩余密度值線性插值得到，滿足關(guān)系式：，（）式中，、和表示四面體體元四個頂點的剩余

9、密度值對于任意形狀的四面體，利用式（）計算地面點的重力異常值，不易得到、和三個方向的積分域，因此也不利于重力異常的計算本文借鑒文獻(xiàn)提出的水平截面法，首先將四面體切分為若干個水平面層，計算每一層所產(chǎn)生的重力異常，然后疊加得到一個四面體所引起的重力異常值與文獻(xiàn)的不同之處在于：四面體四個頂點的密度是不同的，其每一個面層也是變密度的圖為變密度的水平截面被劃分為若干個三角形的示意圖如圖所示，為地面計算點在目標(biāo)地質(zhì)體的一個面層上的投影，則面層被分為個三角形，以為例，計算一個三角形對點產(chǎn)生的（）當(dāng)線段平行于軸（）時，與的組成項變?yōu)?，槡槡；?dāng)線段平行于軸（）時，與的組成項變?yōu)闃殬殻ǎ?、和三點在一條直線上，即

10、，則；（）根據(jù)格林公式求解與時，定義逆時針為正，要求、與、一致，即或時表示順時針，系數(shù)，采用公式（）、（）、（）和（）計算出系數(shù)陣后，即可計算一個三角形質(zhì)面對地面點產(chǎn)生的重力異常，按逆時針方向計算組成一個截面的個三角形對點產(chǎn)生的異常，將其疊加即為水平截面對地面點產(chǎn)生的重力異常假設(shè)四面體被分為了層水平截面，、和分別為相應(yīng)埋深、和的水平截面所產(chǎn)生的重力異常值，則一個四面體對地面點產(chǎn)生的重力異常為綜上所述，最終可建立目標(biāo)體個剩余密度值與個地面重力異常的關(guān)系式，即算例分析目標(biāo)體的正演計算分別采用常規(guī)塊體剖分方法和剖分方法對目標(biāo)體建模，并計算了位于目標(biāo)體正上方范圍為××的空間格網(wǎng)點的

11、重力異常，常規(guī)塊體方法的正演計算公式可以參見文獻(xiàn)，兩種方法的計算結(jié)果如圖所示從圖可以看出：兩種方法的計算結(jié)果相差較小，最大差值為×對目標(biāo)體設(shè)計分塊數(shù)為時，分別采用常規(guī)塊體剖分方法和剖分方法對其建模，正演計算了位于目標(biāo)體正上方處個格網(wǎng)點的重力異常，其差值的比較結(jié)果如表所示從表可以看出：隨著分塊數(shù)的增加，兩種方法的計算結(jié)果并無明顯差異，最大差值在×以內(nèi)，均方誤差在×以內(nèi)因此，對于規(guī)則形狀且連續(xù)密度變化的目標(biāo)體，采用常規(guī)塊體剖分方法和剖分方法均能獲得較好的重力正演結(jié)果目標(biāo)體的正演計算分別采用常規(guī)塊體剖分方法和剖分方法對目標(biāo)體建模，計算了與節(jié)相同區(qū)域的重力異常，計算結(jié)果如

12、圖所示從圖可以看出：對于非規(guī)則形狀變密度的目標(biāo)體，兩種方法的計算結(jié)果差異較大，其差值在（）×范圍之內(nèi)，主要差異體現(xiàn)在雙傾斜臺階組合模型體的公共斜平面上，其原因是常規(guī)塊體方法改變了目標(biāo)體的幾何結(jié)構(gòu)和物性特征，對地質(zhì)體的建模不夠準(zhǔn)確類似于目標(biāo)體，對目標(biāo)體設(shè)計分塊數(shù)為，分別采用常規(guī)塊體方法和剖分方法對目標(biāo)體建模并計算其產(chǎn)生的重力異常，兩種方法的比較結(jié)果見表結(jié)果表明：隨著分塊數(shù)的增加，常規(guī)塊體剖分方法對目標(biāo)體的建模逐漸逼近真實地質(zhì)體，使得兩者計算結(jié)果的差值逐漸減小從表和表可以看出：平均值與均方誤差量級相當(dāng)，說明常規(guī)塊體方法和剖分方法之間存在明顯的系統(tǒng)誤差，并且系統(tǒng)誤差隨著分塊數(shù)的增加而減小這

13、主要是目標(biāo)體的建模誤差和式（）中劃分層數(shù)的取值所導(dǎo)致根據(jù)文獻(xiàn)和，的取值應(yīng)根據(jù)研究區(qū)域及其地質(zhì)條件來確定，本文的取值為目標(biāo)體的反演計算剖分算法不僅可以應(yīng)用于地質(zhì)體的重力正演計算，也可以應(yīng)用于地質(zhì)體的三維物性反演針對地球物理的反演問題，和最早對連續(xù)線性問題做了深入研究；和引入深度加權(quán)函數(shù)，克服了三維反演結(jié)果的“趨附”現(xiàn)象；等將光滑度引入到正則化反演；在重磁反演求解大型線性方程組時，為了不降低計算效率及消耗過多計算機(jī)資源，采用共軛梯度法進(jìn)行三維反演本文基于剖分算法，采用共軛梯度法加密度約束條件對目標(biāo)體進(jìn)行三維物性反演地質(zhì)體剩余密度值與地面重力異常值之間的線性關(guān)系式（）往往是欠定的，不易得到真實物性反

14、演結(jié)果，為此本文增加了密度約束條件定義目標(biāo)函數(shù)的表達(dá)式為這里為觀測值，為實例模型結(jié)點剩余密度值，為該線性方程組的系數(shù)陣，、與分別對應(yīng)于式（）的、與，和分別表示密度約束的下界和上界，在本算例的數(shù)值分別為采用最小二乘法求解式（）時，其矩陣求逆時間較長，而且有可能出現(xiàn)病態(tài)方程，因此本文利用共軛梯度法求解假定該區(qū)域的地質(zhì)條件已知，背景密度為，以節(jié)基于本文方法正演計算的重力異常作為地面觀測值，總共個地面重力異常值，其等值線圖見圖；采用個地面重力異常值反演目標(biāo)體的個剩余密度值，其反演結(jié)果見圖；利用密度反演結(jié)果計算得到的重力異常見圖；圖為圖和圖的重力異常之差，用于評價三維密度反演效果圖的結(jié)果表明，采用共軛梯

15、度法加密度約束條件得到的反演結(jié)果非常接近于真實目標(biāo)地質(zhì)體的密度分布；從圖可以看出：利用反演密度值計算的重力異常（見圖）與模擬重力異常（見圖相差很小，其量級為以上結(jié)果初步驗證了本文的重力建模方法和三維密度反演方法的正確性和有效性結(jié)論與展望針對任意形狀和變密度的三度體，本文提出了基于剖分算法的重力建模方法采用剖分算法將三維目標(biāo)地質(zhì)體分解為若干變密度四面體體元，推導(dǎo)了基于四面體體元的重力正演公式；以規(guī)則形狀且連續(xù)密度變化的地質(zhì)體（如長方體）和非規(guī)則形狀變密度的地質(zhì)體（如傾斜臺階組合模型體）為例，比較分析了常規(guī)塊體方法和剖分算法應(yīng)用于重力正演的計算效果，并采用共軛梯度法加密度約束條件對非規(guī)則形狀變密度的地質(zhì)體進(jìn)行了反演計算模擬結(jié)果表明對于規(guī)則形狀且連續(xù)密度變化的地質(zhì)體，基于常規(guī)塊體方法和剖分算法的重力建模均能獲得滿意的計算結(jié)果；對于非規(guī)則形狀變密度的地質(zhì)體，基于剖分算法的重力建模能獲得更加合理的計算結(jié)果；隨著對地質(zhì)體分塊數(shù)的增加，盡管常規(guī)塊體方法的建模結(jié)果將逐漸逼近真實地質(zhì)體