李懷旭，孫曉巖，陶少輝，夏力，項(xiàng)曙光

（青島科技大學(xué)化工學(xué)院，過程系統(tǒng)工程研究所，山東 青島 266042）

引 言

隨著石油加工領(lǐng)域進(jìn)入分子管理時(shí)代，業(yè)內(nèi)逐漸形成一套完整的油品分析模擬體系，即油品分子表征（重構(gòu)）-組分集總（劃分）-解集總（還原）體系。油品分子表征通過現(xiàn)代化分析儀器進(jìn)行[1-3]，在此基礎(chǔ)上依靠計(jì)算機(jī)輔助方式構(gòu)建分子體系的模型化方法來進(jìn)行分子重構(gòu)[4-7]。油品真實(shí)分子組成在帶來更豐富信息的同時(shí)，也極大增加了煉油過程模擬和優(yōu)化模型的規(guī)模，導(dǎo)致計(jì)算量激增，因而如何在盡量保留油品真實(shí)分子組成信息的同時(shí)，有效減少組分?jǐn)?shù)以提高煉油過程的模擬與優(yōu)化效率，是石油分子管理體系亟待解決的關(guān)鍵問題之一。而集總過程則可以根據(jù)一定規(guī)則將油品分子劃分為少數(shù)幾個(gè)虛擬組分，是減小過程模擬計(jì)算規(guī)模的有效方式[8-9]。如在油藏模擬領(lǐng)域，Montel 等[10]首次提出通過迭代的動(dòng)態(tài)聚類（dynamical clustering,DynC）算法將大量組分通過性質(zhì)的相似性進(jìn)行集總，將150 個(gè)組分用少于10 種虛擬組分替代，并提出了虛擬組分物性的計(jì)算方法?；谏倭刻摂M組分進(jìn)行模擬后，為獲得模擬結(jié)果中的產(chǎn)品真實(shí)組成，Leibovici 等[11]通過閃蒸計(jì)算實(shí)現(xiàn)解集總（delumping）來還原油品組分，封松等[12-13]在其構(gòu)建的分子級(jí)組成模型中利用可逆集總方案的映射規(guī)則倒推產(chǎn)品分子組成。

Montel 等[10]對(duì)真實(shí)組分進(jìn)行集總時(shí)，采用了需要事先指定類別數(shù)和初始聚類中心的迭代動(dòng)態(tài)聚類算法，但針對(duì)特定體系，確定虛擬組分?jǐn)?shù)目需要進(jìn)行大量的分析和比較，且聚類結(jié)果對(duì)于初始中心的選擇極為敏感[14-16]。這些弊端使得動(dòng)態(tài)聚類在實(shí)際應(yīng)用時(shí)難以產(chǎn)生穩(wěn)定的結(jié)果，即使Montel 等提出通過人工選取類別數(shù)量最多的組分為中心可以使算法變得更穩(wěn)定，但這是基于經(jīng)驗(yàn)的，依賴于使用者對(duì)數(shù)據(jù)集的了解程度。對(duì)于過程模擬計(jì)算而言，往往希望得到一套穩(wěn)健高效的集總計(jì)算模型，以適用于不同類型的數(shù)據(jù)。

基于此，本工作對(duì)油品分子的虛擬組分劃分?jǐn)?shù)目（即集總數(shù)目）對(duì)模擬計(jì)算的影響進(jìn)行討論分析，將機(jī)器學(xué)習(xí)領(lǐng)域中一種無須事先指定類別數(shù)和聚類中心的密度峰聚類（density peak clustering,DPC）[17]算法引入油品分子的集總過程，以真實(shí)分子的熱力學(xué)性質(zhì)為聚類變量，從而實(shí)現(xiàn)對(duì)表征汽油分子的集總，同時(shí)以相平衡與精餾計(jì)算對(duì)方法的有效性進(jìn)行驗(yàn)證。

1 模型與方法

1.1 DPC算法原理

DPC 算法基于一個(gè)核心假設(shè)，即認(rèn)為給定數(shù)據(jù)集的中心一定是該數(shù)據(jù)集中密度較大且互相遠(yuǎn)離（相對(duì)距離較大）的一些數(shù)據(jù)點(diǎn)。因而該算法的計(jì)算過程取決于數(shù)據(jù)集的距離矩陣dn×n（n為樣本個(gè)數(shù)）和每個(gè)點(diǎn)的密度ρi。距離是聚類算法中常用的相似性度量標(biāo)準(zhǔn)，采用何種距離計(jì)算公式取決于不同種類的算法和目標(biāo)數(shù)據(jù)集的分布情況。這里采用Montel等[10]提出的距離公式來度量數(shù)據(jù)點(diǎn)i與j間的差異，如式（1）所示。

