改進(jìn)的蟻群算法及其在南嶺地區(qū)花崗巖侵入體探測(cè)中的應(yīng)用

劉 雙，陳國(guó)雄，劉天佑

1.中國(guó)地質(zhì)大學(xué)地球物理與空間信息學(xué)院，武漢 430074

2.中國(guó)地質(zhì)大學(xué)地質(zhì)過(guò)程與礦產(chǎn)國(guó)家重點(diǎn)實(shí)驗(yàn)室，武漢 430074

0 引言

南嶺地區(qū)地質(zhì)構(gòu)造復(fù)雜，礦產(chǎn)資源豐富，是世界上研究燕山期大陸成礦體系和花崗巖成巖成礦理論最典型的地區(qū)之一。研究南嶺地區(qū)花崗巖侵入體的侵位特征，對(duì)研究區(qū)深部控礦規(guī)律、深部成礦模式的認(rèn)識(shí)以及深部成礦預(yù)測(cè)意義重大。而重磁勘探是研究該區(qū)區(qū)域深部構(gòu)造、花崗巖侵位以及成礦預(yù)測(cè)的重要手段，前人做了很多相關(guān)工作。秦葆瑚［1－3］指出，湘南一些大型－特大型內(nèi)生多金屬礦田不但具有明顯的特殊地質(zhì)條件，而且在區(qū)域重、磁場(chǎng)上有明顯的異常反映；并對(duì)水口山、東坡、瑤崗仙、黃沙坪4個(gè)大型多金屬礦田的區(qū)域重磁異常特征進(jìn)行了總結(jié)。李少蓮［4］在研究湘南地區(qū)隱伏構(gòu)造中指出了殼型和幔型花崗巖類的重磁場(chǎng)特征：殼型花崗巖類重力異常規(guī)模較大，航磁異常為正負(fù)伴生；而幔型花崗巖類通常處于重力梯度帶或者重力梯度帶交叉部位，航磁ΔT異常規(guī)模較小。段利華［5］利用湘桂粵贛區(qū)域重力垂向二次導(dǎo)數(shù)、重力區(qū)域場(chǎng)垂向二次導(dǎo)數(shù)以及航磁局部場(chǎng)分析淺部和深部的花崗巖帶信息，推斷出區(qū)內(nèi)的隱伏花崗巖體信息，認(rèn)為地幔凹區(qū)是成礦的有利地區(qū)，礦體的形成與隱伏巖體及其形態(tài)相關(guān)。朱自強(qiáng)等［6］利用神經(jīng)網(wǎng)絡(luò)反推算法，對(duì)千里山巖體及其鄰區(qū)進(jìn)行三維重磁反演，給出了巖體的空間形態(tài)及侵位，并預(yù)測(cè)找礦靶區(qū)。饒家榮等［7］根據(jù)南嶺地區(qū)地質(zhì)地球物理資料，研究了南嶺中段主要深斷裂和隱伏花崗巖帶的分布特征，指出大義山－郴州－大寶山北西向隱伏構(gòu)造帶控制著大型－超大型礦床，并提出以地球物理和地球化學(xué)資料為基礎(chǔ)，創(chuàng)建具有“重低磁高”或重磁變異特征的半隱伏－隱伏花崗巖帶的構(gòu)造巖漿系列礦床成礦模式，指導(dǎo)和開展南嶺地區(qū)深部找礦。

蟻群智能優(yōu)化算法（ant colony optimization，ACO）是繼模擬退火算法、遺傳算法、人工神經(jīng)網(wǎng)絡(luò)等啟發(fā)式隨機(jī)搜索算法以后的又一種啟發(fā)式算法。該算法在控制與決策問(wèn)題中取得了較好的效果，例如旅行商問(wèn)題（TSP）［8－13］，背包問(wèn)題（QAP）［14－20］，指派問(wèn)題（JSP）［21］等經(jīng)典組合優(yōu)化問(wèn)題。組合優(yōu)化問(wèn)題是在離散空間尋求最優(yōu)解，然而，地球物理反問(wèn)題是多維連續(xù)目標(biāo)函數(shù)的優(yōu)化，解空間比組合優(yōu)化問(wèn)題大很多，從而增加了蟻群算法在地球物理反演應(yīng)用中的難度。近年來(lái)，有學(xué)者將蟻群算法應(yīng)用于地震反演和地震解釋，如陳雙全等［22］用蟻群算法進(jìn)行地震數(shù)據(jù)AVO反演，嚴(yán)哲等［23］用蟻群算法進(jìn)行地震剖面斷層識(shí)別和層速度反演；而將ACO應(yīng)用于地球物理資料特別是重磁資料反演的實(shí)例還很少。Yuan等［24］在進(jìn)行地球物理數(shù)據(jù)的群智能最優(yōu)化反演時(shí)，采用連續(xù)參數(shù)離散化的策略，將ACO應(yīng)用于地震記錄及重力異常的參數(shù)化反演，在與擬牛頓法、Levenberg－Marquardt（L－M）、模擬退火算法（SA）及遺傳算法（GA）反演結(jié)果的比較中，ACO顯示出更好的并行性、魯棒性及收斂性。筆者改進(jìn)了傳統(tǒng)組合優(yōu)化問(wèn)題中蟻群算法的映射模式，使得蟻群算法更適用于重磁資料的反演；然后根據(jù)兩條1∶2 000高精度重力剖面數(shù)據(jù)，用蟻群算法探測(cè)千里山巖體和九嶷山巖體的延深、賦存形態(tài)，以期為研究南嶺地區(qū)花崗巖侵入模式和多金屬成礦規(guī)律提供參考。

