賈虎，鄧力琿

（油氣藏地質(zhì)與開發(fā)工程國家重點實驗室 西南石油大學(xué)，成都 610500）

0 引言

中國陸相沉積的碎屑巖和海相沉積的碳酸鹽巖油藏，儲集層非均質(zhì)性較強，且目前許多高含水砂巖油藏歷經(jīng)長期水驅(qū)開發(fā)后剩余油分布雜亂且分散，難以有效認(rèn)識水驅(qū)油藏動用規(guī)律，導(dǎo)致注水開發(fā)效率低下。同時地層在經(jīng)歷長時間流體沖刷后易形成優(yōu)勢通道，導(dǎo)致無效水循環(huán)[1]，降低水驅(qū)效率。為了更好描述水驅(qū)油藏特征，前人通過水驅(qū)特征曲線[2]及井間連通性模型[3]對水驅(qū)流場進(jìn)行了表征。同時，侯玉培等[1,4-5]在流場識別方面做了大量研究，建立了流場評價體系，采取層次分析法對流場進(jìn)行了評價，但層次分析法需要人為設(shè)置評價權(quán)重，因此方案決策存在不確定性。在地質(zhì)模型基礎(chǔ)上，前人通過流線模擬對生產(chǎn)數(shù)據(jù)進(jìn)行擬合并分析了流體運移規(guī)律及流動軌跡[6-10]。與根據(jù)孔滲飽分布確定流場強弱的常規(guī)方法相比，流線模擬可更準(zhǔn)確且直觀地顯示出流體主要流通區(qū)域。此外，可依據(jù)注水對采油的貢獻(xiàn)量對水驅(qū)油藏進(jìn)行注水優(yōu)化、水驅(qū)控制、井網(wǎng)優(yōu)化[11-18]。主流的油藏數(shù)值模擬軟件如Petrel RE雖可準(zhǔn)確計算流線分布并依據(jù)上述方法進(jìn)行方案優(yōu)化，但缺乏量化分析流場的方法，且離散的流線分布無法代表實際流場，也無法區(qū)分水相驅(qū)動能力不同的流動區(qū)域。通過對商業(yè)軟件Petrel RE進(jìn)行二次開發(fā)實現(xiàn)了流線聚類功能，可對水相驅(qū)動能力不同的流線進(jìn)行區(qū)分并加以分析，為水驅(qū)流場識別提供了一種新的研究手段。

聚類方法又被稱為無監(jiān)督學(xué)習(xí)，是人工智能方法的分支之一。不同于通過注采井間配位對流線進(jìn)行歸類的分析方式，本文采用聚類算法[19]對流線依據(jù)其空間位置及流線性質(zhì)進(jìn)行聚類，提取出流場潛在的分布結(jié)構(gòu)。聚類算法以最大化類與類間的差異并最小化類內(nèi)差異為原則對流線進(jìn)行聚類，能確保每一類流線均可最大程度表征其所在的流場，同時區(qū)分出不同類型的流線，該類無監(jiān)督學(xué)習(xí)算法對經(jīng)驗公式與人為評價依賴度較小，因此對不同類型油藏均具有較好的適用性。

本文以某碳酸鹽巖油藏為例，通過流線聚類方法，對水相驅(qū)動能力不同的流線進(jìn)行區(qū)分并進(jìn)一步對同一注采井間流線進(jìn)行細(xì)分，開展流場精細(xì)化描述，找出潛在優(yōu)勢流場，為水驅(qū)流場調(diào)整與提高采收率方案決策提供科學(xué)依據(jù)與技術(shù)支撐。

1 流線聚類

流線模擬以流線表示流體長期沖刷而形成的優(yōu)勢通道，模擬地層流體的流動軌跡。本文通過流線模擬方法表征流場，以油藏流體在一定時間步長內(nèi)均沿流線方向流動為基本假設(shè)，流線軌跡通過流線追蹤算法[20]獲取，在獲取流線軌跡后通過求解沿流線方向的一維質(zhì)量守恒方程，即可得到流體沿流線的飽和度分布，并顯示油藏流體流動規(guī)律。

流線模擬結(jié)果表征了流體性質(zhì)在離散流線軌跡上的分布，因此需對流線信息加以提取才可得到流場分布的具體描述。本文通過流線聚類的方法將流場分為不同流動區(qū)域，識別出不同類型流場分布以及該類流場特性，通過這些信息便可針對性地進(jìn)行后續(xù)油藏開發(fā)方案設(shè)計。流線聚類步驟包括特征參數(shù)提取，聚類分析及聚類評價。