式中，wk為不同性質(zhì)（維度）的權(quán)重值；pik與pjk為數(shù)據(jù)點(diǎn)的性質(zhì)值；N為選取性質(zhì)的個(gè)數(shù)（即數(shù)據(jù)集的維數(shù)）。采用式（1）度量性質(zhì)之間的相似性旨在擴(kuò)大性質(zhì)間的微小差距。對(duì)于數(shù)據(jù)點(diǎn)的密度，本文采用高斯核計(jì)算數(shù)據(jù)點(diǎn)的局部密度[18-19]：

式中，dc為截?cái)嗑嚯x，其物理意義是用于判定數(shù)據(jù)點(diǎn)的疏密程度，定義為將距離矩陣的元素升序排列，其1%～5%處對(duì)應(yīng)的距離。需要指出，dc為該算法唯一需要人工指定的參數(shù)，Rodriguez 等[17]的研究結(jié)果表明，其在一定范圍內(nèi)調(diào)節(jié)有助于更好地發(fā)現(xiàn)中心點(diǎn)。

相對(duì)距離δi定義為數(shù)據(jù)點(diǎn)i與所有密度大于點(diǎn)i的數(shù)據(jù)點(diǎn)之間的最小距離，對(duì)于密度最高的數(shù)據(jù)點(diǎn)，相對(duì)距離表達(dá)式為

根據(jù)各數(shù)據(jù)點(diǎn)的ρ和δ可做出二維平面決策圖，在圖中明顯離群的數(shù)據(jù)點(diǎn)即為類中心點(diǎn)，后文實(shí)例應(yīng)用時(shí)將會(huì)詳細(xì)展示如何尋找中心點(diǎn)。得到中心點(diǎn)后，將其余點(diǎn)分配給比其密度高的最近鄰點(diǎn)所在的類別即可。DPC算法步驟如下所示。

（1） 對(duì)數(shù)據(jù)集進(jìn)行標(biāo)準(zhǔn)化處理，將所有性質(zhì)變量變換到[0,1]區(qū)間；

（2） 通過式（1）計(jì)算距離矩陣dn×n：

（3） 通過式（2）計(jì)算每個(gè)數(shù)據(jù)點(diǎn)的密度；

（4） 通過式（3）、式（4）計(jì)算所有數(shù)據(jù)點(diǎn)的δ；

（5） 做出決策圖并選取類中心；

（6） 按照分配機(jī)制分配其他數(shù)據(jù)點(diǎn)。

由上述步驟可知，DPC 算法應(yīng)用到油品分子的集總劃分時(shí)的優(yōu)點(diǎn)為：（1）無須預(yù)定中心點(diǎn)，這使得在確定虛擬組分?jǐn)?shù)目時(shí)有足夠的參考；（2）算法所得到的類中心是隸屬于數(shù)據(jù)集的真實(shí)樣本，且將每一類看作是數(shù)據(jù)空間中被低密度區(qū)域分割開的稠密區(qū)域，故能檢測(cè)到任意形狀的數(shù)據(jù)集，對(duì)于非球形數(shù)據(jù)集有很好的適用性，一些研究已經(jīng)對(duì)此加以證明[20-25]；（3）DPC 算法流程中無迭代過程，類中心的發(fā)現(xiàn)和其余點(diǎn)的分配都是一次性完成的。這意味著在預(yù)先計(jì)算好距離矩陣后，DPC 算法效率不受數(shù)據(jù)集維度影響，而動(dòng)態(tài)聚類法雖然在處理大量數(shù)據(jù)時(shí)有良好的伸縮性，但需要反復(fù)地迭代計(jì)算距離并重置中心點(diǎn)，且迭代結(jié)果受初始中心的影響很大，每次的聚類結(jié)果也都存在差異[26-28]。

1.2 集總性質(zhì)選擇

Montel 等[10]提出以Peng-Robinson（PR）狀態(tài)方程參數(shù)a、b、m以及分子量四個(gè)參數(shù)作為集總性質(zhì)，其中參數(shù)a、b、m分別通過式(5)～式(7)求得。

可見，參數(shù)a和b均由臨界壓力Pc和臨界溫度Tc求出，m由偏心因子ω求解，四個(gè)參數(shù)除了能很好地反映純物質(zhì)的臨界特征和極性特征外，還能體現(xiàn)各個(gè)組分在熱力學(xué)過程中的性質(zhì)變化，且與進(jìn)行模擬計(jì)算時(shí)選取的熱力學(xué)計(jì)算方法保持一致。此外，Montel 等[10]建議a、b、m以及分子量（MW）在計(jì)算相似度時(shí)的權(quán)重值分別為1.0、1.0、0.1、0.5，這里也采用此權(quán)重值。

1.3 虛擬組分物性估算

對(duì)于集總得到的虛擬組分的物性計(jì)算，Leibovici 等[29]提出了集總系統(tǒng)特性估計(jì)方程對(duì)集總組分信息損失進(jìn)行最小化處理，該方法準(zhǔn)確且與使用的熱力學(xué)模型完全一致。但在集總組分過多的情況下，該方法計(jì)算過程過于復(fù)雜，且本文側(cè)重于研究集總劃分方法而非熱力學(xué)計(jì)算，故在虛擬組分物性計(jì)算方面，本文采用一組簡(jiǎn)單的石油餾分物性估算模型集合來進(jìn)行計(jì)算，具體如下所述。

