楊曉峰,謝 巍,張浪文?

(1.華南理工大學(xué)自動化科學(xué)與工程學(xué)院,廣東廣州 510640;2.廣州市標(biāo)準(zhǔn)化研究院,廣東廣州 510110)

1 引言

現(xiàn)代化學(xué)過程通常由相互連接的操作單元組成,它們通過能量和信息流緊密地集成在一起.先進(jìn)控制算法在應(yīng)用于這些過程時應(yīng)滿足預(yù)定的目標(biāo),達(dá)到預(yù)期的安全,可持續(xù)性和經(jīng)濟(jì)性水平[1].然而,傳統(tǒng)的集中或分散控制架構(gòu)難以適用于復(fù)雜過程系統(tǒng)[2].

近十年來,為提高計算效率,維護(hù)靈活性和容錯性,分布式控制架構(gòu)得到快速發(fā)展[3].分布式控制系統(tǒng)設(shè)計包括兩個關(guān)鍵步驟[4]:1)合理地將整個系統(tǒng)劃分成較小的子系統(tǒng),需要保證子系統(tǒng)的能控性、能觀性和靈活性等;2)針對每個子系統(tǒng)設(shè)計局部控制器及其協(xié)調(diào)算法.在過去的十年里,分布式控制算法(尤其是分布式模型預(yù)測控制)吸引了非常廣泛的研究關(guān)注[5].然而,與其同等重要的子系統(tǒng)劃分問題受到的關(guān)注相對較少[6].

近幾年來,子系統(tǒng)劃分方法取得了一些成果[7–9].文獻(xiàn)[7]提出了一種基于層次聚類的分散控制系統(tǒng)分解方法.文獻(xiàn)[8]提出了一種既考慮耦合強(qiáng)度又考慮結(jié)構(gòu)緊密性的輸入–輸出配對方法.另外,文獻(xiàn)[9]研究了基于無權(quán)重有向圖的社區(qū)發(fā)現(xiàn)子系統(tǒng)劃分算法.在社區(qū)發(fā)現(xiàn)算法中,將一個大規(guī)模系統(tǒng)看作一個網(wǎng)絡(luò),通過將網(wǎng)絡(luò)分解為較小的子網(wǎng)絡(luò),建立面向分布式控制或估計的子系統(tǒng)模型.

考慮到不同狀態(tài)和輸出變量之間的結(jié)構(gòu)緊密性,已有研究嘗試將非線性系統(tǒng)分解為較小的單元,然后設(shè)計分布式狀態(tài)估計[10–11].文獻(xiàn)[10]提出了一種面向分布式狀態(tài)估計的過程網(wǎng)絡(luò)子系統(tǒng)劃分方法.在文獻(xiàn)[11],分布式狀態(tài)估計和控制在一個架構(gòu)內(nèi)被共同考慮.上述面向分布式狀態(tài)估計的子系統(tǒng)劃分方法均只考慮過程系統(tǒng)的物理拓?fù)?而未考慮不同變量之間的連接強(qiáng)度.研究發(fā)現(xiàn),如果不考慮變量間連通性的強(qiáng)弱,可能導(dǎo)致子系統(tǒng)劃分的次優(yōu)或不適當(dāng),特別是對于連接強(qiáng)度差異較大的系統(tǒng)影響更大[12].

基于上述分析,本文提出一種面向分布式狀態(tài)估計的子系統(tǒng)劃分方法,構(gòu)建了加權(quán)有向圖,既考慮結(jié)構(gòu)的緊密性,又考慮系統(tǒng)狀態(tài)/輸出之間的連接強(qiáng)度.基于加權(quán)有向圖,利用社區(qū)發(fā)現(xiàn)算法將所有系統(tǒng)變量劃分為較小的群組,使得劃分得到的子系統(tǒng)內(nèi)部關(guān)聯(lián)較強(qiáng),而子系統(tǒng)之間的耦合強(qiáng)度較弱.

2 問題描述

考慮如下復(fù)雜非線性系統(tǒng)的子系統(tǒng)劃分問題:

其中:x ∈Rnx表示狀態(tài)向量,y ∈Rny表示輸出向量,f和h分別表示非線性系統(tǒng)狀態(tài)和輸出方程.

本研究的目標(biāo)是將系統(tǒng)(1)劃分成多個子系統(tǒng),得到第i個子系統(tǒng)的模型,以設(shè)計分布式狀態(tài)估計器

