基于加速度泰勒展开的动力学方程显式积分方法

文 颖1,2,陶 蕤1

(1.中南大学土木工程学院,长沙 410075;2.重载铁路工程结构教育部重点实验室(中南大学),长沙 410075)

摘 要:该文旨在提出兼顾适用性、可靠性与高效性的结构振动时域积分算法。基于加速度的泰勒展开式,引入截断系数考虑高阶项的影响,提出了具有4阶精度的加速度公式;通过积分并考虑典型时间步初始时刻系统动力平衡条件,建立了位移和速度的单步递推公式,运用终止时刻系统运动方程修正加速度。与多步积分法相比,单步积分法无需记录当前时间步以外时刻响应。稳定性分析表明,临界步长相比中心差分法增加40%。通过线性系统振动响应计算发现,当步长-系统固有周期(荷载周期)比达到0.2时,该文方法的振幅衰减率和周期延长率均小于5%;对于非线性系统,为降低算法阻尼和周期误差的影响,需控制步长周期比小于0.1。

关键词:显式时域积分;泰勒展开;稳定性分析;算法阻尼;周期延长

结构动力学从数学分析角度看归结为求解时空耦合的四维初边值问题。为便于问题的解决,采取状态变量(位移、速度和加速度)在空间和时间域分别离散的策略。以有限单元法为代表的数值方法建立了状态变量的空间描述,受结构外形、边界约束和荷载分布等因素影响,一般采用细密网格满足精度要求。因此,状态变量的时域描述方法决定了结构振动响应计算精度和规模[1]。选取较小时间步长虽能得到准确响应,但以显著增加机时和数据存储量(如以中心差分法(下文称 C.D.法)为代表的多步显式积分法)为代价,计算效率低。另一方面,加大步长必须考虑算法稳定性和精度问题。国内外学者在追求时域积分算法高稳定性、高精度和高效率方面已取得丰硕成果,相关探索仍在持续。

现有系统运动状态时域描述方法主要分为 3类:1)在典型时间步(tt+dt)内进行Lagrange插值,逼近t+dt时刻系统状态变量的真实解;2)围绕t时刻状态变量进行泰勒级数展开,按精度要求进行截断;3)采用差分公式求得t+dt时刻系统状态变量的近似值。根据t+dt时刻状态变量是否完全取决于状态变量的历史值,这3类方法又分为显式算法(如第2和第3类,包括显式差分法[2],常氏法[3],Chen-Ricles(下文称 C-R)法[4]等)和隐式算法(如第 1类和第3类,包括Newmark法[5],Wilson-θ[5],Houbolt法[5],HHT-α[6],Bathe 法[7]等)。隐式方法尽管具有良好稳定性,增大步长后将影响分析精度,且计算量较大(需求解联立方程组,对于非线性问题需迭代求解)[8];显式方法的计算量虽小[9],但稳定性和精度稍差。不少学者提出了改进的显式算法:吕和祥等[10]提出了线性、非线性系统动力分析的 4阶精度精细显式积分法;王进廷和杜修力[11]由2次Lagrange位移插值函数建立具有2阶精度的显式差分法,稳定性略逊于 C.D.法;冉田苒等[12]基于HHT-α法提出一种无条件稳定显式算法,由于沿用HHT-α法框架,仍需求解联立方程组;陈学良和袁一凡[13]将Newmark法递推公式和3阶位移泰勒展开式相结合,建立比C.D.法稳定性要求更苛刻的显式积分法;杜晓琼等[14]基于C-R法提出具有2阶精度无条件稳定的显式算法;Noh和 Bathe[15]提出了波动问题显式积分的子步法,具有较高精度。多数研究关注系统自由振动响应计算存在的振幅耗散与周期变化问题,鲜有研究全面分析不同外荷载频率成分和时变/非线性系统参数条件下步长周期比对算法性能的影响。此外,在提高算法精度的同时应明确与之对应的临界步长[16―17]。因为,增加步长即使保证解不发散,而获得具有较高精度的响应是非常困难的[8]。逐步积分法的实质是通过离散时刻系统运动状态分析逼近状态变量的连续变化,当选择的步长无法描述系统固有周期运动及荷载频率特征时,不可能得到准确结果[18]