1 重磁反演蟻群優(yōu)化算法

1.1 蟻群行為

螞蟻屬于群居昆蟲，個(gè)體行為極其簡(jiǎn)單，而群體行為卻相當(dāng)復(fù)雜。相互協(xié)作的一群螞蟻很容易找到從蟻巢到食物源的最短路徑，而單個(gè)螞蟻則不能。人們通過(guò)大量的研究發(fā)現(xiàn)，螞蟻個(gè)體之間是通過(guò)在其所經(jīng)過(guò)的路上留下一種稱之為“信息素”的物質(zhì)來(lái)進(jìn)行信息傳遞的［8－9］。隨后的螞蟻遇到信息素時(shí)，不僅能檢測(cè)出該物質(zhì)的存在及多少，而且可根據(jù)信息素的濃度來(lái)選擇前進(jìn)方向。同時(shí)，該物質(zhì)隨著時(shí)間的推移會(huì)逐漸揮發(fā)，于是路徑的長(zhǎng)短及該路徑上通過(guò)螞蟻的多少就對(duì)殘留信息素的強(qiáng)度產(chǎn)生不同的影響，反過(guò)來(lái)信息素的強(qiáng)弱又指導(dǎo)著其他螞蟻的行動(dòng)方向。因此，某一路徑上走過(guò)的螞蟻越多，則后來(lái)者選擇該路徑的概率就越大。這就構(gòu)成了螞蟻群體行為表現(xiàn)出的一種信息正反饋現(xiàn)象，螞蟻個(gè)體之間就是通過(guò)這種信息交流搜索最優(yōu)路徑的。

1.2 重磁最優(yōu)化反演目標(biāo)函數(shù)

重磁資料最優(yōu)化反演的目標(biāo)函數(shù)φ由數(shù)據(jù)約束φd和模型約束φm兩部分組成［25－26］，即

式中：λ為正則化參數(shù)；d＝（d1，d2，…，dm）為觀測(cè)數(shù)據(jù)矢量；m＝（m1，m2，…，mn）為模型參數(shù)矢量，m和n分別為觀測(cè)數(shù)據(jù)個(gè)數(shù)和模型參數(shù)個(gè)數(shù)；G為敏感度矩陣；為觀測(cè)數(shù)據(jù)權(quán)重矩陣。對(duì)觀測(cè)數(shù)據(jù)進(jìn)行加權(quán)，若觀測(cè)數(shù)據(jù)含均值為0的高斯白噪聲，且相互獨(dú)立［25－26］，則

式中：σ為觀測(cè)數(shù)據(jù)噪聲的標(biāo)準(zhǔn)差；I為單位矩陣。由于重磁異常幅值隨測(cè)點(diǎn)到場(chǎng)源的距離增大而急劇衰減，若不對(duì)模型進(jìn)行加權(quán)約束，磁性或密度的分布將集中分布于地表，出現(xiàn)“趨膚效應(yīng)”。為抵消重磁異常急劇衰減的影響，需對(duì)模型進(jìn)行深度加權(quán)［25－26］，即

式中：z為單元體深度；z0是與深度有關(guān)的常數(shù)；β為與重磁異常隨深度的衰減速率有關(guān)的常數(shù)。Li和 Oldenburg［25－26］進(jìn)行3D重力數(shù)據(jù)反演時(shí)，設(shè)β＝2.0；進(jìn)行3D磁力數(shù)據(jù)反演時(shí)，設(shè)β＝3.0。

1.3 重磁反演蟻群優(yōu)化算法

蟻群最優(yōu)化算法首先應(yīng)用于求解組合優(yōu)化問(wèn)題（COPs），如旅行商問(wèn)題（TSP）等。然而，對(duì)于重磁資料反演則是多維的連續(xù)的目標(biāo)函數(shù)優(yōu)化問(wèn)題，解空間維數(shù)比COPs大得多，增加了最優(yōu)化問(wèn)題的難度。目標(biāo)函數(shù)優(yōu)化問(wèn)題與TSP的相互對(duì)比如表1所示。

表1 重磁反演中目標(biāo)函數(shù)優(yōu)化問(wèn)題與TSP優(yōu)化問(wèn)題的對(duì)比Table 1 Comparison of optimization to potential field data inversion and TSP

