應用IMOABC空調房間CO2濃度二自由度內模分數階PI控制性能的研究

楊 睿 李紹勇 魏明松 王少波

（1.蘭州理工大學土木工程學院 蘭州 730050；2.廣州市恒盛建設工程有限公司 廣州 510080）

0 引言

伴隨世界經濟快速增長，人們的生活質量得到了明顯改善。同時，對建筑環境的舒適度、空氣品質的關注度也與日俱增，使得中央空調系統在各類建筑環境中得到了廣泛應用。由于定風量空調系統（Constant Air Volume,CAV）具備結構簡單、初始投資低、易操作和送風量穩定等特點，在高大空間和室內人員密度大的建筑場所，例如大型會議室、商場和辦公大廳等，廣泛應用[1]。對于CAV系統中的關鍵設備CAVAHU而言，新風量qoa、回風量qra以及排風量qea經空調工藝設計后，往往保持不變[2]。當CAVAHU實際運行，常出現空調房間內的額定人員超員或不足的狀況，導致qoa不夠或過量，產生室內CO2濃度超標（Cn＞Cn,set）或處理新風所消耗的冷、熱負荷增加的問題，使得室內空氣品質下降或增大CAVAHU運行能耗。

室內CO2濃度變化呈現非線性、多干擾等特性，采用傳統控制方式很難得到一個滿意的控制效果。因此，室內CO2濃度的控制問題，是當前暖通空調自動化領域的一個有研究價值的課題。Kusuda[3]針對傳統的機械通風系統，首次提出基于室內CO2平均濃度的需求控制通風策略，將房間CO2平均濃度作為調節新風量大小的控制指標。Igor 等[4]在通風系統上設計了室內CO2濃度內模控制系統，對輸送到室內的新風量qoa進行調節。結果顯示，qoa能夠根據室內CO2濃度的變化快速做出響應。文獻[5]利用直接反饋線性化設計了基于室內CO2濃度變化的機械通風控制系統。結果表明，該控制系統對于室內CO2濃度控制的速度和穩定性等方面都優于傳統PID 控制系統。

分數階PID（Fractional Order PID Controller,FOPID）控制器比整數階PID 控制器[6]更具有設計自由度廣、自適應好和強魯棒性等特性[7]，在許多工業領域得到了應用，呈現出了更好的控制性能[8-10]。將內模控制與FOPID 控制相結合能使得控制系統獲得更好的設定值跟蹤和抗干擾等性能[11,12]。此外，控制器的參數整定對于控制系統的穩定性和快速性具有顯著意義。目前，許多研究工作都集中在控制器參數整定的單目標優化問題[9-13,14]。但在實際應用中，控制系統所輸出的被控參數的多個性能指標之間往往是相互制約的[15]。因此，解決控制器參數整定的多目標優化問題，顯得十分迫切。文獻[15]提出了一種基于改進多目標極值優化算法的自動穩壓器系統FOPID 控制器設計方法。通過積分絕對誤差，絕對穩態誤差和調節時間三個性能指標，建立多目標優化函數來整定FOPID 控制器參數的最優值。結果表明，該方法使自動穩壓器系統的控制精度和魯棒性方面有明顯提升。針對旋轉式倒立擺系統在穩定階段的控制問題，文獻[16]選擇擺角絕對積分時間誤差和轉角絕對積分時間誤差兩個目標函數，建立多目標優化問題，并設計了一種基于改進多目標粒子群算法來實施FOPID 控制器參數的整定。通過試驗驗證了該算法的有效性，使得倒立擺系統具有良好的穩定性和動態品質。