本文基于t时刻加速度的泰勒展开式,引入截断系数φ考虑高阶项的贡献。通过积分具有4阶精度的加速度近似式并考虑tt+dt时刻系统动力平衡条件,建立了单步显式积分公式。稳定性分析表明,与C.D.法相比,临界步长周期比增加40%。数值算例部分给出了不同阻尼比、荷载频率成分和非线性系统参数条件下单步显式积分法、Newmark-β法、C.D.法和C-R法的振幅衰减率和周期延长率比较,通过列车-桥梁时变系统振动分析算例,得到了响应精度随计算机时的变化规律,说明本文方法具有较高精度和分析效率。

1 加速度的截断泰勒展开式

描述系统运动状态需要确定任意时刻位移、速度和加速度,一般做法是初始化某个状态变量(含待定参数),由微分、积分运算得到其余变量,依据tt+dt时刻系统动力平衡条件更新初始化变量。为了避免微分运算降低计算精度,本文通过积分加速度近似式得到速度和位移。以tt+dt时间步为例,根据泰勒展开定理得到加速度的级数表达式如下:

式中:右上标的希腊数字表示对时间t求导的阶数。取加速度近似式如下:

式中:αβ表示待定系数矩阵;f代表仅依赖τ的列阵(假定t时刻系统运动状态为已知)。令式(1)和式(2)表示的加速度在t+dt时刻相等,可得:

式中:φ称为截断系数,近似考虑o(dt4)的影响,故由式(2)和式(3)确定的加速度具有4阶精度。考虑无阻尼单自由度系统(ω=2π/(rad/s))由初始干扰=0引起的自由振动,建立10个固有周期内峰值响应最大误差与φ的关系(如图1所示)。当dt减小到o(dt4)可以忽略时,φ的取值对精度影响较小;当dt为有限小量时,取φ=0.53能确保不同步长周期比条件下的响应具有较高精度。

图1 截断系数φ取值
Fig.1 Determination of truncation parameterφ

为了确定αβf,首先对t时刻系统运动方程微分3次,得:

式中:mc、和k表示系统质量矩阵、阻尼矩阵和刚度矩阵;p表示外荷载列阵。将式(4)代入式(3),得:

将式(5)代入式(2)并积分,得t+dt时刻速度和位移:

其中:

式中:I表示单位矩阵。需要注意的是,式(6)等号右边均为与t时刻相关的已知项。由式(6)和t+dt时刻系统动力平衡条件,得到修正的加速度列式:

用式(7)修正式(2)的必要性可由下面例子进行说明。仍然考虑前述无阻尼单自由度系统,分别由式(7)和式(2)计算自由振动加速度时程,如图 2所示,修正后的加速度与解析解高度吻合,说明由式(6)和式(7)计算的系统状态变量已具有较高精度,降低了加速度初始化列式近似性带来的影响,提高了整体分析精度。

图2 初始化与修正加速度精度比较(dt/T=0.1)
Fig.2 Comparisons of accuracy between initialized and updated accelerations(dt/T=0.1)

式(6)和式(7)构成了系统状态变量计算的单步显式递推公式。与多步差分法相比,无需调用t时刻以外的状态变量。由式(6)和式(7)可知,本文方法无需求解系数矩阵包含质量、阻尼和刚度特性的大型代数方程组,面对非线性运动方程求解问题也无需迭代,有效提高计算效率。此外,阻尼矩阵列式不局限于瑞雷(Rayleigh)假定,变换阻尼矩阵形式不增加求解负担。然而,式(6)需要计算t时刻荷载的前3阶导数,当荷载无法表示为解析函数或某些时刻荷载导数不存在时,一般使用较高采样频率得到荷载时程的离散值,再通过各阶差分公式得到荷载导数的近似值,如MATLAB软件拥有计算离散数据序列导数的专用函数。当然,计算荷载导数的采样步长与时程分析步长不必一致。

