基于帝國競爭算法和LQG的磁流變阻尼器-建筑結構系統(tǒng)的振動控制

林秀芳， 陳淑梅, 黃國榮(.福州大學 機械工程及自動化學院，福州 3506;2.福建農林大學 金山學院，福州 350002)

提高建筑結構在地震等外部擾動下的抗震性以及減輕動力作用危害是土木結構工程領域面臨的重大課題之一。近些年，作為一種智能的半主動控制裝置，磁流變阻尼器已經在振動領域獲得越來越多的關注[1-3]?？刂扑惴ㄊ前胫鲃涌刂葡到y(tǒng)的核心，它直接決定了控制系統(tǒng)的控制效果。現(xiàn)有的利用MRD進行振動控制的算法主要包括H∞控制、滑模控制[4]、最優(yōu)控制[5]、模糊控制[6]等。其中，LQG控制是一種適用性很強的最優(yōu)控制，其最大的優(yōu)點是可得到狀態(tài)線性反饋的最優(yōu)控制律。又因其性能指標物理意義明確，最優(yōu)解可以解析表達，因此被廣泛運用于振動控制中。

LQG控制算法的設計是在加權矩陣選定的基礎上進行的，其設計效果直接取決于加權矩陣的選取。如果加權矩陣選取不當則無法得到滿意的設計結果。其結果會造成控制裝置不能發(fā)揮最優(yōu)性能和資源的浪費，甚至對結構響應產生不利的影響。因此，加權矩陣的選取是LQG設計中的一個關鍵問題。但LQG控制系統(tǒng)的加權矩陣的設計過程往往依賴于人工經驗和領域專家的先驗知識，而且人工設計的效率很低且不能保證系統(tǒng)達到最優(yōu)。近些年，隨著智能優(yōu)化算法的興起，用其對線性二次型控制系統(tǒng)的加權矩陣進行優(yōu)化選擇已成為熱點。例如，宋春生等[7]利用遺傳算法(GA)優(yōu)化最優(yōu)控制的加權矩陣，實現(xiàn)對復雜雙層磁懸浮精密系統(tǒng)的隔振控制。結果表明優(yōu)化后的控制系統(tǒng)的隔振效果有顯著提高。景益娟[8]在設計大型雷達天線的抗風擾伺服控制系統(tǒng)中利用人工蜂群算法優(yōu)化LQG控制器的加權矩陣。仿真結果表明優(yōu)化過的LQG控制器能夠顯著地提高風擾情況下的跟蹤精度。王磊等[9]利用人工魚群算法求解基于主動質量阻尼器的隨機結構的控制系統(tǒng)中加權矩陣的優(yōu)化問題。

帝國競爭算法(ICA)由Atashpaz-Gargari等[10]于2007年提出，是一種富有生命力和高性能的群智能優(yōu)化算法。該算法與多種智能優(yōu)化算法相比都表現(xiàn)出較高的優(yōu)化效率和魯棒性，因此近年來已被越來越多的學者關注和應用。例如，Kaveh等[11]提出基于ICA的骨架結構的優(yōu)化設計，與其它多種優(yōu)化算法(如粒子群算法、蟻群算法等)比較的結果表明了該算法的優(yōu)越性。Abd-Elazim等[12]利用ICA優(yōu)化多機環(huán)境下的靜止同步補償器，仿真結果表明該算法的性能超過了GA。Shabani等[13]為了實現(xiàn)對載荷頻率的魯棒控制，提出了基于ICA優(yōu)化的PID控制器。仿真結果表明其性能由于基于基于GA優(yōu)化的PID控制器和基于神經網絡優(yōu)化的PID控制器。

作為一種近幾年來新興的智能優(yōu)化算法，雖然ICA已在多領域得到成功應用，但尚未發(fā)現(xiàn)其在LQG控制器優(yōu)化這方面的文獻。本文基于優(yōu)化設計的思想，在選中某種優(yōu)化指標的基礎上，嘗試采用ICA對LQG的加權矩陣進行優(yōu)化選取，從而獲得結構振動控制所需的理想控制力。此外，由于磁流變阻尼器的力不可控，只有電壓可控，所以本文還利用剪枝最優(yōu)控制(COC)算法將理想控制力轉化成MRD所需的控制電流，從而實現(xiàn)基于MRD的建筑結構振動的半主動控制。

1 基于LQG控制的受控系統(tǒng)的數(shù)學模型

