基于能量一致积分的拉索-阻尼器实时混合试验方法

杨 格， 孙红硕， 吴 斌， 潘天林， 王 贞

(1.武汉理工大学 土木工程与建筑学院，武汉 430070;2.天津大学 中国地震局地震工程综合模拟与城乡抗震韧性重点实验室，天津 300350；3.东北电力大学 建筑工程学院，吉林 132012)

斜拉索因其结构柔、质量轻，在风荷载等外力作用下非常容易发生振动。由于拉索振动易导致拉索锚固端发生疲劳破坏等问题，需要采取措施来抑制其大幅度振动，工程中广泛采用的措施是在拉索锚固位置附近安装阻尼器。为了研究阻尼器对拉索的减振效果，一般采用数值模拟[1]或室内试验[2]的方式来获取拉索-阻尼器系统的动力响应。数值模拟应用方便，但现有的阻尼器数值模型均存在一定的简化，不能真实反应阻尼器在实际结构中的性能。室内试验因试验场地规模和技术等问题，很难对较长拉索进行真实工况下的模拟。

结构混合试验方法可以将结构中无法准确模拟的部分作为试验子结构，其他部分采用数值模型，两者在线交互完成对整体结构的模拟。该方法解决了因实验室场地规模限制等因素而无法开展大型结构试验的问题，被广泛应用于结构抗震试验研究中[3-6]。其中，实时混合试验方法能够考虑加载速度的影响，可实现对速度相关型试件力学性能的准确反映。鉴于实时混合试验的这一优点，可以将拉索作为数值子结构，阻尼器作为试验子结构，开展拉索-阻尼器系统的实时混合试验。该方法对研究拉索-阻尼器系统中阻尼器的减振效果有着重要意义。

实时混合试验对数值模型的计算效率有着极高的要求。为了保证计算效率，现有关于拉索-阻尼器实时混合试验的研究大多采用简化的数值模型来模拟拉索。Jung等[7]对拉索-阻尼器系统在不同风速下的响应进行了实时混合试验，其拉索数值模型采用伽辽金方法得到。伽辽金方法是一种简化计算方法，其在抛物线假设的基础上引入一个参数来考虑索垂度以减少计算量；而传统有限元分析是通过自重作用下索的形状来得到垂度，因此伽辽金方法计算精度要比传统有限元分析计算精度差。Duan等[8]建立了基于向量式有限元的拉索模型，阻尼器采用线性黏滞阻尼器数值模型来代替，对拉索-阻尼器系统的一阶模态振动进行了实时混合试验数值仿真。相对于传统有限元方法，基于向量式有限元的拉索模型具有更高的计算效率，但其在进行质点运动方程求解时，采用了显式的中心差分法。中心差分法是一种对线性体系有条件稳定的算法，而拉索-阻尼器系统属于强非线性系统，一方面拉索振动具有较强的几何非线性;另一方面阻尼器也往往具有很强的非线性，因此无法保证中心差分法在拉索-阻尼器系统动力分析中的稳定性。当拉索存在高频振动时，通常需要很小的时间步长来满足算法的计算精度和稳定性要求[9]，仍然存在计算时间可能超过时间步长而导致试验失败的问题。

因此，在拉索-阻尼器系统的实时混合试验中需要对非线性体系无条件稳定的时间积分方法。然而，目前混合试验中常用的无条件稳定时间积分方法主要是针对线性体系，如无条件稳定的显式CR方法(Chen and Ricles)[10-12]、显式Chang方法族[13]、隐式平均加速度法[14]、隐式中点法[15]等。对于非线性体系，Crisfield等[16]采用平均加速度法求解桁架单元的动力方程时，发现结果会出现发散。潘天林[17]通过对具有几何非线性的桁架体系进行分析，发现隐式平均加速法和隐式中点方法不能保证无条件稳定。Kuhl等[18]发现基于桁架单元的能量耗散积分方法也会存在计算结果不稳定的现象。为了实现对非线性体系的无条件稳定，Wu等[19]提出了能量一致积分方法，并将基于梁柱单元的能量一致积分方法应用到足尺钢框架结构的混合试验中。由于能量一致积分方法为隐式算法需要迭代，目前该积分方法方法尚未应用于实时混合试验中。考虑到近年来混合试验为应对速度相关型试件的精细化模拟需求，逐渐趋于有限元化和实时化[20]，因此，将能量一致积分方法应用于实时混合试验中具有重要意义。

