基于自適應(yīng)雜交遺傳算法的CO2地質(zhì)封存的儲層參數(shù)反演研究

郝艷軍，楊頂輝*，程遠鋒

1 清華大學(xué)數(shù)學(xué)科學(xué)系，北京 1000842 中國石化集團勝利油田分公司物探研究院，山東東營 257001

?

基于自適應(yīng)雜交遺傳算法的CO2地質(zhì)封存的儲層參數(shù)反演研究

郝艷軍1，楊頂輝1*，程遠鋒2

1 清華大學(xué)數(shù)學(xué)科學(xué)系，北京 1000842 中國石化集團勝利油田分公司物探研究院，山東東營 257001

二氧化碳地質(zhì)封存是減少溫室氣體排放和減緩溫室效應(yīng)的重要手段.二氧化碳封存的一個重要組成部分是地震監(jiān)測，即用地震的方法監(jiān)測封存后的二氧化碳的分布變化.為了實現(xiàn)這個目標(biāo)，需要建立儲層參數(shù)與地震性質(zhì)之間的關(guān)系(巖石物理模型)和從地震監(jiān)測數(shù)據(jù)中反演獲得儲層流體的飽和度等參數(shù).首先，本文以Biot理論為基礎(chǔ)，結(jié)合多相流模型研究了多個物理參數(shù)(孔隙度、二氧化碳飽和度、溫度和壓力等)對同時含有二氧化碳和水的孔隙介質(zhì)的波速和衰減等屬性的影響.結(jié)果表明：孔隙度和二氧化碳飽和度對巖石的頻散和衰減屬性影響強烈，而溫度和壓力通過孔隙流體性質(zhì)對巖石的波速產(chǎn)生影響.然后，本文基于含多相流的Biot理論，應(yīng)用抗干擾能力強、且具有更好的局部搜索能力和抗早熟能力的自適應(yīng)雜交遺傳算法對實際數(shù)據(jù)進行了反演研究.對巖心實驗數(shù)據(jù)的反演研究表明了算法的有效性，而且表明含多相流的Biot理論能夠很好地解釋水和二氧化碳飽和巖石的波速特征.最后，我們將自適應(yīng)雜交遺傳算法應(yīng)用于實際封存項目的地震監(jiān)測數(shù)據(jù)，獲得了封存后不同時期的二氧化碳飽和度，達到了用地震方法監(jiān)測二氧化碳分布的目的.

二氧化碳地質(zhì)封存；Biot理論；多相流；遺傳算法；儲層參數(shù)反演

1 引言

從工業(yè)革命開始，人類對化石能源的大規(guī)模利用使大量的二氧化碳被排入大氣層，由此導(dǎo)致的碳循環(huán)不平衡引起了嚴(yán)重的溫室效應(yīng).二氧化碳捕獲和封存(CCS)是減少二氧化碳排放和減緩溫室效應(yīng)的有效途徑(Metz et al.,2005;Bachu,2008;郝艷軍和楊頂輝,2012).經(jīng)過二十多年的發(fā)展，眾多的試驗性二氧化碳封存項目證明CCS可以大幅減少從煤電廠和水泥廠等工業(yè)設(shè)施排放的二氧化碳，從而減緩溫室效應(yīng).可行的封存方式主要有三種：地質(zhì)封存，海洋封存和礦化封存，其中地質(zhì)封存具有最大的封存潛力，并且實施起來也最容易.礦化封存的封存量較少，而海洋封存會使海水酸化，對環(huán)境以及生態(tài)造成未知的影響，因此二者在實用性上都遠不如地質(zhì)封存(Metz et al.,2005).完整的地質(zhì)封存包括三方面：從碳源(煤電廠、水泥廠等)將二氧化碳從廢氣中分離并運輸?shù)椒獯娴攸c，將二氧化碳持續(xù)注入地層中，以及注入之后對二氧化碳的分布進行監(jiān)測以保證安全性.適合于封存二氧化碳的地質(zhì)構(gòu)造包括深部咸水層、枯竭油氣田和不可開采的煤層，其中深部咸水層具有最大的潛力，據(jù)估計可以封存多達104Gt二氧化碳(Metz et al.,2005).為了監(jiān)測注入地下的二氧化碳的分布變化，首先需要對封存之后地層中的所有物理和化學(xué)過程有一個全面的認(rèn)識.注入二氧化碳后，儲層流體(水、CO2、油、氣等等)形成復(fù)雜的多相流，在重力以及各種流相的壓力差驅(qū)動下流動，較輕的CO2會上浮并受到蓋層的阻擋而留在儲層中.二氧化碳溶解于水形成碳酸，進而與儲層巖石發(fā)生反應(yīng)形成新的礦物，在一個相當(dāng)長的時期內(nèi)改變儲層的巖性、孔隙度和滲透率等屬性.由于二氧化碳的物理性質(zhì)對于溫度和壓力十分敏感，注入二氧化碳導(dǎo)致的溫度和壓力變化會使整個物理和化學(xué)過程變得更加復(fù)雜.針對影響這些復(fù)雜物理和化學(xué)過程的重要因素(諸如相對滲透率、毛管壓力、化學(xué)反應(yīng)種類和速率、儲層溫度和壓力等)的模擬研究已有很多(Rutqvist et al.,2001；Pruess,2005；Doughty,2007;Gaus et al.,2008；Birkholzer et al.,2009)，這些工作揭示了二氧化碳地質(zhì)封存的物理和化學(xué)過程的特殊性，為儲層地震監(jiān)測打下了基礎(chǔ).