LQG控制算法實際上是輸出反饋控制算法，是由最優(yōu)狀態(tài)反饋器(Linear Quadratic Regulator, LQR)和基于Kalman濾波的狀態(tài)觀測器組成。

1.1 基于LQR控制的受控系統(tǒng)的數(shù)學模型

在一個具有n個自由度的結構上安裝p個MRD，則基于MRD的受控結構的運動方程可以表示為

(1)

(2)

式中，

y(t)=Cz(t)+D1f(t)

(3)

其中，

綜合式(2)和式(3)，則相應于運動方程(1)的控制系統(tǒng)的狀態(tài)方程可以用標準狀態(tài)方程表示為

(4)

由于位移響應是與安全性有關的參數(shù)，加速度響應是與舒適性有關的參數(shù)，為了同時兼顧建筑物的安全性和其中人員的舒適性，同時限制控制力的幅值，控制問題的目標性能函數(shù)定義為

(5)

將目標性能函數(shù)重新改寫為

(6)

式中：Q為狀態(tài)變量的加權矩陣，R為控制變量的加權矩陣，這兩個矩陣中的參數(shù)將由ICA優(yōu)化確定。根據(jù)LQR最優(yōu)控制率，最終得到控制器的反饋為

(7)

式中：G是最優(yōu)反饋增益，P由如下的Riccati方程求出

-PA-ATP+PBR-1BTP-Q=0

(8)

1.2 Kalman濾波器

LQR控制器需要全部的狀態(tài)變量，但在實際控制中，一般只有加速度能夠被簡單地測得。因此，需要從可測量的信號中估計出所有的狀態(tài)變量，這里釆用Kalman濾波器估計系統(tǒng)的狀態(tài)。Kalman濾波器的性能函數(shù)定義如下

(9)

(10)

(11)

借助式(11)求解得到濾波器的增益

(12)

其中，Qe=qe，Re=re×eye(n,n),eye(n,n)為一個n×n的單位矩陣。最終，系統(tǒng)的最優(yōu)控制力可表述為狀態(tài)估計值的反饋

(13)

2 基于ICA優(yōu)化的LQG控制

2.1 ICA的優(yōu)化原理

ICA是在模擬人類社會發(fā)展過程中帝國競爭過程的基礎上提出的一個新穎的元啟發(fā)式算法。該算法通過國家間的競爭來實現(xiàn)優(yōu)化，“國家”為群體中的個體，每個“國家”代表待優(yōu)化問題的一個解，即每一個“國家”是一組加權矩陣Q和R中的待優(yōu)化參數(shù)的集合。該算法的基本流程包括：帝國初始化、帝國的同化與更新、帝國集團間的競爭和帝國的滅亡。

(1) 帝國初始化：在一個N維的優(yōu)化問題中，一個國家即為一個1×N的矩陣，該矩陣定義為

country=[p1,p2,p3,…,pN]

(14)

式中:p1,p2,p3,…pN為待優(yōu)化的變量。每個國家的實力大小通過代價函數(shù)來衡量

cost=f(country)=f(p1,p2,p3,…,pN)

(15)

代價函數(shù)值越小，國家實力越大。首先，設置算法開始時國家的數(shù)量Npop，選取其中實力最強的Nimp個國家作為帝國主義國家，剩下的Ncol個國家作為帝國主義國家的殖民地。然后，為了形成初始帝國集團，每個帝國主義國家根據(jù)實力大小，按照輪盤賭選擇法獲得相應的殖民地。

(2) 帝國的同化與更新：每一個帝國為了壯大實力，會盡可能地同化所占領的殖民地，使殖民地和帝國相融合。這個過程可以表達為殖民地在二維平面內沿著坐標軸向帝國隨機移動，如圖1所示。x表示殖民地向帝國移動的單元，它是一個隨機數(shù)，服從均勻分布x～U(0,β×d)。d表示殖民地與帝國之間的距離，β是大于1的系數(shù)，表示殖民地從兩邊向帝國靠近。但帝國的同化過程并不一定會使殖民地直接向帝國移動，因此引入一個隨機角度θ模擬該運動，θ是一個服從θ～U(-r,r)的均勻分布，r用于調整偏移初始方向的角度。一般取β=2，r=π/4。

當殖民地向帝國移動的過程中，可能出現(xiàn)殖民地的實力大于帝國的實力，在這種條件下，帝國和殖民地角色互換，即殖民地成為帝國，原先的帝國淪落為殖民地。

圖1 殖民地向其帝國移動示意圖

