本發(fā)明涉及磁場(chǎng)梯度張量數(shù)值模擬方法
技術(shù)領(lǐng)域：
，具體的涉及一種任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法。
背景技術(shù)：
：磁法勘探是地球物理勘探中一種重要的技術(shù)手段，在區(qū)域地質(zhì)調(diào)查、普查找礦工作中發(fā)揮著重要作用。該方法通過測(cè)量由巖石、礦石等磁性體的磁化率差異所引起的磁異常，對(duì)測(cè)量數(shù)據(jù)進(jìn)行反演、解釋，進(jìn)而研究地質(zhì)構(gòu)造和礦產(chǎn)資源的分布規(guī)律。隨著傳感器技術(shù)的不斷發(fā)展，磁場(chǎng)梯度張量測(cè)量技術(shù)日漸成熟。相比磁異常和磁場(chǎng)三分量，磁場(chǎng)梯度張量具有更高的分辨率。利用磁場(chǎng)梯度張量數(shù)據(jù)進(jìn)行反演、解釋，是磁法勘探的重要研究方向。磁性體磁場(chǎng)梯度張量的數(shù)值模擬是進(jìn)行反演的基礎(chǔ)。任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量快速、高精度數(shù)值模擬一直是個(gè)難點(diǎn)問題。針對(duì)磁場(chǎng)數(shù)值模擬，眾多國(guó)內(nèi)外學(xué)者進(jìn)行了研究。數(shù)值模擬首先對(duì)研究區(qū)域進(jìn)行剖分，然后根據(jù)剖分方式，采用某種方法計(jì)算磁場(chǎng)。文獻(xiàn)(姚長(zhǎng)利,郝天珧,管志寧,張聿文.重磁遺傳算法三維反演中高速計(jì)算及有效存儲(chǔ)方法技術(shù).地球物理學(xué)報(bào),2003.46(2):252-258.)采用結(jié)構(gòu)化剖分方式，根據(jù)離散化后數(shù)學(xué)問題的特點(diǎn)，提出了“格架分離”技術(shù)和“格架等效計(jì)算方案”，較好解決了計(jì)算效率和計(jì)算精度問題，但對(duì)于大規(guī)模剖分情形，該文獻(xiàn)所給出的數(shù)值模擬方法的計(jì)算效率仍然比較低；文獻(xiàn)(Tontini,F.C.,L.Cocchi,C.Carmisciano.Rapid3-DforwardmodelofpotentialfieldswithapplicationtothePalinuroSeamountmagneticanomaly(southernTyrrhenianSea,Italy).JournalofGeophysicalResearch,2009.114.)采用結(jié)構(gòu)化剖分方法，采用三維傅里葉變換，給出了任意密度或者磁化率分布情形下重磁數(shù)值模擬的波數(shù)域表達(dá)式，借助三維快速傅里葉變換算法，實(shí)現(xiàn)了快速數(shù)值模擬，該方法效率極高，但為克服截?cái)嘈?yīng)，使用該方法前需要對(duì)剖分區(qū)域進(jìn)行擴(kuò)邊，影響了數(shù)值模擬精度；文獻(xiàn)(Wu,L.,Tian,G.High-precisionFourierforwardmodelingofpotentialfields.Geophysics,2014,79(5):G59-G68.)提出了一種重磁數(shù)值模擬的Gauss-FFT方法，該方法有效克服了傳統(tǒng)傅里葉變換方法的截?cái)嘈?yīng)問題，提高了數(shù)值模擬精度，但計(jì)算效率降低，且計(jì)算場(chǎng)源內(nèi)部磁場(chǎng)時(shí)誤差比較大。剖分方式和計(jì)算方法共同決定了數(shù)值模擬的效率和精度。數(shù)值模擬的效率和精度是一對(duì)矛盾體，目前已有的數(shù)值模擬方法存在的最大問題是不能同時(shí)保證計(jì)算效率和精度，無(wú)法滿足大規(guī)模數(shù)據(jù)三維磁化率精細(xì)反演、人機(jī)交互建模和解釋的需求；同時(shí)，目前的方法不能準(zhǔn)確計(jì)算場(chǎng)源內(nèi)部的磁場(chǎng)梯度張量，不能滿足井地聯(lián)合反演的需求。因此，尋找一種計(jì)算效率高、同時(shí)能保證計(jì)算精度、且能精確計(jì)算全空間磁場(chǎng)梯度張量的數(shù)值模擬方法具有重要的現(xiàn)實(shí)意義。技術(shù)實(shí)現(xiàn)要素：本發(fā)明的目的在于提供一種任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，該發(fā)明解決了目前磁場(chǎng)梯度張量數(shù)值模擬方法計(jì)算精度低、計(jì)算時(shí)間長(zhǎng)，不能準(zhǔn)確計(jì)算場(chǎng)源內(nèi)部磁場(chǎng)，無(wú)法滿足大規(guī)模觀測(cè)數(shù)據(jù)精細(xì)反演、人機(jī)交互建模解釋和井地聯(lián)合反演需求的技術(shù)問題。本發(fā)明提供了一種任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，包括以下步驟：步驟S100：建立包含目標(biāo)區(qū)域的棱柱體模型，將棱柱體模型剖分為多個(gè)棱柱體，對(duì)各棱柱體設(shè)置磁化率，并計(jì)算各棱柱體的磁化強(qiáng)度得到組合棱柱體模型；步驟S200：按公式(5)采用二維離散卷積法計(jì)算各層棱柱體的磁場(chǎng)梯度張量其中，表示第r層(r＝1,2,…,L)棱柱體在高度面z0產(chǎn)生的磁場(chǎng)梯度張量，my(ξp,ηq,ζr)為坐標(biāo)為(ξp,ηq,ζr)的棱柱體在x方向的磁化強(qiáng)度分量、my(ξp,ηq,ζr)為坐標(biāo)為(ξp,ηq,ζr)的棱柱體在y方向的磁化強(qiáng)度分量，mz(ξp,ηq,ζr)為坐標(biāo)為(ξp,ηq,ζr)的棱柱體在z方向的磁化強(qiáng)度分量，hxz(xm-ξp,yn-ηq,z0-ζr)為磁化強(qiáng)度XZ分量的加權(quán)系數(shù)、hyz(xm-ξp,yn-ηq,z0-ζr)為磁化強(qiáng)度YZ分量的加權(quán)系數(shù)、hzz(xm-ξp,yn-ηq,z0-ζr)為磁化強(qiáng)度ZZ分量的加權(quán)系數(shù)，M為目標(biāo)區(qū)域x方向棱柱體的剖分個(gè)數(shù)，N為目標(biāo)區(qū)域y方向棱柱體的剖分個(gè)數(shù)；步驟S300：按公式(16)累加各層棱柱體的磁場(chǎng)梯度張量得到組合棱柱體模型的模擬磁場(chǎng)梯度張量其中，L表示目標(biāo)區(qū)域z方向棱柱體剖分個(gè)數(shù)。進(jìn)一步地，加權(quán)系數(shù)按公式(6)～(8)計(jì)算：其中，(xm,yn,z0)為觀測(cè)點(diǎn)坐標(biāo)，z0為常值；x1＝ξp-0.5Δx-xm，x2＝ξp+0.5Δx-xm，y1＝ηq-0.5Δy-yn，y2＝ηq+0.5Δy-yn，z1＝ζr-0.5Δz-z0，z2＝ζr+0.5Δz-z0，μijk＝(-1)i(-1)j(-1)k，i＝1,2，j＝1,2，k＝1,2。進(jìn)一步地，棱柱體模型為規(guī)則棱柱型。進(jìn)一步地，根據(jù)目標(biāo)區(qū)域的磁化率分布將每個(gè)棱柱體的磁化率設(shè)置為常值，并將位于目標(biāo)區(qū)域空氣部分的棱柱體的磁化率設(shè)為零。進(jìn)一步地，棱柱體磁場(chǎng)強(qiáng)度的計(jì)算方法包括以下步驟：步驟S110：根據(jù)地球主磁場(chǎng)模型IGRF，計(jì)算棱柱體中心點(diǎn)的地球主磁場(chǎng)X軸的分量Tx(ξp,ηq,ζr)、Y軸的分量Ty(ξp,ηq,ζr)、Z軸的分量Tz(ξp,ηq,ζr)；步驟S120：按式(1)(2)(3)計(jì)算磁化強(qiáng)度mx(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Tx(ξp,ηq,ζr)(1)my(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Ty(ξp,ηq,ζr)(2)mz(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Tz(ξp,ηq,ζr)(3)其中，(ξp,ηq,ζr)表示目標(biāo)區(qū)域中編號(hào)為(p,q,r)的棱柱體幾何中心坐標(biāo)，χ(ξp,ηq,ζr)表示該棱柱體的磁化率值，Tx(ξp,ηq,ζr)表示(ξp,ηq,ζr)處地球主磁場(chǎng)的X軸上的分量、Ty(ξp,ηq,ζr)表示(ξp,ηq,ζr)處地球主磁場(chǎng)的Y軸上的分量、Tz(ξp,ηq,ζr)表示(ξp,ηq,ζr)處地球主磁場(chǎng)的Z軸上的分量，mx(ξp,ηq,ζr)表示(ξp,ηq,ζr)處磁化強(qiáng)度的X軸上的分量，my(ξp,ηq,ζr)處磁化強(qiáng)度的Y軸上的分量、mz(ξp,ηq,ζr)處磁化強(qiáng)度的Z軸上的分量。本發(fā)明的技術(shù)效果：1、本發(fā)明提供的任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，通過建立包含所有目標(biāo)區(qū)域的規(guī)則棱柱體模型，使得目標(biāo)區(qū)域(包含起伏地形)能完全嵌入在該棱柱體模型中，能很容易地刻畫任意磁化率分布復(fù)雜磁性體以及起伏地形，從而提高后續(xù)模型模擬的計(jì)算精度；而且所構(gòu)建模型剖分方法簡(jiǎn)單、靈活，能提高模擬效率。2、本發(fā)明提供的任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，能夠?qū)崿F(xiàn)磁場(chǎng)梯度張量的快速、高精度數(shù)值模擬，可以滿足大規(guī)模數(shù)據(jù)三維磁化率精細(xì)反演、人機(jī)交互建模和解釋的需求；3、本發(fā)明提供的任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，用于處理大規(guī)模數(shù)值模擬時(shí)，不但計(jì)算效率和計(jì)算精度高，并且所需計(jì)算機(jī)內(nèi)存?。?、本發(fā)明提供的任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，不僅可以應(yīng)用于磁性體感應(yīng)磁化產(chǎn)生的磁場(chǎng)梯度張量數(shù)值模擬，還可以應(yīng)用于磁性體剩余磁化產(chǎn)生的磁場(chǎng)梯度張量數(shù)值模擬；5、本發(fā)明提供的任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法，不但可以精確模擬場(chǎng)源外部的磁場(chǎng)梯度張量，還可以精確模擬場(chǎng)源內(nèi)部的磁場(chǎng)梯度張量，有利于開展井地聯(lián)合反演的研究。具體請(qǐng)參考根據(jù)本發(fā)明的食品泵軸提出的各種實(shí)施例的如下描述，將使得本發(fā)明的上述和其他方面顯而易見。附圖說(shuō)明圖1為本發(fā)明提供任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法流程示意圖；圖2為本發(fā)明提供任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法優(yōu)選實(shí)施例流程圖；圖3為本發(fā)明優(yōu)選實(shí)施例中所處理復(fù)雜磁性體的模型表示圖；圖4為本發(fā)明優(yōu)選實(shí)施例中磁性體剖分后所得棱柱磁性體剖面圖；圖5為本發(fā)明優(yōu)選實(shí)施例的磁場(chǎng)梯度張量計(jì)算值；圖6為本發(fā)明優(yōu)選實(shí)施例的磁場(chǎng)梯度張量理論值；圖7為本發(fā)明優(yōu)選實(shí)施例中x＝0m觀測(cè)線對(duì)應(yīng)理論值與采用本發(fā)明提供方法所得計(jì)算值對(duì)比示意圖。圖例說(shuō)明：L：表示z方向剖分小棱柱體個(gè)數(shù)；M：表示x方向剖分小棱柱體個(gè)數(shù)；N：表示y方向剖分小棱柱體個(gè)數(shù)；χ：表示磁化率。具體實(shí)施方式構(gòu)成本申請(qǐng)的一部分的附圖用來(lái)提供對(duì)本發(fā)明的進(jìn)一步理解，本發(fā)明的示意性實(shí)施例及其說(shuō)明用于解釋本發(fā)明，并不構(gòu)成對(duì)本發(fā)明的不當(dāng)限定。參見圖1，本發(fā)明提供的任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬方法包括以下步驟：步驟S100：建立包含目標(biāo)區(qū)域的棱柱體模型，將棱柱體模型剖分為多個(gè)棱柱體，對(duì)各棱柱體設(shè)置磁化率，并計(jì)算各棱柱體的磁化強(qiáng)度得到組合棱柱體模型；復(fù)雜磁性體模型表示：建立包含目標(biāo)區(qū)域的規(guī)則棱柱體模型，使得目標(biāo)區(qū)域(包含起伏地形)完全嵌入在該棱柱體模型中；將該棱柱體模型劃分成許多小棱柱體，每個(gè)小棱柱體磁化率為常值，不同棱柱體磁化率取值不同，以此刻畫任意磁化率分布復(fù)雜磁性體模型；將位于空氣部分的小棱柱體的磁化率值設(shè)為零，以此刻畫起伏地形。采用這種剖分形式結(jié)合后續(xù)所用步驟能在保證剖分過程簡(jiǎn)單快捷的情況下，保證計(jì)算結(jié)果的精確性。有利于提供計(jì)算效率。磁化強(qiáng)度計(jì)算：可以按現(xiàn)有方法進(jìn)行計(jì)算，優(yōu)選的包括以下步驟：步驟S110：根據(jù)公布的地球主磁場(chǎng)模型IGRF，計(jì)算剖分所得各棱柱體中心點(diǎn)的地球主磁場(chǎng)三分量Tx(ξp,ηq,ζr),Ty(ξp,ηq,ζr),Tz(ξp,ηq,ζr)；步驟S120：計(jì)算磁化強(qiáng)度mx(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Tx(ξp,ηq,ζr)(1)my(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Ty(ξp,ηq,ζr)(2)mz(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Tz(ξp,ηq,ζr)(3)式(1)(2)(3)中，(ξp,ηq,ζr)表示目標(biāo)區(qū)域中編號(hào)為(p,q,r)的棱柱體幾何中心坐標(biāo)；χ(ξp,ηq,ζr)表示該棱柱體的磁化率值；Tx(ξp,ηq,ζr)、Ty(ξp,ηq,ζr)、Tz(ξp,ηq,ζr)分別表示(ξp,ηq,ζr)處地球主磁場(chǎng)的x、y、z分量；mx(ξp,ηq,ζr)、my(ξp,ηq,ζr)、mz(ξp,ηq,ζr)分別表示(ξp,ηq,ζr)處磁化強(qiáng)度的x、y、z分量。按此方法計(jì)算磁化強(qiáng)度，計(jì)算效率得到提高。為了提高棱柱體計(jì)算精度，本發(fā)明提供了組合棱柱體模型磁場(chǎng)梯度張量計(jì)算公式為其中，L表示z方向棱柱體剖分個(gè)數(shù)；M表示x方向棱柱體剖分個(gè)數(shù)；N表示y方向棱柱體剖分個(gè)數(shù)。本發(fā)明提出的公式(4)，為計(jì)算組合棱柱體模型產(chǎn)生的場(chǎng)的精確公式，因而本發(fā)明提供的方法具有較高的計(jì)算精度。同時(shí)配合前述剖分棱柱體組合模型，能實(shí)現(xiàn)對(duì)復(fù)雜形狀磁性體的刻畫準(zhǔn)確，從而提高了場(chǎng)就計(jì)算準(zhǔn)確。步驟S200，按公式(5)采用二維離散卷積法計(jì)算各層棱柱體的磁場(chǎng)梯度張量其中，表示第r層(r＝1,2,…,L)棱柱體組合模型在高度面z0產(chǎn)生的磁場(chǎng)梯度張量，my(ξp,ηq,ζr)、my(ξp,ηq,ζr)、mz(ξp,ηq,ζr)為坐標(biāo)為(ξp,ηq,ζr)的棱柱體在x方向、y方向、z方向的磁化強(qiáng)度分量，hxz(xm-ξp,yn-ηq,z0-ζr)、hyz(xm-ξp,yn-ηq,z0-ζr)、hzz(xm-ξp,yn-ηq,z0-ζr)分別為對(duì)應(yīng)磁化強(qiáng)度三分量的加權(quán)系數(shù)，M表示目標(biāo)區(qū)域x方向棱柱體剖分個(gè)數(shù)，N表示目標(biāo)區(qū)域y方向棱柱體剖分個(gè)數(shù)。優(yōu)選的，加權(quán)系數(shù)計(jì)算：加權(quán)系數(shù)計(jì)算公式為其中，(xm,yn,z0)表示觀測(cè)點(diǎn)坐標(biāo)，z0為常值；其它符號(hào)含義如下x1＝ξp-0.5Δx-xm，x2＝ξp+0.5Δx-xm，y1＝ηq-0.5Δy-yn，y2＝ηq+0.5Δy-yn，z1＝ζr-0.5Δz-z0，z2＝ζr+0.5Δz-z0，μijk＝(-1)i(-1)j(-1)k，i＝1,2，j＝1,2，k＝1,2。采用二維離散卷積計(jì)算方法來(lái)計(jì)算任一層(相對(duì)z方向而言)棱柱體組合模型磁場(chǎng)梯度張量。此處采用二維離散卷積計(jì)算方法配合所提出公式(4)能有效提高該方法的整體處理效率。二維離散卷積計(jì)算方法現(xiàn)有技術(shù)中主要用于處理數(shù)字圖像處理，首次用于處理磁性體模擬計(jì)算。二維離散卷積計(jì)算方法可以按現(xiàn)有方法步驟進(jìn)行，優(yōu)選的，式(5)中包含三個(gè)二維離散卷積，步驟S200中的二維離散卷積快速計(jì)算方法，包括以下步驟：(1)將加權(quán)系數(shù)h(x1-ξp,y1-ηq,z0-ζr)排列成矩陣t，記為式(9)中，矩陣元素ti,j與加權(quán)系數(shù)h(x1-ξp,y1-ηq,z0-ζr)存在關(guān)系ti,j＝h(x1-ξi+1,yj+1-ηq,z0-ζr)(10)(2)將矩陣t補(bǔ)零擴(kuò)展成矩陣將矩陣分成四個(gè)塊矩陣，記為將塊矩陣互換位置，得到矩陣cext(3)將第r層對(duì)應(yīng)加權(quán)系數(shù)h(x1-ξp,y1-ηq,z0-ζr)的磁化強(qiáng)度分量m(ξp,ηq,ζr)，其中p＝1,2,…,M，q＝1,2,…,N，排列成矩陣g，矩陣元素gi,j與密度值存在關(guān)系gi,j＝m(ξi,ηj,ζr)(14)將矩陣g補(bǔ)零擴(kuò)展成矩陣gext式(15)中，0M×N表示M×N零矩陣；(4)計(jì)算式中，fft2()表示二維快速傅里葉變換；(5)計(jì)算式中，“.*”表示對(duì)應(yīng)元素相乘運(yùn)算；(6)計(jì)算式中，ifft2()表示二維快速傅里葉反變換；(7)提取矩陣fext的前M行前N列，構(gòu)成矩陣f，當(dāng)矩陣f中h分別取hxz，hyz，hzz，而m分別對(duì)應(yīng)取mx，my，mz時(shí)，即可得到式(5)中三個(gè)二維離散卷積：的累加結(jié)果。步驟S300，將各層棱柱體組合模型磁場(chǎng)梯度張量進(jìn)行累加，得到整個(gè)組合模型的磁場(chǎng)梯度張量，即本發(fā)明提供方法：(1)通過對(duì)模型剖分形狀和方式進(jìn)行改進(jìn)，給出了精確的棱柱體組合模型磁場(chǎng)梯度張量計(jì)算公式(4)，從而在保證處理效率的同時(shí)提高了計(jì)算精度；(2)本發(fā)明提供方法通過采用計(jì)算效率較高的二維離散卷積算法，進(jìn)一步提高了處理效率?；谏鲜龌A(chǔ)，本發(fā)明提供方法實(shí)現(xiàn)了任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量數(shù)值模擬在計(jì)算效率和計(jì)算精度上的統(tǒng)一。本發(fā)明提供的方法還能用于處理其余分量，以上僅為磁場(chǎng)梯度張量分量Bzz的處理說(shuō)明。參見圖2，下面結(jié)合附圖對(duì)本發(fā)明提供的方法作進(jìn)一步說(shuō)明。1、復(fù)雜磁性體模型表示：首先，建立包含所有目標(biāo)區(qū)域的規(guī)則棱柱體模型，確定棱柱體在x，y，z方向的起始位置，使得目標(biāo)區(qū)域(包含起伏地形)完全嵌入在該棱柱體模型中；其次，根據(jù)實(shí)際問題需求，將棱柱體劃分成許多規(guī)則小棱柱體(如圖3所示)，確定小棱柱體的幾何尺寸Δx，Δy，Δz；最后，根據(jù)目標(biāo)區(qū)域的磁化率分布，對(duì)每個(gè)小棱柱體磁化率進(jìn)行賦值，位于空氣部分的小棱柱體，其磁化率設(shè)為零；2、磁化強(qiáng)度計(jì)算首先，根據(jù)公布的地球主磁場(chǎng)模型IGRF，計(jì)算目標(biāo)區(qū)域的主磁場(chǎng)三分量Tx(ξp,ηq,ζr),Ty(ξp,ηq,ζr),Tz(ξp,ηq,ζr)；其次，計(jì)算磁化強(qiáng)度mx(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Tx(ξp,ηq,ζr)(1)my(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Ty(ξp,ηq,ζr)(2)mz(ξp,ηq,ζr)＝χ(ξp,ηq,ζr)Tz(ξp,ηq,ζr)(3)其中，(ξp,ηq,ζr)表示編號(hào)為(p,q,r)的小棱柱體幾何中心坐標(biāo)；χ(ξp,ηq,ζr)表示該棱柱體的磁化率值；Tx(ξp,ηq,ζr)、Ty(ξp,ηq,ζr)、Tz(ξp,ηq,ζr)分別表示(ξp,ηq,ζr)處地球主磁場(chǎng)的x、y、z分量；mx(ξp,ηq,ζr)、my(ξp,ηq,ζr)、mz(ξp,ηq,ζr)分別表示(ξp,ηq,ζr)處磁化強(qiáng)度的x、y、z分量；3、加權(quán)系數(shù)計(jì)算：加權(quán)系數(shù)計(jì)算公式為其中，(xm,yn,z0)表示觀測(cè)點(diǎn)坐標(biāo)，z0為常值；其它符號(hào)含義如下x1＝ξp-0.5Δx-xm，x2＝ξp+0.5Δx-xm，y1＝ηq-0.5Δy-yn，y2＝ηq+0.5Δy-yn，z1＝ζr-0.5Δz-z0，z2＝ζr+0.5Δz-z0，μijk＝(-1)i(-1)j(-1)k，i＝1,2，j＝1,2，k＝1,24、組合棱柱體模型磁場(chǎng)梯度張量計(jì)算磁場(chǎng)梯度張量計(jì)算公式為其中，L表示z方向棱柱體剖分個(gè)數(shù)；M表示x方向棱柱體剖分個(gè)數(shù)；N表示y方向棱柱體剖分個(gè)數(shù)；首先，采用二維離散卷積快速計(jì)算方法來(lái)計(jì)算一層(相對(duì)z方向而言)棱柱體組合模型磁場(chǎng)梯度張量，其計(jì)算公式為：其中，表示第r層(r＝1,2,…,L)棱柱體組合模型在高度面z0產(chǎn)生的磁場(chǎng)梯度張量；其次，將各層棱柱體組合模型磁場(chǎng)梯度張量進(jìn)行累加，得到整個(gè)組合模型的磁場(chǎng)梯度張量，即其中包含三個(gè)二維離散卷積，以第一個(gè)為例，步驟4中所述的二維離散卷積快速計(jì)算方法，包括以下步驟：(1)將加權(quán)系數(shù)hxz(x1-ξp,y1-ηq,z0-ζr)排列成矩陣t，記為式(9)中，矩陣元素ti,j與加權(quán)系數(shù)hxz(x1-ξp,y1-ηq,z0-ζr)存在關(guān)系ti,j＝h(x1-ξi+1,yj+1-ηq,z0-ζr)(10)(2)將矩陣t補(bǔ)零擴(kuò)展成矩陣將矩陣分成四個(gè)塊矩陣，記為將塊矩陣互換位置，得到矩陣cext(3)將第r層磁化強(qiáng)度mx(ξp,ηq,ζr)(p＝1,2,…,M，q＝1,2,…,N)排列成矩陣g，矩陣元素gi,j與密度值存在關(guān)系gi,j＝m(ξi,ηj,ζr)(14)將矩陣g補(bǔ)零擴(kuò)展成矩陣gext式(15)中，0M×N表示M×N零矩陣；(4)計(jì)算式中，fft2()表示二維快速傅里葉變換；(5)計(jì)算式中，“.*”表示對(duì)應(yīng)元素相乘運(yùn)算；(6)計(jì)算式中，ifft2()表示二維快速傅里葉反變換；(7)提取矩陣fext的前M行前N列，構(gòu)成矩陣f，即為二維離散卷積計(jì)算結(jié)果。參見圖4～7，以下結(jié)合實(shí)例對(duì)本發(fā)明提供方法所得結(jié)果進(jìn)行分析說(shuō)明：為了說(shuō)明本發(fā)明所提出的方法用于計(jì)算任意磁化率分布復(fù)雜磁性體磁場(chǎng)梯度張量時(shí)的效率和精度，設(shè)計(jì)了如下模型算例(圖4所示)：目標(biāo)區(qū)域有一個(gè)棱柱型磁性體，目標(biāo)區(qū)域范圍為：x方向從-1000m到1000m，y方向從-1000m到1000m，z方向從0m到500m(z軸向下為正)；磁性體展布范圍為：x方向從-500m到500m，y方向從-500m到500m，z方向從100m到400m；磁化率為0.01(SI)；目標(biāo)區(qū)域地球主磁場(chǎng)為50000nT，磁偏角為0度，磁傾角為90度。將目標(biāo)區(qū)域剖分成400×400×200個(gè)大小相同的小棱柱體，計(jì)算高度為150m平面(圖4中虛線所示，該平面經(jīng)過場(chǎng)源)上的磁場(chǎng)梯度張量，計(jì)算點(diǎn)個(gè)數(shù)為400×400。本發(fā)明提供方法采用Fortran語(yǔ)言編程實(shí)現(xiàn)，運(yùn)行程序所用的個(gè)人臺(tái)式機(jī)配置為：CPU為i7-2620，主頻為2.7GHz，內(nèi)存為32GB，四核八線程。運(yùn)行所需時(shí)間約為7秒，即可得的模擬結(jié)果，由此可見本發(fā)明提供方法計(jì)算效率高。磁場(chǎng)梯度張量模擬算法計(jì)算值和理論值分別如圖5～6所示，兩者形態(tài)一致，說(shuō)明本發(fā)明提供的模擬方法精度較高。圖7所示為x＝0m觀測(cè)線上計(jì)算值和理論值對(duì)比，兩者吻合度較高。用理論值減去本發(fā)明提供方法所得計(jì)算值，并對(duì)所得差值進(jìn)行統(tǒng)計(jì)，統(tǒng)計(jì)結(jié)果列于表1中。表1本發(fā)明提供方法所得計(jì)算值與磁場(chǎng)梯度張量理論值的誤差統(tǒng)計(jì)表(單位：nT/m)最大值最小值均值均方值理論值1.9-1.91.7*10-160.4誤差6.6*10-13-6.4*10-136.1*10-162.9*10-13可知本發(fā)明提供方法的模擬結(jié)果與理論值誤差較小，精度高。且由于所處理實(shí)例中所處理的觀測(cè)面通過場(chǎng)源內(nèi)部，計(jì)算得結(jié)果即為場(chǎng)源內(nèi)部的場(chǎng)，因而本發(fā)明提供的方法也能精確計(jì)算場(chǎng)源內(nèi)部磁場(chǎng)梯度張量。本領(lǐng)域技術(shù)人員將清楚本發(fā)明的范圍不限制于以上討論的示例，有可能對(duì)其進(jìn)行若干改變和修改，而不脫離所附權(quán)利要求書限定的本發(fā)明的范圍。盡管己經(jīng)在附圖和說(shuō)明書中詳細(xì)圖示和描述了本發(fā)明，但這樣的說(shuō)明和描述僅是說(shuō)明或示意性的，而非限制性的。本發(fā)明并不限于所公開的實(shí)施例。通過對(duì)附圖，說(shuō)明書和權(quán)利要求書的研究，在實(shí)施本發(fā)明時(shí)本領(lǐng)域技術(shù)人員可以理解和實(shí)現(xiàn)所公開的實(shí)施例的變形。在權(quán)利要求書中，術(shù)語(yǔ)“包括”不排除其他步驟或元素，而不定冠詞“一個(gè)”或“一種”不排除多個(gè)。在彼此不同的從屬權(quán)利要求中引用的某些措施的事實(shí)不意味著這些措施的組合不能被有利地使用。權(quán)利要求書中的任何參考標(biāo)記不構(gòu)成對(duì)本發(fā)明的范圍的限制。當(dāng)前第1頁(yè)1 2 3