對(duì)式（1）目標(biāo)函數(shù)進(jìn)行蟻群算法最優(yōu)化時(shí)，首先需將模型參數(shù)離散化［27－29］。在模型參數(shù)mi的取值范圍內(nèi)，將其劃分為Ni個(gè)節(jié)點(diǎn)（i＝1，2，…，n），將每個(gè)模型參數(shù)看成1層，共計(jì)n層，個(gè)節(jié)點(diǎn)，代表TSP中的個(gè)“城市”（圖1）。劃分節(jié)點(diǎn)時(shí)，可以不必進(jìn)行等分，且每一層的劃分?jǐn)?shù)量和參數(shù)范圍也可以不一樣。但是，為了使最優(yōu)化問(wèn)題簡(jiǎn)化，通常在每個(gè)模型參數(shù)的取值范圍內(nèi)進(jìn)行N等分點(diǎn)，且每個(gè)模型參數(shù)的等分?jǐn)?shù)量相同，因此共計(jì)n層n×N個(gè)節(jié)點(diǎn)。N的數(shù)量決定了最優(yōu)化反演的精度：N越大，反演精度越高，但反演的計(jì)算量就越大；反之，N越小，反演的精度就越低，反演計(jì)算的計(jì)算量就越小。

初始時(shí)刻，將M只螞蟻隨機(jī)地放在m1的N個(gè)節(jié)點(diǎn)上。搜索開始后，螞蟻按照轉(zhuǎn)移概率

從（i－1，j）節(jié)點(diǎn)轉(zhuǎn)移至（i，k）節(jié)點(diǎn)（2≤i≤n，1≤j≤N，1≤k≤N）。其中：τ（i，k）（t）為螞蟻在t時(shí)刻節(jié)點(diǎn)（i，k）的信息素殘留總量；η（i，k）（t）為螞蟻在t時(shí)刻節(jié)點(diǎn)（i，k）的啟發(fā)信息；w1、w2分別表示信息素和啟發(fā)信息的權(quán)重。每只螞蟻從第1層m1開始，層層轉(zhuǎn)移，一直達(dá)到最后一層mn。則第s只螞蟻經(jīng)過(guò)的路徑就對(duì)應(yīng)一個(gè)解ms，可計(jì)算出對(duì)應(yīng)的目標(biāo)函數(shù)值φs。所有螞蟻完成遍歷后，更新節(jié)點(diǎn)（i，j）的信息素為

其中：ρ∈（0，1），稱為揮發(fā)度；▽?duì)樱╥，j）（t）為所有螞蟻在第t次遍歷在節(jié)點(diǎn)（i，j）殘留的信息素。揮發(fā)度避免了系統(tǒng)信息素的無(wú)節(jié)制積累，同時(shí)使得蟻群系統(tǒng)有機(jī)會(huì)弱化舊的信息，突出新的信息。單個(gè)節(jié)點(diǎn)的信息數(shù)殘留量為所有螞蟻在該節(jié)點(diǎn)殘留的信息素之和，即

螞蟻信息素殘留量的取值極為重要，影響到反演的收斂速度和優(yōu)化能力。在TSP最優(yōu)化問(wèn)題中，Dorigo等［8－9，30－31］提出 3 種模型，分別是 ant－cycle，ant－density和ant－quantity模型。其中，應(yīng)用最廣泛的ant－cycle模型定義為

其中：Q是一個(gè)常量，用來(lái)表示螞蟻完成一次完整的路徑搜索后釋放的信息素總量，該值對(duì)搜索結(jié)果沒(méi)有影響；Ls是第s只螞蟻的總花費(fèi)，它等于第s只螞蟻經(jīng)過(guò)各段路徑上所需花費(fèi)的總和，即Ls＝式中，li為第i段線路行走單位里程的費(fèi)用為第s只螞蟻經(jīng)過(guò)的第i段線路的路程。antcycle說(shuō)明，如果螞蟻的路徑總花費(fèi)越高，那么其在單位路徑上所釋放的信息素濃度就越低，反之，殘留的信息數(shù)濃度就越高。

蟻群系統(tǒng)的搜索機(jī)理說(shuō)明，如果模型ms對(duì)應(yīng)的目標(biāo)函數(shù)φs越小，那么螞蟻在經(jīng)過(guò)路徑上殘留的信息素就越多，從而指導(dǎo)其他螞蟻的搜索。當(dāng)蟻群系統(tǒng)完成一次遍歷后，更新節(jié)點(diǎn)信息素，開始下一次搜索，蟻群系統(tǒng)重新回到第一層，開始新的搜索，直到超過(guò)εs的螞蟻遍歷在同一條路徑為止（0＜εs≤1，為收斂系數(shù)）。最終，蟻群系統(tǒng)將趨于同一條最短的路徑，對(duì)應(yīng)的最優(yōu)解為m＊。

當(dāng)蟻群系統(tǒng)收斂至同一路徑以后，搜索基本停滯。若此時(shí)反演結(jié)果不能滿足要求，則可以逐步縮小剖分范圍，在解的臨域

內(nèi)進(jìn)行重新剖分，重建解空間，并進(jìn)行新一次的蟻群搜索，直到觀測(cè)數(shù)據(jù)擬合誤差小于ε為止（ε為收斂的擬合誤差）。δ為臨域半徑，

式中：為模型參數(shù)mi的最優(yōu)解；，分別為模型參數(shù)mi的最小值和最大值；κ為與縮小剖分范圍速度的系數(shù)。蟻群算法的流程圖如圖2所示。

2 重磁反演蟻群優(yōu)化算法改進(jìn)