綜上分析，空調房間CO2濃度與輸配到每個空調房間的新風量密切相關。基于保持CAVAHU 輸出的總送風量qsa恒定，同時動態調節新風量qoa和新風比m，確保Cn=Cn,set的目標，本文提出一種室內CO2濃度二自由度內模分數階PI 控制策略，發揮控制器參數少，系統調節靈活，快速跟蹤和魯棒性強等優勢。此外，對于室內CO2濃度二自由度內模分數階PI 控制器參數的整定問題，依據人工蜂群算法，對雇傭蜂階段引入自適應慣性權重，同時對觀察蜂階段施加精英組策略，分別進行非線性遞減和柯西變異的演變，并結合觀察蜂搜索特性，將最小粒子角度[17]引入外部檔案集，獲取相應的Pareto 解集，設計得到IMOABCA。考慮調節時間tc、穩態誤差絕對值ESS和絕對積分時間誤差（Integrated Time Absolute Error,ITAE），選擇min（ITAE,ESS,tc）作為IMOABCA 目標函數，設計室內CO2濃度二自由度內模分數階PI 控制器參數整定算法。仿真結果表明：該參數整定算法可獲取3 個控制器參數和3 個目標函數的Pareto 最優解。從而提升室內CO2濃度調節品質，確保Cn=Cn,set和滿足室內空氣品質需求。同時，二自由度內模分數階PI 控制器作用下的CAVAHU，其輸出的新風量qoa是動態變化的，可適應空調房間內的額定人員超員或不足的實際情況。

1 CAVAHU 作用下的空調房間CO2濃度控制

對于新風+一次回風的空調方式，CAVAHU 輸送的總送風量qsa包括新風量qoa和一次回風量qra，即qsa=qoa+qra，新風比m=qoa/qsa。同時為了維持室內的微正壓，空調工藝的設計使得新風量等于排風量qea，即qoa=qea，且qsa＞qra。由上述可知，CAVAHU輸配的qoa對空調區域的CO2濃度Cn影響顯著。所以綜合考慮CAVAHU 輸出的qsa恒定和空調房間內的額定人員超員或不足的實際狀況，本文提出一種室內CO2濃度二自由度內模分數階PI 控制策略，能夠動態調節新風量qoa，實現Cn=Cn,set的目的。相應的空調工藝測控流程圖，如圖1所示。

圖1 CAVAHU 運行工藝測控流程圖Fig.1 Flowchart of measurement and control for CAVAHU operation process

當CAVAHU 運行空調冬/夏季工況時，設置在一次回風管道上CT 將反映n個空調房間平均CO2濃度測量值信號Cn上傳到CC，與經過CF 的室內CO2濃度設定值Cn,set進行比較，求偏差e=Cn,set-Cn，對其進行內模分數階PI 運算后，輸出3 路控制信號p1、p2和p3，分別到新風電動調節閥、一次回風電動調節閥和排風電動調節閥。這里，p1=p3，p2=pmax+pmin-p1。p1、p2、p3∈[pmin，pmax]，pmin，pmax-CC 輸出的控制信號的下限、上限。動態改變輸送到空調房間的新風量qoa，經過通風置換，實現Cn=Cn,set的目的，適應空調房間內的額定人員超員或不足的實際情況。同時CAVAHU 以變新風比（m＜100%）運行，輸出的總送風量qsa是恒定的；且維持室內的微正壓。這是由于當e＜0（＞0）時，接受控制指令p1和p3，新風閥開度v1和排風閥開度v3是同向增大（減小），即v1=v3=a%（某個開度值），qoa=qea；而接受控制指令p2，回風閥開度v2（100%-a%）是反向減小（增大），即v1+v2=100%，qoa+qra=qsa。

這樣，該室內CO2濃度二自由度內模分數階PI 控制策略，可以實現CAVAHU 對室內CO2濃度的實時調整，既動態調節新風量qoa，實現Cn=Cn,set的目的；又使得CAVAHU 以輸出的qsa是恒定的和變新風比（m＜100%）運行。

2 室內CO2濃度二自由度內模分數階PI 控制系統

2.1 空調房間CO2 濃度的數學模型

分析文獻[5]，可得空調房間CO2濃度質量平衡方程如下：

式中，qoa為新風量，m3/s；Cw和Cn為室外和室內CO2濃度，%；mn為室內CO2發生量，m3/s；V為房間體積，m3。