1.1 特征參數(shù)提取

聚類算法需依據(jù)樣本特征對其進(jìn)行聚類，因此需提取出具有區(qū)分度且便于實際應(yīng)用的流線特征，本文提取的流線特征及公式見表1。其中流線油水體積比由Vow表示，若該值較大，則代表沿該流線方向未被注入水波及到的油較多。流線油水流速比由參數(shù)vow表示，該值越大，則表示沿該流線方向水驅(qū)油能力越強。

表1 流線特征及公式

因為流線在z軸方向跨度不大，故選擇流線中點x，y軸坐標(biāo)表征流線位置。同時，由于流線模擬計算的數(shù)值耗散問題，部分流線可能出現(xiàn)油水兩相流速均接近0的情況，因此無法通過流體流速直接表征流線位置的水相驅(qū)動能力。此外，對于水驅(qū)后期的流場，水相滲透率往往較高，而流體飽和度為線性增加的屬性，因此對驅(qū)動能力區(qū)分效果不明顯，故本文采用油水體積比及油水流速比表征流線位置水相驅(qū)動能力，該方式可對水相飽和度較低而油相飽和度較高的流線進(jìn)行細(xì)分，同時對水相飽和度大于一定程度的流線進(jìn)行一定的粗化，從而最大程度反映不同水相驅(qū)動能力的流場分布。

參數(shù)選取的意義在于將流線依據(jù)其空間位置、驅(qū)替效率、波及效率進(jìn)行聚類，將具有類似性質(zhì)的流線歸為一類，從而識別出具有開發(fā)價值的油藏區(qū)域，為后期注水優(yōu)化、井網(wǎng)層系調(diào)整、深部調(diào)剖等方案決策提供科學(xué)依據(jù)和技術(shù)支撐。

在進(jìn)行特征提取后，需對每條流線的特征參數(shù)進(jìn)行歸一化處理以便之后的聚類分析。

1.2 聚類分析

聚類分析是一種用于尋找數(shù)據(jù)之間內(nèi)在結(jié)構(gòu)的技術(shù)。聚類算法有密度峰值、K-means、譜聚類、層次聚類算法等，其中由于密度峰值聚類算法對數(shù)據(jù)結(jié)構(gòu)區(qū)分能力較強且其結(jié)果不具有隨機性，因此本文采用此法。密度峰值聚類算法[19]的基本原理為，依據(jù)與樣本性質(zhì)類似的樣本數(shù)量的多少作為樣本密度，選取局部區(qū)域內(nèi)密度最大點作為聚類中心，同時確保其與其他密度更大的樣本之間的“距離”相對更大。密度峰值聚類算法流程圖見圖1。

圖1 密度峰值聚類算法流程圖

樣本i的局部密度計算公式為[19]：

設(shè)dc可使樣本平均局部密度為總樣本數(shù)的 1%～2%[19]，dij通過歐幾里得距離公式計算。（1）式計算的結(jié)果是離散值，且易影響算法穩(wěn)定性，因此筆者采用指數(shù)函數(shù)的形式計算局部密度：

同時，需保證樣本與其他比該樣本密度更大的樣本有較大距離，設(shè)iqρ為樣本局部密度降序排列的樣本序號，即，則可定義樣本分離距離iqδ為：

從而得到不同樣本的局部密度及分離距離分布（見圖2）。

聚類中心應(yīng)同時具有較大的局部密度與分離距離值，因此圖 2內(nèi)偏離橫、縱坐標(biāo)軸較遠(yuǎn)處的點可作為聚類中心點，第i個樣本的中心點評價公式為：

圖2 樣本局部密度與分離距離分布圖

若要將總體樣本分為nc類，即可選取依據(jù)γ降序排列在前nc的樣本作為聚類中心點。Rodriguez等[21]通過對人工隨機生成數(shù)據(jù)集分析發(fā)現(xiàn)，若以降序排序后的γ值序號作為橫軸，lgγ值作為縱軸，則lgγ曲線從聚類中心過渡到非聚類中心的曲線趨勢有明顯跳躍，可為聚類中心選取提供參考，因此本文取lgγ值排序曲線趨勢出現(xiàn)跳躍的位置作為總樣本類數(shù)，即聚類數(shù)nc。確定聚類中心后，將非聚類中心樣本與距離最近的聚類中心歸為一類。