臨界溫度和臨界壓力采用Riazi-Daubert 模型計(jì)算，模型關(guān)聯(lián)式如式（8）、式（9）所示[30]。

式中，Tb為正常沸點(diǎn)；Sg為相對(duì)密度。

偏心因子采用Lee-Kesler 模型計(jì)算，關(guān)聯(lián)式如式（10）所示[31]。

式中，xi表示第i個(gè)組分；c為組分總數(shù)。相對(duì)密度（Sg）可以用分子量和沸點(diǎn)通過關(guān)聯(lián)式估算，故這里不再收集每個(gè)組分的密度，其他諸如液相摩爾體積、水溶性、黏度、汽化潛熱等性質(zhì)的計(jì)算模型見文獻(xiàn)[32]。

2 實(shí)例應(yīng)用與討論

2.1 數(shù)據(jù)來源

油品數(shù)據(jù)來源于某催化汽油吸附脫硫工藝體系的287 種真實(shí)組分，其中包括C2～C12烯烴、C3～C12正構(gòu)烷烴、C4～C10異構(gòu)烷烴、C5～C10環(huán)烷烴、C6～C12芳烴以及二十余種含硫烴類化合物。

為方便比較，以a、b、m和MW 四個(gè)參數(shù)為數(shù)據(jù)特征。圖1 顯示了287 個(gè)數(shù)據(jù)點(diǎn)的分布情況，每一個(gè)子圖顯示了兩個(gè)性質(zhì)間的相互關(guān)系。從對(duì)角線的直方圖和其余散點(diǎn)圖可以大致判斷出由油品分子性質(zhì)組成的數(shù)據(jù)集有三個(gè)主要特征：(1)a、b、m和MW 歸一化數(shù)值都大量集中在0.25～0.75 之間；（2）組分兩兩性質(zhì)之間的線性相關(guān)性較差（擬合優(yōu)度在0.79～0.95 間）；(3) 數(shù)據(jù)集在二維平面上呈現(xiàn)出不規(guī)則的條狀或橢圓狀形態(tài)，大部分?jǐn)?shù)據(jù)點(diǎn)處于高度密集且重疊的狀態(tài)，不同類別之間沒有明顯的界限。顯然這對(duì)于一般的聚類方法而言很不友好。

圖1 標(biāo)準(zhǔn)化后的數(shù)據(jù)分布Fig.1 Standardized data distribution

2.2 集總數(shù)目對(duì)計(jì)算結(jié)果的影響

對(duì)于動(dòng)態(tài)聚類而言，集總數(shù)目的選擇會(huì)影響最終的集總效果，因此下文將通過指定壓力（0.1 MPa）和汽化分率以計(jì)算物流溫度的閃蒸過程，討論不同集總數(shù)目對(duì)計(jì)算結(jié)果的影響。

圖2 展示了使用動(dòng)態(tài)聚類的集總方法將287 種真實(shí)組分分別集總為4、6、8、16、21種虛擬組分的閃蒸計(jì)算結(jié)果。顯然，計(jì)算結(jié)果的準(zhǔn)確性與集總數(shù)目密切相關(guān)。當(dāng)集總數(shù)目為4 個(gè)時(shí)（4 components, 4 Comp），虛擬組分的閃蒸曲線的斜率小于真實(shí)組分，因?yàn)檫^少的集總數(shù)目會(huì)使得每個(gè)虛擬組分中包含許多性質(zhì)差異較大的組分，從而使其熱力學(xué)狀態(tài)變化難以描述實(shí)際情況，同樣6 Comp 時(shí)亦如此。再者，計(jì)算誤差也會(huì)因系統(tǒng)中虛擬組分?jǐn)?shù)量太少而擴(kuò)大。隨著集總數(shù)目的增加，真實(shí)組分得到更細(xì)致的分配，其閃蒸溫度曲線的斜率更接近實(shí)際情況，但當(dāng)集總數(shù)目太多時(shí)（21 Comp），許多本該屬于同一類別的組分被強(qiáng)行分配到不同類別，特別是對(duì)于密集分布的數(shù)據(jù)集，這會(huì)導(dǎo)致集總類別之間出現(xiàn)鏈?zhǔn)叫?yīng)，即相鄰的每個(gè)虛擬組分的性質(zhì)都會(huì)產(chǎn)生較大的誤差。

圖2 不同集總數(shù)目在0.1 MPa下的閃蒸計(jì)算結(jié)果Fig.2 Flash calculation results at 0.1 MPa for different lumping

