基于帝国竞争算法和LQG的磁流变阻尼器-建筑结构系统的振动控制

林秀芳， 陈淑梅, 黄国荣(.福州大学 机械工程及自动化学院，福州 3506;2.福建农林大学 金山学院，福州 350002)

提高建筑结构在地震等外部扰动下的抗震性以及减轻动力作用危害是土木结构工程领域面临的重大课题之一。近些年，作为一种智能的半主动控制装置，磁流变阻尼器已经在振动领域获得越来越多的关注[1-3]。控制算法是半主动控制系统的核心，它直接决定了控制系统的控制效果。现有的利用MRD进行振动控制的算法主要包括H∞控制、滑模控制[4]、最优控制[5]、模糊控制[6]等。其中，LQG控制是一种适用性很强的最优控制，其最大的优点是可得到状态线性反馈的最优控制律。又因其性能指标物理意义明确，最优解可以解析表达，因此被广泛运用于振动控制中。

LQG控制算法的设计是在加权矩阵选定的基础上进行的，其设计效果直接取决于加权矩阵的选取。如果加权矩阵选取不当则无法得到满意的设计结果。其结果会造成控制装置不能发挥最优性能和资源的浪费，甚至对结构响应产生不利的影响。因此，加权矩阵的选取是LQG设计中的一个关键问题。但LQG控制系统的加权矩阵的设计过程往往依赖于人工经验和领域专家的先验知识，而且人工设计的效率很低且不能保证系统达到最优。近些年，随着智能优化算法的兴起，用其对线性二次型控制系统的加权矩阵进行优化选择已成为热点。例如，宋春生等[7]利用遗传算法(GA)优化最优控制的加权矩阵，实现对复杂双层磁悬浮精密系统的隔振控制。结果表明优化后的控制系统的隔振效果有显著提高。景益娟[8]在设计大型雷达天线的抗风扰伺服控制系统中利用人工蜂群算法优化LQG控制器的加权矩阵。仿真结果表明优化过的LQG控制器能够显著地提高风扰情况下的跟踪精度。王磊等[9]利用人工鱼群算法求解基于主动质量阻尼器的随机结构的控制系统中加权矩阵的优化问题。

帝国竞争算法(ICA)由Atashpaz-Gargari等[10]于2007年提出，是一种富有生命力和高性能的群智能优化算法。该算法与多种智能优化算法相比都表现出较高的优化效率和鲁棒性，因此近年来已被越来越多的学者关注和应用。例如，Kaveh等[11]提出基于ICA的骨架结构的优化设计，与其它多种优化算法(如粒子群算法、蚁群算法等)比较的结果表明了该算法的优越性。Abd-Elazim等[12]利用ICA优化多机环境下的静止同步补偿器，仿真结果表明该算法的性能超过了GA。Shabani等[13]为了实现对载荷频率的鲁棒控制，提出了基于ICA优化的PID控制器。仿真结果表明其性能由于基于基于GA优化的PID控制器和基于神经网络优化的PID控制器。

作为一种近几年来新兴的智能优化算法，虽然ICA已在多领域得到成功应用，但尚未发现其在LQG控制器优化这方面的文献。本文基于优化设计的思想，在选中某种优化指标的基础上，尝试采用ICA对LQG的加权矩阵进行优化选取，从而获得结构振动控制所需的理想控制力。此外，由于磁流变阻尼器的力不可控，只有电压可控，所以本文还利用剪枝最优控制(COC)算法将理想控制力转化成MRD所需的控制电流，从而实现基于MRD的建筑结构振动的半主动控制。

1 基于LQG控制的受控系统的数学模型

LQG控制算法实际上是输出反馈控制算法，是由最优状态反馈器(Linear Quadratic Regulator, LQR)和基于Kalman滤波的状态观测器组成。

1.1 基于LQR控制的受控系统的数学模型

在一个具有n个自由度的结构上安装p个MRD，则基于MRD的受控结构的运动方程可以表示为

(1)

(2)

式中，

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

(3)

其中，

综合式(2)和式(3)，则相应于运动方程(1)的控制系统的状态方程可以用标准状态方程表示为

(4)