1.3 聚類評價

聚類數(shù)對聚類效果有著極大影響：聚類數(shù)過小將粗化數(shù)據(jù)，丟失有效信息；聚類數(shù)過大將致使聚類結(jié)果無法有效壓縮，造成分析困難。同時，不同聚類算法如密度峰值、K-means、譜聚類、層次聚類算法等結(jié)果往往不盡相同，因此，有必要對聚類結(jié)果進(jìn)行評價，從而選擇合理的聚類算法及聚類數(shù)。Liu Y等[22]在已知樣本真實類別數(shù)的情況下，分析了不同聚類評價算法如Calinski-Harabasz、指數(shù)、Dunn、輪廓系數(shù)、Davies-Bouldin、Xie-Beni、SD、S-Dbw等是否具有識別出真實樣本聚類數(shù)的能力。輪廓系數(shù)與S-Dbw算法同其他聚類評價算法相比，對不同的樣本類型均能有效識別出樣本實際的類別數(shù)，然而S-Dbw算法對離散點較敏感，結(jié)果易出現(xiàn)不連續(xù)情況。因此本文采用輪廓系數(shù)算法[23]，計算樣本的輪廓系數(shù)，根據(jù)已評價聚類后的樣本與同類樣本間平均距離以及樣本與非同類樣本間平均距離，對聚類效果進(jìn)行評價，輪廓系數(shù)值為-1～1，趨近于1時代表聚類效果相對較優(yōu)。

2 流場表征

圖 3為某碳酸鹽巖油藏內(nèi)孔隙度及滲透率分布，由于溶孔溶洞的存在，部分區(qū)域孔隙度及滲透率較高，為流體提供了較好的滲流條件。圖4為用Petrel RE軟件生成的該油藏注水開發(fā) 6年后通過流線模擬得到的流線分布。圖4a為依據(jù)流線起點進(jìn)行分類的結(jié)果，圖4b為油相飽和度在流線上的分布。在流線分布圖中，流線越密集則流動強度越大，由圖 4流線模擬結(jié)果可看出流體流動區(qū)域主要集中于滲透率、孔隙度均較高，即溶孔溶洞較為發(fā)育的區(qū)域，如 I4井至 P4、P6井間區(qū)域，這些區(qū)域內(nèi)滲流條件較好，因此易形成注入水流動的優(yōu)勢通道。

圖3 某碳酸鹽巖油藏孔隙度、滲透率分布

圖4 通過Petrel RE軟件生成的開發(fā)6年后流線分類（a）及流線油相飽和度分布（b）

2.1 聚類參數(shù)及算法

在獲取了流線模擬結(jié)果后，筆者通過自行開發(fā)基于Ocean平臺的插件實現(xiàn)了與Petrel RE的對接，能導(dǎo)出表示流線油藏流體、水相流體以及油相流體流動速率的屬性數(shù)據(jù)，表示流線油相、水相飽和度分布的屬性數(shù)據(jù)以及流線位置數(shù)據(jù)。

插件將Petrel RE流線數(shù)據(jù)導(dǎo)出至txt或excel文件，之后即可通過具有快速開發(fā)能力的Python編程語言進(jìn)行后續(xù)數(shù)據(jù)處理及聚類分析。依據(jù)表 1公式提取出流線特性后，進(jìn)行密度峰值聚類計算，當(dāng)設(shè)定截斷距離為0.01時樣本平均局部密度為1.405 3%。由于流線聚類結(jié)果對聚類數(shù)極為敏感，不同聚類數(shù)下聚類結(jié)果有較大差異（見圖 5），因此有必要通過合理的方法選取聚類數(shù)。

圖5 不同聚類數(shù)下流線聚類結(jié)果

圖6 不同樣本序號下lgγ分布

密度峰值聚類算法通過計算流線樣本局部密度與距離屬性后依據(jù)二者的乘積進(jìn)行排序，并選擇lgγ出現(xiàn)跳躍處序號作為聚類數(shù)。圖 6中紅色虛線處為跳躍較為明顯位置，表示聚類數(shù)為14，但除該點外，樣本序號為 4時同樣出現(xiàn)非常明顯的跳躍，因此需進(jìn)一步通過聚類評價算法對聚類結(jié)果合理性進(jìn)行評價。本文采用聚類評價能力較好的輪廓系數(shù)聚類評價算法計算不同聚類數(shù)下的輪廓系數(shù)得分，并以該得分對水驅(qū)開發(fā)6年后的流線聚類效果進(jìn)行評價（見圖7）。