考慮多個(gè)集總數(shù)目下，不同汽化分率時(shí)的閃蒸溫度的平均絕對(duì)誤差（mean absolute deviation,MAE），得到如圖3 所示的變化趨勢(shì)（對(duì)0 個(gè)集總數(shù)目下的MAE做了縮小處理）。這說明集總結(jié)果的準(zhǔn)確程度并不是隨著集總數(shù)目的增加持續(xù)上升，這中間存在或大或小的波動(dòng)。理論上，集總數(shù)目越多，計(jì)算規(guī)模便越大，復(fù)雜度自然也越大。本例中，集總數(shù)目在10～16 之間時(shí)，MAE 較小，集總數(shù)目持續(xù)上升至21時(shí)，MAE出現(xiàn)激增，而后又出現(xiàn)驟降。

圖3 閃蒸溫度的平均絕對(duì)誤差和計(jì)算復(fù)雜度隨集總數(shù)目增加的變化Fig.3 Variation of MAE and computational complexity of the flash temperature with the increase of the number of lumps

一般地，以誤差的大小衡量集總結(jié)果的準(zhǔn)確性，得到圖4（a）所示的變化趨勢(shì)。Opt（optimal）代表最優(yōu)集總數(shù)目，此時(shí)較少的虛擬組分?jǐn)?shù)量卻有較高的準(zhǔn)確性；Lim（limit）代表極限集總數(shù)目，此時(shí)的集總數(shù)目等于組分個(gè)數(shù)，相當(dāng)于每一個(gè)組分作為一個(gè)集總。小于Opt的部分代表許多性質(zhì)差異大的組合更多地被集總在一起。Opt 到Lim 之間代表最優(yōu)集總數(shù)目下的虛擬組分被拆開，其子組分先是被互相劃分到與其性質(zhì)差異較大的類別中，后又得到了更精細(xì)的組合，此時(shí)的復(fù)雜度也更大。圖4（b）描述了一種簡(jiǎn)化情況來表示Opt-Lim 間虛擬組分的劃分情況，P1,P2,P3 為Opt 中的某三個(gè)虛擬組分，實(shí)線圈起的部分代表原本正確的劃分被拆開，虛線圈起的部分代表被拆開的子組分重新組合為新的虛擬組分（假設(shè)新的虛擬組分比之前劃分的更加精細(xì)）。

圖4 計(jì)算精度與組分分配情況隨集總數(shù)目變化的示意圖Fig.4 Schematic diagram of the variation of calculation accuracy and component distribution with the number of lumps M, S—增加的集總數(shù)目; v, r, z—真實(shí)子組分

此發(fā)現(xiàn)具有一定意義，說明在組分集總過程中，確實(shí)存在一種兼顧效率和準(zhǔn)確性的可能，而DPC 集總方法的作用就在于自適應(yīng)地發(fā)現(xiàn)圖中所示的Opt。需要指出，當(dāng)物性估算模型依賴于摩爾分?jǐn)?shù)權(quán)重時(shí)，這一規(guī)律并不明顯。因?yàn)閷?duì)于含量很少的組分，即使被錯(cuò)誤分配，對(duì)虛擬組分的性質(zhì)影響也不是很大。但對(duì)于含量高的組分（尤其是重組分），差距將會(huì)很明顯。

2.3 DPC集總過程

首先根據(jù)式（1）、式（2）計(jì)算得到距離矩陣d287×287和局部密度ρi，其中dc取1%計(jì)算（經(jīng)試驗(yàn)，本例dc在1%～5%內(nèi)變化對(duì)最終聚類結(jié)果影響很小），再根據(jù)式（3）、式（4）得到δi，得到如圖5（a）所示的平面決策圖。圖中彩色的數(shù)據(jù)點(diǎn)為明顯離群的點(diǎn)，根據(jù)算法假設(shè)，這些數(shù)據(jù)點(diǎn)的共同特征是同時(shí)具有較大的δi和ρi，這意味著它們彼此相互遠(yuǎn)離且擁有各自的稠密區(qū)域，故將此視為數(shù)據(jù)集的10 個(gè)中心點(diǎn)。為了更清晰地展示類中心點(diǎn)與其他數(shù)據(jù)點(diǎn)區(qū)別，將每個(gè)數(shù)據(jù)點(diǎn)的ρ2和δ相乘，記作γ并降序排列，采用Lv 等[33]提出的一種通用閾值計(jì)算方法計(jì)算γ閾值（詳細(xì)算法過程參見該文獻(xiàn)），如圖5（b）所示，從決策圖中獲取的10 個(gè)中心點(diǎn)在圖5（b）中依然得到了突出顯示。就此將油品中含有的287個(gè)烴類分子劃分為10類，即10個(gè)虛擬組分。

圖5 集總數(shù)目及中心組分的決策圖Fig.5 Decision diagram of the number of lumps and center components