由于空調CO2濃度變化是非線性的，利用式（1）對空調房間CO2濃度模型在平衡點（穩定點）建立線性過程，得到相應的傳遞函數為[4]：

式中，a=-(qoa*/V)，b=(Cw-Cn*)/V。

這里，qoa*、Cn*-qoa、Cn的平衡點。兩者的關系如下：

式中，m*-qoa*對應的新風比。

2.2 空調房間CO2 濃度測量變送器傳遞函數

室內CO2濃度測量變送器的作用是實時測量空調房間內的CO2濃度Cn大小，輸送標準信號到CC控制器與Cn,set進行比較，其傳遞函數如下式[18]：

2.3 新風量、排風量與一次回風量電動調節閥傳遞函數

由于新風量qoa可直接影響房間CO2濃度Cn，所以設置新風量電動調節閥，接受輸出的控制指令p1（0～10mA·DC 或4～20mA·DC），通過改變其開度v1來控制qoa大小。對于排風量qea和一次回風量qra，由第1 節可知qea=qoa，qoa+qra=qsa，且p1=p3，p2=pmax+pmin-p1。所以，三個風閥本文均選擇直線流量特性的電動調節閥，其可看作放大環節處理[19]：

2.4 室內CO2濃度二自由度內模分數階PI 控制器傳遞函數

從上文可知，新風量電動調節閥與空調房間所構成的廣義被控對象為一階慣性數學模型：

式中，K=10b/-a；T=b/-a。

分析文獻[11]可以得出，內模分數階PI 控制器的分數特性是通過施加的分數階濾波器體現的，其傳遞函數如下：

式中，T1、α為分數階濾波器的時間常數和階次的非整數項（0＜α＜1）。

從上式可以看出，內模分數階PI 控制器是由一個整數階PI 控制器與一個分數階積分器串聯而成，其中KP=T/KT1，KI=1/KT1。

此外，為了提升參數控制系統對設定值的跟蹤特性，加快響應速度和降低超調量，引入二自由度內模控制理論[20]，構建設定值前饋濾波器，其傳遞函數如下：

式中，T2為一階低通濾波器的時間常數；T1、α同上。

所以，考慮內模分數階PI 控制器和設定值前饋濾波器的優勢，本文提出的室內CO2濃度二自由度內模分數階PI 控制器，其結構如圖2所示。

圖2 室內CO2 濃度二自由度內模分數階PI 控制器結構Fig.2 Block diagram of structure for two-degree of freedom internal model fractional order PI controller for the indoor CO2 concentration

如上所述，該控制器可輸出3 路控制信號p1、p2和p3，分別作用于新風量、一次回風量和排風量電動調節閥，在保證送風量qsa不變的情況下，動態調節新風量qoa，實現Cn=Cn,set的目的。

因此，本文提出的基于IMOABCA 空調房間CO2濃度二自由度內模分數階PI 控制系統，如圖3所示。

圖3 基于IMOABCA 空調房間CO2 濃度二自由度內模分數階PI 控制系統方框圖Fig.3 Block diagram of two degree of freedom internal model fractional order PI control system based on IMOABCA for the concentration of CO2 in air-conditioning room

圖3為二自由度內模分數階PI 控制器作用下的空調房間CO2濃度控制系統方框圖，并采用IMOABCA 對控制器參數進行整定。考量ITAE，ESS和tc，確定min（ITAE,ESStc）作為IMOABCA的目標函數，進而對控制器的三個參數[T1,T2,α]進行整定，獲取相應的Pareto 最優解：[T1*,T2*,α*]和[ITAE*,ESS*tc*]，提升調節質量，使得Cn=Cn,set，滿足室內空氣品質要求；且動態調節新風量qoa，減少CAVAHU 處理新風所消耗的冷、熱負荷。

3 控制器參數整定

3.1 多目標人工蜂群算法改進方法

多目標優化問題是指當在解決一個問題時，需要對多個目標進行最大或最小優化，但這幾個目標都是相互制約，相互矛盾的。因此，在多目標優化過程中會在多個目標問題之間尋求一個平衡值，從而得出這些目標函數的Pareto 解集[17]。