本次對(duì)蟻群算法的改進(jìn)在于目標(biāo)函數(shù)與信息素之間的映射函數(shù)。盡管式（7）ant－cycle模型成功應(yīng)用于 TSP、QAP和JSP等組合優(yōu)化問(wèn)題［8－9，30］，但是它在重磁資料反演中的應(yīng)用效果還不好。重磁反演問(wèn)題與TSP不同，重磁反演的解空間是通過(guò)連續(xù)域參數(shù)離散化得到的，而TSP本身就是由各個(gè)城市組成的離散化解空間；因此，TSP的解空間比重磁反演問(wèn)題的解空間小得多，前者蟻群系統(tǒng)更容易搜索到最優(yōu)解。在TSP中應(yīng)用的ant－cycle模型尋優(yōu)能力不強(qiáng)，沒(méi)有實(shí)時(shí)更好地突出螞蟻個(gè)體（即各個(gè)解）之間的差異，蟻群系統(tǒng)容易搜索到一個(gè)局部最優(yōu)解，導(dǎo)致搜索過(guò)程過(guò)早停滯。由于蟻群系統(tǒng)搜索停滯，陷入局部極值，即使增加迭代次數(shù)或改變反演參數(shù)，也難從實(shí)質(zhì)上提高收斂的速度和反演的精度。

圖1 蟻群算法最優(yōu)化目標(biāo)函數(shù)示意圖Fig.1 Diagram of ACO optimizing objective function

筆者認(rèn)為，當(dāng)蟻群系統(tǒng)完成一次遍歷后，有必要對(duì)此次蟻群系統(tǒng)的搜索結(jié)果進(jìn)行統(tǒng)計(jì)和評(píng)估。根據(jù)本次蟻群系統(tǒng)的搜索結(jié)果，合理計(jì)算信息素殘留量，進(jìn)而引導(dǎo)蟻群的下一步搜索；這樣勢(shì)必會(huì)增加算法的穩(wěn)定性和收斂速度等。對(duì)蟻群系統(tǒng)搜素結(jié)果的統(tǒng)計(jì)，最主要是計(jì)算每只螞蟻對(duì)應(yīng)的目標(biāo)函數(shù)值，所以對(duì)蟻群系統(tǒng)的統(tǒng)計(jì)就是對(duì)所有螞蟻對(duì)應(yīng)的目標(biāo)函數(shù)值的統(tǒng)計(jì)。因此，筆者提出新的映射模型。新的映射模型為Gauss函數(shù)形式，稱為Gauss模型，即

式中：μ（φ）和σ（φ）是整個(gè)蟻群系統(tǒng)目標(biāo)函數(shù)的均值和方差；A為與信息素總量有關(guān)的常數(shù)，該常數(shù)對(duì)搜索結(jié)果沒(méi)有影響。該式對(duì)每次搜索的結(jié)果進(jìn)行數(shù)理統(tǒng)計(jì)，根據(jù)統(tǒng)計(jì)的目標(biāo)函數(shù)分布，計(jì)算對(duì)應(yīng)的信息素增量，突出了螞蟻個(gè)體之間的差異，有利于加快收斂的速度，尋找全局最優(yōu)解。

圖2 蟻群算法最優(yōu)化目標(biāo)函數(shù)流程圖Fig.2 Flowchart of ACO optimizing objective function

圖3是二維矩形截面棱柱體模型磁異常的Gauss模型和ant－cycle模型蟻群算法反演結(jié)果。兩者螞蟻數(shù)量、節(jié)點(diǎn)層數(shù)、揮發(fā)系數(shù)等反演參數(shù)均是相同的。反演時(shí)，將1 000m×500m的介質(zhì)范圍劃分為40×20＝800個(gè)網(wǎng)格單元，即模型參數(shù)的個(gè)數(shù)為800。將模型參數(shù)剖分為0和1兩等分（0代表無(wú)磁性，1代表有磁性），因此共有2 800種可能的解。兩者的螞蟻數(shù)量為200，揮發(fā)系數(shù)為0.7；正則化參數(shù)（λ）也相同，均為103。由于Gauss模型和antcycle模型反演的初始模型是隨機(jī)生成的，如果螞蟻數(shù)量、等分?jǐn)?shù)量相同，則隨機(jī)生成的初始模型也是相同的。雖然兩者控制信息素總量的常數(shù)不同（A＝1，Q＝107），但它們只控制蟻群系統(tǒng)的信息素總量，不會(huì)影響轉(zhuǎn)移概率，所以對(duì)反演的結(jié)果不會(huì)產(chǎn)生影響。由圖3可知，Gauss模型的反演結(jié)果比antcycle模型反演結(jié)果要好。前者磁性單元的聚集程度高，分辨率高，形狀、邊界與理論模型吻合很好；后者磁性單元分布較散，反演的磁性體分布和理論模型差別較大。

圖4和表2是Gauss模型和ant－cycle模型反演的收斂過(guò)程對(duì)比。由圖4、表2可知，Gauss模型的反演效率比ant－cycle模型高。Gauss模型目標(biāo)函數(shù)下降快，收斂穩(wěn)定；而ant－cycle模型目標(biāo)函數(shù)下降緩慢，收斂速度慢。前者迭代97次收斂，耗時(shí)314s，擬合觀測(cè)數(shù)據(jù)的相對(duì)誤差達(dá)到2.0%；而后者迭代1 000次仍沒(méi)收斂，耗時(shí)4 052s，擬合觀測(cè)數(shù)據(jù)的相對(duì)誤差為18.0%，反演的效率比Gauss模型低很多。這充分說(shuō)明改進(jìn)的Gauss模型比傳統(tǒng)的antcycle模型效果好，這說(shuō)明了算法改善的有效性。