(3) 帝國集團間的競爭：出于本性每個帝國總是希望盡可能多地占領殖民地，因此不同帝國間不可避免地產生競爭。每一個帝國集團的實力主要由帝國的實力所決定，但殖民地的實力也對帝國集團有所影響。所以，將衡量一個帝國集團的實力大小的代價函數(shù)歸一化為

(16)

其中，ξ的大小決定了殖民地國家對整個帝國實力的影響程度。Mean表示該帝國所有殖民地代價函數(shù)值的平均值。每個帝國集團都試圖兼并其它帝國集團的殖民地。在帝國競爭過程中，每個集團都有占領最弱集團殖民地的可能性。

(4) 帝國滅亡及算法結束：實力弱的帝國經過帝國競爭后逐漸失去殖民地。當一個帝國失去所有殖民地時，該帝國滅亡。最終只剩一個帝國集團時，算法結束，輸出最優(yōu)解；否則，返回(2)。

2.2 優(yōu)化目標

為了同時兼顧建筑物的安全性和其中人員的舒適性，ICA優(yōu)化的目標函數(shù)定義為

obj=wJ1+(1-w)J2

(17)

其中

(18)

(19)

2.3 優(yōu)化流程

用ICA優(yōu)化LQG加權函數(shù)的參數(shù)的具體步驟如下。

(1) 根據(jù)受控對象的屬性和LQG的控制目標，確定加權矩陣Q和R的結構以及矩陣中待優(yōu)化參數(shù)的數(shù)量。

(2) 確定待優(yōu)化參數(shù)的范圍和優(yōu)化算法的編碼形式，本文采用實值編碼。

(3) 初始化算法和帝國：初始化算法的參數(shù)，如國家的數(shù)量Npop，帝國的數(shù)量Nimp以及迭代的次數(shù)k等。然后，隨機產生初始國家，若國家(解)對應的理想控制力超出MRD的最大量程，令代價函數(shù)值為1，否則，根據(jù)式(18)計算代價函數(shù)值。最后根據(jù)代價函數(shù)值的大小分配帝國和殖民地。

(4) 迭代開始：根據(jù)2.1節(jié)的第(2)步進行每個帝國的同化與更新，然后按照上一步的方法計算每個帝國集團的代價函數(shù)值，最后根據(jù)2.1節(jié)的第(3)步進行帝國間的競爭。

(5) 判斷是否滿足收斂標準：如果滿足，則算法結束，輸出最優(yōu)解，否則，返回(4)。

3 ICA-LQG-COC半主動控制

要實現(xiàn)對MRD的應用，通常需要兩個步驟：首先根據(jù)結構響應計算理想控制力，然后將理想力轉換為MRD所需的輸入信號。本文所提的ICA-LQG-COC半主動控制系統(tǒng)是由ICA-LQG控制器和COC這兩部分組成，如圖2所示。圖中，ICA-LQG控制器根據(jù)建筑物的響應計算出理想控制力，COC通過調整輸入電流使得MRD能夠跟蹤理想控制力。需要注意的是圖2中的MRD在仿真中由動力學模型替代。

圖2 本文所提的ICA-LQG-COC策略

3.1 MRD的動力學模型

本文使用的是Lord公司生產的RD-8040-1型MRD。該阻尼器屬于剪切閥式單出桿型，其結構原理如圖3所示。當活塞桿在缸體中往復運動時，通過控制輸入線圈電流的大小來改變環(huán)形孔周圍的磁場強度，從而改變磁流變液的流變特性，產生不同的阻尼力。對MRD進行性能測試所使用的是INSTRON公司生產的INSTRON1342電液伺服疲勞試驗機，圖4為MRD的測試系統(tǒng)。測試時通過Lord公司生產的RD-300-2-03型電流控制器調節(jié)MRD輸入電流的大小。采用正弦位移激勵對MRD進行拉壓測試，當頻率取2.4 Hz時，幅值分別取3 mm、5 mm、6.5 mm、10 mm；當頻率取1.6 Hz時，幅值分別取2.17 mm、3.5 mm、7.5 mm。測試過程中提供的輸入電流大小在0 A到1 A之間變化，變化間隔為0.2 A。

圖3 磁流變阻尼器的結構原理圖

通過采集MRD在不同激勵性質及控制電流下的力學性能數(shù)據(jù)，分析其動力學特性，從而進行MRD動力學模型的辨識及驗證。Bouc-Wen模型因其精度較高，數(shù)字處理方便簡單而得到較為廣泛的應用[14-15]。因此，本文用該模型表征MRD的動力學特性，模型中的阻尼力f表示為