能量一致积分方法是一种隐式方法，将其应用于实时混合试验时，会遇到迭代导致作动器加载速度波动较大的问题。对于阻尼器试件，加载速度波动较大会导致测得的阻尼器出力严重失真，致使试验失败。为此，本文提出了基于能量一致积分的拉索-阻尼器系统实时混合试验方法，一方面通过Jung等[21]提出的固定迭代次数并对迭代位移进行插值的方式来求解隐式差分方程，实现平滑加载；另一方面为保证试验过程中拉索-阻尼器系统的能量一致，对试验测得的阻尼器出力进行恢复力修正。在不考虑拉索抗弯刚度的情况下，拉索可由若干个桁架单元模拟，本文将基于桁架单元的能量一致积分方法应用于拉索-阻尼器系统的实时混合试验中，对拉索-阻尼器系统进行了一阶模态振动下的实时混合试验数值仿真，验证了方法的可行性。

1 基于桁架单元的能量一致积分方法基本原理

对于桁架单元，连续的动力方程可以表示为

(1)

式中,m、v、r和f分别为质量矩阵、速度向量、节点力向量和外荷载向量。其中节点力向量r可表示为

r=N[-ee1]T

(2)

式中:N为单元的轴力;e1为桁架单元的轴向单位向量

(3)

式中:x=xQ-xP为节点坐标差，xP与xQ分别为桁架两个节点的坐标向量;L为单元的长度。

桁架单元的节点力考虑了几何非线性和材料非线性的影响，其中几何非线性由向量e1考虑，材料非线性则由N考虑。为了实现单元的能量一致，需要对恢复力进行非线性修正，其修正格式为

(4)

ri+β=r[(1-β)ui+βui+1]，

ri+1-β=r[βui+(1-β)ui+1]

(5)

式中：β为修正系数；ui、ui+1分别为第i、第i+1个积分点时刻的位移向量。然后分别对式(2)中两项非线性相关向量进行离散，对e1离散得到

(6)

式中,β1为几何非线性修正系数。对N离散得到

(7)

式中,β2为材料非线性修正系数。能量一致积分方法采用平均加速度法对位移和速度的假设

(8)

(9)

式中：Δt为时间积分步长；vi、vi+1分别为第i、i+1个积分点时刻速度向量；ai、ai+1分别为第i、第i+1个积分点时刻加速度向量。

基于此，式(1)的离散形式可表示为

(10)

式中:fi、fi+1分别为第i、第i+1个积分点时刻外荷载向量;εi、εi+1分别为第i、第i+1个积分点时刻的应变;A0为桁架单元原截面面积;L0为桁架单元原长;σ为工程应力;第三式和第四式分别为关于β1和β2的非线性方程。由于能量一致积分方法为隐式算法，整体节点力方程的求解通过牛顿迭代实现。

2 基于能量一致积分的拉索-阻尼器实时混合试验方法

2.1 基于能量一致积分的拉索-阻尼器实时混合试验方法基本原理

在拉索-阻尼器系统的实时混合试验中，以拉索作为数值子结构，阻尼器作为试验子结构，连续的动力方程可以表示为

(11)

式中:rN为数值子结构恢复力向量；rE为试验子结构恢复力向量，即阻尼器的出力。采用基于桁架单元的能量一致积分方法，式(11)可以离散为

(12)

式中:nN为桁架单元数目;βN1、βN2分别为几何非线性修正系数和材料非线性性系数。本文假定拉索材料为线弹性，故第四式中所有桁架单元的材料非线性修正系数βN2=0。由于测得的试验子结构恢复力rE不是连续的，故修正后的等效恢复力rEβ没有理论解，需根据试验数据结合第五式计算。

(13)

(14)

式中:k=0，1，…，n-1,n为每一积分步固定迭代次数,n=Δt/δt，Δt为积分步长，δt为采样步长；系数m=(k+1)/n；v0、u0分别为初始速度向量和初始位移向量。在计算第一个积分步时，因为没有上一积分点时刻的位移，位移命令可通过式(14)获得。

在纯时滞与总时滞相比非常小的情况下，作动器系统的动力性能可用不包含纯时滞环节的二阶传递函数来代替，传递函数可写为

(15)

式中：s为Laplace变量；ξ为作动器系统的阻尼比；ω为作动器系统的频率。时滞会影响实时混合试验仿真精度[22]，但在阻尼器作为试验子结构的实时混合试验中，时滞对仿真精度影响一般较小，故本文不进行作动器的时滞补偿。