2 算法稳定性分析

为了得到控制积分有界的临界步长,由式(6)和式(7)建立典型时间步状态变量的递推公式:

式中:I表示单位矩阵。根据稳定性的物理概念,经历任意小的干扰后,算法仍能保持t+dt时刻解的有界性,这要求的谱半径ρ()小于 1。为了简化讨论,考虑由初始干扰引起的单自由度系统状态变量递推公式系数矩阵谱半径的计算,故有:

其中:

当忽略系统阻尼(ξ=0)时,可得的非零特征值。

不难得到ρ()=max|λ|与ωdt的关系如图 3所示。相比 C.D.法,无阻尼系统临界步长增加了22.5%,超过部分显式算法的稳定性界限[19]。阻尼加速了系统响应的衰减,对算法稳定性有利,阻尼最多使临界步长增大12.3%。当步长满足ωdt=2.45时,的谱半径最小,系统响应衰减最快。

图3 不同阻尼比条件下递推矩阵谱半径
Fig.3 Spectral radiusρ(A)in cases of distinct damping ratio

3 算法精度分析

本节继续考虑前节无阻尼单自由度系统自由振动响应,揭示振幅衰减率和周期延长率随步长周期比的变化规律。

由图4可知,本文方法的周期延长率可以忽略,即使采用较大步长(dt/T=0.2),计算结果与解析解吻合良好。C.D.法、Newmark-β法和C-R法的计算精度依次下降。各算法振幅衰减率和周期延长率的比较见表1。

表1 各算法关于无阻尼单自由度系统自由振动振幅衰减率和周期延长率的对比
Table 1 Comparison of distinct algorithms in terms of amplitude decay(AD)and period elongation(PE)of undamped free vibrations of single degree-of-freedom system

注:AD:振幅衰减,PE:周期延长(负号表示缩短)。

4 数值算例

下面探讨上节未涉及的因素,如阻尼比、荷载的频谱特征、非线性阻尼和刚度特性以及线性时变阻尼和刚度特性对积分算法性能的影响。

图4 数值积分算法精度对比
Fig.4 Comparison of accuracy of integration algorithms

4.1 不同阻尼条件下单自由度系统自由振动分析

考虑阻尼比取0.05、0.1和0.3时单自由度系统自由振动响应的衰减特征。

1)ξ=0.05,dt/T=0.1

由图5(a)可知,本文结果与解析解完全吻合,C.D.法和Newmark-β法存在周期误差累积,平均每个周期的延长率分别增长-2%和3.5%,C-R法不存在周期误差累积,周期延长率保持5%。Newmark-β法的精度最低,最大峰值响应误差为 5%,其次是C-R法,最大误差为2%,C.D.法的精度与本文方法接近,均不超过1%。

2)ξ=0.1,dt/T=0.2

如图5(b)所示,本文方法仍具有较好的周期保持性和精度,相比无阻尼系统,阻尼加快了 C.D.法和Newmark-β法的周期误差累积速度,周期延长率分别增长-6%和10%,C-R法的周期延长率增加到9%,本文算法的周期延长率为2%。C.D.法的振幅最大误差为10%,Newmark-β法和C-R法的最大峰值响应误差均超过40%。

图5 有阻尼系统数值积分算法性能
Fig.5 Performance of distinct integration algorithms for analyzing damped system

3)ξ=0.3,dt/T=0.32

此时,步长已超过C.D.法的稳定性界限而导致响应发散,如图5(c)所示。本文算法的周期延长率达到 10%,Newmark-β法和 C-R法分别为 31%和-20%。本文方法的振幅衰减率为 2%,Newmark-β法和C-R法的振幅最大计算误差均超过100%且为负,说明算法阻尼与实际阻尼的效应相反[20],前两种阻尼工况也观察到相同的现象。