此外，決策圖中的黑色菱形標(biāo)記也是數(shù)據(jù)集中明顯離群的數(shù)據(jù)點(diǎn)，但與彩色的中心點(diǎn)不同，菱形點(diǎn)只有較大的δ，其ρ卻異常?。ń咏?），這說明此類數(shù)據(jù)點(diǎn)彼此孤立且相互遠(yuǎn)離，在整個(gè)數(shù)據(jù)集中不屬于任何一類，本文將其定義為“離群組分”。對(duì)于其他ρ極小的點(diǎn)（緊靠縱坐標(biāo)軸），因?yàn)槠洇牟⒉煌怀觯什荒茏鳛殡x群組分，這類點(diǎn)通常處于“團(tuán)簇”的邊緣。通過查找，3 個(gè)離群組分分別對(duì)應(yīng)正十二烷、碳二苯并噻吩、1,4-二異丙苯（按δ值降序），其部分性質(zhì)已經(jīng)列于表1，與其他組分比較發(fā)現(xiàn)，它們具有個(gè)別較為突出的性質(zhì)，因此將其視作“離群”是合理的。對(duì)于離群組分，本文的處理方式與算法要求不同（Rodriguez 等[17]認(rèn)為每個(gè)離群點(diǎn)應(yīng)單獨(dú)為一類），本文建議兩種處理方法：若離群組分?jǐn)?shù)量太多，則將其歸為一類，作為一個(gè)虛擬組分；若其數(shù)量較少、含量低，則將其包含在較近的類別中。本例中僅含2個(gè)離群組分，且為重組分，故采取第二種處理方法。

表1 離群組分的部分性質(zhì)Table 1 Partial properties of outlier components

需要指出，本例中含有硫化氫氣體和二氧化硫氣體，其理論上也為離群組分，但算法并未偵測(cè)到兩者，在此有必要做出解釋：DPC 集總過程使用了式（1）來度量組分間的相似性，式（1）雖能擴(kuò)大性質(zhì)間的差異，但有一個(gè)潛在問題，即對(duì)于較小的數(shù)值并不靈敏。如圖6 所示，點(diǎn)Ⅰ、Ⅲ代表兩個(gè)離群組分，點(diǎn)Ⅱ代表某正常組分，且假設(shè)A、B、C為三個(gè)組分的性質(zhì)。其中點(diǎn)Ⅰ的性質(zhì)比點(diǎn)Ⅱ低50%（反向差距），點(diǎn)Ⅲ的性質(zhì)比點(diǎn)Ⅱ高50%（正向差距），若采用歐氏距離計(jì)算空間中的真實(shí)距離，則兩者與點(diǎn)Ⅱ的相似度是一致的，兩者都應(yīng)該計(jì)算出較小的ρ，進(jìn)而得到較大的δ，最后形成孤立點(diǎn)。但當(dāng)采用式（1）時(shí)，反向差距大的組分得不到較大的相對(duì)距離，即使它看起來與正向差距大的組分并無區(qū)別，這便導(dǎo)致了所有類似的組分不能在決策圖中突出顯示。本例中，H2S 和SO2具有很小的分子量，很低的臨界溫度和很大的臨界壓力，故其a、b、m、MW 均很小，算法沒有對(duì)其進(jìn)行挖掘。實(shí)際上，諸如SO2等不凝氣或一些雜質(zhì)氣體往往單獨(dú)從油品分子中劃分出來，所以在檢測(cè)離群組分時(shí)可以不考慮這一類反向差距很大的組分；相反，重組分的偵測(cè)更加重要，同時(shí)也更難偵測(cè)，而式（1）對(duì)于正向差距大的組分敏感性更強(qiáng)，故這種處理方法是合理的。

圖7（a）為最終集總結(jié)果，圖中的不同著色代表不同的類別，即反映不同真實(shí)組分所屬的虛擬組分，對(duì)角線的直方圖匯總了每種虛擬組分的性質(zhì)分布情況。通過觀察圖1 發(fā)現(xiàn)，參數(shù)a、MW 相關(guān)性較其他更高，為更直觀地顯示集總結(jié)果，暫時(shí)忽略性質(zhì)a，僅對(duì)b、m、MW 進(jìn)行三維可視化，如圖7（b）所示。經(jīng)集總得到的虛擬組分的部分物性由式（8）～式（14）計(jì)算，與其中心組分的性質(zhì)一并列于表2。

表2 中心組分與對(duì)應(yīng)的虛擬組分的部分性質(zhì)Table 2 Partial properties of center components and corresponding pseudo-components

圖8展示了每個(gè)虛擬組分中子組分的分子類型和碳原子數(shù)?？梢?，PC4 和PC5 分別主要包含C8異構(gòu)烷、C10異構(gòu)烷，兩者數(shù)量較少，即所謂的“小簇”。關(guān)注數(shù)據(jù)集中的“小簇”極為重要，這在一定程度上能反映算法對(duì)數(shù)據(jù)進(jìn)行挖掘的徹底性。從圖7（b）中又可以發(fā)現(xiàn)兩個(gè)“小簇”都分布在PC6和PC10中，而PC6 包含較多C8～C9的異構(gòu)烷及環(huán)烷，PC10 則主要為C11以上的芳烴。