人工蜂群算法是一種簡單、高效的群體仿生優化算法，具有收斂快、魯棒性強等特點[21]。但由于采用隨機選擇鄰居的策略，導致Akbari 等[23]所設計的MOABCA在雇傭蜂階段的前期全局搜索能力弱，后期局部探索能力不足。因此，加入余弦自適應權重因子w，利用余弦函數在[π/2,π]呈現單調遞減變化，使得雇傭蜂：（1）前期具有較強的全局搜索能力。（2）隨著算法迭代次數的增加，后期也具有充足的局部探索能力。這樣，雇傭蜂在全局搜索任務中的能力得以提升。相應的w及改進后的雇傭蜂搜索公式如下：

式中，it、itmax為算法的當前迭代次數和最大迭代次數；φi,j為[-1.1]之間的隨機數；xi,j為當前蜜源位置；xk,j為隨機蜜源位置；vi,j為新蜜源位置。

當觀察蜂對較好蜜源進行挖掘時，MOABCA往往存在陷入局部最優或早熟的問題。因此，引入柯西算子Cauchy對其步長進行變異，使得觀察蜂：（1）當出現陷入局部最優的趨向時，利用較長的步長跳出局部極值。（2）當在最優解附近挖掘時，利用較小的步長加速收斂。標準的柯西分布的概率密度表達式如下：

與此同時，本文對MOABCA 的外部檔案（External Files,EF）的網格方法應用最小粒子角度[17]，并基于精英組策略，利用觀察蜂“優中尋優”的特點對網格中的Pareto 解集進行優化，進而提升Pareto 解的質量。為了保證算法在搜索過程中的多樣性，加入控制參數MR[22]進行修正。改進后的網格結構（基于兩目標函數優化）和觀察蜂搜索公式分別如圖4和公式（12）所示：

圖4 基于最小粒子角度的網格結構圖Fig.4 Grid structure based on minimal particle angles

式中，xi,j為雇傭蜂發送的優質蜜源；leader為xi,j對應EF 中的最優蜜源；repk,j為EF 中隨機生成的蜜源；Cauchy為柯西算子；MR為控制參數，一般取0.5[22]。

分析圖4，可知觀察蜂從雇傭蜂搜索的蜜源中選擇較優的蜜源，進行了開采。與此同時，通過計算所選的開采蜜源與網格中的蜜源角度，選取EF中與其蜜源角度最小的蜜源，作為最優蜜源進行柯西變異，達到優化EF 中的Pareto 解集目的。

因此，基于改進的網格結構以及搜索公式，本文在MOABCA 的基礎上構建出了IMOABCA。

3.2 室內CO2濃度二自由度內模分數階PI 控制器參數整定算法

本節選擇調節時間tc、穩態誤差絕對值ESS和絕對積分時間誤差ITAE[9]作為IMOABCA 的三目標函數，即F(x)=minJ=min(ITAE,ESS,tc)，來設計室內CO2濃度二自由度內模分數階PI 控制器參數整定算法。基于IMOABCA 改進內容，相應運算流程如圖5所示。

圖5 控制器參數整定算法流程Fig.5 Parameter tuning algorithm flow of the controller

Step1：設置MALAB/Simlink 工具組態中的室內CO2濃度二自由度內模分數階PI 控制器參數變量[T1,T2,α]，并設定T1、T2和α的上限、下限。

Step2：初始化：設置決策空間維數D，種群規模N，雇傭蜂數目Ne=觀察蜂數目Ns=N，閾值limit，外部檔案規模nEF，最大迭代次數itmax以及多目標函數F(x)=minJ。

Step3：輸入單位閉環負反饋室內CO2濃度二自由度內模分數階PI 控制系統的反饋信號e(t)，基于min(ITAE,ESS,tc)，計算每個初始蜜源的適應度值F，完成網格非支配排序。

