基于实时混合试验误差累积规律的时滞影响修正方法

周子豪，李 宁,2，李忠献,2

(1.天津大学建筑工程学院/滨海土木工程结构与安全教育部重点实验室(天津大学)，天津 300350；2.中国地震局地震工程综合模拟与城乡抗震韧性重点实验室(天津大学)，天津 300350)

1969年日本东京大学Motohiko课题组实现了第一个混合试验[1]，开启了混合试验方法的研究历程。经过30余年发展，混合试验结构划分的试验理念得到了广泛认同：将试验对象中本构关系明确、易于模拟的部分划分为数值子结构，而试验对象中难于建模却又对整体结构的响应起到重要作用的部分划分为物理子结构，并在试验室中进行试验以获得其实际的响应。从拟动力混合试验到快速拟动力混合试验再到实时混合试验，混合试验的应用愈发成熟。近十年来，涌现了许多实时混合试验成果。Karavasilis等[2]利用实时混合试验对基于性能抗震设计的装有阻尼器的钢框架进行了性能评估测试。Chen等[3]利用实时混合试验对安装大型磁流变阻尼器的钢框架进行了分析。Shao等[4]对大型结构的振动台和作动器联合实时混合试验方法进行了研究。Sorkhabi等[5]利用实时混合试验对TLD在地震作用下的性能进行了研究。Chen等[6]研发了实时混合试验平台，考虑了作动器、控制架构、数值积分算法以及时滞补偿等，并对大型实时混合试验的多作动器时滞自适应补偿[7]进行了探索。Chae等[8]利用实时混合试验对大型MR阻尼器控制策略进行了研究。Xu等[9]利用实时混合试验对结构主动控制性能进行了研究。Eem等[10]对MRE智能基础隔振的性能进行了研究和评估。Asai等[11]利用实时混合试验对高耸结构智能外伸阻尼系统进行了研究。Saouma等[12]进了非延性混凝土结构的实时混合试验，并将实时混合试验的结果与振动台试验进行了对比。Gao等[13]对多自由度实时混合试验的一般流程进行了研究和试验验证。Cha等[14]对比了多种MR阻尼器半主动控制实时混合试验结果和纯模拟结果。Jiang等[15]针对复杂数值子结构的情况开展了安装MR阻尼器的桥梁实时混合试验。

目前实时混合试验几乎都是在少数几个设备精良的高校实验室中完成的。这些试验室中环境干扰少，试验人员专业知识和试验经验较丰富，能够较好地对试验的过程进行把控，对试验中的不足可以及时发现和修正。还有一个重要的原因是，由于实时混合试验需保证试验过程的实时性，对实时机、控制器和作动器都有较高的要求，造成试验设备的成本高，这在一定程度上也阻碍了实时混合方法的发展和应用[16−17]。此外，随着作动器最大负载的增加，要达到实时的难度也随之增加。目前大部分实时混合实验中物理子结构依然是缩尺结构模型。对大型的结构而言，如桥墩、高层建筑等，体系并不复杂，仅仅是体形较大，受限于加载设备能力，致使实时混合试验难于开展原型试验。近几年，随着实时子结构试验技术和设备的发展，国内外也出现了许多针对大型结构的实时混合试验，例如，Liu等[18]基于实时混合试验探究了钢筋混凝土薄壁高墩桥在地震作用下的结构响应以及提出了该结构失效模型预测方法；Wu等[19]进行了足尺钢框架的界面模型更新的实时混合试验；Mei等[20]进了高墩混凝土桥梁的基于模型更新的实时混合试验。

为了保证实时混合试验的实时性，问题主要集中在如何有效地对试验中的时滞影响进行消除。在实时混合试验中，相对于数值计算和数据传输的时滞，作动器时滞是占主导地位的。为了消除时滞的影响，国内外的学者进行了非常多的探索。其中代表性的有Horiuchi等[21]提出的多项式外插时滞补偿方法、Chen[22]提出的一阶反补偿时滞补偿算法以及Phillips等[23]提出的基于模型的时滞补偿方法。笔者所在的科研团队也在时滞补偿算法的开发上做出了自己的贡献。例如，将在线时滞预估[24]分别与多项式外插时滞补偿[21]和一阶反补偿[22]相结合，提出了二阶在线迭代时滞补偿算法[25]以及自适应时滞补偿算法[26]；引入模型误差提出了二阶的基于模型的时滞补偿方法[27]。但是，目前时滞补偿算法依旧无法消除时滞，特别是对于高频信号，补偿的效果仍然不理想。