圖7 287個(gè)組分的集總結(jié)果分布及結(jié)果可視化Fig.7 Lumped results of 287 components and visualization of results

表面上，不管從數(shù)據(jù)集分布形態(tài)還是具體組分來看，似乎PC4 應(yīng)隸屬于PC6。但實(shí)際發(fā)現(xiàn)，PC4 中10 個(gè)組分的平均偏心因子為0.38，而PC6 中39 個(gè)組分的平均偏心因子僅為0.30，且有近半數(shù)的組分的偏心因子小于0.30，這說明PC6 中組分的極性與PC4 有明顯差距，將此團(tuán)數(shù)據(jù)集分成兩類是符合實(shí)際的。對(duì)于PC5 而言，結(jié)合圖8（a）和圖8（b）可以很直觀地發(fā)現(xiàn)其與PC10 中的組分性質(zhì)差距更大，PC10 以C11以上的芳烴為主，其極性較大，而PC5 以C10以下的異構(gòu)烷烴為主，其極性要小于芳烴。再者，通過比較發(fā)現(xiàn)，PC5 中組分的平均臨界壓力為2.17 MPa，平均臨界溫度為337℃，而PC10 為2.52 MPa、415℃，兩者存在明顯差距。而DPC 集總方法恰好可以有效地對(duì)PC4 和PC5 進(jìn)行挖掘，這得益于該算法擅長(zhǎng)處理具有重疊性的非球形數(shù)據(jù)集。

圖8 虛擬組分中子組分的烴類分子構(gòu)成和碳原子數(shù)分布Fig.8 Hydrocarbon composition and carbon number distribution of the subcomponents of pseudo-component

2.4 相平衡計(jì)算驗(yàn)證

本文重在考察虛擬組分在汽液相平衡計(jì)算中的表現(xiàn)。Leibovici 等[11]在其體系已證實(shí)了基于動(dòng)態(tài)聚類的集總方法在閃蒸計(jì)算中的可靠性。需要指出，基于動(dòng)態(tài)聚類的油品集總方案需要預(yù)先指定集總數(shù)目，為便于比較，將基于DPC 算法集總得到的10 個(gè)虛擬組分作為其預(yù)先設(shè)定的集總數(shù)目。下面展示了DPC 集總方法和動(dòng)態(tài)聚類方法與真實(shí)組分的對(duì)比。

很明顯，DPC 方法取得了預(yù)期效果。圖9 中，其閃蒸溫度曲線與真實(shí)組分的貼合程度較高；而采用動(dòng)態(tài)聚類的集總方法即使預(yù)先給定了一個(gè)近乎合理的集總數(shù)目，其效果與DPC 方法相比也有一些差距，尤其是在接近露點(diǎn)溫度時(shí)誤差較大，原因是動(dòng)態(tài)聚類無法捕捉不規(guī)則或具有層疊性的數(shù)據(jù)，故部分?jǐn)?shù)據(jù)點(diǎn)被劃分到了與自身性質(zhì)差異較大的類別中，進(jìn)而在宏觀上引起了計(jì)算結(jié)果的差異。從圖9得出，DPC 集總方法所得到的閃蒸溫度的MAE 較動(dòng)態(tài)聚類方法而言下降36.0%，計(jì)算精度有明顯提升。

圖9 不同集總方法在0.1 MPa下的閃蒸溫度曲線對(duì)比Fig.9 Comparison of flash temperature curves at 0.1 MPa for different lumped methods

圖10為不同集總方法在0.1 MPa下的閃蒸熱負(fù)荷的絕對(duì)誤差隨汽化分率的變化。與真實(shí)情況相比，DPC 方法得到的各點(diǎn)誤差均在較小范圍之內(nèi)，而基于動(dòng)態(tài)聚類的方法只有在汽化分率0.7～0.8 之間才能達(dá)到與其持平的效果。從圖10得出，DPC 集總方法所得到的閃蒸熱負(fù)荷的MAE 較動(dòng)態(tài)聚類方法而言下降64.2%，計(jì)算精度有明顯提升。熱負(fù)荷是對(duì)虛擬組分的熱容、焓值等性質(zhì)的綜合考察，故此結(jié)果在一定程度上說明了DPC 集總方法的可靠性。

圖10 不同集總方法下閃蒸熱負(fù)荷的絕對(duì)誤差（與真實(shí)組分下的熱負(fù)荷相比）比較Fig.10 Comparison of the absolute error of the flash heat duty(compared to the heat load of the real component) of different lumped methods

2.5 中心組分和虛擬組分的比較