Step4：進入雇傭蜂階段，根據公式（10）生成新蜜源vi,j，計算目標函數值F，并與舊蜜源進行支配關系判斷：若新蜜源支配舊蜜源，則替代舊蜜源，成為非劣解，雇傭蜂階段結束。非劣解被加入EF 中進行網格非支配排列，并刪除網格中被支配蜜源。

Step5：利用輪盤賭方式，在雇傭蜂階段獲取的蜜源F中，為觀察蜂選擇開采對象。通過計算開采對象與EF 中的蜜源之間的角度，選取角度最小的蜜源作為最優蜜源。將開采對象和最優蜜源同時發送給觀察蜂。

Step6：進入觀察蜂階段，根據公式（12）生成新蜜源vi,j，計算目標函數值F，并與收到的開采對象、最優蜜源進行支配關系判斷：若新蜜源支配開采對象或最優蜜源，則替代它，成為非劣解，觀察蜂階段結束。非劣解被加入EF 中進行網格非支配排列，并刪除被支配蜜源。

Step7：進入偵察蜂階段，若某蜜源的開采次數超過了閾值limit，則隨機生成新蜜源替代舊蜜源，并計算目標函數值F。反之，則繼續保留舊蜜源。

Step8：若滿足it=itmax，輸出控制器參數最優值[T1*,T2*,α*]，ITAE、tc和ESS的Pareto 最優解minJ*以及相應的閉環負反饋控制系統動態響應曲線。否則，返回Step4。

需要說明的是，借助MATLAB 軟件，將上述的控制器參數整定流程進行編程，保存為一個獨立程序，命名為IMOABCA.m。該IMOABCA 是與后面章節中的控制系統組態Simulink 模型同步運行的。

4 數值模擬

4.1 室內CO2濃度二自由度內模分數階PI 控制系統的組態與仿真

對于蘭州市一個大型辦公區域，其容積V=1200m3，額定人員90 人，人員的活動強度為輕度，人均CO2發生率為0.0173m3/h。根據GB/T 18883-2022 標準，室內CO2濃度不應超過0.1%[25]，故本文選擇室內CO2濃度設定值Cn,set=0.1%，相應的CAVAHU 總送風量qsa=4.1m3/s。對該CAVAHU配置本文設計的室內CO2濃度二自由度內模分數階PI控制系統。基于圖3和MATLAB中的Simulink工具，對該控制系統進行組態，相應的Simulink模型如圖6所示。

圖6 基于IMOABCA 的室內CO2 濃度二自由度內模分數階PI 控制系統組態圖Fig.6 Configuration diagram of two-degree of freedom internal model fractional order PI control system for the indoor CO2 concentration based on IMOABCA

根據空調工藝要求[25]，Cn,set=0.1%，假定室內CO2濃度初始值Cn,0=0.11%，IMOABCA.m 中的相關參數設置如下：N=100，nEF=3，控制器參數整定的數目D=3，itmax=200，limit=300，w∈[0,1]，T1∈[90,110]，T2∈[10,12]，α∈[0.4,0.6]，仿真時間設為500s。

在MATLAB 中的Command Window 界面和Simulink 環境，同步運行IMOABCA.m 和圖6中的Simulink 模型，可獲得室內CO2濃度的動態響應Cn(t)，min(ITAE,ESS,tc)Pareto 前沿變化和3 組控制器參數Pareto 最優值[T1*,T2*,α*]，分別如圖7，圖8和表1所示。

圖7 三組控制器參數Pareto 最優解下的室內CO2 濃度動態響應Fig.7 Dynamic responses of indoor CO2 concentration under three Pareto optimal solutions of parameters of the controller

圖8 min（ITAE，ESS，tc）Pareto 前沿變化Fig.8 Pareto front variation for min(ITAE,ESS,tc)

分析圖7可知，基于IMOABCA 得到的控制器三組最優參數均能使得Cn=Cn,set，說明IMOABCA 對控制器的參數整定是可行的。從局部放大圖可見，輸出三組Cn(t)在超調量，調節時間表現出了不同的差異，也與表1中的三個目標函數值所表征的相呼應。