考虑到时滞普遍存在于实时混合试验中，针对由于实验设备性能限制或结构负载较大而导致的时滞补偿不理想等情况，本文提出了一种基于双显式数值积分算法误差累积规律的试验后时滞修正方法。可在试验后对结果进行修正以消除时滞影响，降低试验成本。首先，对双显示数值积分算法的误差累积规律进行了探索，指出在实时混合试验过程中，时滞导致的误差是按照一定的规律进行累积传播的；其次，根据得到的误差累积规律提出了一种试验后的时滞补偿方法，在试验后对试验结果进行修正，以达到减小时滞影响的目的；最后对方法的有效性通过4种不同的物理子结构工况进行了验证。

1 数值积分算法的误差累积规律

数值积分算法在实时混合试验中是关键部分。物理子结构的反馈进入数值子结构之后，通过数值积分算法得到下一步的位移指令，进而发送到控制器驱动作动器运动。由于时滞的存在，在试验过程中作动器的位置相对于位移指令滞后，导致物理子结构反馈滞后，此时积分得到的结果也就必然与理想无时滞状态存在误差。因此，无时滞理想状态下与有时滞状态下数值子结构计算的结果是不同的。无时滞、有时滞单自由度体系运动方程可写为式(1)和式(2)。

双显式数值积分算法的通用表达式为式(4)、式(5)所示，速度和位移均为显式递推得出。

当数值积分算法分别在无时滞理想情况和有时滞情况下，可分别改写为式(6)、式(7)及式(8)、式(9)。

式(8)减去式(6)，式(9)减去式(7)，得式(10)、式(11)。

式(10)、式(11)即为双显式数值积分的位移和速度误差累积公式，其中：Ed、Ev和Ea分别表示位移误差、速度误差和加速度误差。可以看到，其表达式与积分算法的形式是一致的，而且位移和速度的误差累积公式是显式的。将运动方程式(2)减去式(1)，可以得到式(12)、式(13)。

综上所述，式(14)、式(15)及式(16)组成的方程组：

综合上述的推导，将整体的修正过程总结如图1所示。

图1 试验后修正方法框图Fig.1 Flow chart of post-correction method

Gao等[28]对时滞试验结果的影响进行了分析，其时滞微分方程列式如下：

令x=Asinωt，忽略Gao等[28]引入的幅值误差，只考虑时滞，则xm=Asinω(t−δt)。其中：M、C、K分别为质量、阻尼和刚度；下角标N、P分别代表数值和物理子结构；x¨g为地震动加速度；xm为作动器量测位移；A和ω分别为结构振动的幅值和频率；δt为时滞。假设δt很小，将xm进行泰勒展开得式(22)～式(23)：

将其代入时滞微分方程得：

观察式(25)，发现物理子结构的刚度向运动方程中引入了负阻尼CKp=−KPδt，向系统中引入了能量，这与Horiuchi等[21]得到的结论一致。从式(25)还可以发现，物理子结构的质量向运动方程中引入了正阻尼CMp=MPω2δt，耗散了系统的能量；物理子结构的阻尼向运动方程中引入了负质量MCp=−CPδt，向系统中引入了能量。

因此，当存在时滞时，物理子结构的质量会增加系统的能耗、降低系统的能量；而物理子结构的刚度和阻尼则会向系统中引入能量。故此，当物理子结构质量耗散的能量小于刚度和阻尼引入的能量时，实时混合试验的结构响应要大于无时滞的结构响应。因此，在物理子结构的质量相对较小的情况下，可以在有时滞的实时混合试验结果中找到无时滞条件下的结构响应。

至此得到应用本文提出方法的第一个限制条件：物理子结构质量配比相对较小，并使得物理子结构质量耗散的能量小于物理子结构刚度和阻尼引入的能量。这一点在实际的应用中很容易满足，其原因在于质量带来的惯性力可以比较准确的进行测量，且质量在试验过程中一般是不变的，其物理规律是比较清晰的。故在实时混合试验中一般不会将质量占总结构比例非常大的构件作为物理子结构[29]。