4.2 任意频率荷载作用下系统振动分析

1)采用文献[21]中两自由度系统参数,计算其强迫振动响应。系统运动方程可表示为:

其中:

式中,ω2表示系统最高阶模态的圆频率。取dt/T1=0.075(T1表示系统卓越周期),各算法性能比较见图6。C-R法、C.D.法和Newmark-β法描述高频振动的能力较差,它们要达到本文算法的精度应降低步长为目前步长的1/6、1/3和1/2。与单自由度系统自由振动分析的结论一致,C.D.法的精度好于Newmark-β法和C-R法。

图6 不同频率荷载组合下算法性能比较(dt/T1=0.075)
Fig.6 Comparison of algorithm performance in the case of mixed loads having distinct frequencies(dt/T1=0.075)

2)将第1种工况的外荷载改成:

式中,ω1表示系统基频。系统将出现共振而使响应逐步发散,如图7所示。即使增大积分步长1倍,C-R法仍表现出较高精度,C.D.法和Newmark-β法则出现“人工”拍振现象,说明两者只能近似描述荷载频率成分(当荷载频率与系统基频接近,但不相等时便出现“拍振”)。Newmark-β法和 C.D.法达到本文算法精度应降低步长为当前步长的 1/5和1/3。

图7 单频率同步荷载作用下算法性能比较(dt/T1=0.15)
Fig.7 Comparison of algorithm performance in the case of synchronous loads having single frequency(dt/T1=0.15)

4.3 非线性刚度和阻尼条件下系统振动分析

1)考虑去掉小位移假定的单摆自由振动问题,其运动方程为:

式中,l是摆长,取ω=1。运用泰勒展开表示sinθ,可得:

式中,带下划线项反映了几何非线性效应对单摆刚度的修正。由此看出,角位移θ的近似值将产生“虚拟”几何刚度,改变单摆固有特性。误差累积是非线性运动方程积分算法面临的挑战。整理式(13)得:

各算法均基于式(14)进行求解。假定初始条件为即单摆摆臂从水平静止位置自由落下。图8(a)和图8(b)给出了考虑和忽略非线性刚度效应条件下各算法计算的振动响应,可以发现非线性刚度效应对各算法振幅衰减率的影响较小,但周期延长率却呈现不同影响规律,本文方法、C.D.法和C-R法的周期延长率略有增大,其中C.D.法和C-R法的周期误差由线性分析时的缩短转变为延长,表明非线性刚度效应削弱了系统刚度特性,有延长系统周期的趋势;Newmark-β法的周期误差反而减小,说明非线性刚度效应产生的周期缩短抵消了线性分析的周期延长效应。

图8 单摆非线性运动方程积分(dt/T=0.08)
Fig.8 Integration of nonlinear equation of motion for a pendulum(dt/T=0.08)

2)考虑如下强非线性耗散系统[10],其运动方程为:

整理式(15)得:

式中:ξ=0.5,ω=1。式(15)包含非线性阻尼和刚度。图9(a)表明本文方法与精确解(dt/T=0.05时的结果)吻合良好,C-R法的周期缩短率初始值达到 7%并以10%的速度递增,Newmark-β法和C.D.法的周期缩短率分别以 5%和 2%的速度递增。C-R法、Newmark-β法和C.D.法达到本文算法精度应分别降低步长为目前的1/50、1/5和1/5。倘若忽略非线性刚度效应,如图9(b)所示,C-R法、Newmark-β法和C.D.法的周期误差均下降,说明相同步长条件下非线性阻尼和刚度效应降低计算精度。

图9 强非线性耗散系统数值积分(dt/T=0.05)
Fig.9 Numerical integration of equation for a strong nonlinear dissipative system(dt/T=0.05)

4.4 列车-桥梁线性时变系统垂向振动分析