在二氧化碳封存過程中和封存后，需要長時期的監(jiān)測以確保封存的安全性，避免泄露對周圍人類和環(huán)境的危害，因此監(jiān)測方法至關(guān)重要.經(jīng)過長時期的實踐，地震方法已被證明是最有效且便宜的對二氧化碳分布進行監(jiān)測的手段，且已經(jīng)在很多封存項目中被成功運用(Michael et al.,2010).在美國德克薩斯州Frio構(gòu)造進行的實驗性CO2封存項目中，時移VSP方法和井間地震層析成像方法對CO2的監(jiān)測結(jié)果與數(shù)值模擬所預(yù)測的結(jié)果相符，顯示出大量CO2因相對滲透率降低而被封存在儲層中(Hovorka et al.,2006).挪威Sleipner的CO2封存項目是4D地震成功實施的典范，由于CO2被咸水層中的薄泥巖層所阻擋，不同時期的地震剖面上顯示出非常明顯的反射振幅增強以及時間下推效應(yīng)(Chadwick et al.,2009).In Salah氣田的CO2封存項目綜合使用了地震方法、InSAR以及地球化學(xué)方法等來監(jiān)測封存后儲層的變化，而4D地震方法為項目提供了高質(zhì)量的儲層剖面圖像(Mathieson et al.,2011).Zhang等(2013)應(yīng)用全波形反演的方法對德國Ketzin封存項目的時移地震數(shù)據(jù)進行了反演研究，反演的CO2分布結(jié)果與數(shù)值模擬結(jié)果相吻合.Bergmann等(2016)總結(jié)了在德國Ketzin地區(qū)的CO2實驗性封存項目的地球物理監(jiān)測結(jié)果.在這個歐洲首個陸上CO2封存項目中，研究者綜合使用了各種地球物理方法來監(jiān)測地下CO2的運移，例如4D地震調(diào)查、VSP方法、井間調(diào)查、MSP方法以及電阻率層析成像等方法.其中，4D地震技術(shù)被證明是最為有效的方法，為觀測井位置選取以及封存量估計提供了重要依據(jù)(Bergmann et al.,2016).

二氧化碳注入之后，儲層的各種屬性，如流體組分、飽和度、壓力、溫度、孔隙度、滲透率、巖石骨架和巖性等，都會經(jīng)歷一定程度的變化.這些變化都會反映在儲層的地震屬性上，所以我們需要發(fā)展能夠描述這些地質(zhì)屬性變化對地震屬性影響的理論，即巖石物理模型.巖石物理模型是連接儲層中流體組分和飽和度等屬性與地震屬性的橋梁，也是地震監(jiān)測成功的關(guān)鍵.適合于封存二氧化碳的儲層屬于典型的含流體孔隙介質(zhì)，針對這種雙相(流-固耦合)介質(zhì)中的波傳播理論研究已有半個多世紀(jì)的歷史.在大量的研究工作中，Gassmann(1951)和Biot(1956a,b)的工作，即被后人稱為Biot-Gassmann理論的工作，奠定了含流體多孔隙雙相介質(zhì)理論的基礎(chǔ).Biot通過引入地震波導(dǎo)致的流體宏觀流動(Biot流)，建立了雙相介質(zhì)中地震波傳播的Biot(1956a,b)理論，并預(yù)測了慢P波的存在.模型的有效性需要實驗的驗證，近幾年涌現(xiàn)出了很多測量含二氧化碳的巖石儲層中的波速和衰減實驗研究(Lei and Xue,2009;Siggins et al.,2010;Shi et al.,2011;Lebedev et al.,2013;Nakagawa et al.,2013;Mikhaltsevitch et al.,2014;Zhang et al.,2015).目前多數(shù)研究更偏向于測量超聲波段(0.2～1 MHz)的數(shù)據(jù)，而低頻的波速和衰減測量則存在技術(shù)上的困難(Zhao et al.,2013).直接使用超聲波段測量的數(shù)據(jù)對地震數(shù)據(jù)進行校正是有問題的，因為波速存在頻散的現(xiàn)象.就我們所知，Mikhaltsevitch等(2014)針對含有二氧化碳和水的砂巖儲層進行的低頻(0.1～100 Hz)頻散和品質(zhì)因子的測量是我們目前已知唯一的地震頻段實驗.

地震監(jiān)測的目的是獲得封存之后的隨時間變化的二氧化碳分布(MacBeth et al.,2006;Azuma et al.,2011;Quei?er and Singh,2013)，而儲層參數(shù)反演可以達到這個目的.巖石物理模型提供了儲層參數(shù)與地震數(shù)據(jù)之間的關(guān)系，因此應(yīng)用儲層參數(shù)反演方法可以從地震數(shù)據(jù)中反演得到儲層參數(shù)，尤其是CO2飽和度分布.根據(jù)選取的巖石物理模型的不同，儲層參數(shù)反演可能是有很強非線性的計算過程，一般的局部搜索算法例如共軛梯度法、最速下降法和單純形法等可能很難得到全局最優(yōu)解.在這種情況下，隨機算法可以發(fā)揮更大的作用.遺傳算法作為一種隨機算法已經(jīng)被用到很多反演問題中，在地球物理領(lǐng)域也有廣泛的應(yīng)用(楊磊等，2014；Fang and Yang,2015).但是諸如遺傳算法等隨機算法也面臨著容易早熟和易陷入局部最優(yōu)解等問題.為了解決這些問題，F(xiàn)ang 和 Yang(2015)提出了一種自適應(yīng)雜交遺傳算法，通過引入罰函數(shù)來自適應(yīng)地調(diào)整個體的適應(yīng)值和交叉概率，然后利用模擬退火的方法對接受概率進行評價.通過這些綜合的改進，有效提高了遺傳算法的局部搜索能力和預(yù)防早熟的能力，從而在進行反演的時候可以有效提高反演成功率和計算效率.

本文以Biot理論為基礎(chǔ)，結(jié)合多相流模型來分析CO2和水飽和巖石中的多個物理參數(shù)對波速和衰減等屬性的影響，并應(yīng)用雜交遺傳算法對Mikhaltsevitch等(2014)的巖石物理數(shù)據(jù)進行參數(shù)反演研究，驗證了模型的有效性.最后，應(yīng)用雜交遺傳算法對實際的二氧化碳封存項目的儲層地震監(jiān)測數(shù)據(jù)(Quei?er and Singh,2013)進行儲層參數(shù)反演，獲得了不同時期的二氧化碳飽和度分布，實現(xiàn)了地震監(jiān)測的目的.

2 含多相流的Biot理論及其分析

在本節(jié)中，我們將詳細闡述含多相流的Biot理論，并分析多個物理參數(shù)：CO2飽和度、孔隙度、溫度、孔隙壓力和有效壓力對于快P波和S波的速度和衰減的影響.

2.1 含多相流的Biot理論

部分應(yīng)力形式的三維Biot方程為(Biot,1956a,b;楊磊,2014)

(1)

其中u和U分別是固體和流體的位移向量，λdry=Kdry-2μdry/3和μdry是巖石骨架的拉梅系數(shù)，Kdry是巖石骨架的體積模量，φ和k分別是巖石的孔隙度和滲透率，η是流體的黏度.α=1-Kdry/K0稱為Biot系數(shù)，其中K0是巖石礦物的體積模量.F=[1/Kf+1/(φ Q)]-1，其中Kf是流體的體積模量，參數(shù)Q的表達式為[(1-α)/Kdry-(1-α+φ)/K0]-1.密度項可以表達為ρ11=(1-φ)ρs+ρa，ρ22=φ ρf+ρa，ρ12=-ρa，其中ρs和ρf分別是巖石礦物和流體的密度，ρa是固流耦合密度.將諧波解代入Biot方程中，通過分析可以得到兩種P波和S波的頻散和衰減表達式，詳細表達式見附錄.