图2 搜索方式示意图Fig.2 Schematic diagram of search method

3 模拟验证

本文采用4种典型工况对所提出方法进行模拟验证：物理子结构分别为线性刚度、线性阻尼、非线性刚度及非线性阻尼构件。其中非线性刚度构件采用的是非线性弹簧[30]以及Bouc-Wen模型[31]，而非线性阻尼构件则是采用MR阻尼器模型[32]。

首先，分析了时滞对试验结果的影响。对于单自由度体系，取数值子结构的MN、CN、KN分别为7 kg、36 (N·s)/m、100 N/m；物理子结构的MP、CP、KP分别为3 kg、4(N·s)/m、900 N/m；引入如图3所示的时滞量，根据目前作动器的性能，其时滞数值的取值在15/1024 s～25/1024 s波动；地震激励为El-Centro波。有时滞与无时滞的模拟结果对比如图4所示。有时滞情况得到的试验结果与无时滞情况得到的试验结果误差如图5所示。

图3 人为引入时滞Fig.3 Artificially introduced time delay signal

图4 有/无时滞试验结果对比Fig.4 Test results with/without time delay

图5 有时滞与无时滞试验结果误差Fig.5 Errorsbetween testing resultswith/without delay

3.1 线性刚度、阻尼构件

表1 数值与物理子结构参数Table 1 Parametersof numerical and physical substructures

从图6～图7可以看出，在时滞相同、物理子结构为线性刚度的情况下，随着物理子结构刚度占总体刚度比例的增加，有时滞与无时滞试验结果的位移、速度和加速度误差呈近似指数增长。且通过上述的模拟结果可以发现，在相同的时滞条件下，随着物理子结构刚度占比的增加，有时滞与无时滞情况下的准确值相差越来越大。上述模拟的结果均为理想模拟得到，在实际的试验过程中，由于环境噪声、系统扰动等因素的干扰下误差势必会进一步增加。

图7 线性阻尼物理子结构修正前后误差对比Fig.7 RMScomparison for linear damping physical substructure tests before and after correction

当物理子结构为线性阻尼的情况下，随着物理子结构阻尼占总体阻尼比例的增加，有时滞与无时滞的试验结果表明，位移、速度和加速度误差同物理子结构为线性刚度时都呈现出增大趋势。当物理子结构为线性刚度时，物理子结构的刚度占总体结构的比例为20%时，其RMS值已经达到了400%，试验结果已不具备任何参考价值；而当物理子结构为线性阻尼时，物理子结构的阻尼占总体结构的比例为5%～80%，但是其RMS却远远小于物理子结构为线性刚度的情况。而且，物理子结构为线性阻尼时，其占总体阻尼比例的增加对位移、速度和加速度的影响是依次递增的，并不象物理子结构为线性刚度时，物理子结构刚度占比的增加对位移、速度和加速度误差增速的影响是基本一致的。

3.2 非线性刚度构件

对于非线性的刚度物理子结构，我们采用了两种不同的模型：1)高阶非线性弹簧模型[30]，如式(27)所示；2) Bouc-Wen模型[31]如式(28)、式(29)所示。模型中物理子结构和数值子结构参数的取值见表2及表3。高阶弹簧及Bouc-Wen模型的位移-力曲线如图8所示，其中Bouc-Wen模型采用6种不同幅值的正弦激励，A代表正弦激励的幅值。

图8 高阶弹簧及Bouc-Wen模型特性Fig.8 Characteristicsof high-order spring and Bouc-Wen model

表2高阶弹簧及Bouc-Wen模型参数Table2 Parametersof high-order spring and Bouc-Wen model

表3数值子结构参数Table 3 Parametersof numerical substructure

为了验证算法的普适性，以及保证实时混合试验的稳定性，本文将物理子结构的反馈力乘以一个系数之后再反馈到数值子结构中，分别称为高阶弹簧系数及Bouc-Wen系数。随着系数不断增大，物理子结构的刚度占总体刚度的比例不断增加。图9为有时滞情况下高阶弹簧及Bouc-Wen模型修正前后的误差(RMS值)对比。通过图9可以发现，类似于线性刚度情况，时滞对于试验结果的影响随着物理子结构的刚度占比增加而增加。而且，本文方法在物理子结构为高阶弹簧时表现出了与纯刚度情况下一样的准确性，这是因为物理子结构为高阶弹簧时，根据位移就可以准确地对应其反馈力。其修正前、后的位移、速度和加速度误差对比如图10～图12所示。

