本發(fā)明涉及磁力勘探和地球物理勘探領(lǐng)域，尤其涉及一種磁化率反演方法及裝置。

背景技術(shù)：
磁測(cè)數(shù)據(jù)在地球物理勘探領(lǐng)域中有著廣泛的應(yīng)用。例如，在礦物勘探中，磁方法可以用于決定體參數(shù)；在石油勘探中，磁方法可以用來映射地下的沉積特性和缺點(diǎn)，用來控制沉積盆地的沉積背景。物理參數(shù)的反演，例如計(jì)算磁化率和磁化強(qiáng)度，是使用磁場(chǎng)數(shù)據(jù)的主要科學(xué)問題。實(shí)際應(yīng)用中，能夠反演目標(biāo)地理區(qū)域的磁化率，得到該目標(biāo)地理區(qū)域的磁化率模型，通過該磁化率模型能夠反映該目標(biāo)地理區(qū)域內(nèi)各個(gè)地理位置的磁化率，進(jìn)而將磁化率符合預(yù)設(shè)要求的某地理位置選擇出來，作為進(jìn)一步的研究對(duì)象。相關(guān)技術(shù)中通常基于總磁場(chǎng)強(qiáng)度(TMI)數(shù)據(jù)反演求解磁化率。然而，由于總磁場(chǎng)強(qiáng)度(TMI)數(shù)據(jù)包含的磁測(cè)數(shù)據(jù)不夠全面，因此導(dǎo)致相關(guān)技術(shù)中反演得到的磁化率不夠精確。