由于巖石中同時含有兩種流體(水和二氧化碳)，為了使用Biot理論計算波傳播性質(zhì)需要將兩種流體等效為一種流體，將兩組流體參數(shù)整合為一組流體參數(shù).混和流體的等效流體密度等于兩種流體相對于飽和度的加權(quán)平均，即ρf=Slρl+Sgρg.其中S代表飽和度，下腳標(biāo)l和g分別代表液相和氣相，在我們所考慮的問題中水是液相，二氧化碳是氣相(雖然在深部咸水層中二氧化碳一般處于超臨界流體狀態(tài)，但是為了簡便仍用氣相來稱呼).對于混合流體的黏度，我們采用Carcione等(2006)中的關(guān)系式來計算：η=ηg(ηl/ηg)Sl，其中ηl和ηg分別是液相和氣相的黏度.混合流體的體積模量在計算波度和衰減時是一個很關(guān)鍵的參數(shù).通常來講，兩種流體混合之后的等效模量處于Ruess平均和Voigt平均之間(Mavko et al.,2003)，公式為

(2)

(3)

由Brie公式計算的模量介于Ruess平均和Voigt平均之間.研究表明在計算波速時使用Brie公式計算的混合流體模量能夠得到和實際數(shù)據(jù)更匹配的結(jié)果(Lumley,2010)，但是也要注意其中的參數(shù)e的選擇.在實際中，巖石的彈性性質(zhì)是受到多方面因素控制的，例如巖石的壓實、沉積歷史、黏土含量、巖性和孔隙形狀，孔隙流體的黏度、密度、浸潤性、流體種類和比例，以及環(huán)境的溫度和壓力等等.因此Biot理論的適用性要取決于儲層巖石本身的物性，在應(yīng)用時需要根據(jù)儲層條件謹(jǐn)慎選擇合適的巖石物理理論.

2.2 CO2飽和度對波頻散和衰減的影響

為了分析CO2飽和度對波頻散和衰減的影響，我們根據(jù)Carcione等(2006)中二氧化碳封存的地震模擬的一個算例選取如下一組參數(shù)：巖石干骨架體積模量為Kdry=1.37 GPa，剪切模量為μdry=0.85 GPa；巖石礦物的模量為K0=40 GPa，密度為ρs=2600 kg·m-3.孔隙度為φ=0.3，滲透率為k=10-12m2.水的密度為ρl=995.72 kg·m-3，聲速為vl=1549.4 m·s-1，黏度為ηl=6.3019×10-4Pa·s；CO2的密度為ρg=665.37 kg·m-3，聲速為vg=278.0 m·s-1，黏度為ηg=5.1743×10-5Pa·s.固流耦合密度是ρa=200 kg·m-3.在公式(3)中，選取指數(shù)e=2.0.巖石骨架的這組參數(shù)對應(yīng)于高孔隙度高滲透率的未固結(jié)砂巖.

圖1和圖2給出了當(dāng)CO2的飽和度分別為0.1，0.3，0.65和0.8時，快P波和S波的頻散和衰減曲線.其中，圖1a和1b分別是當(dāng)飽和度取4個不同值時快P波波速和逆品質(zhì)因子隨頻率變化的曲線.從圖1a中可見，在某一固定頻率下，快P波速度隨著CO2飽和度的增加而減小，這說明快P波速度對巖石整體的可壓縮性較為敏感，而當(dāng)CO2飽和度增加時，巖石整體的可壓縮性減??；當(dāng)CO2飽和度為定值時，隨著頻率的增加，快P波的速度略有增加.4條曲線對應(yīng)的最大頻散分別為(按CO2飽和度從小到大排列)：28.9 m·s-1，12.4 m·s-1，2.0 m·s-1，12.0 m·s-1.從圖1b中可見，逆品質(zhì)因子在較高頻率時(>104Hz)隨著CO2的飽和度增加而減小，但是頻率較低時沒有表現(xiàn)出類似的關(guān)系.逆品質(zhì)因子的峰值所對應(yīng)的頻率隨著CO2飽和度的增加而降低，這可以通過對附錄A中的臨界頻率的表達式進行分析而得知；逆品質(zhì)因子的峰值隨著CO2飽和度的增加先減小，后增加.對比圖1a中4條曲線的最大頻散值可以發(fā)現(xiàn)，頻散越大，逆品質(zhì)因子的峰值越高.

圖1 當(dāng)CO2飽和度取不同值時快P波的(a)頻散和(b)衰減曲線Fig.1 (a)Dispersion and(b)attenuation curves of fast P-waves at different CO2 saturation values

圖2 當(dāng)CO2飽和度取不同值時S波的(a)頻散和(b)衰減曲線Fig.2 (a)Dispersion and(b)attenuation curves of S-waves at different CO2 saturation values

圖2a和2b分別是當(dāng)CO2飽和度取不同值時的S波波速和逆品質(zhì)因子隨頻率變化的曲線.從圖2a可見，在某一固定頻率下，S波速度隨著CO2飽和度的增加而增加，說明S波速度對巖石整體的密度更加敏感，對于CO2飽和度的增加而引起的巖石可壓縮性變化不敏感；當(dāng)CO2飽和度為定值時，S波的速度隨著頻率的增加.4條曲線對應(yīng)的最大頻散分別為(按CO2飽和度從小到大排列)：27.4 m·s-1，24.9 m·s-1，20.7 m·s-1，19.0 m·s-1.從圖2b可見，較高頻率下(>104Hz)，逆品質(zhì)因子隨著CO2的飽和度增加而減小，在較低頻率時則呈現(xiàn)相反的趨勢.逆品質(zhì)因子的峰值所對應(yīng)的頻率隨著CO2飽和度的增加而降低，這和快P波的表現(xiàn)相同；逆品質(zhì)因子的最大值隨著CO2飽和度的增加而減小.與圖2a中4條曲線的最大頻散值對比可以發(fā)現(xiàn)，頻散越大，逆品質(zhì)因子的峰值越高.

2.3 飽和度和孔隙度對波速的影響