考虑 6辆编组列车匀速(v=27.78 m/s)通过单跨简支梁,车辆采用文献[22]中两轴模型模拟,车辆动力学参数如下:车体质量Mv=54 t,转动惯量Iv=13800 t·m2;悬挂参数kv=41350 kN/m,cv=0.0;轮对质量Mw=6.24 t;轮间距d1=17.5 m,相邻车辆轮间距d2=2.5 m。桥梁参数与文献[22]基本一致,只是本文考虑桥梁阻尼比ξ=0.05。为了突出桥梁自振周期对其动力性能的影响,采用模态综合法建立桥梁模型。为了简化分析模型,忽略轨道结构的参振作用,但是计入轨道不平顺的影响,假定轨面粗糙度由短波长单一谐波型不平顺描述,其表达式为:

式中,L表示桥梁跨度。分别运用本文提出的单步显式积分法和Newmark-β法计算桥梁跨中位移时程。如图10所示(N表示选取的模态阶数,T0表示桥梁卓越周期),以 Newmark-β法(N=10,dt/T0=0.001)的结果为准确解,说明本文方法(N=3,dt/T0=0.04)具有良好精度(峰值响应相对误差小于1%),相同条件下,Newmark-β法的峰值响应相对误差达到7%,受步长限制,高频响应难以准确模拟。

在dt/T0=0.04的基础上减小计算步长dt,直到计算结果与准确解完全吻合,记录响应峰值和计算机时。图 11表明,同等精度条件下,本文方法所需计算机时只有Newmark-β法的1/2~2/3。

图10 列车-桥梁线性时变系统垂向振动响应
Fig.10 Vertical vibration response of linear time-varying vehicle-bridge interactive system

图11 桥梁峰值响应精度与计算耗时关系
Fig.11 Relations between the accuracy of peak response and computer time

5 结论

结构动力学方程时域积分算法发展目标是追求计算精度和执行效率的高度平衡。显式积分算法因较高的精度和执行效率引起了广泛关注。本文基于加速度泰勒展开式,建立了具有4阶精度的加速度计算公式;通过积分加速度及考虑tt+dt时刻系统动力平衡条件,得到了系统状态变量递推公式。随后开展了算法稳定性以及不同阻尼比、荷载频率组成和时变/非线性阻尼与刚度条件下精度分析,对算法性能进行了比较。本文提出的单步显式积分法具有如下特点:

(1)与多步显式积分法相比,无需存储t时刻以外的状态变量。系统阻尼无需满足Rayleigh假定,使用任意形式的阻尼矩阵不影响计算效率。

(2)相比C.D.法,控制算法稳定性的临界步长在无阻尼条件下增长22.5%,临界步长随阻尼比的增加而增加,最多增长40%。

(3)步长周期比达到0.2时,无阻尼单自由度系统自由振动振幅衰减率和周期延长率为0.5%。有阻尼系统的振幅计算误差小于 2%,周期延长率随阻尼比的增加而增大,最大值为Newmark-β法的1/3。

(4)对具有任意频率组成荷载的适应性较好,步长-荷载周期比达到0.2,振幅计算误差和周期延长率均小于 5%。非线性阻尼对精度影响较小,非线性刚度的误差敏感性使周期延长率增加,需控制步长周期比小于0.1。

参考文献:

[1]Argyris J H,Mlenjnek H P.Dynamics of structures [M].New York: North Holland Publishing House,1991.

[2]Thompson W T,Dahleh M D.Theory of vibration with applications,fifth edition [M].New York: Pearson Education Press,1997.

[3]Chang S Y.An explicit method with improved stability[J].International Journal for Numerical Methods in Engineering,2009,77(8): 1100-1120.

[4]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.

[5]Bathe K J.Finite element procedures [M].New York:Prentice Hall Publishing Company,1996.

