剖分算法建模分析

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

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

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

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

5、平面離散點(diǎn)的有效作用范圍，首先定義了二維平面上的圖，年由將圖演化出了更易于分析應(yīng)用的三角網(wǎng)圖和三角網(wǎng)為目前普遍接受和廣泛采用的分析研究區(qū)域離散數(shù)據(jù)的有力工具，已應(yīng)用到石油勘探、地質(zhì)、礦業(yè)、城市規(guī)劃和環(huán)境監(jiān)測(cè)等領(lǐng)域?qū)崿F(xiàn)了網(wǎng)生成算法，對(duì)邊界恢復(fù)算法進(jìn)行改進(jìn)，陳曉勇、李清泉等研究了四面體格網(wǎng)結(jié)構(gòu)（）模型的生成算法剖分算法的剖分結(jié)果為四面體體元，四個(gè)頂點(diǎn)為空間離散點(diǎn)，包含坐標(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)具有如下特點(diǎn)：生成的四面體之

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

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

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

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

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

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

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

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

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

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