現(xiàn)在我們來分析飽和度和孔隙度兩個參數(shù)對快P波和S波速度的影響.頻率固定為50Hz，其余參數(shù)和2.2節(jié)中相同.圖3a和3b分別展示了當(dāng)CO2飽和度為0.1、0.3、0.5、0.7和0.9時快P波和S波的波速隨孔隙度的變化.從圖中可見，當(dāng)孔隙度固定時，隨著CO2飽和度的增加，快P波的速度減小，而S波的速度增加，這與圖1a和圖2a一致.當(dāng)CO2飽和度小于等于0.7時，快P波速度隨著孔隙度的增加而減?。坏钱?dāng)CO2飽和度為0.9時，快P波速度隨著孔隙度的增加先減小，后增加.對于S波并沒有這樣的表現(xiàn)，即若CO2飽和度固定，S波速度隨著孔隙度的增加而增加.綜合圖1、2和3可知，P波和S波波速受到巖石可壓縮性和密度的共同作用，而孔隙度和流體飽和度的變化對于巖石可壓縮性和密度的影響是互相耦合的，所以二者對波速的影響是比較復(fù)雜的.

2.4 溫度和孔隙壓力對波速的影響

由于二氧化碳的各種物理參數(shù)：密度、聲速和黏度，對溫度和壓力比較敏感，所以含有CO2的巖石的波速受溫度和壓力的影響比較大.為了模擬溫度和孔隙壓力對快P波和S波速度的影響，需要獲得水和CO2在不同溫度和壓力下的各種物理參數(shù)值.一般來講，需要通過狀態(tài)方程來計算某一種流體的物理參數(shù)，不過這是一個很耗時和復(fù)雜的過程，所以我們使用NIST Chemistry WebBook的數(shù)據(jù)庫(Linstrom and Mallard,2001)來獲得不同溫度和壓力下水和CO2的物理參數(shù).設(shè)定二氧化碳的飽和度為0.5，頻率為50 Hz，除了水和CO2的物理參數(shù)通過數(shù)據(jù)庫獲得外，其他參數(shù)和2.2節(jié)中相同.我們分別計算了溫度T為303.15 K，308.15 K，313.15 K，318.15 K，323.15 K，孔隙壓力從7 MPa增加到15 MPa時的快P波速度和S波速度.圖4a和4b分別給出了快P波速度和S波速度.圖4a顯示，在較低壓力下快P波速度隨溫度的增加而增加；在較高壓力下，快P波速度隨溫度的增加而減小.對于同一溫度，快P波速度隨著孔隙壓力的增加先減小，后增加.對于S波，速度隨溫度和壓力變化的規(guī)律比較清晰.對于同一壓力，S波速度隨著溫度的增加而增加；對于同一溫度，速度隨著壓力的增加而減小.從圖4a和4b中可見當(dāng)T=303.15 K，壓力在7至8 MPa之間時，波速變化最為劇烈，而且顯示出隨壓力的非單調(diào)性變化.這是因為CO2的物性在其臨界點(304.13 K，7.38 MPa)附近變化很大，密度等物理量隨壓力和溫度是呈現(xiàn)非單調(diào)性變化的，而在其他區(qū)域變化則較為平緩(Linstrom and Mallard,2001).從P波和S波不同的表現(xiàn)來看，P波速度對于CO2的聲速和密度變化都比較敏感，所以在CO2的臨界點處會發(fā)生非單調(diào)的變化；S波則只對CO2的密度變化比較敏感，所以其變化趨勢則比較單一.

圖3 不同CO2飽和度下(a)快P波波速和(b)S波波速隨孔隙度的變化Fig.3 (a)Fast P-wave and(b)S-wave velocity as functions of porosity at CO2 saturations

圖4 不同溫度下(a)快P波波速和(b)S波波速隨孔隙壓力的變化Fig.4 (a)Fast P-wave velocity and(b)S-wave velocity as functions of pore pressure at different temperatures

3 應(yīng)用自適應(yīng)雜交遺傳算法對實測數(shù)據(jù)進行反演

遺傳算法的基本原理是將優(yōu)化問題的解編碼轉(zhuǎn)化為遺傳空間的基因，初始化生成一定數(shù)量的種群，計算它們的適應(yīng)值，然后通過選擇算子、交叉算子和變異算子生成新一代的個體，如此迭代直到最優(yōu)解(適應(yīng)值最大).傳統(tǒng)的二進制編碼的遺傳算法的解精度不高、存儲量大、計算效率低，容易產(chǎn)生“海明懸崖”問題，而實數(shù)編碼的遺傳算法能夠解決這些問題.此外，傳統(tǒng)遺傳算法還具有易早熟、容易陷入局部最優(yōu)解等問題.為此，F(xiàn)ang 和 Yang(2015)提出了一種自適應(yīng)的雜交遺傳算法，引入罰函數(shù)調(diào)整個體的適應(yīng)值和交叉概率，以保護種群多樣性和避免早熟現(xiàn)象，同時利用模擬退火的方法對接受概率進行評價.對合成數(shù)據(jù)和實測數(shù)據(jù)的反演研究表明，這種自適應(yīng)雜交遺傳算法具有更好的局部搜索能力，更高的精度，更快的收斂率和計算速度.我們將應(yīng)用此算法對實測數(shù)據(jù)進行反演研究.

3.1 對巖心數(shù)據(jù)的反演

我們選取Mikhaltsevitch等(2014)的數(shù)據(jù)進行反演研究.針對澳大利亞Donnybrook砂巖樣本，Mikhaltsevitch等(2014)測量了超聲頻段的干巖石、水飽和巖石以及水和CO2混合飽和巖石的波速隨有效壓力變化的數(shù)據(jù)，以及三種巖石的低頻(0.1～100 Hz)下的波頻散和衰減的數(shù)據(jù).

首先考察巖石干骨架的波速隨有效壓力的變化關(guān)系.有效壓力是圍壓和孔隙壓力的差值，Eberhart-Phillips等(1989)指出，砂巖的波速是有效壓力的指數(shù)函數(shù).相應(yīng)地，巖石的體積模量和剪切模量也是有效壓力的指數(shù)函數(shù).因此我們采用如下的指數(shù)函數(shù)關(guān)系，公式為

(5)

用P波和S波速度隨有效壓力變化的數(shù)據(jù)來分別反演其中的參數(shù)：a1，b1，c1，a2，b2，c2.定義目標(biāo)函數(shù)為

(6)

