COMPARISONS OF MODEL-BASED EXPLICIT INTEGRATION ALGORITHMS IN REAL-TIME SUBSTRUCTURE TESTING
-
摘要: 实时子结构试验技术中的一个关键问题是求解数值子结构的动力响应,而这一过程可选取合适的数值积分算法来实现。针对目前已建立的三种基于模型的显式积分算法(Chang法、CR法和实时子结构RST法),对比分析了各算法在线性系统和非线性系统中的数值特性。结果表明:三种算法在线性系统和具有刚度软化特性的非线性系统中均是无条件稳定的,在具有刚度硬化特性的非线性系统中变为有条件稳定。当结构阻尼比为零时,三种算法均无数值阻尼,且周期延长率结果完全一致,并随Ω的增加而增大;当结构阻尼比不为零时,三种算法均存在数值阻尼,且CR法的数值阻尼绝对值最小,而RST法的周期延长率最小。两个算例表明,RST法和Chang法的精度要优于CR法,但因Chang法是半显式的,因此RST法更适于实时子结构试验中数值子结构的仿真计算。Abstract: In real-time substructure testing, a key issue is the numerical simulation of the numerical substructure, which can be done by applying appropriate integration algorithms to solve the step-by-step structural equations of motion. Three well-developed model-based explicit integration algorithms, that is, Chang, CR and RST methods, are discussed and compared in terms of the numerical properties in linear and non-linear systems. Stability analysis indicates that the three methods are unconditionally stable for linear and instantaneous stiffness softening systems while they are only conditionally stable for instantaneous stiffness hardening systems. The three methods show exactly the same characteristics for zero-damping systems as no numerical damping is produced and the period error (PE) increases with the increase of Ω. However, the three methods introduce numerical damping into the simulation for damped systems, and the numerical damping caused by the CR method is the smallest at a certain value of Ω. Meanwhile, the RST method has a smaller PE than the other two methods for damped systems. Two examples indicate that the RST method and Chang method have a better accuracy than the CR method in solving dynamic problems of both linear and non-liner systems. Since the velocity of the Chang method is in an implicit form, the RST is more advantageous in real-time substructure testing.
-
实时子结构试验作为近年来新兴的试验技术得到了广泛关注,该技术的核心思想是基于结构动力分析中的子结构概念,将整体结构划分为若干个子结构,其中难以模拟或者具有强非线性的构件作为试验子结构进行物理模型试验,而其余部分作为数值子结构进行数值仿真计算,通过子结构之间相互作用传递,获得整体结构的动力反应[1-2]。在求解数值子结构的运动方程时,通常需要选定合适的数值积分算法以保证实时子结构试验顺利开展,因此研究数值积分算法的精度和稳定性就显得尤为重要。
实时子结构试验的数值积分算法有显式和隐式之分,当第i+1步的位移和速度均可由第i步及其以前的解表达时称之为显式算法,反之为隐式算法[3]。显式算法的最大缺点在于条件稳定性,特别是多自由度体系子结构试验的条件稳定性对时间步长要求十分严格;隐式算法的主要优点是无条件稳定和能量耗散性能高,但缺点是需要迭代[4]。近年来,不少学者都致力于开展无条件稳定显式积分算法的研究,提出了几种基于模型的积分算法,即算法的积分参数与结构模型有关。2002年,台湾学者Chang[5-6]最先提出了基于模型的无条件稳定半显式积分算法Chang法,其位移表达式为显式,但速度表达式因含有下一步的加速度而为隐式;2008年,Chen和Ricles[7]利用离散传递函数的概念构造了一种新算法−CR法,实现了位移和速度的显式表达;随后,部分学者[8-9]保留CR法的位移和速度表达式,在积分参数的表达式中引入控制参数,扩大了CR法的应用范围;Kolay和Ricles[10]则利用离散控制理论提出了具有可控数值阻尼的无条件稳定显式KR-α法,该算法可以为高阶模态提供数值阻尼,而忽略低阶模态的数值阻尼,但计算过程十分复杂;Tang和Lou[11-12]在CR法的基础上进一步改进,提出了一种实时子结构试验RST法,并验证了该算法在线性系统中的无条件稳定性,但是缺乏算法在非线性系统中稳定性和精度的研究。本文对比分析Chang法(CEM)、CR法(CRM)及RST法在线性系统中的精度问题及在非线性系统中的数值特性,以寻求更适用于实时子结构试验的数值积分算法。
1 基于模型的显式积分算法
1.1 单自由度体系
当数值子结构和试验子结构均为单自由度(SDOF)体系时,结构运动微分方程可以写成离散形式为:
mN¨ui+1+cN˙ui+1+kNui+1+RE(ui+1,˙ui+1)=fi+1 (1) 式中:
mN 、cN 、kN 分别为数值子结构的质量、阻尼和刚度;¨ui+1 、˙ui+1 、ui+1 分别为数值子结构在第i+1步的加速度、速度和位移;RE(ui+1,˙ui+1) 为试验子结构在第i+1步的反作用力;fi+1 为第i+1步的外部激励。将位移和速度表达式统一写成:{ui+1=ui+α1Δt˙ui+α2Δt2¨ui˙ui+1=˙ui+α3Δt¨ui+α4Δt¨ui+1 (2) 式中:
Δt 为时间步长;基于模型的显式积分算法中系数α1 、α2 、α3 、α4 分别为:{α1=4+4ξΩ4+4ξΩ+Ω2,α2=24+4ξΩ+Ω2,α3=12,α4=12(CEM)α1=1,α2=44+4ξΩ+Ω2,α3=44+4ξΩ+Ω2,α4=0(CRM)α1=44+4ξΩ+Ω2,α2=4−2ξΩ−8ξ24+4ξΩ+Ω2,α3=1,α4=0(RST) (3) 式中:
Ω=Δtωn=Δt√(kN+kE)/m ;ξN=cN/(2mωn) ;ξE=cE/(2mωn) ;ξ=ξN+ξE 。1.2 多自由度体系
多自由度(MDOF)体系的结构运动微分方程可以写成矩阵的形式:
M⋅¨U(t)+C⋅˙U(t)+K⋅U(t)=F(t) (4) 式 中:M、C和K分别为多自由度体系的质量矩阵、阻尼矩阵和刚度矩阵;
¨U 、˙U 和U 分别为加速度向量、速度向量和位移向量;F 为力向量。为简单起见,这里假定矩阵C为比例阻尼矩阵。对于一个MDOF体系,各自由度之间的动力反应是相互耦合的,因此需要通过转换到模态坐标系进行解耦,考虑到模态的正交性,将式(2)和式(4)改写为:
M∗⋅¨Y(t)+C∗⋅˙Y(t)+K∗⋅Y(t)=F∗(t) (5) {Yi+1=Yi+α∗1Δt˙Yi+α∗2Δt2¨Yi˙Yi+1=˙Yi+α∗3Δt2¨Yi (6) 式(5)和式(6)均代表N个解耦的联立方程,其中
M∗=ΦTMΦ 、C∗=ΦTCΦ 、K∗=ΦTKΦ 分别为模态质量矩阵、模态阻尼矩阵和模态刚度矩阵,均为对角阵;F∗(t)=ΦT⋅F(t) 为模态力向量;Φ 为N×N阶振型矩阵,Φ= [ϕ1ϕ2⋯ϕN ] (i=1,2,⋯,N )为第i阶振型向量,N为模态数;Y(t) 为模态坐标向量,其与几何坐标向量U(t) 的关系为U(t)=Φ⋅Y(t) ;α∗1=Φ−1α1Φ 、α∗2=Φ−1α2Φ 、α∗3=Φ−1α3Φ 为对角积分参数矩阵。MDOF体系在模态坐标系下的动力反应是解耦的,因此对应于每阶模态的运动方程求解都可以采用同SDOF体系完全相同的步骤进行,然后再利用两个坐标系下各矩阵关系,最终得到几何坐标系下的解。
2 积分算法在线性系统中的数值特性
研究积分算法的数值特性主要从稳定性和精度两个方面着手。首先将式(1)~式(3)写成递推矩阵的形式:
Xi+1=AXi+LFi+1 (7) 式中:
Xi=[ui,Δt˙ui,Δt2¨ui]T ;L=[0,0,Δt2/m]T ;A为放大矩阵。以SDOF无阻尼系统为模型考核算法稳定性是有意义的,文献[4]中已经验证了各算法在线性系统中的无条件稳定性,因此不再赘述。此处主要讨论算法的精度问题,通常采用数值损耗和飘移来衡量,这两项指标是由算法的数值阻尼引起的,具体表现在计算中会出现振幅衰减(或增长)和周期延长(或缩短)。当算法无条件稳定时,
λ1,2 是一对共轭复根,记作:λ1,2=σ±iε=exp(−ˉξ⋅ˉΩ±iˉΩd) (8) 式中:
i=√−1 ;σ=A1 ;ε=√A2−A21 ;ˉΩd=ˉωdΔt= arctan(|ε/σ|) ;ˉΩ=ˉωnΔt=ˉΩd/√1−ˉξ2 ;ˉξ=−ln(σ2+ ε2)/2ˉΩ 为数值阻尼比。定义数值阻尼比与结构阻尼比的差值为算法阻尼比KSI,周期延长率PE为:
KSI=ˉξ−ξ,PE=ˉT−TT=ωˉω−1 (9) 在线性系统中,三种算法具有相似的数值特性,图1给出了三种算法在结构阻尼比为0.1时,算法阻尼比随Ω的变化曲线,可见当结构阻尼比为0时,各算法的数值阻尼比始终为0,说明此时三种算法均无数值阻尼;当结构阻尼比不为零时,算法阻尼比的绝对值均随Ω的增加而增大,CR法和RST法的算法阻尼比重合,且绝对值小于Chang法。图2显示,三种算法的周期延长率在线性条件下保持一致,且其值随Ω的增加而增大,随ξ的增大而减小。
3 积分算法在非线性系统中的数值特性
为了描述积分算法在非线性系统中的数值特性,需要定义结构在第i+1步结束时的瞬时刚度为
ki+1 ,其与结构初始刚度k0 的比值称为瞬时非线性度δi+1 。当δi+1=1 时,表示第i+1步结束时的瞬时刚度等于初始刚度;0<δi+1<1 时,表示第i+1步结束时发生了瞬时刚度软化;而δi+1>1 时,表示第i+1步结束时发生了瞬时刚度强化。由稳定条件和放大矩阵的特征值可以得到各算法在特定瞬时非线性度δi+1 下的稳定性,结果可知三种算法对于刚度软化系统是无条件稳定的,对于刚度硬化系统有条件稳定,其中各算法在δi+1=2.0 时谱半径随Ω的变化曲线如图3所示。图3中可见,在同一ξ条件下Chang法的稳定界限最大,其次是RST法;Chang法和RST法的稳定界限随ξ的增加而增大,CR法的稳定界面与ξ无关,始终保持为2。图4和图5分别给出了在ξ=0.1且δ分别等于0.5和2.0条件下,三种算法的算法阻尼比和周期延长率随Ω的变化规律。4 算例分析
4.1 线性体系
本算例为一个三层隔振体系,其中主体结构为剪切型框架(数值子结构),底层为隔振层并设有一个阻尼器(试验子结构)。数值子结构的各层质量
m=160×103kg ,刚度kN=360×106N/m ,阻尼cN=0N⋅s/m ;试验子结构的刚度kE=9× 106 N/m,阻尼cE=30×105N⋅s/m 。初始条件为u0 =1 cm,˙u0=0cm/s 。使用三种算法对结构进行自由振动下的结反应计算,积分时间间隔Δt=0.02 s,Newmark显式算法(Δt=0.001 s)的计算结果作为准确解用于对比。图6给出了结构作自由振动时顶层位移、速度和加速度的时程反应曲线。图6中可见,RST法和Chang法的精度优于CR法;当t=0时刻,由各数值积分算法得到的速度反应均小于原体系的速度,说明数值积分算法会改变原体系初始条件;在0.1 s内由各数值积分算法得到的顶层加速度反应与准确解的误差较大;0.1 s后与准确解吻合较好。
4.2 非线性体系
本算例为一个具有刚度软化特性的两层剪切型框架结构,将底层框架作为试验子结构,而二层框架作为数值子结构。各层质量分别为
mE=3× 103kg ,mN=2.5×102kg ;各层初始刚度分别为kE,0=107N/m ,kN,0=104N/m ;第i层的刚度表达式为ki=ki,0(1−√|Δui|) ,其中ki,0 为第i层的初始刚度,Δui 为第i层的层间位移;各层阻尼分别为cE=cN=0N⋅s/m 。计算可得前两阶初始频率分别为ω1=6.32rad/s ,ω2=57.76rad/s 。由三种算法求解结构的动力反应,积分时间步长Δt=0.02;外荷载为sin(2t) ,初始条件为u0=0cm ,˙u0=0cm/s 。Newmark显式算法(Δt=0.001 s)的计算结果作为准确解用于对比。图7给出了三种积分算法计算的结构顶层反应时程曲线,可以得到与精度分析相同的结论。5 结论
数值积分算法是求解结构运动方程的有效方法,也是实现实时子结构试验的关键技术之一。本文选取三种基于模型的显式积分算法:Chang法、CR法和RST法,对比分析了各算法在线性系统和非线性系统中的稳定性和精度。稳定性分析可知,三种算法在线性系统和具有刚度软化特性的非线性系统中均是无条件稳定的,在具有刚度硬化特性的非线性系统中变为有条件稳定,且Chang法和RST法的稳定界限随Ω的增加而增大,CR法的稳定界限始终为2。精度分析可知,当结构阻尼比为零时,三种算法均无数值阻尼;当结构阻尼比不为零时,三种算法均存在数值阻尼,并且在同一个Ω下CR法的数值阻尼比的绝对值最小。三种算法的周期延长率均随Ω的增加而增大,并且在同一个Ω下RST法的周期延长率最小。两个算例表明,RST法和Chang法的精度要优于CR法,但因Chang法是半显式的,因此RST法更适于实时子结构试验中数值子结构的仿真计算。
-
-
[1] 陈永盛, 王涛, 张令心, 等. 结构混合试验新进展: 16WCEE相关论文综述[J]. 地震工程与工程振动, 2017, 37(3): 136 − 143. Chen Yongsheng, Wang Tao, Zhang Lingxin, et al. Advance in hybrid simulation of structures: Review of the papers in 16th World Conference on Earthquake Engineering [J]. Earthquake Engineering and Engineering Dynamics, 2017, 37(3): 136 − 143. (in Chinese)
[2] 王进廷, 徐艳杰, 金峰. 实时耦联动力试验理论与实践[M]. 北京: 中国建筑工业出版社, 2014. Wang Jinting, Xu Yanjie, Jin Feng. Real-time dynamic hybrid testing theory and practice [M]. Beijing: China Architecture & Building Press, 2014. (in Chinese)
[3] Mahin S. and Shing P. Pseudo dynamic method for seismic testing [J]. Journal of Engineering Mechanics, 1985, 111(7): 1482 − 1503.
[4] 唐玉. 实时子结构振动台试验技术中若干问题的研究[D]. 上海: 同济大学, 2017. Tang Yu. Theory and numerical analysis on key issues of real-time substructure testing [D]. Shanghai: Tongji University, 2017. (in Chinese)
[5] Chang S Y. Explicit pseudodynamic algorithm with unconditional stability [J]. Journal of Engineering Mechanics, 2002, 128(9): 935 − 947. doi: 10.1061/(ASCE)0733-9399(2002)128:9(935)
[6] Chang S Y. Improved explicit method for structural dynamics [J]. Journal of Engineering Mechanics, 2007, 133(7): 748 − 760. doi: 10.1061/(ASCE)0733-9399(2007)133:7(748)
[7] Chen C, Ricles J M. Development of direct integration algorithms for structural dynamics using discrete control theory [J]. Journal of Engineering Mechanics, 2008, 134(8): 676 − 683. doi: 10.1061/(ASCE)0733-9399(2008)134:8(676)
[8] Gui Y, Wang J T, Jin F, et al. Development of a family of explicit algorithms for structural dynamics with unconditional stability [J]. Nonlinear Dynamics, 2014, 77(4): 1157 − 1170. doi: 10.1007/s11071-014-1368-3
[9] 傅博. 使用基于模型的积分算法的子结构振动台试验方法[D]. 上海: 同济大学, 2017. Fu Bo. Substructure shaking table testing method using model-based integration algorithm [D]. Shanghai: Tongji University, 2017. (in Chinese)
[10] Kolay C, Ricles J M. Development of a family of unconditionally stable explicit direct integration algorithms with controllable numerical energy dissipation [J]. Earthquake Engineering & Structural Dynamics, 2014, 43(9): 1361 − 1380.
[11] Tang Y, Lou M L. New unconditionally stable explicit integration algorithm for real-time hybrid testing [J]. Journal of Engineering Mechanics, 2017, 143(7): 1 − 15.
[12] Tang Y, Lou M L. Closure to “new unconditionally stable explicit integration algorithm for a real-time hybrid testing by Yu Tang and Menglin Lou” [J]. Journal of Engineering Mechanics, 2018, 144(10): 1 − 11.
-
期刊类型引用(11)
1. 康厚军,钱登宇,苏潇阳,张晓宇,丛云跃. 刚构体系斜拉桥塔柱计算长度系数求解方法研究. 工程力学. 2025(01): 44-52 . 本站查看
2. 张晓笛,段冰,吴健,王金昌,杨仲轩,龚晓南,徐荣桥. 混凝土芯水泥土复合桩竖向承载特性分析方法. 岩土力学. 2024(01): 173-183 . 百度学术
3. 江杰,柴文成,张探,龙逸航. 基于T-P模型的桩-桩水平振动相互作用研究. 工程力学. 2023(01): 87-99 . 本站查看
4. 张晓笛,王金昌,杨仲轩,龚晓南,徐荣桥. 基于状态空间法的阶梯型变截面水平受荷桩分析方法. 岩土工程学报. 2023(09): 1944-1952 . 百度学术
5. 黄申,翟恩地,许成顺,孙毅龙. 考虑附加抗力影响的单桩水平受力分析方法. 工程力学. 2022(06): 156-168 . 本站查看
6. 尹平保,王翱,赵衡,杨铠波,赵明华. 斜坡段桥梁基桩受力与变形分析的传递矩阵法. 湖南大学学报(自然科学版). 2022(11): 216-224 . 百度学术
7. 王洋,龚维明,竺明星,戴国亮,万志辉. 砂土中水平循环荷载下单桩式海上风电体系自振频率偏移规律. 中南大学学报(自然科学版). 2022(11): 4372-4380 . 百度学术
8. 张小玲,赵景玖,孙毅龙,许成顺. 基于圆孔扩张理论的桩基水平承载力计算方法. 工程力学. 2021(02): 232-241+256 . 本站查看
9. 王董平,赵成光. 较大水平荷载桩基在储煤棚设计中的应用. 煤炭工程. 2021(01): 33-37 . 百度学术
10. 付臣志,江杰,欧孝夺,龙逸航,张探. 砂土地基中刚性单桩水平极限承载力改进计算方法. 中南大学学报(自然科学版). 2021(10): 3668-3679 . 百度学术
11. 邢康宇,吴文兵,张凯顺,刘浩. 基于改进应变楔模型的大直径桩基水平承载力分析方法. 安全与环境工程. 2020(03): 200-207 . 百度学术
其他类型引用(5)