前文已述，通過DPC 算法得到聚類中心隸屬于真實(shí)組分，其在一定程度上能代表子組分的性質(zhì)，故有必要討論中心組分是否可以表征油品分子的真實(shí)性質(zhì)。圖11、圖12展示了中心組分和虛擬組分分別表征油品時(shí)的結(jié)果對(duì)比。由圖可見，在汽化分率小于0.1 時(shí)，中心組分代替真實(shí)組分的效果尚好，但隨著汽化分率的增大，其閃蒸溫度總是低于真實(shí)情況，這說明中心組分中有過多的輕組分（相比于虛擬組分而言較輕）。從圖12 的PT 包絡(luò)圖來看，臨界條件的低估又說明中心組分的極性和沸點(diǎn)難以反映真實(shí)情況。表3中列出了每個(gè)虛擬組分與中心組分的沸點(diǎn)差距，可見PC6 的中心組分的沸點(diǎn)比虛擬組分低28.33 K，PC8 的中心組分的沸點(diǎn)比虛擬組分低32.60 K。

圖11 中心組分與虛擬組分在不同汽化分率下的閃蒸溫度對(duì)比（0.1 MPa）Fig.11 Comparison of the flash temperature of the center component and the pseudo-component at different vaporization fractions at 0.1 MPa

圖12 中心組分與虛擬組分的PT相包絡(luò)圖對(duì)比（1 bar=0.1 MPa）Fig.12 Comparison of the PT phase envelopes of center component and pseudo-component

表3 各虛擬組分與其中心組分的沸點(diǎn)差距Table 3 Boiling point gap between each pseudocomponent and its center component

雖然中心組分目前不能完全反映虛擬組分的性質(zhì)，但分析其中的原因尚有價(jià)值。實(shí)際上，通過算法得到的虛擬組分，可看作是空間中一塊稠密程度從中心向周圍不規(guī)則發(fā)散的區(qū)域。如圖13（a）所示，以b和m為特征性質(zhì)，將PC8 中子組分的局部密度ρ映射到氣泡圖中，圖中面積最大的點(diǎn)為中心組分，可見，ρ較大的組分全部靠近中心組分，這些組分的性質(zhì)與中心組分更相似；而ρ較小的組分距離中心組分較遠(yuǎn)，其與中心組分的相似度也就較差（相比于本類別中的ρ較大的組分而言）。更直觀地，圖13（b）的密度分布熱力圖展示了中心組分對(duì)其子組分的描述能力大小。熱力圖的核心部分（黑色區(qū)域）為中心組分，它對(duì)其周圍的密度接近自己的點(diǎn)能很好反映；但對(duì)稍微遠(yuǎn)離的點(diǎn)（綠色及藍(lán)色區(qū)域），其代表性是有限的。對(duì)于包含子組分較多的虛擬組分，由于當(dāng)前數(shù)據(jù)分布的不規(guī)則性，其中心組分很難描述到每一個(gè)組分的性質(zhì)。而用虛擬組分來表征真實(shí)油品時(shí)，實(shí)際上等同于用圖中低密度的點(diǎn)去修正中心點(diǎn)，只是修正過程涉及了其區(qū)域中的全部點(diǎn)。

圖13 PC8中每個(gè)真實(shí)組分的局部密度與所在空間位置的可視化Fig.13 Display of local density and spatial location of each real component in PC8

2.6 精餾模擬驗(yàn)證

下面使用ASPEN 進(jìn)一步驗(yàn)證集總方法的性能。使用本文方法（DPC）和動(dòng)態(tài)聚類法（DynC）分別得到的10個(gè)虛擬組分進(jìn)行常壓精餾模擬，并與真實(shí)組分的結(jié)果進(jìn)行對(duì)比。此外，還將動(dòng)態(tài)聚類法得到的兩種不恰當(dāng)集總的結(jié)果也加入比較。進(jìn)料條件及塔參數(shù)設(shè)置如圖14所示，模擬結(jié)果的對(duì)比展示于圖15 中。顯然，使用虛擬組分替代真實(shí)油品進(jìn)行精餾模擬時(shí)，方法面臨了嚴(yán)峻挑戰(zhàn)：隨著D/F（塔頂采出量與進(jìn)料量之比）的增大，4 個(gè)集總結(jié)果的塔頂采出溫度均被高估。其中DynC（10 Comp）方法的塔頂溫度MAE 為5.03℃，塔釜溫度MAE 為4.34℃，相比之下DPC 方法取得了更可令人接受的結(jié)果，其塔頂溫度的MAE 為1.95℃，精度提高61.2%；塔釜溫度的MAE 為2.01℃，精度提高53.7%，這為下一步的“組分還原”奠定了相對(duì)可靠的基礎(chǔ)。

圖14 進(jìn)料條件及精餾塔參數(shù)設(shè)置Fig.14 Feeding conditions and distillation column parameter settings

圖15 不同集總結(jié)果下的精餾計(jì)算結(jié)果Fig.15 Calculation results of distillation of different lumped results

2.7 模擬計(jì)算效率