一般情況下標(biāo)準(zhǔn)的遺傳算法求解的是一個極大值問題，因此通過下面的公式將(6)式轉(zhuǎn)化為極大值函數(shù)(楊磊等,2014)，求解相應(yīng)的優(yōu)化問題，公式為

(7)

種群個體數(shù)量設(shè)置為40，初始交叉概率為0.85，初始變異概率為0.1，新一代個體保留上一代的10%的精英個體，限制每個個體被選中的次數(shù)為4次.待反演參數(shù)的取值范圍分別是：0

圖5 巖石干骨架的P波和S波速度隨有效壓力變化的實驗數(shù)據(jù)與從反演結(jié)果計算的曲線的比較Fig.5 Comparison of experimental data of dry rock P-wave and S-wave velocity changing with effective pressure and calculated data using inversion results

Vpdry=4188.755-1419.161exp(-0.0791Peff),

Vsdry=2630.520-1050.705exp(-0.0760Peff).

(8)

P波數(shù)據(jù)和S波數(shù)據(jù)反演結(jié)果的適應(yīng)值分別為0.9990，0.9985.反演結(jié)果和實測數(shù)據(jù)的比較如圖5所示，最大的相對誤差為1.69%，說明反演結(jié)果和實驗數(shù)據(jù)吻合得非常好.

為了比較自適應(yīng)遺傳算法與傳統(tǒng)遺傳算法的反演效果，我們同時采用傳統(tǒng)遺傳算法對Mikhaltsevitch等(2014)的干巖石隨有效壓力變化的數(shù)據(jù)進行了反演，并在圖6中展示了兩種算法的目標(biāo)函數(shù)(公式(6))收斂過程.在圖6中，橫坐標(biāo)表示迭代次數(shù)，縱坐標(biāo)表示目標(biāo)函數(shù)值；圖6a和6b分別給出了反演P波和S波速度數(shù)據(jù)時的目標(biāo)函數(shù)收斂過程.從圖6中可見，隨著迭代步數(shù)的增加，目標(biāo)函數(shù)值趨于穩(wěn)定.在P波和S波數(shù)據(jù)的反演過程中，自適應(yīng)雜交遺傳算法的目標(biāo)函數(shù)值一直小于傳統(tǒng)遺傳算法的目標(biāo)函數(shù)值，在趨于穩(wěn)定之后，自適應(yīng)雜交遺傳算法得到了比傳統(tǒng)遺傳算法更精確的反演結(jié)果.這表明，自適應(yīng)雜交遺傳算法比傳統(tǒng)遺傳算法精度更高，收斂速度更快.

下面應(yīng)用多相流Biot理論對水和二氧化碳飽和巖石的P波和S波隨有效壓力變化數(shù)據(jù)進行反演.根據(jù)Mikhaltsevitch等(2014)，實驗時巖心的溫度為315.15 K，孔隙壓力為10 MPa，實驗進行時有效壓力從10 MPa增加到60 MPa.巖石孔隙度為φ=0.1154，滲透率k=0.28 mD，巖石礦物的密度為ρs=2593.3 kg·m-3，礦物模量為Ks=31.7 GPa，固流耦合密度使用ρa=(α-1)φρf計算，曲折度為1.7，波的頻率為0.5 MHz.巖石干骨架的體積模量可以由(8)式和巖石干骨架的密度計算得到.在實驗的溫度和孔隙壓力條件下，水的密度為995.72 kg·m-3，聲速為1549.4 m·s-1，黏度為6.3019×10-4Pa·s；CO2的密度為582.0 kg·m-3，聲速為248.7 m·s-1，黏度為4.3×10-5Pa·s(Linstrom and Mallard,2001).二氧化碳的飽和度為0.6，水飽和度為0.4，由此計算混合流體的密度為747.5 kg·m-3，黏度為1.2585×10-4Pa·s.給定上述參數(shù)，我們來反演混合流體的體積模量.類似于(6)，定義目標(biāo)函數(shù)為

圖6 目標(biāo)函數(shù)的收斂過程(a) P波的反演過程；(b) S波的反演過程.Fig.6 Convergence process of the objective function(a) Inversion process of P-waves；(b) Inversion process of S-wave.

minE(Kf)=

(9)

適應(yīng)值函數(shù)與(7)式的構(gòu)造方式相同.種群個體數(shù)設(shè)置為60，其余參數(shù)和前例相同.待反演的混和流體體積模量的取值范圍很自然地規(guī)定為CO2和水的體積模量形成的區(qū)間，即36 MPa≤Kf≤2.39 GPa.使用自適應(yīng)雜交遺傳算法得到的反演結(jié)果為Kf=695.12 MPa，適應(yīng)值為0.9974.由公式(2)計算的Ruess平均和Voigt平均分別為59.4 MPa和977.8 MPa，而反演結(jié)果處于二者之間，這表明了反演結(jié)果的合理性.圖7比較了用反演結(jié)果計算的波速和實驗數(shù)據(jù)，其中橫坐標(biāo)表示有效壓力，縱坐標(biāo)表示波速.從圖7可見，由反演結(jié)果計算的波速能夠很好地擬合實驗數(shù)據(jù)，通過計算可得P波速度的最大相對誤差為0.56%，S波速度的最大相對誤差為3.06%.這樣的結(jié)果說明含多相流的Biot理論能夠很好地解釋水和CO2飽和的Donnybrook砂巖的速度數(shù)據(jù).

Mikhaltsevitch等(2014)還測量了同一巖石樣本(40%水和60%CO2)的低頻頻散數(shù)據(jù).我們將上面反演得到的混和流體體積模量結(jié)果應(yīng)用于解釋這組低頻頻散數(shù)據(jù).圖8對比了快P波和S波的低頻頻散實驗數(shù)據(jù)和用反演結(jié)果計算的低頻頻散數(shù)據(jù)，其中星形表示實驗數(shù)據(jù)，圓圈表示由反演結(jié)果計算的結(jié)果.圖8表明，使用反演結(jié)果計算的快P波波速和實驗數(shù)據(jù)擬合得非常好(最大相對誤差約為1.85%)，而S波的計算結(jié)果與實驗數(shù)據(jù)的偏差比快P波要大一些(最大相對誤差約為3.45%)，但仍然處在實驗數(shù)據(jù)的誤差范圍以內(nèi).這說明反演的混和流體體積模量結(jié)果和多相流Biot理論能夠很好地描述同時含有水和二氧化碳的Donnybrook砂巖的低頻頻散數(shù)據(jù).