3 理論試算

筆者設(shè)計(jì)兩組磁性體模型來(lái)檢驗(yàn)ACO算法的效果，分別是單個(gè)二度矩形截面棱柱體模型和組合二度板狀體組成的向斜模型。兩組磁性體單元的磁化強(qiáng)度M＝100A/m，磁化傾角I＝45°，測(cè)線方位角為0°。地面觀測(cè)資料不含噪聲，剖面長(zhǎng)度為1 000 m，點(diǎn)距為20m，觀測(cè)點(diǎn)數(shù)為51。矩形截面棱柱體模型ACO反演的正則化參數(shù)為103，向斜模型ACO反演的正則化參數(shù)為5×102。

3.1 單個(gè)模型

圖3 Gauss模型和ant－cycle模型反演結(jié)果對(duì)比Fig.3 Comparisons of inversion results using Gauss system and ant－cycle system

圖5是矩形截面棱柱體模型ACO反演的迭代收斂過(guò)程和反演結(jié)果。初始時(shí)刻，剖面的磁性分布是隨機(jī)的、任意的（圖5a）；隨著蟻群系統(tǒng)的搜索，磁性單元的分布逐漸收攏，并與真實(shí)模型逐漸接近（圖5b、c、d）；當(dāng)經(jīng)過(guò)97次搜索后，收斂停止，反演終止（圖5e）。反演結(jié)果與真實(shí)模型非常接近，反演效果很好。

圖4 Gauss模型和ant－cycl模型的收斂過(guò)程對(duì)比Fig.4 Comparisons of convergence processes using Gauss system and ant－cycle system

圖5 二維矩形截面棱柱體模型蟻群算法反演結(jié)果Fig.5 ACO inversion results for synthetic two－dimensional rectangular prism model

表2 Gauss模型和ant－cycle模型反演結(jié)果對(duì)比Table 2 Comparison of inversion results of Gauss system and ant－cycle system

3.2 組合模型

圖6是板狀體組成的向斜模型ACO反演的迭代收斂過(guò)程和反演結(jié)果。同單個(gè)矩形截面棱柱體模型的反演類似，初始時(shí)刻，剖面的磁性分布是隨機(jī)的、任意的（圖6a）；隨著蟻群系統(tǒng)的搜索，磁性單元的分布逐漸收攏，并與真實(shí)模型逐漸接近（圖6b、c、d）；當(dāng)經(jīng)過(guò)191次搜索后，收斂停止，反演終止（圖6e）。反演結(jié)果與真實(shí)模型非常接近，反演效果很好。

4 南嶺花崗巖侵入體構(gòu)造探測(cè)

4.1 地質(zhì)與地球物理概況

圖6 二維向斜模型體模型蟻群算法反演結(jié)果Fig.6 ACO inversion results for synthetic two－dimensional syncline model

南嶺地區(qū)在構(gòu)造上橫跨華夏板塊與揚(yáng)子板塊，位于華南陸塊構(gòu)造巖漿活動(dòng)帶，歷經(jīng)多期次構(gòu)造－巖漿－成礦作用，區(qū)內(nèi)廣泛發(fā)育著各種花崗巖和其他巖漿巖，是我國(guó)重要的金屬礦成礦帶之一，鎢、錫、銅、鉛、鋅等許多內(nèi)生金屬礦床與花崗巖有密切關(guān)系?；◢弾r密度較低，侵入到沉積地層中表現(xiàn)為負(fù)異常特征。根據(jù)前人所做的物性統(tǒng)計(jì)工作［7］得知，區(qū)內(nèi)巖石密度分布有如下特征：1）碳酸鹽巖密度大于砂頁(yè)巖密度，灰?guī)r平均密度為2.70g/cm3，頁(yè)巖平均密度為2.53g/cm3；大理巖化灰?guī)r、白云巖密度值為2.73～2.75g/cm3；矽卡巖平均值密度為3.27 g/cm3。2）淺源重熔型花崗巖密度常見值為2.6g/cm3，比一般圍巖低，剩余密度為－0.05～－0.10g/cm3；而深源同熔型花崗巖其密度常見值為2.72g/cm3，與圍巖無(wú)明顯的密度差。根據(jù)巖石密度資料統(tǒng)計(jì)，研究區(qū)花崗巖的平均密度為2.64g/cm3，古生界到中生界平均密度為2.69～2.73g/cm3，因此取花崗巖剩余密度為－0.05～－0.09g/cm3。

4.2 千里山巖體