圖4 磁流變阻尼器測試系統(tǒng)

(20)

(21)

式中:x是活塞的相對位移；c0和k0代表高速時的黏性阻尼和剛度；x0是補償器氣體彈簧的初始位移;α是滯環(huán)算子的比例因子;γ,β,A和n是調整滯環(huán)的參數(shù)，其中n通常取2。

為了簡化參數(shù)識別的工作，需要讓每個待識別參數(shù)與電流產生關系。在對實驗數(shù)據(jù)進行參數(shù)擬合的過程中，為了盡可能減小預測阻尼力和實驗阻尼力之間的誤差，待識別參數(shù)與電流的關系表征如下

c0=c01+c02i

(22)

k0=k01+k02i

(23)

x0=x01+x02i

(24)

α=α1+α2i+α3i2

(25)

γ=γ1e(i/γ2)+γ3

(26)

β=β1e(i/β2)+β3

(27)

A=A1+A2i+A3i2

(28)

由式(22)～式(28)可見，總共有18個參數(shù)需要識別。使用MATLAB中的遺傳算法工具箱，取激勵頻率1.6 Hz、激勵幅值7.5 mm的實驗數(shù)據(jù)進行參數(shù)識別，得到的Bouc-Wen模型的具體參數(shù)數(shù)值見表1所示。在此基礎上，圖5顯示了預測阻尼力和實驗阻尼力的對比結果，從圖中可看出模型預測結果與實驗結果非常吻合，說明識別到的Bouc-Wen模型可以準確地表征該MRD的動力學特性。

3.2 剪枝最優(yōu)控制(Clipped-optimal control, COC)

MRD這種控制裝置具有很強的非線性，其阻尼力不能直接控制。為此，本文采用Spencer等提出的COC方法將理想控制力轉化為MRD的控制電流，控制電流i表達如下

表1 Bouc-Wen動力學模型的參數(shù)值

(a) 力-位移圖

(b) 力-速度圖

Fig.5 Comparison of the Bouc-Wen model predicted results and experimental results (1.6 Hz, ±7.5 mm)

i=ImaxH{(fc-f)f}

(29)

式中：Imax是MRD的最大電流，H{}是階躍函數(shù)，fc是理想控制力，f是MRD的控制力。由上式可知，COC的核心控制率是：當f=fc時，控制電流值維持原值；當f

4 數(shù)值計算及分析

本節(jié)以一個三層的建筑物結構為例，分別用三種地震波對其進行激勵，為驗證ICA-LQG-COC方法的有效性，將其與被動控制(0 A)、被動控制(2 A)以及LQG-COC半主動控制進行比較。

4.1 仿真條件

基于地震波激勵的安裝了MRD的框架結構如圖6所示，MRD安裝在在第一層和地面之間。

圖6 基于磁流變阻尼器減震的三層結構

該結構的質量M和剛度矩陣K分別為

阻尼矩陣的計算公式如下

C=αM+βK

(30)

式中:α和β分別為質量阻尼系數(shù)和剛度阻尼系數(shù)。這兩個系數(shù)可通過振型阻尼比計算得到，即：

(31)

式中:ωi和ωj分別為結構的第i和第j階固有頻率；ξi和ξj為相應于第i和第j振型的阻尼比。通常情況下取i=1，j=2。最終計算得到阻尼矩陣為

根據(jù)式(5)和(6)定義的目標性能函數(shù)，LQG控制器的反饋增益G的加權矩陣Q和R表達如下

R=r

其中，q1、q2和r是加權矩陣的待優(yōu)化參數(shù)，這三個參數(shù)的優(yōu)化范圍如下

q1∈[1,105],q2∈[10-5,1],r∈[10-7,10-2]

由于濾波器增益Ke的加權矩陣參數(shù)較易選取，且取值對結果影響較小，故將qe和re分別取為1和10-7。

ICA的具體參數(shù)為：國家的數(shù)量Npop=30、帝國的數(shù)量Nimp=3、迭代次數(shù)=100、β=2、r=π/4、ξ=0.1。令式(17)中的權重w=0.5。算法運行次數(shù)=5，取5次中的最優(yōu)解。至于如何更好地選擇ICA的參數(shù)才得到更優(yōu)化的解，超出了本文的研究范疇。此外，LQG-COC算法中加權矩陣的參數(shù)分別是：q1=10,q2=0.01,r=10-6。

4.2 仿真結果和討論