其中:x(i)∈Rnx(i)表示第i個子系統(tǒng)的狀態(tài)向量,y(i)∈Rny(i)表示第i個子系統(tǒng)的測量輸出,X(i)表示所有與子系統(tǒng)i直接關(guān)聯(lián)的子系統(tǒng)狀態(tài),i=1,···,p,p是劃分得到子系統(tǒng)的個數(shù).在本設(shè)計中,假設(shè)p是已知的且數(shù)量少于測量輸出的個數(shù),即pny.

由于本文子系統(tǒng)的劃分是面向分布式狀態(tài)估計設(shè)計,對于劃分的子系統(tǒng),需要進(jìn)行能觀性的判斷.為了檢驗非線性系統(tǒng)(1)的能觀性,可以判斷該系統(tǒng)的能觀矩陣是否滿秩[13],能觀矩陣通過如下方式構(gòu)建:

其中:Lfh 表示方程h 對方程f 的Lie導(dǎo)數(shù),其定義為表示h對f的r階Lie導(dǎo)數(shù),定義為如果對于x ∈X,有rank(Q(x))=nx,那么,系統(tǒng)(1)在范圍X內(nèi)是局部可觀的[13],否則系統(tǒng)(1)是不可觀的.

社區(qū)發(fā)現(xiàn)算法是一種將大規(guī)模復(fù)雜網(wǎng)絡(luò)劃分成多個子網(wǎng)絡(luò)的有效工具[14],通過由文獻(xiàn)[15]定義的模塊度來評價一個網(wǎng)絡(luò)劃分的好壞.考慮一個帶有N個節(jié)點的網(wǎng)絡(luò),其不加權(quán)的鄰接矩陣A是一個N ×N矩陣,定義元素Aij為

通常,不考慮節(jié)點自身的回路,也就是說鄰接矩陣的主對角線上的元素為0.用模塊度來評價劃分好壞,模塊度Q的范圍為0到1,社區(qū)發(fā)現(xiàn)算法通過最大化模塊度來找到最優(yōu)群組(子系統(tǒng))結(jié)構(gòu),使得劃分的子系統(tǒng)內(nèi)部的連接數(shù)最多,而群組之間的連接數(shù)最少.一般模塊度Q在0.3到0.7之間時被認(rèn)為是較好的劃分方法,否則表示該系統(tǒng)不適合再進(jìn)行劃分.

然而,實際應(yīng)用的子系統(tǒng)劃分中,僅考慮劃分得到的子系統(tǒng)內(nèi)部連接邊數(shù)最多是不夠的,還需要考慮連接邊的連接權(quán)重.本文將對系統(tǒng)構(gòu)造帶權(quán)重的有向圖,并利用模塊度最大化方法對系統(tǒng)進(jìn)行劃分,以保證各個劃分子系統(tǒng)之間的連接權(quán)重最小.

3 基于加權(quán)有向圖的社區(qū)發(fā)現(xiàn)子系統(tǒng)劃分

本節(jié)將提出一種系統(tǒng)化的子系統(tǒng)劃分方法,考慮將非系統(tǒng)(1)劃分成多個子系統(tǒng),并設(shè)計滾動時域估計算法對子系統(tǒng)狀態(tài)進(jìn)行分布式狀態(tài)估計.

本文方法的主要結(jié)構(gòu)如圖1所示.首先,將非線性系統(tǒng)(1)用加權(quán)有向圖進(jìn)行描述,其中系統(tǒng)狀態(tài)和測量輸出變量看成網(wǎng)絡(luò)中的節(jié)點,這些節(jié)點通過加權(quán)連接邊進(jìn)行連接,連接邊上的權(quán)重反映了節(jié)點之間的連接強(qiáng)度.然后,基于加權(quán)有向圖,利用社區(qū)發(fā)現(xiàn)算法找到最大模塊度,將系統(tǒng)劃分成多個子系統(tǒng)群組.最后,通過檢測劃分子系統(tǒng)的能觀性,獲得子系統(tǒng)劃分結(jié)果.

3.1 加權(quán)連接邊的定義