表1 控制器參數及其三目標函數的Pareto 解Table 1 Pareto solution of parameters of the controller and three objective functions

分析圖8，容易觀察到3 個目標函數值所對應的不同控制器參數的Pareto 解之間的聯系，曲面上的控制器參數Pareto 解的分布位置對應了不同的Cn(t)控制特征。

分析表1，可知控制器參數的三組Pareto 最優解所對應的三個目標函數值之間是相互制約或矛盾的。例如，第1 組中的ITAE 最大，但tc和ESS最小；第3 組中的ITAE 最小，但ESS和tc最大。由此說明IMOABCA 可以使得室內CO2濃度二自由度內模分數階PI 控制器獲得側重點不同的Cn(t)控制效果，以滿足不同控制需求。

相應的新風量qoa、一次回風量qra變化和新風比m變化，如圖9所示。

圖9 三組Pareto 控制器參數最優解下的響應變化Fig.9 Response changes under three Pareto optimal solutions of the controller

分析圖9的（a）和（b），控制器的3 組Pareto參數最優解均能夠使控制器根據室內CO2濃度的變化快速做出反應，發送控制信號作用到新風量、一次回風量和排風量電動調節閥完成新風量的調節，達到Cn=Cn,set的目的。圖（c）反映了控制器的3 組Pareto 參數最優解的新風比m變化過程，計算得到3 組平均新風比均約為14%，滿足m≧10%[25]空調工藝要求。

4.2 與現有算法的比對

此外，為了說明本文IMOABCA 的優越性，分別將IMOABCA 和MOABCA 作用于室內CO2濃度二自由度內模分數階PI 控制系統對控制器參數進行整定。可獲得其控制器參數Pareto 最優解[T1*,T2*,α*]及其三個目標函數值如表2所示。

表2 控制器參數及其相應的三目標函數的Pareto 解Table 2 Pareto solutions of parameters of controllers and the corresponding three objective functions

相應的室內CO2濃度的動態響應Cn(t)，如圖10所示。

分析表2，從min(ITAE,ESS,tc)Pareto 最優值不難看出，相比MOABCA，IMOABCA 對本文所設計的控制系統表現出了更好的收斂性和控制精度。同時，從圖10可以看出，IMOABCA 作用下的控制系統具有更快的響應速度和更小的超調量，且能夠很快的達到穩定狀態，說明IMOABCA 在解決多目標優化問題上更具有優越性，可以使系統呈現較強的控制性能，也表明IMOABCA 的設計是可行、有效的。

圖10 兩算法下的室內CO2 濃度動態響應Fig.10 Dynamic responses of indoor CO2 concentration under two algorithms

5 結論

針對CAVAHU 作用下的空調房間CO2濃度的控制問題，本文提出了一種基于IMOABCA 的二自由度內模分數階PI 控制策略。依據人工蜂群算法，對雇傭蜂階段引入自適應慣性權重，同時對觀察蜂階段施加精英組策略，分別進行非線性遞減和柯西變異的演變，并結合觀察蜂搜索特性，將最小粒子角度引入外部檔案集，獲取相應的Pareto 解集，設計了IMOABCA，進而對室內CO2濃度二自由度內模分數階PI 控制器的3 個參數[T1,T2,α]進行整定。結果表明：（1）該IMOABCA 是可行的，可獲取控制器的3 個參數相應的最優解[T1*,T2*,α*]，從而提升室內CO2濃度的調節品質，實現Cn=Cn,set的目的。（2）該IMOABCA 可以使得控制器獲得側重點不同的室內CO2濃度控制效果，以滿足不同控制需求。（3）與現有算法相比，IMOABCA對本文所設計的系統呈現出了明顯的優越性，表明了IMOABCA 的設計是可行、有效的。此外，CAVAHU 是以變新風比m、保持總送風量qsa恒定的方式運行的，為CAVAHU 運行、管理提供了新方式。