對該框架結構采用El-Centro地震波(加速度峰值為4.946 m/s2)進行激勵，采樣頻率為0.02 s，激勵時間為20 s。由于建筑物是一個縮尺的模型，按照結構相似度原則，把地震波加速度加速5倍作為該模型的激勵。經過ICA優(yōu)化確定的參數(shù)分別為q1=73 618,q2=0.604 3,r=4.380 3×10-7。

除了第一層以外，ICA-LQG-COC方法對其它樓層的控制效果都優(yōu)于被動控制(2 A)，但是被動控制(2 A)所需的最大控制力是所有控制方法中最大的。這說明即使將阻尼器的能量調到最大，不一定會得到更好的減震效果。從另一個角度而言，ICA-LQG-COC方法不僅控制效果好，而且更節(jié)能。這里補充一點，在優(yōu)化過程中應滿足的約束條件是，優(yōu)化得到的參數(shù)對應的理想控制力不能大于被動控制(2 A)時的最大控制力(1 876 N)。另外，雖然被動控制(0 A)也能產生一定的減震效果，但其產生的所有峰值減小率是所有方法中最小的。通過對這幾種方法的比較，得出的結論是：控制效果從好到差的排序依次為ICA-LQG-COC、被動控制(2 A)、LQG-COC以及被動控制(0 A)。

表2El-Centro地震波激勵下的響應峰值的比較

Tab.2PeakresponsesunderElCentroearthquakewithdifferentcontrolstrategies

響應無控被動0A被動2ALQG-COCICA-LQG-COCxi/cm0.2940.5610.7410.2590.4880.6310.1100.2950.4200.2200.3900.4880.2070.2440.2650.2940.2590.1100.2200.207di/cm0.2700.2290.2120.2120.1030.1810.1460.1250.1360.071x··i/(cm·s-2)83410011642627914132441081311335898961115566479563F/N33918768021542

為了更直觀地說明本文方法的減震效果，圖7和圖8分別比較了無控時和ICA-LQG-COC控制時第三層的位移時程響應和加速度時程響應。

圖7 ICA-LQG-COC控制和無控時的第三層位移響應的比較

Fig.7 Time histories of the 3rd floor displacement responses with ICA-LQG-COC strategies and without control

對這個框架結構分別采用加速度峰值調為3.5 m/s2的Tianjin地震波和Taft地震波激勵，利用上述各種方法進行控制，得到各樓層的響應峰值分別如圖9和圖10所示。類似于El-Centro波一樣，這兩個例子中的ICA-LQG的控制效果都不如被動控制(2 A)但都優(yōu)于被動控制(0 A)。此外，由圖9可見，除了第一層的位移和加速度，ICA-LQG-COC對其它樓層的所有波峰響應的控制效果都優(yōu)于被動控制(2 A)。由圖10可見，ICA-LQG-COC對最大位移和加速度響應的控制優(yōu)于被動控制(2 A)，只有對第一層的位移響應的控制不如被動控制(2 A)?？傊?，相比其它幾種控制方法，ICA-LQG-COC展示出最優(yōu)的綜合控制性能。

圖8 ICA-LQG-COC控制和無控時的第三層加速度響應的比較

Fig.8 Time histories of the 3rd floor acceleration responses with ICA-LQG-COC strategies and without control

(a) 位移峰值

(b) 加速度峰值

圖9 Tianjin波激勵時受控和無控下的各樓層波峰響應

Fig.9 Peak responses of each floor with different control strategies and without control under Tianjin earthquake

(a) 位移峰值

(b) 加速度峰值

圖10 Taft波激勵時受控和無控時的各樓層波峰響應

Fig.10 Peak responses of each floor with different control strategies and without control under Taft earthquake

5 結 論

本文設計了一種ICA和LQG控制方法相結合的新智能控制方法，并將其與MRD控制方法中的COC相結合，形成了利用MRD進行結構減震的ICA-LQG-COC方法。得到了如下結論：

(1) MRD的拉壓動態(tài)試驗表明：利用依賴電流的Bouc-Wen模型表征MRD的動力學模型，可以比較精確地描述該阻尼器的動力學特性。

(2) 以安裝了MRD的框架結構為例，仿真結果表明：相比LQG-COC半主動控制方法和其它兩種被動控制方法，本文提出的ICA-LQG-COC半主動控制方法在降低樓層的響應峰值時展現(xiàn)出更好的綜合性能。雖然優(yōu)化LQG加矩陣的目標是同時最大限度地降低最大位移和最大加速度，但ICA-LQG-COC方法還很有效地降低了其它樓層的位移和加速度響應，甚至層間位移響應。