[6]Hilbert H M,Hughes T J R,Taylor R L.Improved numerical dissipation for time integration algorithms in structural dynamics [J].Earthquake Engineering & Structural Dynamics,1977,5(3): 283-292.

[7]Bathe K J,Noh G.Insight into an implicit time integration scheme for structural dynamics [J].Computers & Structures,2012,98(1): 1-6.

[8]Chopra A K.Dynamics of structures [M].5th ed.New York: Prentice Hall Publishing Company,2016.

[9]宋佳,许成顺,杜修力,等.基于精细时程积分的 u-p格式饱和两相介质动力问题的显-显式时域算法[J].工程力学,2017,34(11): 9-17.Song Jia,Xu Chengshun,Du Xiuli,et al.A time-domain explicit-explicit algorithm based on the precise time-integration for solving the dynamic problems of fluid-saturated porous media in u-p form [J].Engineering Mechanics,2017,34(11): 9-17.(in Chinese)

[10]吕和祥,蔡志勤,裘春航.非线性动力学问题的一个显式精细积分算法 [J].应用力学学报,2001,18(2):34-40.Lu Hexiang,Cai Zhiqin,Qiu Chunhang.An explicit precise integration algorithm for nonlinear dynamics problems [J].Chines Journal of Applied Mechanics,2001,18(2): 34-40.(in Chinese)

[11]王进廷,杜修力.有阻尼体系动力分析的一种显式差分法 [J].工程力学,2002,19(3): 109-111.Wang Jinting,Du Xiuli.An explicit difference method for dynamic analysis of a structure system with damping[J].Engineering Mechanics,2002,19(3): 109-111.(in Chinese)

[12]冉田苒,王涛,周惠蒙,等.一种无条件稳定的显式数值积分算法[J].防灾减灾工程学报,2016,36(1): 50-54.Ran Tianran,Wang Tao,Zhou Huimeng,et al.An explicit algorithm with unconditional stability for pseudo dynamic testing [J].Journal of Disaster Prevention and Mitigation Engineering,2016,36(1): 50-54.(in Chinese)

[13]陈学良,袁一凡.求解振动方程的一种显式积分格式及其精度与稳定性[J].地震工程与工程振动,2002,22(3): 9-14.Chen Xueliang,Yuan Yifan.Explicit integration formula for vibration equation and its accuracy and stability [J].Earthquake Engineering and Engineering Vibration,2002,22(3): 9-14.(in Chinese)

[14]杜晓琼,杨迪维,赵永亮.一种无条件稳定的结构动力学显式算法[J].力学学报,2015,47(2): 310-319.Du Xiaoqiong,Yang Diwei,Zhao Yongliang.An unconditional stable explicit algorithm for structural dynamics [J].Chinese Journal of Theoretical and Applied Mechanics,2015,47(2): 310-319.(in Chinese)

[15]Noh G,Bathe K J.An explicit time integration scheme for the analysis of wave propagations [J].Computers & Structures,2013,129(1): 178-193.

[16]袁驷,袁全,闫维明,等.运动方程自适应步长求解的一个新进展——基于EEP超收敛计算的线性有限元法[J].工程力学,2018,35(2): 13-20.Yuan Si,Yuan Quan,Yan Weiming,et al.New development of solution of equations of motion with adaptive time-step size-Linear FEM based on EEP super convergence technique [J].Engineering Mechanics,2018,35(2): 13-20.(in Chinese)

[17]Rossi D F,Ferreira W G,Mansur W J,et al.A review of automatic time-stepping strategies on numerical time integration for structural dynamic analysis [J].Engineering Structures,2014,80(1): 118-136.

[18]曾庆元,周智辉,文颖.结构动力学讲义 [M].北京:人民交通出版社,2015.Zeng Qinyuan,Zhou Zhihui,Wen Ying.Lectures on dynamics of structures [M].Beijing: China Communication Press,2015.(in Chinese)