圖7 水和二氧化碳飽和巖石的快P波和S波波速隨有效壓力變化的實驗數(shù)據(jù)和反演結(jié)果的比較Fig.7 Comparison of experimental data of fast P-wave and S-wave velocity changing with effective pressure of the rock saturated with water and CO2 and calculated data using inversion results

圖8 水和二氧化碳飽和巖石的快P波和S波低頻頻散數(shù)據(jù)和反演結(jié)果的比較Fig.8 Comparison of low-frequency dispersion data of fast P-wave and S-wave of the rock saturated with water and CO2 calculated data using inversion results

3.2 對二氧化碳地震監(jiān)測數(shù)據(jù)的反演

挪威Sleipner氣田的二氧化碳封存項目是世界上第一個大規(guī)模商業(yè)性的封存項目.從1996年開始，二氧化碳被注入到深800～1100 m的Utsira深部咸水層中.為了封存的安全性，四維地震技術(shù)被用來監(jiān)測儲層中二氧化碳的分布變化.封存之前的地震調(diào)查于1994年進行，之后的1999、2001、2002、2004、2006和2008年分別進行了三維地震調(diào)查.對四維地震數(shù)據(jù)進行研究獲得的結(jié)果清晰表明，注入的二氧化碳對深部咸水層的地震剖面有十分明顯的影響.Sleipner封存項目是成功運用四維地震技術(shù)的典型例子(Quei?er and Singh,2013).

Quei?er 和 Singh(2013)應(yīng)用非線性全波形反演方法對Sleipner氣田的二氧化碳封存區(qū)域的一個2D垂直剖面的地震數(shù)據(jù)進行了反演研究，獲得了1994、1999和2006年的P波速度剖面圖.反演結(jié)果顯示深部咸水層P波波速剖面有顯著的變化，標(biāo)志著封存的二氧化碳分布的變化.我們應(yīng)用含多相流的Biot方程和自適應(yīng)雜交遺傳算法對Quei?er 和 Singh(2013)中的深部咸水層P波速度剖面進行儲層參數(shù)反演，以得到1999和2006年二氧化碳飽和度的分布.深部咸水層的砂巖儲層中夾雜了數(shù)層泥巖，對CO2上浮形成障礙，有利于封存.根據(jù)Quei?er和Singh(2013)，Utsira砂巖的平均孔隙度φ=0.35，滲透率k=2000 mD，礦物模量Ks=40 GPa，干骨架體積模量Kdry=1.37 GPa，剪切模量μdry=0.85 GPa.泥巖的平均孔隙度φ=0.18，滲透率k=6 mD，礦物模量Ks=20 GPa，干骨架體積模量Kdry=4.7 GPa，剪切模量μdry=0.99 GPa.巖石礦物的密度ρs=2600 kg·m-3.儲層的溫度從頂部的29 ℃增加到底部的39.3 ℃，孔隙壓力從頂部的7.9 MPa增加到底部的10.8 MPa，從地震監(jiān)測數(shù)據(jù)中得到地震波的主頻為38 Hz.

根據(jù)Quei?er 和 Singh(2013)可知儲層由砂巖和泥巖共同組成，所以巖石的等效模量可以通過砂巖和泥巖的模量相對于各自體積分?jǐn)?shù)的Backus平均來計算(Mavko et al.,2003).另外由于泥巖的滲透率非常低(只有6 mD)，所以流體替換主要在砂巖中發(fā)生，故反演CO2飽和度只在砂巖中進行.我們的反演按如下過程進行：利用Backus平均，根據(jù)1994年的速度剖面逐點反演砂巖的體積分?jǐn)?shù)；然后對1999年和2006年的速度剖面逐點反演CO2飽和度，在計算P波速度時利用了之前反演得到的砂巖體積分?jǐn)?shù)以及Backus平均.反演得到的1999年和2006年的飽和度結(jié)果如圖9所示.從圖9可以看出，CO2的分布明顯形成層狀結(jié)構(gòu)，體現(xiàn)出Utsira深部咸水層中砂巖的層狀結(jié)構(gòu)；通過與Quei?er和Singh(2013)中的速度剖面圖進行比較，可以發(fā)現(xiàn)含CO2飽和度高的區(qū)域與1999、2006年的速度剖面相對于1994年變化較大的區(qū)域相吻合，說明了反演結(jié)果的正確性；對比1999年和2006年的CO2分布，可以發(fā)現(xiàn)有明顯變化的區(qū)域也對應(yīng)于速度剖面發(fā)生變化的區(qū)域，這說明反演算法能夠正確反演CO2飽和度分布的變化.通過對Sleipner實際封存項目的地震監(jiān)測數(shù)據(jù)進行反演研究，表明了含多相流的Biot理論和自適應(yīng)雜交遺傳算法在二氧化碳封存問題上的有效性.需要注意的是，Biot理論有其適用范圍，只有當(dāng)儲層巖石的頻散和衰減屬性符合Biot理論的預(yù)測時，才可以應(yīng)用Biot理論進行儲層參數(shù)反演.在實際問題中，需要根據(jù)儲層巖石本身的屬性選擇合適的巖石物理模型.

圖9 (a)1999年和(b)2006年的CO2飽和度反演結(jié)果Fig.9 Inversion results of CO2 saturation of(a)year 1999 and (b)year 2006

4 討論和結(jié)論