由于位移响应是与安全性有关的参数，加速度响应是与舒适性有关的参数，为了同时兼顾建筑物的安全性和其中人员的舒适性，同时限制控制力的幅值，控制问题的目标性能函数定义为

(5)

将目标性能函数重新改写为

(6)

式中：Q为状态变量的加权矩阵，R为控制变量的加权矩阵，这两个矩阵中的参数将由ICA优化确定。根据LQR最优控制率，最终得到控制器的反馈为

(7)

式中：G是最优反馈增益，P由如下的Riccati方程求出

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

(8)

1.2 Kalman滤波器

LQR控制器需要全部的状态变量，但在实际控制中，一般只有加速度能够被简单地测得。因此，需要从可测量的信号中估计出所有的状态变量，这里釆用Kalman滤波器估计系统的状态。Kalman滤波器的性能函数定义如下

(9)

(10)

(11)

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

(12)

其中，Qe=qe，Re=re×eye(n,n),eye(n,n)为一个n×n的单位矩阵。最终，系统的最优控制力可表述为状态估计值的反馈

(13)

2 基于ICA优化的LQG控制

2.1 ICA的优化原理

ICA是在模拟人类社会发展过程中帝国竞争过程的基础上提出的一个新颖的元启发式算法。该算法通过国家间的竞争来实现优化，“国家”为群体中的个体，每个“国家”代表待优化问题的一个解，即每一个“国家”是一组加权矩阵Q和R中的待优化参数的集合。该算法的基本流程包括：帝国初始化、帝国的同化与更新、帝国集团间的竞争和帝国的灭亡。

(1) 帝国初始化：在一个N维的优化问题中，一个国家即为一个1×N的矩阵，该矩阵定义为

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

(14)

式中:p1,p2,p3,…pN为待优化的变量。每个国家的实力大小通过代价函数来衡量

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

(15)

代价函数值越小，国家实力越大。首先，设置算法开始时国家的数量Npop，选取其中实力最强的Nimp个国家作为帝国主义国家，剩下的Ncol个国家作为帝国主义国家的殖民地。然后，为了形成初始帝国集团，每个帝国主义国家根据实力大小，按照轮盘赌选择法获得相应的殖民地。

(2) 帝国的同化与更新：每一个帝国为了壮大实力，会尽可能地同化所占领的殖民地，使殖民地和帝国相融合。这个过程可以表达为殖民地在二维平面内沿着坐标轴向帝国随机移动，如图1所示。x表示殖民地向帝国移动的单元，它是一个随机数，服从均匀分布x～U(0,β×d)。d表示殖民地与帝国之间的距离，β是大于1的系数，表示殖民地从两边向帝国靠近。但帝国的同化过程并不一定会使殖民地直接向帝国移动，因此引入一个随机角度θ模拟该运动，θ是一个服从θ～U(-r,r)的均匀分布，r用于调整偏移初始方向的角度。一般取β=2，r=π/4。

当殖民地向帝国移动的过程中，可能出现殖民地的实力大于帝国的实力，在这种条件下，帝国和殖民地角色互换，即殖民地成为帝国，原先的帝国沦落为殖民地。

图1 殖民地向其帝国移动示意图

(3) 帝国集团间的竞争：出于本性每个帝国总是希望尽可能多地占领殖民地，因此不同帝国间不可避免地产生竞争。每一个帝国集团的实力主要由帝国的实力所决定，但殖民地的实力也对帝国集团有所影响。所以，将衡量一个帝国集团的实力大小的代价函数归一化为

(16)

其中，ξ的大小决定了殖民地国家对整个帝国实力的影响程度。Mean表示该帝国所有殖民地代价函数值的平均值。每个帝国集团都试图兼并其它帝国集团的殖民地。在帝国竞争过程中，每个集团都有占领最弱集团殖民地的可能性。

(4) 帝国灭亡及算法结束：实力弱的帝国经过帝国竞争后逐渐失去殖民地。当一个帝国失去所有殖民地时，该帝国灭亡。最终只剩一个帝国集团时，算法结束，输出最优解；否则，返回(2)。

2.2 优化目标

为了同时兼顾建筑物的安全性和其中人员的舒适性，ICA优化的目标函数定义为