[19]刘祥庆,刘晶波,丁桦.时域逐步积分算法稳定性与精度的对比分析 [J].岩石力学与工程学报,2007,26增刊1): 3001-3008.Liu Xiangqing,Liu Jingbo,Ding Hua.Comparative analysis of stabilization and accuracy of step-by-step time integration methods [J].Chinese Journal of Rock Mechanics and Engineering,2007,26(Suppl1): 3001-3008.(in Chinese)

[20]王进廷,张楚汉,金峰.有阻尼动力方程显示积分方法的精度研究 [J].工程力学,2006,23(3): 1-5.Wang Jinting,Zhang Chuhan,Jin Feng.On the accuracy of several explicit integration schemes for dynamic equation with damping [J].Engineering Mechanics,2006,23(3): 1-5.(in Chinese)

[21]杨超,肖守讷,阳光武,等.一类非耗散的显式时间积分方法[J].振动工程学报,2015,28(3): 441-449.Yang Chao,Xiao Shoune,Yang Guangwu,et al.Non-dissipative explicit time integration methods of the same class [J].Journal of Vibration Engineering,2015,28(3): 441-449.(in Chinese)

[22]Yang Y B,Wu Y S.A versatile element for analyzing vehicle-bridge interaction response [J].Engineering Structures,2001,23(5): 452-469.

AN EXPLICIT TIME-DOMAIN INTEGRATION SCHEME FOR SOLVING EQUATIONS OF MOTION IN STRUCTURAL DYNAMICS BASED ON A TRUNCATED TAYLOR EXPANSION OF ACCELERATION

WEN Ying1,2,TAO Rui1

(1.School of Civil Engineering,Central South University,Changsha 410075,China;2.The Key Laboratory of Engineering Structures of Heavy Haul Railway(Central South University),Ministry of Education,Changsha 410075,China)

Abstract:The aim of this paper is to present a novel time integration algorithm with a high level of balance among applicability and reliability and computational efficiency for the dynamic analysis of structures.A formula for approximating acceleration with a forth-order degree of accuracy has been developed,based on the Taylor expansion approach.In applying the Taylor expansion method,a truncation parameter is defined to consider the contributions of high-order terms upon the accuracy of predicted results.Through an integration of the obtained acceleration and considering the dynamic equilibrium condition at the initial state of a typical time step,a single-step equation for computing displacement and velocity at the end state is correspondingly developed.A revised acceleration can be obtained from the calculated displacement and velocity through the equations of motion at the end state.In this regard,as compared with the multiple-step integration scheme,it is not required for the present method to temporarily record the state variables of previous steps.From the results of stability analysis,the maximum step length to period ratio within which the obtained responses remain bounded has been increased by 40% in comparison to the central difference method.By carrying out a series of numerical analyses for the purpose of demonstration,it is generally observed from the natural and forced vibration investigations for linear systems that the computational amplitude decay and period elongation were less than 5% even if the ratio between the time step length and system inherent period/load period mounts to 0.2.However,to reduce the effects of amplitude decay and period distortion for the time integration of nonlinear systems,the magnitude of the above mentioned ratio should generally be restricted below 0.1.

Key words:explicit time integration; Taylor expansion; stability analysis; algorithm damping; period elongation

作者简介:陶蕤(1993―),男,湖北人,硕士生,从事车-桥系统垂向振动研究(E-mail: taorui_2014@163.com).

通讯作者:文颖(1981―),男,湖南人,副教授,博士,重载铁路工程结构教育部重点实验室副主任,从事桥梁稳定极限承载力及车桥系统振动稳定性研究(E-mail: ywen_ce@csu.edu.cn).

基金项目:国家自然科学基金高铁联合基金重点项目(U1534206);湖南省科技计划项目(2014FJ6036)

收稿日期:2017-08-30;修改日期:2018-03-19

文章编号:1000-4750(2018)11-0026-09

doi:10.6052/j.issn.1000-4750.2017.08.0661

文献标志码:A

中图分类号:TU311.3