本文以Biot理論為基礎(chǔ)，結(jié)合多相流模型，研究了孔隙度、CO2飽和度、溫度和壓力等參數(shù)對同時含有水和二氧化碳的孔隙介質(zhì)中快P波和S波波速和衰減等屬性的影響.數(shù)值模擬結(jié)果表明，孔隙度和二氧化碳飽和度對巖石的頻散和衰減屬性影響強烈，而溫度和壓力的變化通過孔隙流體性質(zhì)對巖石的波速產(chǎn)生影響.對于孔隙儲層，二氧化碳的飽和度是決定波速最直接和最重要的因素之一.數(shù)值實驗表明，當(dāng)CO2飽和度發(fā)生變化，波的各種屬性也會發(fā)生顯著的變化，正因為如此才使地震監(jiān)測成為可能.同時，由于二氧化碳的物理性質(zhì)(密度、聲速、黏度等)對溫度和壓力很敏感，所以這兩個物理參數(shù)也是決定波速的重要因素.另一方面，我們應(yīng)用含多相流的Biot理論和自適應(yīng)雜交遺傳算法對巖心實驗數(shù)據(jù)和低頻頻散數(shù)據(jù)進行了反演.結(jié)果顯示，巖石干骨架的波速以及水和二氧化碳飽和巖石的波速隨有效壓力變化數(shù)據(jù)的反演結(jié)果和實驗數(shù)據(jù)擬合得很好、頻散數(shù)據(jù)的反演結(jié)果和實驗數(shù)據(jù)非常接近，可很好地用于解釋水和二氧化碳飽和的Donnybrook砂巖的波速隨有效壓力的變化.此外，我們將含多相流的Biot理論和自適應(yīng)雜交遺傳算法用于反演Sleipner二氧化碳封存項目的四維地震數(shù)據(jù)，獲得了1999年和2006年的CO2飽和度分布.CO2的飽和度反演結(jié)果和地震監(jiān)測數(shù)據(jù)十分吻合，CO2飽和度高的區(qū)域很好地對應(yīng)于波速發(fā)生較大變化的區(qū)域，同時反演結(jié)果也反映出CO2飽和度隨時間變化的現(xiàn)象.總之，這些反演研究表明了含多相流的Biot理論和自適應(yīng)雜交遺傳算法的正確性和有效性，能夠用于地震監(jiān)測數(shù)據(jù)的儲層參數(shù)反演.這對于二氧化碳封存的四維地震監(jiān)測具有重要意義，可以實現(xiàn)封存后二氧化碳分布的動態(tài)監(jiān)測.

附錄 Biot理論的頻散和衰減關(guān)系

Biot理論的兩種P波(快P波和慢P波)的波速和逆品質(zhì)因子為(楊磊，2014)

其中：

其中A，B，C分別為

B=

其中Mdry=Kdry+4μdry/3為縱波模量，ρ1=(1-φ)ρs，ρ2=φ ρf,ω為角頻率，ωc=η φ/kρf為特征頻率.其他參數(shù)與第2節(jié)中定義相同.

S波的波速和逆品質(zhì)因子為

其中：

Azuma H,Konishi C,Nobuoka D,et al.2011.Quantitative CO2saturation estimation from time lapse sonic logs by consideration of uniform and patchy saturation.Energy Procedia,4:3472-3477.

Bachu S.2008.CO2storage in geological media:role,means,status and barriers to deployment.Progress in Energy and Combustion Science,34(2):254-273.

Bergmann P,Diersch M,G?tz J,et al.2016.Review on geophysical monitoring of CO2injection at Ketzin,Germany.Journal of Petroleum Science and Engineering,139:112-136.

Biot M A.1956a.Theory of propagation of elastic waves in a fluid-saturated porous solid:I.Low-frequency range.J.Acoust.Soc.Am.,28(2):168-178.

Biot M A.1956b.Theory of propagation of elastic waves in a fluid-saturated porous solid:II.Higher-frequency range.J.Acoust.Soc.Am.,28(2):179-191.

Birkholzer J T,Zhou Q L,Tsang C F.2009.Large-scale impact of CO2storage in deep saline aquifers:a sensitivity study on pressure response in stratified systems.International Journal of Greenhouse Gas Control,3(2):181-194.

Brie A,Pampuri F,Marsala A F,et al.1995.Shear sonic interpretation in gas-bearing sands.∥SPE Annual Technical Conference and Exhibition.Dallas,Texas:SPE,701-710.

Carcione J M,Picotti S,Gei D,et al.2006.Physics and seismic modeling for monitoring CO2storage.Pure and Applied Geophysics,163(1):175-207.

Chadwick R A,Noy D,Arts R,et al.2009.Latest time-lapse seismic data from Sleipner yield new insights into CO2plume development.Energy Procedia,1(1):2103-2110.

Doughty C.2007.Modeling geologic storage of carbon dioxide:comparison of non-hysteretic and hysteretic characteristic curves.Energy Conversion and Management,48(6):1768-1781.

Eberhart-Phillips D,Han D H,Zoback M D.1989.Empirical relationships among seismic velocity,effective pressure,porosity,and clay content in sandstone.Geophysics,54(1):82-89.

Fang Z L,Yang D H.2015.Inversion of reservoir porosity,saturation,and permeability based on a robust hybrid genetic algorithm.Geophysics,80(5):R265-R280.

Gassmann F.1951.Elastic waves through a packing of spheres.Geophysics,16(4):673-685.

Gaus I,Audigane P,André L,et al.2008.Geochemical and solute transport modelling for CO2storage,what to expect from it? International Journal of Greenhouse Gas Control,2(4):605-625.

Hao Y J,Yang D H.2012.Research progress of carbon dioxide capture and geological sequestration problem and seismic monitoring research.Progress in Geophys.(in Chinese),27(6):2369-2383,doi:10.6038/j.issn.1004-2903.2012.06.012.

Hovorka S D,Benson S M,Doughty C,et al.2006.Measuring permanence of CO2storage in saline formations:the Frio experiment.Environmental Geosciences,13(2):105-121.

Lebedev M,Pervukhina M,Mikhaltsevitch,V,et al.2013.An experimental study of acoustic responses on the injection of supercritical CO2into sandstones from the Otway Basin.Geophysics,78(4):D293-D306.

Lei X L,Xue Z Q.2009.Ultrasonic velocity and attenuation during CO2injection into water-saturated porous sandstone:Measurements using difference seismic tomography.Physics of the Earth and Planetary Interiors,176(3-4):224-234.

Linstrom P J,Mallard W G.2001.NIST Chemistry webbook;NIST standard reference database No.69.http:∥webbook.nist.gov/chemistry/.

Lumley D.2010.4D seismic monitoring of CO2sequestration.The Leading Edge,29(2):150-155.

MacBeth C,Floricich M,Soldo J.2006.Going quantitative with 4D seismic analysis.Geophysical Prospecting,54(3):303-317.

Mathieson A,Midgely J,Wright I,et al.2011.In Salah CO2storage JIP:CO2sequestration monitoring and verification technologies applied at Krechba,Algeria.Energy Procedia,4:3596-3603.

Mavko G,Mukerji T,Dvorkin J.2003.The Rock Physics Handbook:Tools for Seismic Analysis of Porous Media.Cambridge:Cambridge University Press.Metz B,Davidson O,De Coninck H C,et al.2005.IPCC special report on carbon dioxide capture and storage. Prepared by Working Group III of the Intergovernmental Panel on Climate Change.Cambridge,United Kingdom and New York,NY,USA:Cambridge University Press.

Michael K,Golab A,Shulakova V,et al.2010.Geological storage of CO2in saline aquifers—a review of the experience from existing storage operations.International Journal of Greenhouse Gas Control,4(4):659-667.