千里山巖體位于茶陵－郴州逆沖地殼斷裂帶和大義山－郴州－大寶山斷裂帶交匯地帶，該巖體礦集區(qū)是南嶺花崗巖有色金屬成礦最為典型的區(qū)域之一，其周圍形成的柿竹園、紅旗嶺、東坡以及金船塘等著名的大型－超大型多金屬礦床［32］（圖7），歷來(lái)是地質(zhì)工作者研究的重點(diǎn)。在1∶20萬(wàn)重力異常圖上，千里山巖體呈顯著負(fù)異常特征，并與香花嶺、騎田嶺及諸廣山巖體負(fù)異常相連，平行茶陵－郴州梯級(jí)帶形成北東向負(fù)異常帶，可能與深部隱伏的巖漿遷移通道有關(guān)。采集了千里山1∶2 000高精度重力剖面數(shù)據(jù)，除去淺表干擾，異常曲線如圖8a所示。剖面南北向，布格重力異常以出露巖體為中心呈明顯負(fù)異常特征，幅值達(dá)－26×10－5m/s2，異常形態(tài)南陡北緩。

對(duì)千里山巖體重力剖面進(jìn)行蟻群算法反演時(shí)，將剖面剖分為40×20個(gè)矩形截面棱柱體單元，蟻群系統(tǒng)的數(shù)量為200，信息素?fù)]發(fā)系數(shù)為0.7，根據(jù)該區(qū)物性統(tǒng)計(jì)結(jié)果取該區(qū)花崗巖剩余密度為－0.09 g/cm3，蟻群系統(tǒng)搜索30次后收斂。圖8b是經(jīng)過(guò)千里山巖體重力剖面的蟻群算法探測(cè)結(jié)果。首先，蟻群算法反演快速、收斂穩(wěn)定，較好地?cái)M合了重力觀測(cè)數(shù)據(jù)（圖8a）。蟻群算法圈定的花崗巖巖體邊界清晰，分布集中，揭示千里山花崗巖巖體深部規(guī)模較大，下底面延深達(dá)到16km。這說(shuō)明出露巖體僅為該巖基的頂角，尤其在出露巖體的北部、毗鄰茶陵－郴州斷裂帶附近有一定的隱伏規(guī)模，推測(cè)為深部巖漿遷移通道，并形成千里山－騎田嶺－香花嶺構(gòu)造巖漿帶。

圖7 千里山巖體區(qū)域地質(zhì)圖［32］及重力測(cè)點(diǎn)位置圖Fig.7 Regional geological map of Qianlishan rocks and positions of gravity observation points

4.3 九嶷山巖體

九嶷山巖體位于我國(guó)著名的有色金屬成礦區(qū)——湘南地區(qū)西南部、湘粵交界處，是南嶺地區(qū)重要的錫成礦遠(yuǎn)景區(qū)之一。區(qū)內(nèi)地質(zhì)構(gòu)造復(fù)雜，斷裂構(gòu)造發(fā)育，巖漿活動(dòng)頻繁，成礦地質(zhì)條件優(yōu)越［32］（圖9）。在1∶20萬(wàn)重力異常圖上，九嶷山巖體周圍呈顯著負(fù)異常，并且負(fù)異常在北部超出出露巖體。采集了九嶷山1∶2 000高精度重力剖面數(shù)據(jù)（圖10a），剖面南北向，布格重力異常在出露巖體區(qū)域呈明顯負(fù)異常特征，幅值達(dá)－24×10－5m/s2，且負(fù)異常在出露巖體北部有一定規(guī)模延伸。

對(duì)九嶷山巖體重力剖面進(jìn)行蟻群算法反演時(shí)，將剖面剖分為40×20個(gè)矩形截面棱柱體單元，蟻群系統(tǒng)的數(shù)量為200，信息素?fù)]發(fā)系數(shù)為0.7，根據(jù)該區(qū)物性統(tǒng)計(jì)結(jié)果取該區(qū)花崗巖剩余密度為－0.09 g/cm3，蟻群系統(tǒng)搜索31次后收斂。圖10b是經(jīng)過(guò)九嶷山巖體重力剖面的蟻群算法探測(cè)結(jié)果。首先，蟻群算法反演快速、收斂穩(wěn)定，較好地?cái)M合了重力觀測(cè)數(shù)據(jù)（圖10a）。蟻群算法圈定的巖體邊界清晰，分布集中，且與區(qū)內(nèi)出露巖體吻合。同時(shí)，大地電磁是揭示花崗巖體形態(tài)的重要手段，本次蟻群算法反演的結(jié)果與該剖面大地電磁測(cè)深反演結(jié)果（圖11）相比，所得到的巖體侵入形態(tài)接近。蟻群算法和大地電磁測(cè)深揭示：該剖面上花崗巖侵入體為3個(gè)巖體的組合，從南往北，該剖面上花崗巖年代分別為侏羅紀(jì)、三疊紀(jì)及侏羅紀(jì)；最大延深達(dá)到22km，對(duì)應(yīng)三疊紀(jì)巖體；而北部為隱伏巖體。反演結(jié)果揭示九嶷山復(fù)式巖體深部侵位形態(tài)較復(fù)雜。

圖8 千里山巖體重力異常（Δg）蟻群算法反演結(jié)果Fig.8 ACO inversion results of Qianlishan rocks of gravity anomalies（Δg）

圖9 九嶷山巖體區(qū)域地質(zhì)圖［32］及重力和MT測(cè)點(diǎn)位置圖Fig.9 Regional geological map of Jiuyishan rocks and positions of gravity and MT observation points

圖10 九嶷山巖體重力異常蟻群算法反演結(jié)果Fig.10 ACO inversion results of Jiuyishan rocks of gravity anomalies