圖 7為依據(jù)輪廓系數(shù)算法得出的流線聚類結(jié)果在不同聚類數(shù)下的得分，得分越接近 1代表聚類效果越好。其中紅色圓點表示聚類數(shù)為 14，可觀察到該點處輪廓系數(shù)得分最高，因此合理聚類數(shù)為14。以水驅(qū)開發(fā)6年后油藏流線分布為例，在取聚類數(shù)均為14的情況下，對比密度峰值聚類算法與K-means、層次聚類、譜聚類算法間的聚類結(jié)果。由圖8可見，4種聚類算法計算結(jié)果存在差異。

圖7 不同聚類數(shù)時密度峰值聚類的輪廓系數(shù)得分

圖8 不同聚類算法結(jié)果對比圖

密度峰值、K-means、層次聚類、譜聚類4種聚類算法計算結(jié)果對應(yīng)的輪廓系數(shù)分別為0.851 914，0.853 310，0.850 220，0.816 990。雖然K-means算法輪廓系數(shù)最高，但K-means算法計算結(jié)果具有隨機性，不利于對流場進(jìn)行確定性分析。其次為密度峰值聚類、層次聚類以及譜聚類。層次聚類算法輪廓系數(shù)略低于密度峰值算法，但該算法受異常點影響較大，穩(wěn)定性較差。此外，譜聚類算法的輪廓系數(shù)以及對流線的區(qū)分效果較其他算法均較低，因此本文選擇密度峰值聚類算法作為流線聚類的主要算法。

2.2 聚類結(jié)果

油水流速比與油水體積比可表示流線方向水驅(qū)油能力及流線方向未被注入水波及到的油相體積，因此本文以油水流速比與油水體積比的平均數(shù)作為水驅(qū)評價系數(shù)。該系數(shù)越大，則代表水相驅(qū)替效果越好，按該系數(shù)由大至小將流線聚類結(jié)果從1到14排序，結(jié)果見圖9。與Petrel RE流線模擬結(jié)果（見圖4）相比，本文的流線聚類結(jié)果對流場進(jìn)行了有效區(qū)分，有利于找出具有開發(fā)潛力的區(qū)域，如I3至P6井間的第2類流線；還可識別出影響注水效率的無效注水循環(huán)區(qū)域，如I4至P6井間的第14類流線以及I6至P5井間的第11類流線；同時可對屬于相同注采井間的流線進(jìn)行細(xì)分，將其分為具有不同驅(qū)動能力的邊流與主流：I4至P6井間，第5類流線為邊流，第14類為主流，I4至P4井間，第6、第7類為邊流，第10、第13類為主流。這些現(xiàn)象進(jìn)一步說明流線聚類對不同水相驅(qū)動能力流線具有較好區(qū)分效果。不同類別流線的流場性質(zhì)見表2。

圖9 開發(fā)6年流線聚類圖

表2 不同流線類別的流場性質(zhì)

由表 2可見，不同類流線的油水流速比與油水體積比相差較明顯，說明流線聚類算法對不同類別的流場具有很好的區(qū)分能力。同時，通過該流線聚類分布即可有針對性地分析水相驅(qū)動效果較差區(qū)域，如第10至14類流線區(qū)域，這些區(qū)域往往位于注采井間連線處（見圖9），易發(fā)生水竄，影響注水波及系數(shù)。因此應(yīng)選擇合理方式調(diào)整注水制度，或是對該區(qū)域進(jìn)行深部調(diào)剖，調(diào)整流體流向。此外，流線聚類方法可依據(jù)水驅(qū)評價系數(shù)識別出注采井間不同類型的流線，如I4井至P4及P6井區(qū)域。該區(qū)域因位于溶洞發(fā)育區(qū)域，因此流量較大，流線分布較密集，第5、第6、第7、第10、第13、第14類流線均位于該區(qū)域，其中水相驅(qū)動能力較低的區(qū)域為第10、第13、第14類流線，其水驅(qū)評價系數(shù)分別為 0.155，0.065，0.010，說明該流線區(qū)域內(nèi)主要流通流體為水相，注采井間形成了優(yōu)勢通道并降低了注入水波及效率，而該區(qū)域內(nèi)第5、第6、第 7類流線水驅(qū)效果較好，水驅(qū)評價系數(shù)分別為1.150，0.710，0.700，說明在同一注采井間存在水相驅(qū)動能力相差較大的流場分布。因此僅對該井進(jìn)行注水優(yōu)化無法進(jìn)一步改善I4井的驅(qū)油效果，可考慮對I4井注入調(diào)剖劑，對具有較大水流量區(qū)域進(jìn)行分流，增大邊流處注水流量以提高波及效率。流線聚類方法同時也可用于尋找油藏內(nèi)具有開發(fā)潛力的區(qū)域，如本次聚類結(jié)果中第1、第2、第3類流線，由于其注采井間無裂縫或溶洞將注采井直接相連，流體流動相對緩慢，因此水驅(qū)過程相對平緩，是具有較好開發(fā)潛力的區(qū)域。