技術(shù)實(shí)現(xiàn)要素：
有鑒于此，本發(fā)明提供了一種磁化率反演方法及裝置，能夠反演得到較為精確的磁化率，緩解相關(guān)技術(shù)中反演得到的磁化率不夠精確的問題。第一方面，本發(fā)明實(shí)施例提供了一種磁化率反演方法，包括：獲取全張量磁梯度數(shù)據(jù)；根據(jù)所述全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型；基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率。結(jié)合第一方面，本發(fā)明實(shí)施例提供了第一方面第一種可能的實(shí)施方式，其中，根據(jù)所述全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型，通過以下公式實(shí)現(xiàn)：Jα(m)＝ρ2(Lm,d)+αΩ(m)Sm＝WmWz其中，Jα(m)表示Tikhonov正則化模型，ρ(Lm,d)表示定義在數(shù)據(jù)域的函數(shù)，Ω(m)表示定義在參數(shù)域的函數(shù)，α表示正則化參數(shù)，L表示離散化緊算子，m表示磁化率向量，d表示所述全張量磁梯度數(shù)據(jù)，Sd表示作用于數(shù)據(jù)的尺度算子，Sm表示作用于模型的尺度算子，k表示所述離散化緊算子L的列號(hào)，i表示所述離散化緊算子L的行號(hào)，diag表示對(duì)角化，Wm表示正則化模型的先驗(yàn)約束，Wz表示作用于模型深度的先驗(yàn)約束，表示大于0的常數(shù)。結(jié)合第一方面，本發(fā)明實(shí)施例提供了第一方面第二種可能的實(shí)施方式，其中，采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：通過以下公式對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率：mk+1＝mk+τkhk,τk＝argτminJα(mk+τhk),其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较颍琂表示所述正則化模型，α表示正則化參數(shù)，表示所述正則化模型的梯度，βFR表示共軛方向參數(shù)變量。結(jié)合第一方面，本發(fā)明實(shí)施例提供了第一方面第三種可能的實(shí)施方式，其中，采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：采用預(yù)條件共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argminJα(mk+τhk),hk＝Pg(m0),ifk＝0,其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示所述正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示所述正則化模型的梯度，表示預(yù)條件的共軛方向參數(shù)變量。結(jié)合第一方面，本發(fā)明實(shí)施例提供了第一方面第四種可能的實(shí)施方式，其中，采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：采用混合共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，混合共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argτminJα(mk+τhk),其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示所述正則化模型，α表示正則化參數(shù)，表示所述正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型共軛方向參數(shù)變量，βDY表示DY共軛方向參數(shù)變量，g表示所述正則化模型的梯度。結(jié)合第一方面，本發(fā)明實(shí)施例提供了第一方面第五種可能的實(shí)施方式，其中，采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：采用預(yù)條件混合共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件混合共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argminJα(mk+τhk),hk＝Pg(m0),ifk＝0,其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较颍琂表示所述正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示所述正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型共軛方向參數(shù)變量，βDY表示DY型共軛方向參數(shù)變量。結(jié)合第一方面上述實(shí)施方式，本發(fā)明實(shí)施例提供了第一方面第六種可能的實(shí)施方式，其中，在基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解之前，所述方法還包括：通過凸集投影技術(shù)限定所述正則化模型的有界域，其中，所述凸集投影技術(shù)通過以下公式實(shí)現(xiàn)：mk+1＝PΠ(mk+τkhk)其中，m表示磁化率向量，k表示迭代次數(shù)，PΠ表示凸集投影算子，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?。第二方面，本發(fā)明實(shí)施例提供了一種磁化率反演裝置，包括：數(shù)據(jù)獲取模塊，用于獲取全張量磁梯度數(shù)據(jù)；模型建立模塊，用于根據(jù)所述全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型；磁化率反演模塊，用于基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率。結(jié)合第二方面，本發(fā)明實(shí)施例提供了第二方面第一種可能的實(shí)施方式，其中，所述磁化率反演模塊用于：采用預(yù)條件混合共軛梯度算法對(duì)所述正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件混合共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argminJα(mk+τhk),hk＝Pg(m0),ifk＝0,其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较颍琂表示所述正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示所述正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型共軛方向參數(shù)變量，βDY表示DY型共軛方向參數(shù)變量。結(jié)合第二方面上述實(shí)施方式，本發(fā)明實(shí)施例提供了第二方面第二種可能的實(shí)施方式，其中，所述裝置還包括：有界域限定模塊，用于通過凸集投影技術(shù)限定所述正則化模型的有界域，其中，所述凸集投影技術(shù)通過以下公式實(shí)現(xiàn)：mk+1＝PΠ(mk+τkhk)其中，m表示磁化率向量，k表示迭代次數(shù)，PΠ表示凸集投影算子，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?。本發(fā)明實(shí)施例中的磁化率反演方法及裝置，首先獲取全張量磁梯度數(shù)據(jù)，然后根據(jù)全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型，最后基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率。由于全張量磁梯度數(shù)據(jù)與基于總磁場(chǎng)強(qiáng)度(TMI)數(shù)據(jù)相比，包含的物理信息更加全面，因此通過本實(shí)施例中的磁化率反演方法及裝置，基于全張量磁梯度數(shù)據(jù)求解磁化率，能夠反演得到較精確的磁化率，緩解相關(guān)技術(shù)中反演得到的磁化率不夠精確的問題。附圖說明為了更清楚地說明本發(fā)明實(shí)施例的技術(shù)方案，下面將對(duì)實(shí)施例中所需要使用的附圖作簡(jiǎn)單地介紹，應(yīng)當(dāng)理解，以下附圖僅示出了本發(fā)明的某些實(shí)施例，因此不應(yīng)被看作是對(duì)范圍的限定，對(duì)于本領(lǐng)域普通技術(shù)人員來講，在不付出創(chuàng)造性勞動(dòng)的前提下，還可以根據(jù)這些附圖獲得其他相關(guān)的附圖。圖1示出本發(fā)明實(shí)施例所提供的磁化率反演方法的流程示意圖；圖2示出本發(fā)明實(shí)施例所提供的利用CPU和GPU并行協(xié)同求解磁化率的方法的流程示意圖；圖3示出本發(fā)明實(shí)施例所提供的磁化率反演裝置的結(jié)構(gòu)示意圖。具體實(shí)施方式為使本發(fā)明實(shí)施例的目的、技術(shù)方案和優(yōu)點(diǎn)更加清楚，下面將結(jié)合本發(fā)明實(shí)施例中附圖，對(duì)本發(fā)明實(shí)施例中的技術(shù)方案進(jìn)行清楚、完整地描述，顯然，所描述的實(shí)施例僅僅是本發(fā)明一部分實(shí)施例，而不是全部的實(shí)施例。通常在此處附圖中描述和示出的本發(fā)明實(shí)施例的組件可以以各種不同的配置來布置和設(shè)計(jì)。因此，以下對(duì)在附圖中提供的本發(fā)明的實(shí)施例的詳細(xì)描述并非旨在限制要求保護(hù)的本發(fā)明的范圍，而是僅僅表示本發(fā)明的選定實(shí)施例?；诒景l(fā)明的實(shí)施例，本領(lǐng)域技術(shù)人員在沒有做出創(chuàng)造性勞動(dòng)的前提下所獲得的所有其他實(shí)施例，都屬于本發(fā)明保護(hù)的范圍?？紤]到相關(guān)技術(shù)中通?；诳偞艌?chǎng)強(qiáng)度(TMI)數(shù)據(jù)反演求解磁化率，由于總磁場(chǎng)強(qiáng)度(TMI)數(shù)據(jù)包含的磁測(cè)數(shù)據(jù)不夠全面，因此導(dǎo)致相關(guān)技術(shù)中反演得到的磁化率不夠精確，本發(fā)明提供了一種磁化率反演方法及裝置，下面結(jié)合實(shí)施例進(jìn)行具體描述。如圖1所示，本發(fā)明實(shí)施例提供了一種磁化率反演方法，該方法包括以下步驟：步驟S102，獲取全張量磁梯度數(shù)據(jù)；步驟S104，根據(jù)全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型；步驟S106，基于CPU(CentralProcessingUnit，中央處理器)和GPU(GraphicsProcessingUnit，圖形處理器)協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率。本發(fā)明實(shí)施例中的磁化率反演方法，首先獲取全張量磁梯度數(shù)據(jù)，然后根據(jù)全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型，最后基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率。由于全張量磁梯度數(shù)據(jù)與基于總磁場(chǎng)強(qiáng)度(TMI)數(shù)據(jù)相比，包含的物理信息更加全面，因此通過本實(shí)施例中的磁化率反演方法，基于全張量磁梯度數(shù)據(jù)求解磁化率，能夠反演得到精確的磁化率，緩解相關(guān)技術(shù)中反演得到的磁化率不夠精確的問題。本領(lǐng)域技術(shù)人員能夠理解，對(duì)正則化模型進(jìn)行求解，反演磁化率，實(shí)際解決的就是正則化模型最小化的問題，滿足正則化模型的值最小的磁化率，就是最終反演得到的磁化率。步驟S102中，能夠通過相關(guān)技術(shù)中的設(shè)備、器材以及多種方法獲取到目標(biāo)地理區(qū)域的全張量磁梯度數(shù)據(jù)，這里不再贅述。考慮到介質(zhì)的磁感應(yīng)強(qiáng)度與磁化率的關(guān)系公式，本實(shí)施例中，將磁化率反演問題歸結(jié)為正則化模型的求解問題。上述步驟S104中，根據(jù)全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型，通過以下公式(1)至公式(7)實(shí)現(xiàn)：Jα(m)＝ρ2(Lm,d)+αΩ(m)(1)Sm＝WmWz(5)公式(1)至公式(7)中，Jα(m)表示Tikhonov正則化模型，ρ(Lm,d)表示定義在數(shù)據(jù)域的函數(shù)，Ω(m)表示定義在參數(shù)域的函數(shù)，α表示正則化參數(shù)，L表示離散化緊算子，m表示磁化率向量，d表示全張量磁梯度數(shù)據(jù)，Sd表示作用于數(shù)據(jù)的尺度算子，Sm表示作用于模型的尺度算子，k表示離散化緊算子L的列號(hào)，i表示離散化緊算子L的行號(hào)，diag(.)表示對(duì)角化，Wm表示正則化模型的先驗(yàn)約束，Wz表示作用于模型深度的先驗(yàn)約束，由函數(shù)w(z)在深度方向離散生成，表示大于0的常數(shù)，w(z)表示基于深度的權(quán)重函數(shù)，z表示埋深，η和z0均表示非負(fù)常量。上述公式(1)至公式(7)中，正則化參數(shù)α>0，優(yōu)選α＝0.001。離散化緊算子L是M×N的矩陣，磁化率向量m的長(zhǎng)度為N，全張量磁梯度數(shù)據(jù)d的長(zhǎng)度為M。在確定作用于模型的尺度算子Sm的過程中，表示將向量m的每個(gè)元素插入公式中得到的對(duì)角形式的新向量，其中表示大于0的較小的常數(shù)，優(yōu)選取1.0×10-10?；谏疃鹊臋?quán)重函數(shù)中，η優(yōu)選2，z0優(yōu)選0?；谏疃鹊臋?quán)重函數(shù)可以使∫V(w(z)m(x,y,z))2dv達(dá)到最小?；谏疃鹊臋?quán)重函數(shù)對(duì)于磁化率反演至關(guān)重要，能夠避免磁化率分布集中在表面的情況。在確定作用于數(shù)據(jù)的尺度算子Sd的過程中，表示將離散化緊算子L的每一行2范數(shù)平方后進(jìn)行對(duì)角化形成的新矩陣，這種選擇方式可以減小反演計(jì)算對(duì)觀測(cè)數(shù)據(jù)(即測(cè)量到的全張量磁梯度數(shù)據(jù)d)的依賴性，得到一個(gè)穩(wěn)定的算法。本領(lǐng)域技術(shù)人員能夠理解，由于Wm和Wz表示正則化模型的先驗(yàn)約束，因此通過上述步驟S104得到的正則化模型與先驗(yàn)約束正則化相關(guān)。本實(shí)施例中，通過建立正則化模型，將物理問題轉(zhuǎn)變成數(shù)學(xué)問題，只需要對(duì)正則化模型的最小化問題進(jìn)行求解，就能夠反演得到磁化率。一種優(yōu)選的實(shí)施方式中，上述步驟S106中，采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：通過以下公式(8)對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率：公式(8)中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示正則化模型，α表示正則化參數(shù)，表示正則化模型的梯度，βFR表示共軛方向參數(shù)變量。需要說明的是，βFR同樣表示Flectcher-Reeves公式。上述公式(8)示出了采用基礎(chǔ)共軛梯度算法求解正則化模型，反演得到磁化率的過程。上述公式(8)中，正則化模型的梯度可以通過公式得到，其中，g(m)表示正則化模型的梯度，α表示正則化參數(shù)，m表示磁化率向量，ρ表示定義在數(shù)據(jù)域的函數(shù)，Ω表示定義在參數(shù)域的函數(shù)?？紤]到通過基礎(chǔ)共軛梯度算法求解正則化模型，由于全張量磁梯度數(shù)據(jù)包含的數(shù)據(jù)量大，導(dǎo)致計(jì)算量大，可能存在算法收斂速度過慢的問題，另一個(gè)優(yōu)選的實(shí)施例中，上述步驟S106中，采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：采用預(yù)條件共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件共軛梯度算法通過以下公式(9)實(shí)現(xiàn)：公式(9)中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示正則化模型的梯度，表示預(yù)條件的共軛方向參數(shù)變量。具體地，在基礎(chǔ)共軛梯度算法的技術(shù)公式上增加預(yù)條件技術(shù)，定義一個(gè)新的變量ν使m＝Sv，能夠由公式(8)推導(dǎo)出公式(9)，將公式(9)稱為預(yù)條件共軛梯度算法，其中P是預(yù)條件矩陣并有P＝SST。對(duì)于適當(dāng)?shù)恼齽t化參數(shù)α>0，可以采用矩陣的對(duì)角形式來確定矩陣P。通過如公式(9)所示的預(yù)條件共軛梯度算法求解正則化模型，具有算法收斂速度快的效果，能夠提高求解速度?？紤]到通過基礎(chǔ)共軛梯度算法求解正則化模型，可能存在步長(zhǎng)小、以至于引起欠定反演問題收斂所需的代次數(shù)多、導(dǎo)致求解時(shí)間過長(zhǎng)的問題，另一個(gè)優(yōu)選的實(shí)施例中，上述步驟S106中，采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：采用混合共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，混合共軛梯度算法通過以下公式(10)至公式(13)實(shí)現(xiàn)：公式(10)至公式(13)中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示正則化模型，α表示正則化參數(shù)，表示正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型(Hestenes-Stiefel)共軛方向參數(shù)變量，βDY表示DY型(Dai-Yuan)共軛方向參數(shù)變量，g表示正則化模型的梯度。通過如公式(10)至公式(13)所示的混合共軛梯度算法，能夠在求解正則化模型的過程中避免小步長(zhǎng)的出現(xiàn)，尤其防止病態(tài)問題中，小步長(zhǎng)可能會(huì)引起誤差累計(jì)并導(dǎo)致迭代步數(shù)增加才能收斂的問題，從而加快求解速度。在另一個(gè)優(yōu)選的實(shí)施方式中，為了加快共軛梯度算法的收斂速度，并減少求解過程中的迭代次數(shù)，上述步驟S106中，采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率，包括：采用預(yù)條件混合共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件混合共軛梯度算法通過以下公式(14)至公式(17)實(shí)現(xiàn)：公式(14)至公式(17)中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型(Hestenes-Stiefel)共軛方向參數(shù)變量，βDY表示DY型(Dai-Yuan)共軛方向參數(shù)變量。通過如公式(14)至公式(17)所示的預(yù)條件混合共軛梯度算法求解正則化模型，能夠加快共軛梯度算法的收斂速度，并減少求解過程中的迭代次數(shù)，達(dá)到快速反演磁化率的效果。本發(fā)明實(shí)施例提供的方法流程中，為了使反演得到的磁化率有界，保證反演得到的磁化率的值符合物理意義，得到高效的共軛梯度算法，在基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解之前，還可以包括以下動(dòng)作：通過凸集投影技術(shù)限定正則化模型的有界域，其中，凸集投影技術(shù)通過以下公式(18)至公式(19)實(shí)現(xiàn)：mk+1＝PΠ(mk+τkhk)(18)公式(18)至公式(19)中，m表示磁化率向量，k表示迭代次數(shù)，PΠ表示凸集投影算子，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较颉＞唧w地，對(duì)于計(jì)算磁化率問題，假設(shè)磁化率m是有界的。在本實(shí)施例中，只考慮鐵磁性材料，模型磁化率值為正。定義正則化模型Jα(m)＝ρ2(Lm,d)+αΩ(m)的可行域?yàn)閷?duì)于迭代方法，一旦出現(xiàn)一個(gè)迭代點(diǎn)超出可行域Π，就會(huì)通過凸集投影算子PΠ投影回可行域中。在有界域Π的凸集投影算子可以通過等式(PΠx)(t)＝χΠ(t)x(t)決定，其中，χΩ(t)是有界域Π的特征方程。凸集投影算子PΠ作用于在Π域外全為0的所有函數(shù)構(gòu)成的子空間。PΠ(x)的第i個(gè)構(gòu)成為假設(shè)當(dāng)前循環(huán)中mk可用，則下次循環(huán)可以通過公式mk+1＝PΠ(mk+τkhk)獲得，其中，hk是搜索方向，τk是由前面提到的非確定線性搜索計(jì)算得到的步長(zhǎng)。通過如公式(18)至公式(19)所示的凸集投影技術(shù)限定正則化模型的有界域，能夠增加正則化模型的解的物理約束，采用基于投影到凸集的梯度投影算法來解決正則化問題，能夠生成高效的共軛梯度算法。另外，本實(shí)施例中，為了保證正則化模型Jα(m)＝ρ2(Lm,d)+αΩ(m)是遞減的，在基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解之前，還可以包括以下動(dòng)作：確定共軛梯度算法的迭代步長(zhǎng)的搜索條件，將搜索條件作為求解正則化模型的約束條件；其中，搜索條件通過以下公式(20)至公式(21)確定：公式(20)至公式(21)中，表示目標(biāo)函數(shù)，τ表示步長(zhǎng)，k表示迭代次數(shù)，γ1和γ2均表示正常數(shù)，表示的一階導(dǎo)數(shù)，J表示正則化模型，α表示正則化參數(shù)，m表示磁化率向量，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?。需要說明的是，γ1和γ2是兩個(gè)正常數(shù)并滿足γ1<γ2<1。優(yōu)選γ1＝0.4，γ2＝0.6。由于φ(τ)是下有界的并因此，由中值定理可以得出，存在一個(gè)參數(shù)τ使公式成立。本實(shí)施例中，基于CPU和GPU協(xié)同并行的方式求解磁化率，利用GPU處理大規(guī)模數(shù)值計(jì)算問題，能夠加快運(yùn)算速度。利用本實(shí)施例中的方法，將TMI數(shù)據(jù)、2D和3D的全張量磁梯度數(shù)據(jù)都用于進(jìn)行對(duì)比，結(jié)合大量實(shí)驗(yàn)表明全張量磁梯度數(shù)據(jù)可以得到更高精度的反演結(jié)果，并且利用CPU\GPU并行方式可以顯著提高計(jì)算效率。本實(shí)施例中，上述步驟S106中，基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率，包括以下過程1)、2)、3)、4)：1)給出初始磁化率向量m0，預(yù)條件矩陣P，輸入磁梯度張量數(shù)據(jù)d，設(shè)置參數(shù)γ1、γ2、α，停止迭代求解的閾值ε>0，并設(shè)置迭代索引k:＝0；2)計(jì)算gk＝g(mk)；如果||gk||≤ε，停止計(jì)算；否則，運(yùn)行步驟3)；其中，ε表示停止迭代求解的閾值；3)計(jì)算mk+1＝mk+τkhk，其中τk、hk由公式(8)、公式(9)、公式(10)或公式(14)決定；4)設(shè)置k:＝k+1，轉(zhuǎn)到過程2)。能夠理解，上述過程1)2)3)4)中，迭代求解的停止標(biāo)準(zhǔn)由公式||gk||≤ε決定，當(dāng)滿足公式||gk||≤ε時(shí)，確定求解得到磁化率，迭代停止。本實(shí)施例中的方法在實(shí)際執(zhí)行時(shí)，可以將正則化模型的建立、正則化模型的求解過程以程序的形式輸入計(jì)算機(jī)，操作人員只需在計(jì)算機(jī)中輸入初始磁化率向量m0，預(yù)條件矩陣P，磁梯度張量數(shù)據(jù)d，設(shè)置γ1＝0.4，γ2＝0.6，α＝0.001，ε>0，并設(shè)置迭代索引k＝0，然后由計(jì)算機(jī)自動(dòng)執(zhí)行上述方法，求解得到磁化率。本發(fā)明實(shí)施例還提供了一種如圖2所示的利用CPU和GPU并行協(xié)同求解磁化率的方法，如圖2所示，該方法包括以下步驟：步驟S202，在CPU中初始化算法：給出初始磁化率向量m0，預(yù)條件矩陣P，輸入磁梯度張量數(shù)據(jù)d，設(shè)置γ1＝0.4，γ2＝0.6，α＝0.001，ε>0，并設(shè)置迭代索引k＝0；步驟S204，將初始化得到的初始磁化率向量，預(yù)條件矩陣，磁梯度張量數(shù)據(jù)復(fù)制到GPU全局存儲(chǔ)器中；步驟S206，將初始化數(shù)據(jù)從GPU全局存儲(chǔ)器復(fù)制到GPU共享存儲(chǔ)器中；步驟S208，開始共軛梯度算法，在GPU端進(jìn)行算法中大規(guī)模矩陣向量的計(jì)算、更新操作；步驟S210，將需要進(jìn)行邏輯運(yùn)算的標(biāo)量數(shù)據(jù)傳回CPU；步驟S212，判斷是否達(dá)到迭代停止標(biāo)準(zhǔn)(如達(dá)到規(guī)定的迭代次數(shù)或誤差精度)，若達(dá)到停止標(biāo)準(zhǔn)轉(zhuǎn)到步驟S214，否則轉(zhuǎn)到步驟S208，繼續(xù)迭代操作；步驟S214，將計(jì)算得到的最優(yōu)化結(jié)果復(fù)制到CPU中，并輸出結(jié)果。與相關(guān)技術(shù)相比，本發(fā)明實(shí)施例具有以下優(yōu)點(diǎn)：使用全張量磁梯度數(shù)據(jù)來反演感興趣的參數(shù)是相對(duì)較新的研究，關(guān)鍵原因之一是全張量磁梯度數(shù)據(jù)很難獲得。本發(fā)明實(shí)施例中，將Tikhonov正則化模型應(yīng)用于反演物理參數(shù)。為了實(shí)現(xiàn)這個(gè)目的，考慮對(duì)模型和數(shù)據(jù)進(jìn)行先驗(yàn)約束。在求解過程中，發(fā)展了一種快速收斂的預(yù)條件的凸集投影混合共軛梯度算法，并通過CPU與GPU聯(lián)合計(jì)算來提高算法過程中大規(guī)模矩陣計(jì)算效率。二維和三維的擬合數(shù)據(jù)實(shí)例證明采用全張量磁梯度數(shù)據(jù)反演的磁化率結(jié)果相對(duì)由TMI數(shù)據(jù)計(jì)算所得的結(jié)果精確度更高。對(duì)應(yīng)上述提供的磁化率反演方法，本發(fā)明實(shí)施例還提供了一種磁化率反演裝置，該磁化率反演裝置用于執(zhí)行上述磁化率反演方法，因此適用于上述方法的描述同樣適用于該磁化率反演裝置。如圖3所示，本實(shí)施例中的磁化率反演裝置包括：數(shù)據(jù)獲取模塊31，用于獲取全張量磁梯度數(shù)據(jù)；模型建立模塊32，用于根據(jù)全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型；磁化率反演模塊33，用于基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率。本發(fā)明實(shí)施例中的磁化率反演裝置，首先獲取全張量磁梯度數(shù)據(jù)，然后根據(jù)全張量磁梯度數(shù)據(jù)建立Tikhonov正則化模型，最后基于CPU和GPU協(xié)同并行方式采用共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率。由于全張量磁梯度數(shù)據(jù)與基于總磁場(chǎng)強(qiáng)度(TMI)數(shù)據(jù)相比，包含的物理信息更加全面，因此通過本實(shí)施例中的磁化率反演裝置，基于全張量磁梯度數(shù)據(jù)求解磁化率，能夠反演得到精確的磁化率，緩解相關(guān)技術(shù)中反演得到的磁化率不夠精確的問題。考慮到介質(zhì)的磁感應(yīng)強(qiáng)度與磁化率的關(guān)系公式，上述裝置中，模型建立模塊32用于：通過以下公式建立Tikhonov正則化模型：Jα(m)＝ρ2(Lm,d)+αΩ(m)Sm＝WmWz其中，Jα(m)表示Tikhonov正則化模型，ρ(Lm,d)表示定義在數(shù)據(jù)域的函數(shù)，Ω(m)表示定義在參數(shù)域的函數(shù)，α表示正則化參數(shù)，L表示離散化緊算子，m表示磁化率向量，d表示全張量磁梯度數(shù)據(jù)，Sd表示作用于數(shù)據(jù)的尺度算子，Sm表示作用于模型的尺度算子，k表示離散化緊算子L的列號(hào)，i表示離散化緊算子L的行號(hào)，diag(.)表示對(duì)角化，Wm表示正則化模型的先驗(yàn)約束，Wz表示作用于模型深度的先驗(yàn)約束，由函數(shù)w(z)在深度方向離散生成，表示大于0的常數(shù)，w(z)表示基于深度的權(quán)重函數(shù)，z表示埋深，η和z0均表示非負(fù)常量。一種優(yōu)選的實(shí)施方式中，上述裝置中，磁化率反演模塊33用于：通過以下公式對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率：mk+1＝mk+τkhk,τk＝argτminJα(mk+τhk),其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较颍琂表示正則化模型，α表示正則化參數(shù)，表示正則化模型的梯度，βFR表示共軛方向參數(shù)變量。考慮到通過基礎(chǔ)共軛梯度算法求解正則化模型，由于全張量磁梯度數(shù)據(jù)包含的數(shù)據(jù)量大，導(dǎo)致計(jì)算量大，可能存在算法收斂速度過慢的問題，另一個(gè)優(yōu)選的實(shí)施例中，磁化率反演模塊33用于：采用預(yù)條件共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argminJα(mk+τhk),hk＝Pg(m0),ifk＝0,其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示正則化模型的梯度，表示預(yù)條件的共軛方向參數(shù)變量。通過預(yù)條件共軛梯度算法求解正則化模型，具有算法收斂速度快的效果，能夠提高求解速度?？紤]到通過基礎(chǔ)共軛梯度算法求解正則化模型，可能存在步長(zhǎng)小、以至于引起欠定反演問題收斂所需的代次數(shù)多、導(dǎo)致求解時(shí)間過長(zhǎng)的問題，另一個(gè)優(yōu)選的實(shí)施例中，磁化率反演模塊33用于：采用混合共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，混合共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argτminJα(mk+τhk),其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较颍琂表示正則化模型，α表示正則化參數(shù)，表示正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型(Hestenes-Stiefel)共軛方向參數(shù)變量，βDY表示DY型(Dai-Yuan)共軛方向參數(shù)變量，g表示正則化模型的梯度。通過混合共軛梯度算法，能夠在求解正則化模型的過程中避免小步長(zhǎng)的出現(xiàn)，尤其防止病態(tài)問題中，小步長(zhǎng)可能會(huì)引起誤差累計(jì)并導(dǎo)致迭代步數(shù)增加才能收斂的問題，從而加快求解速度。在另一個(gè)優(yōu)選的實(shí)施方式中，為了加快共軛梯度算法的收斂速度，并減少求解過程中的迭代次數(shù)，磁化率反演模塊33用于：采用預(yù)條件混合共軛梯度算法對(duì)正則化模型進(jìn)行迭代求解，反演得到磁化率；其中，預(yù)條件混合共軛梯度算法通過以下公式實(shí)現(xiàn)：mk+1＝mk+τkhk,τk＝argminJα(mk+τhk),hk＝Pg(m0),ifk＝0,其中，m表示磁化率向量，k表示迭代次數(shù)，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?，J表示正則化模型，α表示正則化參數(shù)，P表示預(yù)條件矩陣，g表示正則化模型的梯度，βhybrid表示共軛方向混合參數(shù)變量，βHS表示HS型(Hestenes-Stiefel)共軛方向參數(shù)變量，βDY表示DY型(Dai-Yuan)共軛方向參數(shù)變量。通過預(yù)條件混合共軛梯度算法求解正則化模型，能夠加快共軛梯度算法的收斂速度，并減少求解過程中的迭代次數(shù)，達(dá)到快速反演磁化率的效果。為了使反演得到的磁化率有界，保證反演得到的磁化率的值符合物理意義，得到高效的共軛梯度算法，上述裝置還包括：有界域限定模塊，用于通過凸集投影技術(shù)限定正則化模型的有界域，其中，凸集投影技術(shù)通過以下公式實(shí)現(xiàn)：mk+1＝PΠ(mk+τkhk)其中，m表示磁化率向量，k表示迭代次數(shù)，PΠ表示凸集投影算子，τ表示步長(zhǎng)，h表示在負(fù)梯度方向?yàn)槌跏挤较虻乃阉鞣较?。通過凸集投影技術(shù)限定正則化模型的有界域，能夠增加正則化模型的解的物理約束，采用基于投影到凸集的梯度投影算法來解決正則化問題，能夠生成高效的共軛梯度算法。與相關(guān)技術(shù)相比，本發(fā)明實(shí)施例具有以下優(yōu)點(diǎn)：使用全張量磁梯度數(shù)據(jù)來反演感興趣的參數(shù)是相對(duì)較新的研究，關(guān)鍵原因之一是全張量磁梯度數(shù)據(jù)很難獲得。本發(fā)明實(shí)施例中，將Tikhonov正則化模型應(yīng)用于反演物理參數(shù)。為了實(shí)現(xiàn)這個(gè)目的，考慮對(duì)模型和數(shù)據(jù)進(jìn)行先驗(yàn)約束。在求解過程中，發(fā)展了一種快速收斂的預(yù)條件的凸集投影共軛梯度算法，并通過CPU與GPU聯(lián)合計(jì)算來提高算法過程中大規(guī)模矩陣計(jì)算效率。二維和三維的擬合數(shù)據(jù)實(shí)例證明采用全張量磁梯度數(shù)據(jù)反演的磁化率結(jié)果相對(duì)由TMI數(shù)據(jù)計(jì)算所得的結(jié)果精確度更高。本發(fā)明實(shí)施例所提供的磁化率反演裝置可以為設(shè)備上的特定硬件或者安裝于設(shè)備上的軟件或固件等。本發(fā)明實(shí)施例所提供的裝置，其實(shí)現(xiàn)原理及產(chǎn)生的技術(shù)效果和前述方法實(shí)施例相同，為簡(jiǎn)要描述，裝置實(shí)施例部分未提及之處，可參考前述方法實(shí)施例中相應(yīng)內(nèi)容。所屬領(lǐng)域的技術(shù)人員可以清楚地了解到，為描述的方便和簡(jiǎn)潔，前述描述的系統(tǒng)、裝置和單元的具體工作過程，均可以參考上述方法實(shí)施例中的對(duì)應(yīng)過程，在此不再贅述。在本發(fā)明所提供的實(shí)施例中，應(yīng)該理解到，所揭露裝置和方法，可以通過其它的方式實(shí)現(xiàn)。以上所描述的裝置實(shí)施例僅僅是示意性的，例如，所述單元的劃分，僅僅為一種邏輯功能劃分，實(shí)際實(shí)現(xiàn)時(shí)可以有另外的劃分方式，又例如，多個(gè)單元或組件可以結(jié)合或者可以集成到另一個(gè)系統(tǒng)，或一些特征可以忽略，或不執(zhí)行。另一點(diǎn)，所顯示或討論的相互之間的耦合或直接耦合或通信連接可以是通過一些通信接口，裝置或單元的間接耦合或通信連接，可以是電性，機(jī)械或其它的形式。所述作為分離部件說明的單元可以是或者也可以不是物理上分開的，作為單元顯示的部件可以是或者也可以不是物理單元，即可以位于一個(gè)地方，或者也可以分布到多個(gè)網(wǎng)絡(luò)單元上。可以根據(jù)實(shí)際的需要選擇其中的部分或者全部單元來實(shí)現(xiàn)本實(shí)施例方案的目的。另外，在本發(fā)明提供的實(shí)施例中的各功能單元可以集成在一個(gè)處理單元中，也可以是各個(gè)單元單獨(dú)物理存在，也可以兩個(gè)或兩個(gè)以上單元集成在一個(gè)單元中。所述功能如果以軟件功能單元的形式實(shí)現(xiàn)并作為獨(dú)立的產(chǎn)品銷售或使用時(shí)，可以存儲(chǔ)在一個(gè)計(jì)算機(jī)可讀取存儲(chǔ)介質(zhì)中?；谶@樣的理解，本發(fā)明的技術(shù)方案本質(zhì)上或者說對(duì)現(xiàn)有技術(shù)做出貢獻(xiàn)的部分或者該技術(shù)方案的部分可以以軟件產(chǎn)品的形式體現(xiàn)出來，該計(jì)算機(jī)軟件產(chǎn)品存儲(chǔ)在一個(gè)存儲(chǔ)介質(zhì)中，包括若干指令用以使得一臺(tái)計(jì)算機(jī)設(shè)備(可以是個(gè)人計(jì)算機(jī)，服務(wù)器，或者網(wǎng)絡(luò)設(shè)備等)執(zhí)行本發(fā)明各個(gè)實(shí)施例所述方法的全部或部分步驟。而前述的存儲(chǔ)介質(zhì)包括：U盤、移動(dòng)硬盤、只讀存儲(chǔ)器(ROM，Read-OnlyMemory)、隨機(jī)存取存儲(chǔ)器(RAM，RandomAccessMemory)、磁碟或者光盤等各種可以存儲(chǔ)程序代碼的介質(zhì)。應(yīng)注意到：相似的標(biāo)號(hào)和字母在下面的附圖中表示類似項(xiàng)，因此，一旦某一項(xiàng)在一個(gè)附圖中被定義，則在隨后的附圖中不需要對(duì)其進(jìn)行進(jìn)一步定義和解釋，此外，術(shù)語“第一”、“第二”、“第三”等僅用于區(qū)分描述，而不能理解為指示或暗示相對(duì)重要性。最后應(yīng)說明的是：以上所述實(shí)施例，僅為本發(fā)明的具體實(shí)施方式，用以說明本發(fā)明的技術(shù)方案，而非對(duì)其限制，本發(fā)明的保護(hù)范圍并不局限于此，盡管參照前述實(shí)施例對(duì)本發(fā)明進(jìn)行了詳細(xì)的說明，本領(lǐng)域的普通技術(shù)人員應(yīng)當(dāng)理解：任何熟悉本技術(shù)領(lǐng)域的技術(shù)人員在本發(fā)明揭露的技術(shù)范圍內(nèi)，其依然可以對(duì)前述實(shí)施例所記載的技術(shù)方案進(jìn)行修改或可輕易想到變化，或者對(duì)其中部分技術(shù)特征進(jìn)行等同替換；而這些修改、變化或者替換，并不使相應(yīng)技術(shù)方案的本質(zhì)脫離本發(fā)明實(shí)施例技術(shù)方案的精神和范圍。都應(yīng)涵蓋在本發(fā)明的保護(hù)范圍之內(nèi)。因此，本發(fā)明的保護(hù)范圍應(yīng)所述以權(quán)利要求的保護(hù)范圍為準(zhǔn)。