Mikhaltsevitch V,Lebedev M,Gurevich B.2014.Measurements of the elastic and anelastic properties of sandstone flooded with supercritical CO2.Geophysical Prospecting,62(6):1266-1277.

Nakagawa S,Kneafsey T J,Daley T M.2013.Laboratory seismic monitoring of supercritical CO2flooding in sandstone cores using the Split Hopkinson Resonant Bar technique with concurrent x-ray Computed Tomography imaging.Geophysical Prospecting,61(2):254-269.

Pruess K.2005.Numerical studies of fluid leakage from a geologic disposal reservoir for CO2show self-limiting feedback between fluid flow and heat transfer.Geophysical Research Letters,32(4):L14404.

Quei?er M,Singh S C.2013.Full waveform inversion in the time lapse mode applied to CO2storage at Sleipner.Geophysical prospecting,61(3):537-555.

Rutqvist J,B?rgesson L,Chijimatsu M,et al.2001.Thermohydromechanics of partially saturated geological media:governing equations and formulation of four finite element models.International Journal of Rock Mechanics and Mining Sciences,38(1):105-127.

Shi J Q,Xue Z Q,Durucan S.2011.Supercritical CO2core flooding and imbibition in Tako sandstone—Influence of sub-core scale heterogeneity.International Journal of Greenhouse Gas Control,5(1):75-87.

Siggins A F,Lwin M,Wisman P.2010.Laboratory calibration of the seismo-acoustic response of CO2saturated sandstones.International Journal of Greenhouse Gas Control,4(6):920-927.

Yang L.2014.Wave dispersion and attenuation in elastic and viscoelastic media with multiphase flow and its applications [Ph.D.Thesis].Beijing:Tsinghua University.

Yang L,Yang D H,Hao Y J,et al.2014.A study on inversion of reservoir parameters under coupling interaction of multiple physical mechanisms.Chinese J.Geophys.(in Chinese),57(8):2678-2686,doi:10.6038/cjg20140826.

Zhang F J,Juhlin C,Ivandic M,et al.2013.Application of seismic full waveform inversion to monitor CO2injection:Modelling and a real data example from the Ketzin site,Germany.Geophysical Prospecting,61(S1):284-299.

Zhang Y,Nishizawa O,Kiyama T,et al.2015.Saturation-path dependency of P-wave velocity and attenuation in sandstone saturated with CO2and brine revealed by simultaneous measurements of waveforms and X-ray computed tomography images.Geophysics,80(4):D403-D415.

Zhao J G,Tang G Y,Deng J X,et al.2013.Determination of rock acoustic properties at low frequency:A differential acoustical resonance spectroscopy device and its estimation technique.Geophysical Research Letters,40(12):2975-2982.

附中文參考文獻

郝艷軍,楊頂輝.2012.二氧化碳地質(zhì)封存問題和地震監(jiān)測研究進展.地球物理學(xué)進展,27(6):2369-2383,doi:10.6038/j.issn.1004-2903.2012.06.012.

楊磊.2014.含多相流彈性及粘彈性介質(zhì)中波頻散與衰減分析及其應(yīng)用[博士論文].北京:清華大學(xué).

楊磊,楊頂輝,郝艷軍等.2014.多種物理機制耦合作用下的儲層介質(zhì)參數(shù)反演研究.地球物理學(xué)報,57(8):2678-2686,doi:10.6038/cjg20140826.

(本文編輯 張正峰)

Reservoir parameter inversion of CO2geological sequestration based on the self-adaptive hybrid genetic algorithm

HAO Yan-Jun1，YANG Ding-Hui1*，Cheng Yuan-Feng2

1 Department of Mathematical Sciences,Tsinghua University,Beijing 100084,China2 Shengli Geophysical Research Institute of SINOPEC,Shandong Dongying 257001,China

Carbon dioxide geological sequestration is an important approach to reduce carbon emission and to mitigate global warming.An important part of CO2sequestration is seismic monitoring,which is to monitor carbon dioxide distribution variation using the seismic method after sequestration.To achieve this goal,we need to establish the relationship between reservoir parameters and seismic properties (rock physics model),and to obtain parameters like reservoir fluid saturation through inversion of seismic monitoring data.First,based on the Biot model and multiphase model,we investigate the effects of several physical parameters (porosity,CO2saturation,temperature and pressure,et al.) on seismic properties such as wave velocity and attenuation of carbon dioxide and water saturated porous media to obtain the regular understanding.The results show that porosity and CO2saturation have a huge impact on the properties of wave dispersion and attenuation,and temperature and pressure can affect the rock velocity through the properties of porous fluids.Next,based on Biot theory with multiphase flow,we apply the self-adaptive hybrid genetic algorithm,which has stronger anti-interference capacity and better ability of local search and anti-interference,to perform inversion of actual data.The inversion of core experimental data indicates the validity of the algorithm,and shows that Biot theory with multiphase flow can explain the wave velocity characters of CO2and water saturated rock.At last,we apply the self-adaptive hybrid genetic algorithm to the seismic monitoring data of actual sequestration project,and obtain the CO2saturation distribution at different periods after sequestration and achieved the purpose of using the seismic method to monitor carbon dioxide distribution.

Carbon dioxide geological sequestration；Biot theory；Multiphase flow；Genetic algorithm；Reservoir parameter inversion

郝艷軍，楊頂輝，程遠鋒.2016.基于自適應(yīng)雜交遺傳算法的CO2地質(zhì)封存的儲層參數(shù)反演研究.地球物理學(xué)報，59(11):4234-4245,

10.6038/cjg20161125.

Hao Y J，Yang D H，Cheng Y F.2016.Reservoir parameter inversion of CO2geological sequestration based on the self-adaptive hybrid genetic algorithm.Chinese J.Geophys.(in Chinese),59(11):4234-4245,doi:10.6038/cjg20161125.

國家基金重大項目(41390452)的資助.

郝艷軍,男,1986年生,清華大學(xué)數(shù)學(xué)科學(xué)系博士研究生,主要研究方向為二氧化碳地質(zhì)封存的數(shù)學(xué)物理模型及其求解.E-mail:haoyanjun@tsinghua.org.cn

*通訊作者 楊頂輝,教授,主要從事計算地球物理、孔隙介質(zhì)波傳播理論等研究.E-mail:dhyang@math.tsinghua.edu.cn

10.6038/cjg20161125

P631

2016-02-11,2016-09-12收修定稿