圖11 九嶷山巖體2D大地電磁測(cè)深反演結(jié)果Fig.11 2DMT inversion results of Jiuyishan rocks

5 結(jié)論

蟻群智能算法在啟發(fā)式隨機(jī)搜索中展示出優(yōu)異的性能，它是基于信息素的正反饋機(jī)制，具有較好的魯棒性、并行性和可移植性。蟻群智能算法是繼蒙特卡洛算法（Monte Carlo）、遺傳算法（GA）、模擬退火（SA）和人工神經(jīng)網(wǎng)絡(luò)（ANN）之后的一種可以廣泛應(yīng)用于地球物理反演的完全非線性反演方法，在地球物理反演中有較好的應(yīng)用前景。

改進(jìn)的蟻群算法反演效果良好，其中九嶷山巖體的反演結(jié)果與大地電磁測(cè)深所揭示的巖體侵入形態(tài)吻合。蟻群系統(tǒng)探測(cè)千里山巖體延深達(dá)到16 km，北部具有一定隱伏規(guī)模，揭示了千里山－騎田嶺－香花嶺隱伏巖漿構(gòu)造帶的存在；九嶷山巖體延深達(dá)到22km，該剖面上分為3個(gè)巖體，北部為隱伏巖體。蟻群系統(tǒng)較好控制了千里山巖體和九嶷山巖體的延深、形態(tài)，為研究南嶺地區(qū)花崗巖侵入模式和多金屬成礦規(guī)律提供了參考。

（References）：

［1］秦葆瑚.湘南四個(gè)金屬礦田的區(qū)域重磁異常特征［J］.物探與化探，1984，8（1）：34－40.Qin Baohu.Characteristics of Regional Gravity and Magnetic Anomalies of Four Major Polymetal Ore Fields in Southern Hunan，China［J］.Geophysical ＆Geochemical Exploration，1984，8（1）：34－40.

［2］秦葆瑚.湘南區(qū)域重磁異常的地質(zhì)解釋及其在成礦預(yù)測(cè)中的應(yīng)用［J］.湖南地質(zhì)，1984，3（2）：1－14.Qin Baohu.The Interpretations of Regional Gravity and Magnetic Anomalies and Their Applications of Minerogenetic Prognoses in Southern Hunan［J］.Hunan Geology，1984，3（2）：1－14.

［3］秦葆瑚.南嶺區(qū)域重磁異常的地質(zhì)解釋［J］.湖南地質(zhì)，1987，1（1）：1－15.Qin Baohu.A Geological Interpretation on the Regional Gravity and Magnetic Anomalies in Nanling Area［J］.Hunan Geology，1987，1（1）：1－15.

［4］李少蓮.湘南地區(qū)隱伏地質(zhì)構(gòu)造與礦產(chǎn)預(yù)測(cè)［J］.礦產(chǎn)與地質(zhì)，1990，4（3）：11－17.Li Shaolian.Concealed Geological Structure and Prediction for Mineral Resource in South Hunan［J］.Mineral Resources and Geology，1990，4（3）：11－17.

［5］段利華.湘桂粵贛區(qū)域性位場(chǎng)資料的數(shù)據(jù)處理及其應(yīng)用［J］.礦產(chǎn)與地質(zhì)，1991，5（4）：318－324.Duan Lihua.Regional Potential Field Data Processing and Their Application in the Area of Hunan，Guangxi，Guangdong and Jiangxi Province［J］.Mineral Resources and Geology，1991，5（4）：318－324.

［6］朱自強(qiáng)，黃國(guó)祥.柿竹園礦田重磁資料的三維反演及成礦預(yù)測(cè)［J］.礦產(chǎn)與地質(zhì)，1996，10（1）：66－72.Zhu Ziqiang，Huang Guoxiang.3DInversion for Potential Field Data and Metallogenic Prognosis in Shizhuyuan Deposits［J］.Mineral Resources and Geology，1996，10（1）：66－72.

［7］饒家榮，金小燕，曾春芳.南嶺中段北段深部構(gòu)造：巖漿巖控礦規(guī)律及找礦方向［J］.國(guó)土資源導(dǎo)刊，2006，3（3）：31－36.Rao Jiarong，Jin Xiaoyan，Zeng Chunfang.Deep Structure of Middle and Northern Section：Magmatic Rocks Ore－Controlling Regularity and Prospecting Direction［J］.Land ＆Resources Herald，2006，3（3）：31－36.

［8］Dorigo M，Maniezzo V，Colorni A.The Ant System：An Autocatalytic Optimizing Process［R］.Milano：Politecnico di Milano，1991.

［9］Dorigo M，Maniezzo V，Colorni A.Ant System：Optimization by a Colony of Cooperating Agents［J］.IEEE Trans on Systems，Man，and Cybernetics：Part B：Cybernetics，1996，26（1）：29－41.

［10］Gambardella L M，Dorigo M.Solving Symmetric and Asymmetric TSPs by Ant Colonies［C］//Proceedings of the 3rd IEEE International Conference on Evolutionary Computation.Nagoya：IEEE Neural Network Council（NNC）and Society of Instrument and Control Engineers（SICE），1996：622－627.