本文將利用加權(quán)有向圖對系統(tǒng)(1)進(jìn)行描述,具體來說,將所有的狀態(tài)和測量輸出變量當(dāng)做有向圖的節(jié)點,這些節(jié)點通過有向邊進(jìn)行連接.令fi表示向量方程f的第i(i=1,···,nx)行,hj表示向量方程h的第j(j=1,···,ny)行,xi表示狀態(tài)向量x的第i(i=1,···,nx)行,yj(j=1,···,ny)表示輸出變量y的第j行.

在已有的無權(quán)重有向圖中,不同變量(節(jié)點)之間的連接強(qiáng)度未被考慮.為獲得權(quán)重,對非線性系統(tǒng)在工作點xs求一階偏導(dǎo)得到有向圖中連接邊的敏感度

圖1 本文子系統(tǒng)劃分結(jié)構(gòu)Fig.1 Flowchart of the proposed subsystem decomposition

定義如下連接邊的權(quán)重:

狀態(tài)變量xi到另一個狀態(tài)變量xl連接邊的權(quán)重定義為(l,i=1,···,nx)

狀態(tài)變量xl到測量輸出變量yj連接邊的權(quán)重定義為(l=1,···,nx,j=1,···,ny)

其中α是一個介于[0 1]的參數(shù),用于調(diào)節(jié)連接邊的權(quán)重在子系統(tǒng)劃分過程中的重要性,當(dāng)α=0時,所有連接邊的權(quán)重為1,即等效于已有的無權(quán)重有向圖;當(dāng)α增加時,連接邊的權(quán)重將在子系統(tǒng)劃分中更多地被考慮.當(dāng)S(xl,yj)=0或S(xi,xl)=0時,連接邊的權(quán)重為無窮,即兩個節(jié)點之間沒有直接連接.

由于不同的α值將影響子系統(tǒng)的劃分結(jié)果,如何選擇合適的參數(shù)α對子系統(tǒng)劃分相當(dāng)重要.總體而言,當(dāng)有向圖中所有連接邊的敏感度差異很小時,在子系統(tǒng)劃分時應(yīng)重點考慮連接邊的數(shù)量,此時選擇較小的α即可.相反,當(dāng)所有連接邊的敏感度差異很大時,應(yīng)更多地考慮連接邊的權(quán)重對子系統(tǒng)劃分的影響[16].

3.2 最短路徑

為實施社區(qū)發(fā)現(xiàn)算法,需要構(gòu)造加權(quán)鄰接矩陣,因此首先應(yīng)找到有向圖中的一個節(jié)點到另一個節(jié)點的最短路徑.令P表示一個狀態(tài)變量到另一個狀態(tài)變量或輸出變量的路徑,在這條路徑P中的所有連接邊表示為e ∈P.路徑L(P)的長度可以計算得到

其中:l=1,···,nx, j=1,···,ny,Plj,Plj表示從狀態(tài)xl到輸出yj所有路徑和其中一條路徑.本文利用Dijkstra’s算法[17]找到任意兩個節(jié)點的最短路徑,其計算復(fù)雜度為O(E log V),其中V 和E分別是節(jié)點的數(shù)量和連接的邊數(shù).

3.3 構(gòu)建鄰接矩陣

本節(jié)將構(gòu)建非線性系統(tǒng)(1)的鄰接矩陣.對于系統(tǒng)(1)所形成的有向圖,總共有na(其中na=nx+ny)個頂點.令ca為狀態(tài)和測量輸出向量,定義為ca=[x1,···,xnx,y1,···,yny].

那么,所構(gòu)建的加權(quán)鄰接矩陣Aw將被用于社區(qū)發(fā)現(xiàn)算法,以對系統(tǒng)(1)進(jìn)行劃分.上述定義的加權(quán)鄰接矩陣(15)是無權(quán)重鄰接矩陣[14]的擴(kuò)展,其主要區(qū)別是:a)加權(quán)鄰接矩陣考慮了兩個節(jié)點之間連接的權(quán)重,而文獻(xiàn)[14]只考慮了兩個節(jié)點之間的連接度;b)本文方法中,考慮了兩個節(jié)點之間的最短路徑,因此,鄰接矩陣中的0元素表征了兩個節(jié)點之間沒有連接.

3.4 尋找最大模塊度