obj=wJ1+(1-w)J2

(17)

其中

(18)

(19)

2.3 优化流程

用ICA优化LQG加权函数的参数的具体步骤如下。

(1) 根据受控对象的属性和LQG的控制目标，确定加权矩阵Q和R的结构以及矩阵中待优化参数的数量。

(2) 确定待优化参数的范围和优化算法的编码形式，本文采用实值编码。

(3) 初始化算法和帝国：初始化算法的参数，如国家的数量Npop，帝国的数量Nimp以及迭代的次数k等。然后，随机产生初始国家，若国家(解)对应的理想控制力超出MRD的最大量程，令代价函数值为1，否则，根据式(18)计算代价函数值。最后根据代价函数值的大小分配帝国和殖民地。

(4) 迭代开始：根据2.1节的第(2)步进行每个帝国的同化与更新，然后按照上一步的方法计算每个帝国集团的代价函数值，最后根据2.1节的第(3)步进行帝国间的竞争。

(5) 判断是否满足收敛标准：如果满足，则算法结束，输出最优解，否则，返回(4)。

3 ICA-LQG-COC半主动控制

要实现对MRD的应用，通常需要两个步骤：首先根据结构响应计算理想控制力，然后将理想力转换为MRD所需的输入信号。本文所提的ICA-LQG-COC半主动控制系统是由ICA-LQG控制器和COC这两部分组成，如图2所示。图中，ICA-LQG控制器根据建筑物的响应计算出理想控制力，COC通过调整输入电流使得MRD能够跟踪理想控制力。需要注意的是图2中的MRD在仿真中由动力学模型替代。

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

3.1 MRD的动力学模型

本文使用的是Lord公司生产的RD-8040-1型MRD。该阻尼器属于剪切阀式单出杆型，其结构原理如图3所示。当活塞杆在缸体中往复运动时，通过控制输入线圈电流的大小来改变环形孔周围的磁场强度，从而改变磁流变液的流变特性，产生不同的阻尼力。对MRD进行性能测试所使用的是INSTRON公司生产的INSTRON1342电液伺服疲劳试验机，图4为MRD的测试系统。测试时通过Lord公司生产的RD-300-2-03型电流控制器调节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在不同激励性质及控制电流下的力学性能数据，分析其动力学特性，从而进行MRD动力学模型的辨识及验证。Bouc-Wen模型因其精度较高，数字处理方便简单而得到较为广泛的应用[14-15]。因此，本文用该模型表征MRD的动力学特性，模型中的阻尼力f表示为

图4 磁流变阻尼器测试系统

(20)

(21)

式中:x是活塞的相对位移；c0和k0代表高速时的黏性阻尼和刚度；x0是补偿器气体弹簧的初始位移;α是滞环算子的比例因子;γ,β,A和n是调整滞环的参数，其中n通常取2。

为了简化参数识别的工作，需要让每个待识别参数与电流产生关系。在对实验数据进行参数拟合的过程中，为了尽可能减小预测阻尼力和实验阻尼力之间的误差，待识别参数与电流的关系表征如下

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个参数需要识别。使用MATLAB中的遗传算法工具箱，取激励频率1.6 Hz、激励幅值7.5 mm的实验数据进行参数识别，得到的Bouc-Wen模型的具体参数数值见表1所示。在此基础上，图5显示了预测阻尼力和实验阻尼力的对比结果，从图中可看出模型预测结果与实验结果非常吻合，说明识别到的Bouc-Wen模型可以准确地表征该MRD的动力学特性。

3.2 剪枝最优控制(Clipped-optimal control, COC)

MRD这种控制装置具有很强的非线性，其阻尼力不能直接控制。为此，本文采用Spencer等提出的COC方法将理想控制力转化为MRD的控制电流，控制电流i表达如下

表1 Bouc-Wen动力学模型的参数值