(3) 利用MRD進行振動控制時不是控制力越大，控制效果就越好。當采用被動控制(2 A)方法時，雖然阻尼力始終大于ICA-LQG-COC方法，但控制效果不如后者。

(4) 本文提出的LQG的優(yōu)化方法簡單易行，能有效地對其加權矩陣的參數(shù)進行優(yōu)化，這種ICA優(yōu)化LQG的方法可以拓展到控制領域的其它應用背景下，具有一定的工程應用價值。另外，ICA的參數(shù)對其代價函數(shù)值有怎樣的影響、如何優(yōu)選這些參數(shù)以及對ICA進行算法改進將是我們今后研究的一個方向。

[1] DYKE S J, SPENCER JR B F, SAIN M K, et al. Modeling and control of magnetorheological dampers for seismic response reduction[J]. Smart Materials and Structures, 1996,5(5): 565-584.

[2] ZEMP R, DE LA LLERA J C, SALDIAS H, et al. Development of a long-stroke MR damper for a building with tuned masses[J]. Smart Materials and Structures, 2016, 25(10): 105006.

[3] FALLAH A Y, TAGHIKHANY T. Robust semi-active control for uncertain structures and smart dampers[J]. Smart Materials and Structures, 2014, 23(9): 095040.

[4] EL-KHOURY O, KIM C, SHAFIEEZADEH A, et al. Mitigation of the seismic response of multi-span bridges using MR dampers: experimental study of a new SMC-based controller[J/OL]. Journal of Vibration and Control, 2016,1077546316633540.

[5] 宗路航. 磁流變阻尼器的動力學模型及其在車輛懸架中的應用研究[D]. 合肥：中國科學技術大學，2013.

[6] LIN X F, CHEN S M. A modified shuffled frog-leaping algorithm-based fuzzy controller for magnetorheological damper-building system[J]. International Journal of Computer Applications in Technology, 2016,53: 279-289.

[7] 宋春生，于傳超，張錦光，等. 基于遺傳算法的復雜雙層磁懸浮精密隔振系統(tǒng)LQR控制研究[J]. 振動與沖擊，2016, 35(16): 99-105.

SONG Chunsheng, YU Chuanchao, ZHANG Jinguang, et al. LQR control of a complex two-stage magnetic suspension active precise isolation system based on the genetic algorithm[J]. Journal of Vibration and Shock, 2016, 35(16): 99-105.

[8] 景益娟.大型雷達天線伺服控制方法及抗風擾特性研究[D].西安：西安電子科技大學，2013.

[9] 王磊，譚平，李森萍. 基于人工魚群算法的隨機結構AMD控制系統(tǒng)LQR權矩陣優(yōu)化[J]. 振動與沖擊，2016, 35(8): 154-158.

WANG Lei, TAN Ping, LI Senping, et al. Optimal analysis of weight matrices of LQR algorithm for stochastic structure-AMD system based on artificial fish algorithm[J]. Journal of Vibration and Shock, 2016, 35(8): 154-158.

[10] ATASHPAZ-GARGARI E, LUCAS C. Imperialist competitive algorithm: an algorithm for optimization inspired by imperialistic competition[C]//2007 IEEE Congress on Evolutionary Computation. 2007: 4661-4667.

[11] KAVEH A, TALATAHARI S. Optimum design of skeletal structures using imperialist competitive algorithm[J]. Computers & Structures, 2010, 88(21): 1220-1229.

[12] ABD-ELAZIM S M, ALI E S. Imperialist competitive algorithm for optimal STATCOM design in a multimachine power system[J]. International Journal of Electrical Power & Energy Systems, 2016, 76: 136-146.

[13] SHABANI H, VAHIDI B, EBRAHIMPOUR M. A robust PID controller based on imperialist competitive algorithm for load-frequency control of power systems[J]. ISA Transactions, 2013, 52(1): 88-95.

[14] SPENCER B F, DYKE S J, SAIN M K, et al. Phenomenological model for magnetorheological dampers[J]. Journal of Engineering Mechanics, 1997, 123(3): 230-238.

[15] ATABAY E, OZKOL I. Application of a magnetorheological damper modeled using the current-dependent Bouc-Wen model for shimmy suppression in a torsional nose landing gear with and without freeplay[J]. Journal of Vibration and Control, 2014,20(11):1622-1644.