［11］Dorigo M，Gambardella L M.Ant Colony System：A Cooperative Learning Approach to the Traveling Salesman Problem［J］.IEEE Transactions on Evolutionary Computation，1997，1（1）：53－66.

［12］Dorigo M，Gambardella L M.Ant Colonies for the Travelling Salesman Problem［J］.BioSystems，1997，43（2）：73－82.

［13］Gambardella L M，Dorigo M.Ant－Q：A Reinforcement Learning Approach to the Traveling Salesman Problem［C］//Proceedings of the 12th International Conference on Machine Learning.Tahoe：The International Machine Learning Society，1995：252－260.

［14］Stutzle T，Hoos H.Improvements on the Ant System：MAX－MIN Ant System［C］//Proceedings of the International Conference on Artificial Neural Networks and Genetic Algorithms.Wien：Springer Verlag，1997：245－249.

［15］Stutzle T，Hoos H.MAX－MIN Ant System and Local Search for the Traveling Salesman Problem［C］//Proceedings of the 4th International Conference on Evolutionary Computation.Indianapolis：IEEE Neural Network Council（NNC）and Evolutionary Programming Society（EPS），1997：308－313.

［16］Maniezzo V，Colorni A.The Ant System Applied to the Quadratic Assignment Problem［J］.IEEE Transactions on Knowledge and Data Engineering，1999，11（5）：769－778.

［17］Gambardella L M，Taillard E，Dorigo M.Ant Colonies for the QAP［R］.Lugano：Istituto Dalle Molle di Studi Sull’Intelligenza Artificiale（IDSIA），1997.

［18］Gambardella L M，Taillard E，Dorigo M.Ant Colonies for the Quadratic Assignment Problem［J］.Journal of the Operational Research Society，1999，50（2）：167－176.

［19］Stützle T，Hoos H.MAX－MIN Ant System and Local Search for Combinatorial Optimization Problems［C］//Stefan V，Ibrahim H O，Catherine R.Meta－Heuristics：Advances and Trends in Local Search Paradigms for Optimization.Norwell：Kluwer Academic Publishers，1997：137－154.

［20］Maniezzo V.Exact and Approximate Nondeterministic Tree－Search Procedures for the Quadratic Assignment Problem［J］.Informs Journal on Computing，1999，11（4）：358－369.

［21］Colorni A，Dorigo M，Maniezzo V，et al.Ant System for Job－Shop Scheduling［J］.Belgian Journal of Operations Research，Statistics and Computer Science，1994，34（1）：39－53.

［22］陳雙全，王尚旭，季敏，等.地震波阻抗反演的蟻群算法實(shí)現(xiàn)［J］.石油物探，2005，44（6）：551－553.Chen Shuangquan，Wang Shangxu，Ji Min，et al，The Ant Colony Algorithm for the Seismic Impedance Inversion［J］.Geophysical Prospecting for Petroleum，2005，44（6）：551－553.

［23］嚴(yán)哲，顧漢明，趙小鵬.基于蟻群算法的非線性AVO反演［J］.石油地球物理勘探，2009，44（6）：700－702.Yan Zhe，Gu Hanming，Zhao Xiaopeng.Non－Linear AVO Inversion Based on Ant Colony Algorithm［J］.Oil Geophysical Prospecting，2009，44（6）：700－702.

［24］Yuan S Y，Wang S X，Tian N.Swarm Intelligence Optimization and Its Application in Geophysical Data Inversion［J］.Applied Geophysics，2009，6（2）：166－174.

［25］Li Y G，Oldenburg D W.3－D Inversion of Magnetic Data［J］.Geophysics，1996，61（2）：394－408.

［26］Li Y G，Oldenburg D W.3－D Inversion of Gravity Data［J］.Geophysics，1998，63（1）：109－119.

［27］Dréo J，Siarry P.A New Ant Colony Algorithm Using the Heterarchical Concept Aimed at Optimization of Multiminima Continuous Functions［C］//Dorigo M，Caro G D，Sampels M.Ant Algorithms.Berlin：Springer Berlin Heidelberg，2002：216－221.

［28］Dréo J，Siarry P.Continuous Interacting Ant Colony Algorithm Based on Dense Heterarchy［J］.Future Generation Computer Systems，2004，20（5）：841－856.

［29］Socha K，Dorigo M.Ant Colony Optimization for Continuous Domains［J］.European Journal of Operational Research，2008，185（3）：1155－1173.

［30］Dorigo M.Optimization，Learning and Natural Algorithms［D］.Milano：Politecnico di Milano，1992.

［31］Dorigo M，Stützle T.The Ant Colony Optimization Metaheuristic： Algorithms， Applications， and Advances［C］//Glover F， Kochenberger G A.Handbook of Metaheuristics.Norwell： Kluwer Academic Publishers，2003：250－285.

［32］於崇文，彭年.南嶺地區(qū)區(qū)域成礦分帶性：復(fù)雜成礦系統(tǒng)中的時(shí)－空同步化［M］.北京：地質(zhì)出版社，2009.Yu Chongwen，Peng Nian.Regional Ore Zonality of the Nanling Region：Spatio－Temporal Synchronization in Complex Ore－Forming Systems［M］.Beijing：Geological Publishing House，2009.