(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{}是阶跃函数，fc是理想控制力，f是MRD的控制力。由上式可知，COC的核心控制率是：当f=fc时，控制电流值维持原值；当f

4 数值计算及分析

本节以一个三层的建筑物结构为例，分别用三种地震波对其进行激励，为验证ICA-LQG-COC方法的有效性，将其与被动控制(0 A)、被动控制(2 A)以及LQG-COC半主动控制进行比较。

4.1 仿真条件

基于地震波激励的安装了MRD的框架结构如图6所示，MRD安装在在第一层和地面之间。

图6 基于磁流变阻尼器减震的三层结构

该结构的质量M和刚度矩阵K分别为

阻尼矩阵的计算公式如下

C=αM+βK

(30)

式中:α和β分别为质量阻尼系数和刚度阻尼系数。这两个系数可通过振型阻尼比计算得到，即：

(31)

式中:ωi和ωj分别为结构的第i和第j阶固有频率；ξi和ξj为相应于第i和第j振型的阻尼比。通常情况下取i=1，j=2。最终计算得到阻尼矩阵为

根据式(5)和(6)定义的目标性能函数，LQG控制器的反馈增益G的加权矩阵Q和R表达如下

R=r

其中，q1、q2和r是加权矩阵的待优化参数，这三个参数的优化范围如下

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

由于滤波器增益Ke的加权矩阵参数较易选取，且取值对结果影响较小，故将qe和re分别取为1和10-7。

ICA的具体参数为：国家的数量Npop=30、帝国的数量Nimp=3、迭代次数=100、β=2、r=π/4、ξ=0.1。令式(17)中的权重w=0.5。算法运行次数=5，取5次中的最优解。至于如何更好地选择ICA的参数才得到更优化的解，超出了本文的研究范畴。此外，LQG-COC算法中加权矩阵的参数分别是：q1=10,q2=0.01,r=10-6。

4.2 仿真结果和讨论

对该框架结构采用El-Centro地震波(加速度峰值为4.946 m/s2)进行激励，采样频率为0.02 s，激励时间为20 s。由于建筑物是一个缩尺的模型，按照结构相似度原则，把地震波加速度加速5倍作为该模型的激励。经过ICA优化确定的参数分别为q1=73 618,q2=0.604 3,r=4.380 3×10-7。

除了第一层以外，ICA-LQG-COC方法对其它楼层的控制效果都优于被动控制(2 A)，但是被动控制(2 A)所需的最大控制力是所有控制方法中最大的。这说明即使将阻尼器的能量调到最大，不一定会得到更好的减震效果。从另一个角度而言，ICA-LQG-COC方法不仅控制效果好，而且更节能。这里补充一点，在优化过程中应满足的约束条件是，优化得到的参数对应的理想控制力不能大于被动控制(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)但都优于被动控制(0 A)。此外，由图9可见，除了第一层的位移和加速度，ICA-LQG-COC对其它楼层的所有波峰响应的控制效果都优于被动控制(2 A)。由图10可见，ICA-LQG-COC对最大位移和加速度响应的控制优于被动控制(2 A)，只有对第一层的位移响应的控制不如被动控制(2 A)。总之，相比其它几种控制方法，ICA-LQG-COC展示出最优的综合控制性能。

图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的拉压动态试验表明：利用依赖电流的Bouc-Wen模型表征MRD的动力学模型，可以比较精确地描述该阻尼器的动力学特性。

(2) 以安装了MRD的框架结构为例，仿真结果表明：相比LQG-COC半主动控制方法和其它两种被动控制方法，本文提出的ICA-LQG-COC半主动控制方法在降低楼层的响应峰值时展现出更好的综合性能。虽然优化LQG加矩阵的目标是同时最大限度地降低最大位移和最大加速度，但ICA-LQG-COC方法还很有效地降低了其它楼层的位移和加速度响应，甚至层间位移响应。

(3) 利用MRD进行振动控制时不是控制力越大，控制效果就越好。当采用被动控制(2 A)方法时，虽然阻尼力始终大于ICA-LQG-COC方法，但控制效果不如后者。

(4) 本文提出的LQG的优化方法简单易行，能有效地对其加权矩阵的参数进行优化，这种ICA优化LQG的方法可以拓展到控制领域的其它应用背景下，具有一定的工程应用价值。另外，ICA的参数对其代价函数值有怎样的影响、如何优选这些参数以及对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] 宋春生，于传超，张锦光，等. 基于遗传算法的复杂双层磁悬浮精密隔振系统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控制系统LQR权矩阵优化[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.