图9 高阶弹簧及Bouc-Wen模型修正前后误差对比Fig.9 Error comparison of high-order spring and Bouc-Wen model before and after correction

图10 修正前后误差对比位移时程图(系数0.8)Fig.10 Displacement comparison before and after correction(physical substructure:Bouc-Wen model,coefficient:0.8)

图12 修正前后误差对比加速度时程图(系数0.8)Fig.12 Acceleration comparison beforeand after correction(physical substructure:Bouc-Wen model,coefficient:0.8)

图11 修正前后误差对比速度时程(系数0.8)Fig.11 Velocity comparison before and after correction(physical substructure:Bouc-Wen model,coefficient:0.8)

但物理子结构为Bouc-Wen模型时，本文所提方法的误差却随着系数的增加而不断的增加。这主要是由于：1)Bouc-Wen模型的反馈力是由位移和速度共同决定，本文提出方法的适用条件之一是位移相关力或者速度相关力占据反馈力的主导地位，而Bouc-Wen模型反馈力由位移和速度共同决定的特性决定了根据位移对其反馈力进行确定会造成误差；2)Bouc-Wen模型的滞回特性是造成误差的另一个重要原因，本文提出的方法是根据位移和速度其中之一进行反馈力的修正，然而在位移曲线的拐点处，Bouc-Wen模型由于滞回的特性，即使位移大小一致，但在拐点前后所对应的反馈力是不一致的，也会造成修正反馈力的误差。但还是可以非常明显地看到本文提出的算法对试验结果的修正效果尚可。

3.3 非线性阻尼构件

非线性阻尼构建的模型采用的是Spencer提出的磁流变阻尼器模型[32]，其表达式如式(30)～式(32)所示，其模型如图13所示，其在正弦激励下的响应如图14所示，其中I为输入磁流变阻尼器的电流。

图13 磁流变阻尼器模型Fig.13 Magnetorheological damper model

图14 磁流变阻尼器性能Fig.14 Performanceof magnetorheological damper

式中，α、c0、k0、k1、c1、γ、β、n和A为阻尼器模型的参数，如表4所示[34]。

表4 磁流变阻尼器模型参数Table 4 Parameters of magnetorheological damper

当物理子结构为磁流变阻尼器时，同样将其反馈力乘以一个系数，以此来调节物理子结构在整体结构中的占比。结果表明(如图15所示)，本文修正方法可以较好地降低时滞引起的试验误差。由磁流变阻尼器的模型可知，模型中也存在Bouc-Wen模型，是位移和速度共同决定的。由于本文提出的方法不能同时根据位移或速度来完全准确的对反馈力进行修正，故修正结果的误差随着物理子结构占比的增加而不断增加，但总体而言其误差还是较大幅度地降低了。

图15 修正前后误差对比Fig.15 Error comparison before and after correction

4 结论

实时混合试验中，在线补偿方法存在以下2个主要的问题：第一，在线时滞补偿本质上是一种预测方法，使得这种方式无法完全消除时滞，目前的时滞补偿方法对于低频信号能够较好的补偿，但是对高频信号效果却不理想，所以实时混合试验中，时滞普遍存在，并对试验结果产生影响；第二，实时混合试验的实时加载，对加载设备和控制系统要求过高，对大型构件，对加载能力要求较高，现有试验条件很难或者需要付出巨大的经济成本才能达到实时要求。

针对上述问题，基于数值子结构积分算法误差累计规律提出了一种试验后时滞修正方法，并进行线性刚度、阻尼、非线性刚度阻尼验证，结论如下：

(1)该方法可同时对试验的位移、速度和加速度结果进行修正；

(2)试验后时滞修正方法的提出，允许在试验的过程中存在一定的时滞，等于放宽了对于加载设备和控制系统的要求，一些性能稍显不足的作动器也可以在保证试验稳定的情况下完成实时混合试验，然后通过本文提出的方法对结果进行修正得到较为准确的结果。

(3)对于试验后发现时滞补偿效果并不理想的试验，也可以通过本文提出的方法对结果进行修正，从而避免必须重复试验的窘境，可以节约试验成本。