采用密度峰值聚類算法對不同注水開發(fā)時刻流線進(jìn)行聚類的結(jié)果見圖10。

圖10反映了流線聚類變化的歷史過程，對應(yīng)時刻的流場聚類性質(zhì)見表3。表3中不同類別流線的油水流動比及油水體積比表明，流線聚類方法對不同時刻的流場均可進(jìn)行較好的區(qū)分，不同類流線間油水流速比與油水體積比差異較大，區(qū)分效果均較好。

圖10 不同注水開發(fā)時刻流線聚類分布

表3 不同注水時間下不同類別流線性質(zhì)

3 結(jié)論

本文為水驅(qū)流場識別提供了一種新的研究方法，以流線聚類的方式將不同類型流場進(jìn)行區(qū)分，同時依據(jù)密度峰值聚類算法及輪廓系數(shù)聚類評價算法確定合理聚類數(shù)，確保聚類結(jié)果準(zhǔn)確性。密度峰值聚類方法結(jié)果較穩(wěn)定，且對流場區(qū)分效果較好，可作為流線聚類主要算法。

油水流速比及油水體積比可表征流場特征，依據(jù)聚類結(jié)果可找出潛在優(yōu)勢流場分布，有效識別油藏中無效注水循環(huán)通道以及具有開發(fā)潛力的區(qū)域，在本文流線聚類結(jié)果中，不同類流場間的油水流速比及油水體積比差異較大，表明區(qū)分效果較好；同時，流線聚類方法可將同一注采井間的流線細(xì)分為主流及支流，進(jìn)一步為注水優(yōu)化、井網(wǎng)層系調(diào)整、深部調(diào)剖等方案提供科學(xué)決策和技術(shù)支撐。流線聚類方法在不同注水時刻均可對不同類型流場進(jìn)行較好的區(qū)分，表明流線聚類算法對流場評價具有較好的適用性。

致謝：感謝斯倫貝謝公司在中國石油工程設(shè)計大賽期間為西南石油大學(xué)提供的Petrel RE軟件及Ocean平臺。

符號注釋：

a——流線；b——位于流線上的點；dc——截斷距離，無因次；dij——樣本i，j間的距離，無因次；ijq qd——樣本qi至樣本qj的距離，無因次；i，j——樣本序號；qi，qj——依據(jù)局部密度降序排列的第i，j個樣本；nc——聚類數(shù)；na——流線a上的總點數(shù)；np——樣本總量；Soab，Swab——第a條流線上b點處油相、水相飽和度，%；Vow——油水體積比，無因次；voab，vwab——第a條流線上b點處油相、水相流動速率，m3/d；vow——油水流速比，無因次；xab，yab——第a條流線上b點處x，y軸坐標(biāo)；Xa，Ya——a流線中點位置x，y軸坐標(biāo)；γ——聚類中心評價系數(shù)，無因次；qδ——樣本的分離距離，無因次；ρ——樣本密度，無因次。

參考文獻(xiàn)：

[1]侯玉培.整裝油藏流場重整提高采收率的數(shù)值模擬研究[D].東營:中國石油大學(xué)(華東), 2011.HOU Yupei.Numerical simulation research of uncompartmentalized oilfields by flow field reforming to improve oil recovery[D].Dongying: China University of Petroleum, 2011.

[2]王繼強, 石成方, 紀(jì)淑紅, 等.特高含水期新型水驅(qū)特征曲線[J].石油勘探與開發(fā), 2017, 44(6): 955-960.WANG Jiqiang, SHI Chengfang, JI Shuhong, et al.New water drive characteristic curves at ultra-high water cut stage[J].Petroleum Exploration and Development, 2017, 44(6): 955-960.