模塊度Qw既考慮了靈敏度,又考慮了邊的數(shù)量.子系統(tǒng)劃分問題等價于通過尋找最大化的模塊化Qw進(jìn)行社區(qū)結(jié)構(gòu)發(fā)現(xiàn).首先,需要進(jìn)行初步的群組配置,群組結(jié)構(gòu)的初始化基于如下假設(shè):i)每個群組至少有一個測量輸出;ii)應(yīng)將直接影響yi的系統(tǒng)狀態(tài)分配給屬于yi同一個群組.初始化群組結(jié)構(gòu)步驟如下:

步驟1定義關(guān)于復(fù)雜系統(tǒng)的群組結(jié)構(gòu)sw(0)=[cx1(0) ··· cxnx(0) cy1(0) ··· cyny(0)],其中cη(0)表示每個節(jié)點xi(i=1,···,nx)或yj(j=1,···,ny)的群體標(biāo)簽.應(yīng)用快速二分法(如文獻(xiàn)[18])尋找最大模塊度Qw時,在sw中的第i個變量被分配到第i個群體,i=1,···,na,即(0)=[1 ··· na].

步驟2對于每個狀態(tài)變量xi,i=1,···,nx,找到所有和狀態(tài)變量xi關(guān)聯(lián)的測量輸出,在ˉsw(0)中將相應(yīng)的元素標(biāo)為相同編號.

步驟3對于每個測量輸出yj,j=1,···,ny,找到所有和輸出方程hj(·)關(guān)聯(lián)的狀態(tài)變量,并將這些狀態(tài)標(biāo)量標(biāo)為和yj相同的編號.

步驟4所得到的分組(0)即可作為實施社區(qū)發(fā)現(xiàn)算法的初始群組.

得到的初始群組結(jié)構(gòu)將用于尋找子系統(tǒng)分解中模塊度的最大值.在所提出的方法中,子系統(tǒng)(即式(2)中的p)的數(shù)量應(yīng)預(yù)先設(shè)定.在本文中,快速二分法用于最大化模塊度[18],其步驟如下:

步驟1對每個節(jié)點i實施如下算法:

步驟1.1對于節(jié)點i的每個鄰居節(jié)點j,通過將節(jié)點i從當(dāng)前群組移動到節(jié)點j所在的群組,計算模塊度的變化值?Qw.

步驟1.2找到上一步中最大的?Qw>0,將頂點i放入到相應(yīng)的群組中(即將節(jié)點i編入最大?Qw相應(yīng)的分組),得到新的分組.

步驟2令nc(k+1)作為節(jié)點集結(jié)后的群組數(shù)量,如果nc(k+1)

3.5 分布式滾動時域估計設(shè)計

為對子系統(tǒng)的劃分結(jié)果進(jìn)行分析,本文將對所得到的子系統(tǒng)應(yīng)用于分布式滾動時域估計算法,并說明考慮加權(quán)有向圖的必要性.本文提出分布式滾動時域估計算法,局部估計器通過共享的通信網(wǎng)絡(luò)進(jìn)行信息交互.子系統(tǒng)i的局部估計器為求解如下優(yōu)化問題:

4 實例應(yīng)用

本節(jié)考慮將所提出的子系統(tǒng)劃分及分布式滾動時域估計算法應(yīng)用于反應(yīng)–分離過程系統(tǒng)[19],該過程包括兩個連續(xù)攪拌反應(yīng)釜和一個快速分離釜(見圖2).第一個連續(xù)攪拌反應(yīng)釜包含反應(yīng)物J,其給定蒸汽為F10,然后產(chǎn)出產(chǎn)品K.然后,產(chǎn)品K在第2個反應(yīng)釜中轉(zhuǎn)換成最終產(chǎn)品L.

圖2 反應(yīng)–分離過程系統(tǒng)Fig.2 Reactor-separator process

其中:xAi,xBi表示產(chǎn)品J,K在容器i(i=1,2,3)中的質(zhì)量分?jǐn)?shù);xC3表示產(chǎn)品L在反應(yīng)釜3中的質(zhì)量分?jǐn)?shù);Ti容器i(i=1,2,3)中液體的溫度;T10,T20表示輸入到容器1和2的蒸汽溫度;F1,F2表示物料從容器1和2流出的流速;F10,F20表示物料從容器1和2流出的穩(wěn)態(tài)流速;Fr,Fp分別表示產(chǎn)品的循環(huán)和廢棄流速;Vi表示容器i(i=1,2,3)的體積;E1和E2表示反應(yīng)過程1和2的催化劑;k1和k2表示反應(yīng)過程1和2的指前因子值;?H1和?H2表示反應(yīng)過程1和2的熱量;αA,αB,αC表示物質(zhì)J,K,L的相對揮發(fā)率;Q1,Q2,Q3分別表示容器i(i=1,2,3)的熱量輸入;Cp熱容量;R表示氣體常數(shù);ρ表示溶液濃度.