基于能量一致积分的拉索-阻尼器系统实时混合试验实施过程,如图1所示。第i+1个积分步的第k次迭代实施过程可简述为：

图1 拉索-阻尼器系统实时混合试验实施过程示意

重复上述步骤，直至达到预定的固定迭代次数n，完成第i+1个积分步的模拟，并依此实现0～t时间内的仿真。

2.2 恢复力修正模块计算方法

为保证与系统能量一致，需要通过恢复力修正模块对测得的试验子结构恢复力进行修正，计算等效恢复力rEβ。在拉索-阻尼器系统中，第i+1个积分步中前j次迭代阻尼器出力所做的实际总功为

(16)

(17)

(18)

通过式(18)可以看出，分子是第i+1个积分步中前j次迭代的阻尼器出力所做总功的近似值，分母是第i+1个积分步中前j次迭代的总位移增量。通过式(18)对测得的阻尼器出力进行修正获得等效恢复力，实现整个过程系统能量一致。

3 拉索-阻尼器系统实时混合模拟

3.1 模型的建立

表1 J26号拉索结构参数

图2 拉索-阻尼器系统模型示意图

3.2 拉索模型的几何非线性和正确性验证

对拉索模型跨中处施加一个周期的正弦位移，幅值取1 m。通过静力分析可得拉索跨中节点恢复力-位移曲线如图3所示。可以看出节点恢复力与位移并非呈线性关系，验证了拉索模型的几何非线性。

图3 拉索跨中节点恢复力-位移曲线

记录模拟得到的拉索1/4跨处位移时程，最大振动幅值约为2.5 mm。通过图3可以看出，在小位移下拉索模型的恢复力-位移曲线基本处于线性阶段，此时可忽略几何非线性的影响。随后对该位移时程进行快速傅里叶变换，可得拉索前三阶模态频率如图4所示。与Duan等的研究中基于ANSYS软件分析得到的拉索前3阶模态频率对比，可看出两者所得前3阶模态频率基本一致，验证了本文拉索模型的正确性,如表2所示。

图4 位移功率谱密度函数

表2 拉索模态频率对比

3.3 实时混合试验的数值仿真

本文假定作动器系统频率和阻尼比分别为ω=3.14 rad/s、ξ=0.8[24],以式(15)所示传递函数来模拟实际作动器系统。在积分步长Δt=0.05 s下，取固定迭代次数n=15，对拉索-阻尼器系统进行实时混合试验数值仿真。为了尽可能真实地再现实际试验的情况，阻尼器出力采用迭代点时刻的实际速度来计算，在第i+1个积分步第k个迭代点时刻速度为

(19)

(20)

将频率为拉索第1阶模态频率的正弦荷载300×sin(2π×0.5t)N施加于拉索的23个结点上，荷载方向均垂直于拉索，激励100 s后释放荷载。在阻尼器作用下，振动自由衰减，可获得前150 s的拉索跨中位移时程，其仿真结果如图5和图6所示。从图5、图6可以看出，本文所提实时混合试验方法的仿真结果与整体模拟的数值解吻合较好。计算得到拉索跨中位移的均方根误差为1.96%，满足精度要求。

图5 拉索跨中位移时程曲线对比(1阶模态)

通过图6可以看出，本文方法计算幅值比数值解偏大。因此，本文方法的计算结果偏于保守，有利于保证结构的安全性。

在计算效率方面，该实时混合试验数值仿真时长为150 s，计算总耗时为9.45 s，每一积分步平均计算耗时0.003 15 s。可以看出本文所提方法计算效率较高，满足实时要求。其中，MATLAB软件版本为MATLAB2017b，计算机配置如表3所示。

表3 计算机配置

为了观察基于能量一致积分的拉索-阻尼器系统实时混合试验方法的迭代效果，取第二个积分步发给作动器模型的位移命令，如图7所示。可以看出发给作动器的位移命令时程曲线较光滑，因此本文所提方法可保证对阻尼器试件的平滑加载。

该仿真中阻尼器出力与速度关系曲线如图8所示。阻尼器出力与速度关系曲线满足本文阻尼器力学模型。

4 结 论

本文提出了基于能量一致积分的实时混合试验方法，并应用于拉索-阻尼器系统的实时混合试验仿真中，主要结论如下：

(1) 基于桁架单元的能量一致积分方法可考虑拉索的几何非线性，实现对拉索-阻尼器系统的振动过程数值模拟。

(2) 基于能量一致积分的实时混合试验方法可以实现速度相关型试件的光滑加载，可应用于拉索-阻尼器系统的实时混合试验中。