[3]趙輝, 康志江, 孫海濤, 等.水驅(qū)開發(fā)多層油藏井間連通性反演模型[J].石油勘探與開發(fā), 2016, 43(1): 99-106.ZHAO Hui, KANG Zhijiang, SUN Haitao, et al.An interwell connectivity inversion model for waterflooded multilayer reservoirs[J].Petroleum Exploration and Development, 2016, 43(1):99-106.

[4]呂琦.基于流線的油藏數(shù)值模擬研究[D].東營: 中國石油大學(xué)(華東), 2009.LYU Qi.The research of numerical simulation of reservoir by streamline[D].Dongying: China University of Petroleum, 2009.

[5]張喬良, 姜瑞忠, 姜平, 等.油藏流場評價體系的建立及應(yīng)用[J].大慶石油地質(zhì)與開發(fā), 2014, 33(3): 86-89.ZHANG Qiaoliang, JIANG Ruizhong, JIANG Ping, et al.Establishment and application of oil reservoir flow-field evaluation system[J].Petroleum Geology & Oilfield Development in Daqing,2014, 33(3): 86-89.

[6]ANWAR A R, MOHAMED B, MUSTAFA A, et al.An effective streamline history matching technique[R].SPE 188073, 2017.

[7]MILLIKEN W J, EMANUEL A S, CHAKRAVARTY A.Applications of 3D streamline simulation to assist history matching[J].SPE Reservoir Evaluation & Engineering, 2001, 4(6): 502-508.

[8]CHENG H, OYERINDE A S, DATTA-GUPTA A, et al.Compressible streamlines and three-phase history matching[R].SPE 99465, 2007.

[9]SINGH A P, MAUCEC M, KNABE S P.History matching using streamline trajectories[R].SPE 172146, 2014.

[10]BADDOURAH M A, ALTURKI A A, AL-ZAWAWI A S, et al.Streamlines guided assisted history matching methodology: A new approach[R].SPE 178007, 2015.

[11]THIELE M R, BATYCKY R.Water injection optimization using a streamline-based workflow[R].SPE 84080, 2003.

[12]SPYROU C E, ROSA A P L, KHATANIAR S K, et al.An approach to alternative waterflood designs and operations using streamline simulation: Application to an oil field in the North German Basin[R].SPE 183879, 2017.

[13]BUKHAMSEEN N Y, ALNAJEM A A, SAFFAR A H, et al.An injection optimization decision-making tool using streamline based fuzzy logic workflow[R].SPE 186021, 2017.

[14]BATYCKY R P, THIELE M R.Technology update: Mature flood surveillance using streamlines[J].Journal of Petroleum Technology,2016, 68(5): 22-25.

[15]LI J, LEI Z, LI S, et al.Optimizing water flood performance to improve injector efficiency in fractured low-permeability reservoirs using streamline simulation[R].SPE 182779, 2016.

[16]M?YNER O, KROGSTAD S, LIE K A.The application of flow diagnostics for reservoir management[R].SPE 171557, 2015.

[17]MAQUI A F, ZHAI X, NEGREIRA A S, et al.A comprehensive workflow for near real time waterflood management and production optimization using reduced-physics and data-driven technologies[R].SPE 185614, 2017.

[18]AHONSI W, LAGER A, MUHAIRI L S A, et al.Tracer and streamline modelling for pattern validation and improved reservoir modelling[R].SPE 181541, 2016.

[19]RODRIGUEZ A, LAIO A.Clustering by fast search and find of density peaks[J].Science, 2014, 344(6191): 1492-1496.

[20]MATRINGE S F, GERRITSEN M G.On accurate tracing of streamlines[R].SPE 89920, 2004.

[21]RODRIGUEZ A, LAIO A.Supplementary materials for clustering by fast search and find of density[DB/OL].[2017-12-20].http://science.sciencemag.org/content/suppl/2014/06/25/344.6191.14 92.DC1.

[22]LIU Y, LI Z, XIONG H, et al.Understanding of internal clustering validation measures[R].Sydney: IEEE, 2010.

[23]ROUSSEEUW P.Silhouettes: A graphical aid to the interpretation and validation of cluster analysis[J].Journal of Computational and Applied Mathematics, 1987, 20(1): 53-65.