相關(guān)模型參數(shù)參見表1.

表1 模型參數(shù)Table 1 Model parameters

該反應(yīng)–分離過程的其中一個穩(wěn)定工作點為

為實現(xiàn)利用測量輸出Ti(i=1,2,3)對系統(tǒng)狀態(tài)xAi,xBi,Ti(i=1,2,3)的估計,首先對子系統(tǒng)進(jìn)行劃分,并和已有的子系統(tǒng)劃分方法(文獻(xiàn)[14])進(jìn)行比較,得到如表2所示的劃分結(jié)果.

表2 子系統(tǒng)劃分結(jié)果Table 2 Decompositions for the CSTR process

在應(yīng)用分布式估計算法時,如果進(jìn)行實時信息交互,分布式估計器的通信成本較高,因此希望子估計器之間可以在較少通信交互的同時獲得較好的估計效果.如圖3所示,應(yīng)用分布式滾動時域估計算法時(通訊間隔n=1),基于子系統(tǒng)劃分2的估計效果比基于子系統(tǒng)劃分1獲得了更好的效果,本文提出的基于加權(quán)有向圖的子系統(tǒng)劃分方法比未考慮連接權(quán)重的分解算法具有更小的估計誤差.

圖3 當(dāng)n=1時,實際狀態(tài)軌跡(藍(lán)色實線)、基于分解1的狀態(tài)估計(綠色虛線)和基于分解2的狀態(tài)估計(紅色點劃線)Fig.3 The trajectories for n=1 of the actual states(blue solid lines),estimated states based on Decomposition 1(green dashed-dot lines)and Decomposition 2(red dashed lines)

圖4 當(dāng)n=3時,實際狀態(tài)軌跡(藍(lán)色實線)、基于分解1的狀態(tài)估計(綠色虛線)和基于分解2的狀態(tài)估計(紅色點劃線)Fig.4 The trajectories for n=3 of the actual states(blue solid lines),estimated states based on Decomposition 1(green dashed-dot lines)and Decomposition 2(red dashed lines).

考慮不同的通信間隔n=1,2,3,對分布式滾動時域估計算法進(jìn)行測試.如表3和圖4(通信間隔n=3)所示,隨著通信間隔的增加,基于分解1的狀態(tài)估計效果受到較大的影響.而對于基于本文方法的分解2,由于子系統(tǒng)之間的關(guān)聯(lián)較弱,因此較大的通信間隔對估計性能的影響并不大.其現(xiàn)實意義是在保證估計效果的同時,可以減輕分布式狀態(tài)估計設(shè)計過程中子系統(tǒng)通信的次數(shù),降低通信成本.

表3 不同通信間隔對估計性能的影響Table 3 Summary of the estimation performance

5 總結(jié)

本文提出一種基于加權(quán)有向圖的社區(qū)發(fā)現(xiàn)子系統(tǒng)劃分方法,并應(yīng)用于分布式狀態(tài)估計設(shè)計中.基于耦合關(guān)聯(lián)的強(qiáng)度研究了一種系統(tǒng)化的復(fù)雜大系統(tǒng)分解方法,根據(jù)系統(tǒng)的動態(tài)模型自動構(gòu)建系統(tǒng)的加權(quán)有向圖,進(jìn)而獲得系統(tǒng)的劃分,得到的子系統(tǒng)內(nèi)部的關(guān)聯(lián)遠(yuǎn)強(qiáng)于子系統(tǒng)之間的關(guān)聯(lián).這類分解方能夠為在復(fù)雜系統(tǒng)的分布式狀態(tài)估計和控制奠定基礎(chǔ),同時也對研究復(fù)雜網(wǎng)絡(luò)的特性具有重要作用.未來可擴(kuò)展到帶系統(tǒng)輸入的復(fù)雜系統(tǒng)劃分方法,及其對分布式控制策略下的性能提升研究,通過分布式狀態(tài)估計和控制,對該子系統(tǒng)劃分方法進(jìn)行系統(tǒng)的評估.