當(dāng)組分?jǐn)?shù)目增加時(shí)，精餾塔的待求方程組規(guī)模也隨之增大，故采用少數(shù)虛擬組分代替大量真實(shí)組分進(jìn)行模擬計(jì)算和調(diào)和設(shè)計(jì)時(shí)，在計(jì)算效率上無疑是有明顯進(jìn)步的。下面以ASPEN 為模擬計(jì)算工具，使用微型計(jì)算機(jī)對(duì)10 個(gè)虛擬組分和287 個(gè)真實(shí)組分的精餾計(jì)算時(shí)間(CPU 使用時(shí)間)進(jìn)行比較（針對(duì)第2.6節(jié)中的流程），計(jì)算機(jī)的硬件配置如下。

通過MATLAB 與ASPEN 的通信功能獲取到兩種情況的計(jì)算時(shí)間(含程序的打開/關(guān)閉時(shí)間，每次計(jì)算前均保證CPU 占用一致)如下：287 個(gè)真實(shí)組分為47.856 s，10 個(gè)虛擬組分為9.671 s。顯然兩者存在顯著差距。此外，若過程模擬涉及包含精餾塔的撕裂流收斂，或針對(duì)精餾塔進(jìn)行優(yōu)化，則計(jì)算效率的差距將更大。

3 結(jié) 論

將基于密度峰的聚類算法思想引入油品分子集總過程，實(shí)現(xiàn)了集總數(shù)目（中心組分）的自動(dòng)偵測(cè)和脫硫汽油組成的可靠劃分，并以某裝置的脫硫汽油實(shí)例進(jìn)行驗(yàn)證，得到以下結(jié)論。

（1）在油品分子集總過程中，存在一個(gè)兼顧效率和準(zhǔn)確性的集總數(shù)目，而本研究提出的方法借助決策圖即可簡(jiǎn)便地確定這個(gè)集總數(shù)目和少數(shù)離群組分，這在一定程度上克服了估計(jì)集總數(shù)目時(shí)的主觀性和盲目性。此外，該方法可對(duì)油品性質(zhì)數(shù)據(jù)集中分布不規(guī)則的數(shù)據(jù)或“小簇”進(jìn)行有效挖掘，在一定程度上避免了傳統(tǒng)劃分方法造成的分配錯(cuò)誤，從而使油品分子得到更精確的劃分。

（2）中心組分代替真實(shí)油品的能力是有限的；而虛擬組分可綜合其子組分（真實(shí)組分）的性質(zhì)，更好地反映真實(shí)情況。但中心組分作為真實(shí)組分，其熱力學(xué)計(jì)算相比于虛擬組分而言是嚴(yán)格的，關(guān)于中心組分描述真實(shí)油品的能力是否可以優(yōu)化，還有待進(jìn)一步討論。

（3）實(shí)例驗(yàn)證結(jié)果表明，DPC 集總方法得到的虛擬組分可以滿足相平衡計(jì)算的要求，該方法對(duì)于將來構(gòu)建完整的石油組分分子水平劃分具有較大潛力。

本工作在集總性質(zhì)的選取方面，現(xiàn)階段使用了固定的熱力學(xué)性質(zhì)組合及性質(zhì)權(quán)重，并未考慮是否有更優(yōu)的選擇?？紤]到不同煉油工藝過程可能對(duì)所需的性質(zhì)組合和權(quán)重值要求不同，故之后的工作重點(diǎn)在于，通過優(yōu)化性質(zhì)組合提高集總結(jié)果對(duì)油品的可描述性以及不同需求下的最優(yōu)權(quán)重設(shè)計(jì)。

符 號(hào) 說 明

c——虛擬組分中的子組分總數(shù)

dc——截?cái)嗑嚯x

dij——數(shù)據(jù)點(diǎn)i與數(shù)據(jù)點(diǎn)j的相互距離

MAE——平均絕對(duì)誤差

MW——分子量

N——集總性質(zhì)個(gè)數(shù)

n——油品真實(shí)組分總數(shù)

Pb——正常沸點(diǎn)下的壓力，MPa

Pbr——正常沸點(diǎn)下對(duì)比飽和蒸氣壓

Pc——純組分臨界壓力，MPa

pik——數(shù)據(jù)點(diǎn)i的第k個(gè)維度（性質(zhì)）的值

R——?dú)怏w常數(shù)，8.314 J/(mol·K)

Sg——相對(duì)密度

Tb——沸點(diǎn)，K

Tbr——對(duì)比正常沸點(diǎn)

Tc——純組分臨界溫度，K

wk——維度（性質(zhì)）k的權(quán)重值

γi——ρi與δi的乘積

δi——第i個(gè)數(shù)據(jù)點(diǎn)的相對(duì)距離

ρi——第i個(gè)數(shù)據(jù)點(diǎn)的局部密度

ω——偏心因子

下角標(biāo)

i,j——組分的樣本序號(hào)

k——性質(zhì)的樣本序號(hào)