CN115453873A - 非线性结构动力学系统瞬态响应无条件稳定时间积分方法 - Google Patents
非线性结构动力学系统瞬态响应无条件稳定时间积分方法 Download PDFInfo
- Publication number
- CN115453873A CN115453873A CN202211126469.0A CN202211126469A CN115453873A CN 115453873 A CN115453873 A CN 115453873A CN 202211126469 A CN202211126469 A CN 202211126469A CN 115453873 A CN115453873 A CN 115453873A
- Authority
- CN
- China
- Prior art keywords
- nonlinear
- time
- term
- equation
- parameter
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 131
- 230000004044 response Effects 0.000 title claims abstract description 50
- 230000010354 integration Effects 0.000 title claims abstract description 36
- 230000001052 transient effect Effects 0.000 title claims abstract description 32
- 238000013016 damping Methods 0.000 claims abstract description 39
- 230000008878 coupling Effects 0.000 claims abstract description 30
- 238000010168 coupling process Methods 0.000 claims abstract description 30
- 238000005859 coupling reaction Methods 0.000 claims abstract description 30
- 238000004458 analytical method Methods 0.000 claims abstract description 21
- 108010074506 Transfer Factor Proteins 0.000 claims abstract description 7
- 238000005516 engineering process Methods 0.000 claims abstract description 7
- 238000005457 optimization Methods 0.000 claims description 32
- 238000006073 displacement reaction Methods 0.000 claims description 25
- 230000001133 acceleration Effects 0.000 claims description 18
- 239000011159 matrix material Substances 0.000 claims description 15
- 238000004422 calculation algorithm Methods 0.000 claims description 8
- 239000013598 vector Substances 0.000 claims description 8
- 238000012546 transfer Methods 0.000 claims description 5
- 230000005540 biological transmission Effects 0.000 claims description 4
- 230000005284 excitation Effects 0.000 claims description 4
- 150000001875 compounds Chemical class 0.000 claims description 2
- 230000000087 stabilizing effect Effects 0.000 claims 1
- 230000008569 process Effects 0.000 abstract description 4
- 238000000556 factor analysis Methods 0.000 abstract description 2
- 238000004364 calculation method Methods 0.000 description 10
- 238000012937 correction Methods 0.000 description 9
- 239000006185 dispersion Substances 0.000 description 8
- 238000004134 energy conservation Methods 0.000 description 8
- 230000008901 benefit Effects 0.000 description 6
- 238000004088 simulation Methods 0.000 description 6
- 230000001419 dependent effect Effects 0.000 description 4
- 238000005183 dynamical system Methods 0.000 description 4
- 239000000463 material Substances 0.000 description 3
- 230000008859 change Effects 0.000 description 2
- 238000010276 construction Methods 0.000 description 2
- 238000013461 design Methods 0.000 description 2
- 238000005381 potential energy Methods 0.000 description 2
- 230000009471 action Effects 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 230000003111 delayed effect Effects 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000000126 substance Substances 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G05—CONTROLLING; REGULATING
- G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
- G05B13/00—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion
- G05B13/02—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion electric
- G05B13/04—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion electric involving the use of models or simulators
- G05B13/042—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion electric involving the use of models or simulators in which a parameter or coefficient is automatically adjusted to optimise the performance
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Health & Medical Sciences (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Computation (AREA)
- Medical Informatics (AREA)
- Software Systems (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Automation & Control Theory (AREA)
- Complex Calculations (AREA)
Abstract
本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法,属于非线性结构动力学系统响应分析领域。本发明实现方法如下:建立包含非线性几何项、非线性阻尼项以及它们的耦合项的结构动力学方程;结合多分步技术将时间单元离散成三个分步,并设计各分步状态变量更新公式;基于BN稳定性理论、传递因子分析方法、局部截断误差分析方法,以效率最大化、二阶精度、可控耗散和无条件稳定性为目标设计分步公式中的参数;建立三分步时间积分方法在包含非线性几何项、非线性阻尼项以及它们的耦合项在结构动力学系统中的求解流程,实现对包含非线性几何项、非线性阻尼项以及它们的耦合项的结构动力学系统瞬态响应的高效、精确且稳定地时域分析。
Description
技术领域
本发明属于非线性结构动力学系统响应分析领域,尤其涉及一种包含非线性几何项、非线性阻尼项及其耦合项的结构动力学系统的具有二阶精度、可调数值耗散以及无条件稳定性的时间积分方法。
背景方法
随着新技术、新材料、新结构等前沿成果在航空航天、土木工程、机器人等领域的广泛应用,由此带来的材料非线性、几何非线性和接触非线性因素给结构动力学分析带来了严峻挑战。目前,最为常用的结构动力学分析方案是“空间离散”与“时间离散”相结合,其也被广泛应用在商业软件CAE(Computer Aided Engineering)中。
灵活有效的有限单元法(Finite Element Method,FEM)是目前较为流行的空间离散工具,其作用是将由偏微分方程控制的无穷自由度连续系统转变为由常微分方程控制的有限自由度离散系统。完成空间离散后,下一个环节就是时间离散,基于差分思想建立的时间积分方法(Time integration method,TIM)是一种强大的时间离散工具。其求解思路是将待求的时间域离散成一系列连续的时间区间,人为假设状态变量在一个时间区间内的变化规律,从上一时刻已知状态变量信息求解当前时刻未知信息。
Newmark家族方法是时间积分方法中的经典工作,其成员包括著名的梯形法则(Trapezoidal Rule,TR)和中心差分方法(Central Difference Method,CDM)等。对于线性系统,梯形法则是无条件稳定的,但对于简单的非线性问题如刚性单摆转动问题,梯形法则出现了发散现象。由此引发当前时间积分方法领域的热点和难点问题:如何提高时间积分方法在求解非线性结构动力学系统时的稳定性,以及能否为非线性结构动力学系统设计出无条件稳定的时间积分方法。为了解决上述问题,学者们发展了耗散型方法、保能量方法以及结构依赖型方法。
大量的数值结果表明通过引入数值耗散可以有效地推迟时间积分方法在非线性结构动力学系统求解中失稳的时刻。由于Newmark方法的耗散格式仅具有一阶精度,因此学者们通过加权Newmark方法的平衡方程构造出了同时具有二阶精度和数值耗散的参数类方法,比如广义-alpha方法、Wood-Bossak-Zienkiewicz-alpha方法、Hiber-Hughes-Taylor-alpha方法等。但是,由于平衡方程在时间结点上没有得到满足,导致该类方法的加速度只具有一阶精度。为了设计出具有较高低频精度的耗散型时间积分方法,线性多步法和复合方法被提出,这两类方法均有效地提高了精度。但是,线性多步法和复合方法并不能从根本上改变方法的稳定性。总之,耗散型方法仍然存在失稳问题,并且其低频精度比非耗散型方法低。
为了能够给非线性结构动力学系统提供稳定的数值仿真结果,基于能量有界准则的保能量方法被提出。能量有界准则由Belytschko提出,其定义是:当前时刻的动能和势能之和小于等于上一时刻的动能、势能和外力功之和。因此,保能量方法是通过强制系统能量守恒来实现无条件稳定性。受构造方法的限制,现有的保能量方法大都只能处理含有非线性几何项的结构动力学问题,比如Krenk方法和EMM等,而只有少部分如ECM能够处理非线性几何项和阻尼项独立存在的非线性结构动力学问题。目前,能够处理含有非线性几何项和阻尼项耦合的结构动力学系统的保能量方法尚未被构造出来。综上,保能量方法虽然在求解非线性结构动力学系统时不会发散,但其适用对象有限,且具有额外计算量(如离散能量函数修正运算等)致使计算效率降低。
对大型复杂结构进行动力学响应分析时,精度、效率和稳定性三者之间往往难以平衡。为了解决这个问题,2002年Chang提出了结构依赖型时间积分方法,这类方法的算法参数与结构初始特性和时间步长大小紧密相关。为了进一步提高Chang方法的数值性能,包括精度、稳定性和耗散特性,一系列性能更加优异的结构依赖型时间积分方法被提出,比如KR方法、CR方法、Du-Yang-Zhao方法等。理论分析和数值实验结果表明,该类方法对线性系统是无条件稳定的,同时对非线性软刚度系统具有理想的稳定性。目前,受构造方式限制,这类方法局限于具有对角质量矩阵、正定对称刚度和阻尼矩阵的线性系统和具有幂形式非线性刚度的结构动力学系统。总之,结构依赖型方法虽然有效地平衡了精度、效率和稳定性,但其在求解非线性结构动力学系统时仍有失稳的可能,且适用对象有限。
因此,针对包含非线性几何项、非线性阻尼项及其耦合项的结构动力学系统,设计可高精度、高效率、无条件稳定地分析瞬态响应的时间积分方法是现阶段亟待突破的研究内容。其将为非线性结构动力学问题提供性能优越的新解法,且可以促进时间积分方法朝着同时改进稳定性、精度和效率方向发展,服务于国家对复杂结构动力学分析高性能算法的需求,具有重要的理论意义和应用价值。
发明内容
本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法主要目的是:利用有限单元方法建立包含非线性几何项、非线性阻尼项以及它们的耦合项的结构动力学方程;结合多分步技术将一个时间单元离散成三个分步,并设计各分步状态变量(包括位移、速度和加速度)的更新公式;基于BN稳定性理论、传递因子分析方法、局部截断误差分析方法,以效率最大化、二阶精度、可控耗散和无条件稳定性为目标设计分步公式中的参数;在此基础上,建立三分步时间积分方法在包含非线性几何项、非线性阻尼项以及它们的耦合项在结构动力学系统中的求解流程,实现对包含非线性几何项、非线性阻尼项以及它们的耦合项的结构动力学系统瞬态响应的高效、精确且稳定地时域分析。
本发明的目的是通过如下技术方法实现。
本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法,通过利用有限单元方法建立包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学方程;将一个时间单元划分成三个分步,基于差分格式建立状态变量在一个时间单元内的演化规律,依次使用n点后向插值格式(n=1,2,3),建立便于数值性能调控的多参数优化时间步进方程,由于调控参数增多,从而提高参数调节的可控性和精度。所述数值性能包括效率、精度、耗散、稳定性。基于构建的多参数优化时间步进方程,推导便于数值性能调控的三个分步的有效刚度矩阵,以效率最大化为目标建立具有统一系数的参数约束关系,利用建立参数约束关系将所述三个分步的有效刚度矩阵变换为具有统一形式的有效刚度矩阵,提高非线性系统的结构动力学响应时域分析效率。基于构建的多参数优化时间步进方程,针对包含非线性几何项、非线性阻尼项以及它们的耦合项的结构动力学系统,以二阶精度为目标建立局部截断误差为O(Δt3)的参数约束关系,以数值耗散可控为目标建立高频段传递因子可调的参数约束关系,以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系,从而提高非线性系统的结构动力学响应时域分析精度、耗散性和稳定性。将得到的参数约束关系代入建立的多参数优化时间步进方程,再将所述多参数优化时间步进方程代入构建的包含非线性几何项、非线性阻尼项及其耦合项的结构动力学控制方程,实现无条件稳定地分析包含非线性几何项、非线性阻尼项及其耦合项的结构动力学系统的瞬态响应。
本发明公开的非线性结构系统瞬态响应无条件稳定时间积分方法,包括以下步骤:
步骤1:基于有限单元方法FEM建立空间离散的包含非线性几何项、非线性阻尼项及其耦合项的结构动力学控制方程。
基于有限单元方法FEM(Finite Element Method)建立空间离散的包含非线性几何项、非线性阻尼项及其耦合项的结构动力学控制方程,如下所示:
步骤2:基于差分格式建立状态变量在一个时间单元内的演化规律,将一个时间单元划分成三个分步,依次使用n点后向插值格式(n=1,2,3),建立便于数值性能调控的多参数优化时间步进方程,由于调控参数增多,提高参数调节的可控性和精度。所述数值性能包括效率、精度、耗散、稳定性。
基于差分格式,建立状态变量,包括位移、速度和加速度,在一个时间单元[t,t+Δt]内随时间t的演化规律。将一个时间单元[t,t+Δt]划分成三个分步,即[t,t+c1Δt]、[t+c1Δt,t+c2Δt]以及[t+c2Δt,t+c3Δt]。第一个分步[t,t+c1Δt]的步进公式为
其中,xt+c1Δt、和代表t+c1Δt时刻的位移、速度和加速度;xt和代表t时刻的位移和速度;Δt代表时间步长大小;c1为第一个分步引入的算法参数。第二个分步[t+c1Δt,t+c2Δt]的步进公式为
其中,xt+c3Δt、和代表t+c3Δt时刻的位移、速度和加速度;c3、β和η为第三个分步引入的算法参数。时间单元终点时刻t+Δt信息通过对三个已知时刻信息t+c1Δt、t+c2Δt以及t+c3Δt进行加权得到,对应的步进方程为
其中,b1和b2为待定加权参数。
步骤3:基于步骤2构建的多参数优化时间步进方程,推导便于数值性能调控的三个分步的有效刚度矩阵,以效率最大化为目标建立具有统一系数的参数约束关系,利用建立参数约束关系将所述三个分步的有效刚度矩阵变换为具有统一形式的有效刚度矩阵,提高非线性系统的结构动力学响应时域分析效率。根据具有统一形式的有效刚度矩阵得到统一系数的参数约束关系。
从公式(7)至(9)中发现,当参数满足公式(10)中的关系时,三个分步的有效刚度矩阵完全相同。此时,三个分步使用同一个有效刚度矩阵,进而降低计算量。
c1=c2α=c3η (10)
另外,要求时间点t+c1Δt和t+c3Δt关于时间单元中点t+0.5Δt对称,且对应的加权参数相同。
由此可得,
c1+c3=1 (11)
b1=1-b1-b2 (12)
步骤4:基于步骤2构建的多参数优化时间步进方程,针对包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学系统,以二阶精度为目标建立局部截断误差为O(Δt3)的参数约束关系,以数值耗散可控为目标建立高频段传递因子可调的参数约束关系,以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系,从而提高非线性系统的结构动力学响应时域分析的精度、耗散性和稳定性。
步骤4.1:针对包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学系统,推导便于数值性能调控的一个时间单元内的Butcher表。首先,将公式(2)至(5)等价地写为基于一阶方程的形式,如下所示
步骤4.2:在步骤4.1的基础上,建立的基于一阶方程的多参数优化时间步进方程推导其在一个时间单元内的传递方程,如下所示
zt+Δt=A(τ)zt,τ=λΔt (18)
其中,λ代表特征根。式中,传递因子A的表达式为
步骤4.3:在步骤4.1的基础上,基于局部截断误差确定多参数优化时间步进方程的精度阶次,局部截断误差的定义为
σ=A(τ)-Aexact(τ),Aexact(τ)=exp(τ) (20)
根据二阶精度,要求
A(0)=A(1)(0)=A(2)(0)=1 (21)
根据式(21)求解得到如下满足二阶精度要求的时间步进方程的参数关系,即
步骤4.4:在步骤4.1的基础上,以数值耗散可控为目标建立高频段传递因子可调的参数约束关系,如下所示,
对于非耗散格式ρ∞=1,根据式(23)建立如下参数关系,
对于耗散格式0≤ρ∞<1,根据式(23)建立如下参数关系
步骤4.5:在步骤4.1的基础上,以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系。根据BN稳定性定义,建立如下不等式约束条件,即
S1=b1(2c1-b1)≥0 (26)
S2=b2(2c2α-b2)≥0 (27)
S3=(1-b1-b2)[2c3η-(1-b1-b2)]≥0 (28)
S5=S1S3-(1-b1-b2)2[c3(1-β-η)-b1]2≥0 (29)
S6=S2S3-(1-b1-b2)2(c3β-b2)2≥0 (30)
上述不等式约束条件,即公式(26)至公式(31),能够保证在任意时间步长下稳定地计算包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学系统。根据公式(26)至(31),可建立β和c3之间的关系。对于非耗散格式ρ∞=1,根据公式(26)至公式(31),建立如下参数关系,
对于耗散格式0≤ρ∞<1,根据公式(26)—公式(31),建立如下参数关系,
式中,c3=f(ρ∞)为关于高频段耗散量的已知函数。
步骤5:根据步骤3得到具有统一系数的参数约束关系,根据步骤4得到满足二阶精度、可控数值耗散和无条件稳定性的参数约束关系,将步骤3和步骤4得到的参数约束关系代入步骤2建立的多参数优化时间步进方程,将所述多参数优化时间步进方程结合Newton迭代技术,代入步骤1构建的包含非线性几何项、非线性阻尼项及其耦合项的结构动力学控制方程,实现高效、精确、且无条件稳定地计算包含非线性几何项、非线性阻尼项及其耦合项的结构动力学系统瞬态响应。
还包括步骤6:根据步骤1至步骤5得到的适用于包含非线性几何项、非线性阻尼项及其耦合项的结构动力学系统的无条件稳定时间积分方法,能够在结构动力学工程应用中广泛应用于非线性结构动力学系统的瞬态响应分析,改善并预测非线性结构动力学结构和性能,解决相关工程技术问题。
有益效果:
1、本发明公开的一种非线性结构动力学系统瞬态响应的无条件稳定时间积分方法,基于差分格式建立状态变量在一个时间单元内的演化规律,将一个时间单元划分成三个分步,依次使用n点后向插值格式(n=1,2,3),建立便于数值性能(包括效率、精度、耗散性和稳定性)调控的多参数优化时间步进方程,由于调控参数增多,提高参数调节的可控性和精度。
2、本发明公开的一种非线性结构动力学系统瞬态响应的无条件稳定时间积分方法,基于步骤2构建的多参数优化时间步进方程,推导便于数值性能调控的三个分步的有效刚度矩阵,以效率最大化为目标建立具有统一系数的参数约束关系,利用建立参数约束关系将所述三个分步的有效刚度矩阵变换为具有统一形式的有效刚度矩阵,提高非线性结构动力学系统瞬态响应的时域分析效率。
3、本发明公开的一种非线性结构动力学系统瞬态响应的无条件稳定时间积分方法,基于步骤2构建的多参数优化时间步进方程,针对包含非线性几何项、非线性阻尼项以及它们耦合项的非线性结构动力学系统,以二阶精度为目标建立满足局部截断误差为O(Δt3)的参数约束条件;以耗散可控为目标建立高频段传递因子可调的参数约束关系;以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系,提高包含非线性几何项、非线性阻尼项以及它们耦合项的非线性结构动力学系统瞬态响应的时域分析的精度、耗散性和稳定性。
4、本发明公开的一种非线性结构动力学系统瞬态响应的无条件稳定时间积分方法,相较于已有的耗散型方法、保能量方法以及结构依赖型方法,计算量小、精度高、稳定性强,方便使用。即使在缺乏专业知识背景的情况下也能进行操作,能够在结构动力学工程应用中广泛应用于非线性结构动力学系统的瞬态响应分析,改善并预测非线性结构动力学结构和性能,解决相关工程技术问题。
附图说明:
图1为本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法的计算流程图;
图2为考虑非线性格林应变的悬臂梁冲击模型;
图3为广义-alpha方法、ρ∞-Bathe方法和本发明方法在悬臂梁自由端沿z方向位移和速度随时间演化规律(初始速度为20m/s);
图4为广义-alpha方法、ρ∞-Bathe方法和本发明方法在悬臂梁自由端沿z方向位移和速度随时间演化规律(初始速度为40m/s);
图5为弹塑性悬臂板冲击模型;
图6为广义-alpha方法、ρ∞-Bathe方法和本发明方法在悬臂板自由端沿x方向位移随时间演化规律;
图7为广义-alpha方法、ρ∞-Bathe方法和本发明方法在悬臂板自由端沿z方向位移随时间演化规律。
具体实施方式
为了更好地说明本发明的目的和优点,下面结合附图和实例对发明内容做进一步说明。其中,自始至终相同或类似的符号表示相同或类似功能。
实施例1:
如图1所示,本实施例公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法,应用于预测考虑非线性格林应变的悬臂梁冲击模型的动力学响应,具体实现步骤如下:
步骤1:以考虑非线性格林应变的悬臂梁冲击模型为例,如图2所示,来说明本发明方法的稳定性和效率优势。基于有限单元法,利用10个非线性双节点梁单元对考虑非线性格林应变悬臂梁冲击模型进行空间离散,进而可得如下非线性结构动力学方程
其中,v0为在自由端施加的初始速度。
步骤2:在步骤1的基础上,将待求时间域[0,0.5]离散成n个连续的时间单元,即[t,t+Δt]、[t+Δt,t+2Δt]、…、[t+(n-1)Δt,t+nΔt],其中Δt代表一个时间单元的尺寸。基于差分格式,建立状态变量,包括位移、速度和加速度,在一个是时间单元[t,t+Δt]内随时间t的演化规律。将一个时间单元[t,t+Δt]划分成三个分步,即[t,t+c1Δt]、[t+c1Δt,t+c2Δt]以及[t+c2Δt,t+c3Δt]。
第一个分步[t,t+c1Δt]的步进公式为
第二个分步[t+c1Δt,t+c2Δt]的步进公式为
第三个分步[t+c2Δt,t+c3Δt]的步进公式为
时间单元终点时刻t+Δt信息通过对三个已知时刻信息t+c1Δt、t+c2Δt以及t+c3Δt进行加权得到,对应的步进方程为
步骤3:根据具有统一形式的有效刚度矩阵得到统一系数的参数约束关系,如下所示
c1=c2α=c3η (39)
另外,要求时间点t+c1Δt和t+c3Δt关于时间单元中点t+0.5Δt对称,且对应的加权参数相同。由此可得,
c1+c3=1 (40)
b1=1-b1-b2 (41)
步骤4:基于步骤2构建的多参数优化时间步进方程,推导便于数值性能调控的一个时间单元内的Butcher表,如下所示,
推导便于数值性能调控的一个时间单元内的传递因子公式,如下所示
基于局部截断误差确定多参数优化时间步进方程的精度阶次,局部截断误差的定义为
σ=A(τ)-Aexact(τ),Aexact(τ)=exp(τ) (44)
根据二阶精度,要求
A(0)=A(1)(0)=A(2)(0)=1 (45)
根据式(21)求解得到如下满足二阶精度要求的时间步进方程的参数关系,即
以耗散可控为目标建立高频段传递因子可调的参数约束关系,如下所示
对于非耗散格式ρ∞=1,根据式(47)可建立如下参数关系,
对于耗散格式0≤ρ∞<1,根据式(47)可建立如下参数关系
以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系,如下所示
S1=b1(2c1-b1)≥0 (50)
S2=b2(2c2α-b2)≥0 (51)
S3=(1-b1-b2)[2c3η-(1-b1-b2)]≥0 (52)
S5=S1S3-(1-b1-b2)2[c3(1-β-η)-b1]2≥0 (53)
S6=S2S3-(1-b1-b2)2(c3β-b2)2≥0 (54)
根据公式(51)—(55),可建立β和c3之间的关系。对于非耗散格式ρ∞=1,整理可得
对于耗散格式0≤ρ∞<1,整理可得
式中,c3=f(ρ∞)为关于高频段耗散量的已知函数。
考虑非线性格林应变的悬臂梁冲击模型是保守系统。因此,根据动力学特性,非耗散格式ρ∞=1被采用,由公式(39)、(40)、(41)、(46)、(48)和(56)可确定其余参数为
步骤5:将步骤3和步骤4得到的参数约束关系代入步骤2建立的多参数优化时间步进方程,结合Newton迭代技术和所述多参数优化时间步进方程由初始时刻响应依次递推求解步骤1建立的考虑非线性格林应变的悬臂梁冲击模型在每个时间离散点处的状态变量信息。
在第一个分步[t,t+c1Δt]中,首先假定
其中,上标“0”表示迭代次数。将公式(59)代入第一个分步的步进方程(35),可得
进而可得修正的加速度为
在第二个分步[t+c1Δt,t+c2Δt]中,首先假定
其中,上标“0”表示迭代次数。将公式(63)代入第二个分步的步进方程(36),可得
进而可得修正的加速度为
在第三个分步[t+c2Δt,t+c3Δt]中,首先假定
其中,上标“0”表示迭代次数。将公式(67)代入第三个分步的步进方程(37),可得
进而可得修正的加速度为
图3和图4分别表示v0=20m/s和40m/s时自由端沿y方向的位移和速度,随着v0的增大,该动力学模型能够激发出更多的高频模式从而导致时间积分方法更容易失稳。由数值模拟结果可以得出:本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法在整个仿真过程中都给出了稳定且精确的结果,而常用的ρ∞-Bathe方法在早期仿真阶段得到了发散的结果。此外,表1将已有方法和本发明方法的计算效率进行了比较,可以发现本发明方法的计算成本最低。由实施例1可以说明,本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法比传统方法具有更强的稳定性和更高的计算效率。
表1广义-alpha方法、ρ∞-Bathe方法和本发明方法的CPU和迭代次数(容许误差ε=10e-12)
实施例2:
如图1所示,本实施例公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法,应用于预测弹塑性悬臂板在集中冲击力作用下的动力学响应,具体实现步骤如下:
步骤1:以弹塑性悬臂板为例,如图5所示,来说明本发明方法的精度优势。基于有限单元法,利用1296个考虑J2双线性次弹塑性材料板单元对悬臂板模型进行空间离散,进而可得如下非线性结构动力学方程
步骤2:在步骤1的基础上,将待求时间域[0,1]离散成n个连续的时间单元,即[t,t+Δt]、[t+Δt,t+2Δt]、…、[t+(n-1)Δt,t+nΔt],其中Δt代表一个时间单元的尺寸。基于差分格式,建立状态变量,包括位移、速度和加速度,在一个是时间单元[t,t+Δt]内随时间t的演化规律。将一个时间单元[t,t+Δt]划分成三个分步,即[t,t+c1Δt]、[t+c1Δt,t+c2Δt]以及[t+c2Δt,t+c3Δt]。
第一个分步[t,t+c1Δt]的步进公式为
第二个分步[t+c1Δt,t+c2Δt]的步进公式为
第三个分步[t+c2Δt,t+c3Δt]的步进公式为
时间单元终点时刻t+Δt信息通过对三个已知时刻信息t+c1Δt、t+c2Δt以及t+c3Δt进行加权得到,对应的步进方程为
步骤3:根据具有统一形式的有效刚度矩阵得到统一系数的参数约束关系,如下所示
c1=c2α=c3η (76)
另外,要求时间点t+c1Δt和t+c3Δt关于时间单元中点t+0.5Δt对称,且对应的加权参数相同。由此可得,
c1+c3=1 (77)
b1=1-b1-b2 (78)
步骤4:基于步骤2构建的多参数优化时间步进方程,推导便于数值性能调控的一个时间单元内的Butcher表,如下所示,
推导便于数值性能调控的一个时间单元内的传递因子公式,如下所示
基于局部截断误差确定多参数优化时间步进方程的精度阶次,局部截断误差的定义为
σ=A(τ)-Aexact(τ),Aexact(τ)=exp(τ) (81)
根据二阶精度,要求
A(0)=A(1)(0)=A(2)(0)=1 (82)
根据式(82)求解得到如下满足二阶精度要求的时间步进方程的参数关系,即
以耗散可控为优化目标建立高频段传递因子可调的参数约束关系,如下所示
对于非耗散格式ρ∞=1,根据式(84)可建立如下参数关系,
对于耗散格式0≤ρ∞<1,根据式(84)可建立如下参数关系
以无条件稳定性为优化目标建立满足BN稳定性条件的参数约束关系,如下所示
S1=b1(2c1-b1)≥0 (87)
S2=b2(2c2α-b2)≥0 (88)
S3=(1-b1-b2)[2c3η-(1-b1-b2)]≥0 (89)
S5=S1S3-(1-b1-b2)2[c3(1-β-η)-b1]2≥0 (90)
S6=S2S3-(1-b1-b2)2(c3β-b2)2≥0 (91)
根据公式(87)—(92),可建立β和c3之间的关系。对于非耗散格式ρ∞=1,整理可得
对于耗散格式0≤ρ∞<1,整理可得
式中,c3=f(ρ∞)为关于高频段耗散量的已知函数。
考虑非线性格林应变的悬臂梁冲击模型是保守系统。因此,根据动力学特性,非耗散格式ρ∞=1被采用,由公式(76)、(77)、(78)、(83)、(85)和(93)可确定其余参数为
步骤5:将步骤3和步骤4得到的参数约束关系代入步骤2建立的多参数优化时间步进方程,结合Newton迭代技术和所述多参数优化时间步进方程由初始时刻响应依次递推求解步骤1建立的考虑非线性格林应变的悬臂梁冲击模型在每个时间离散点处的状态变量信息。
在第一个分步[t,t+c1Δt]中,首先假定
其中,上标“0”表示迭代次数。将公式(96)代入第一个分步的步进方程(72),可得
进而可得修正的加速度为
在第二个分步[t+c1Δt,t+c2Δt]中,首先假定
其中,上标“0”表示迭代次数。将公式(100)代入第二个分步的步进方程(73),可得
进而可得修正的加速度为
在第三个分步[t+c2Δt,t+c3Δt]中,首先假定
其中,上标“0”表示迭代次数。将公式(104)代入第三个分步的步进方程(74),可得
进而可得修正的加速度为
图6和图7分别绘制了A点(自由端中点)在x和z方向的位移,由数值模拟结果可以得出:本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法在整个仿真过程中都给出了最高的计算精度。由实施例2表明,本发明公开的非线性结构动力学系统瞬态响应无条件稳定时间积分方法比传统方法具有更高的计算精度。
以上所述的具体描述,对发明的目的、技术方案和有益效果进行了进一步详细说明,所应理解的是,以上所述仅为本发明的具体实施例而已,并不用于限定本发明的保护范围,凡在本发明的精神和原则之内,所做的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。
Claims (5)
1.非线性结构动力学系统瞬态响应无条件稳定时间积分方法,其特征在于:包括以下步骤,
步骤1:基于有限单元方法FEM建立空间离散的包含非线性几何项、非线性阻尼项及其耦合项的结构动力学控制方程;
步骤2:基于差分格式建立状态变量在一个时间单元内的演化规律,将一个时间单元划分成三个分步,依次使用n点后向插值格式(n=1,2,3),建立便于数值性能调控的多参数优化时间步进方程,由于调控参数增多,提高参数调节的可控性和精度;所述数值性能包括效率、精度、耗散、稳定性;
步骤3:基于步骤2构建的多参数优化时间步进方程,推导便于数值性能调控的三个分步的有效刚度矩阵,以效率最大化为目标建立具有统一系数的参数约束关系,利用建立参数约束关系将所述三个分步的有效刚度矩阵变换为具有统一形式的有效刚度矩阵,提高非线性系统的结构动力学响应时域分析效率;根据具有统一形式的有效刚度矩阵得到统一系数的参数约束关系;
步骤4:基于步骤2构建的多参数优化时间步进方程,针对包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学系统,以二阶精度为目标建立局部截断误差为O(Δt3)的参数约束关系,以数值耗散可控为目标建立高频段传递因子可调的参数约束关系,以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系,从而提高非线性系统的结构动力学响应时域分析的精度、耗散性和稳定性;
步骤5:根据步骤3得到具有统一系数的参数约束关系,根据步骤4得到满足二阶精度、可控数值耗散和无条件稳定性的参数约束关系,将步骤3和步骤4得到的参数约束关系代入步骤2建立的多参数优化时间步进方程,将所述多参数优化时间步进方程结合Newton迭代技术,代入步骤1构建的包含非线性几何项、非线性阻尼项及其耦合项的结构动力学控制方程,实现高效、精确、且无条件稳定地计算包含非线性几何项、非线性阻尼项及其耦合项的结构动力学系统瞬态响应。
3.如权利要求2所述的非线性结构动力学系统瞬态响应无条件稳定时间积分方法,其特征在于:步骤2实现方法为,
基于差分格式,建立状态变量,包括位移、速度和加速度,在一个是时间单元[t,t+Δt]内随时间t的演化规律;将一个时间单元[t,t+Δt]划分成三个分步,即[t,t+c1Δt]、[t+c1Δt,t+c2Δt]以及[t+c2Δt,t+c3Δt];第一个分步[t,t+c1Δt]的步进公式为
其中,xt+c1Δt、和代表t+c1Δt时刻的位移、速度和加速度;xt和代表t时刻的位移和速度;Δt代表时间步长大小;c1为第一个分步引入的算法参数;第二个分步[t+c1Δt,t+c2Δt]的步进公式为
其中,xt+c3Δt、和代表t+c3Δt时刻的位移、速度和加速度;c3、β和η为第三个分步引入的算法参数;时间单元终点时刻t+Δt信息通过对三个已知时刻信息t+c1Δt、t+c2Δt以及t+c3Δt进行加权得到,对应的步进方程为
其中,b1和b2为待定加权参数。
5.如权利要求4所述的非线性结构动力学系统瞬态响应无条件稳定时间积分方法,其特征在于:步骤4实现方法为,
步骤4.1:针对包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学系统,推导便于数值性能调控的一个时间单元内的Butcher表;首先,将公式(2)至(5)等价地写为基于一阶方程的形式,如下所示
步骤4.2:在步骤4.1的基础上,建立的基于一阶方程的多参数优化时间步进方程推导其在一个时间单元内的传递方程,如下所示
zt+Δt=A(τ)zt,τ=λΔt (18)
其中,λ代表特征根;式中,传递因子A的表达式为
步骤4.3:在步骤4.1的基础上,基于局部截断误差确定多参数优化时间步进方程的精度阶次,局部截断误差的定义为
σ=A(τ)-Aexact(τ),Aexact(τ)=exp(τ) (20)
根据二阶精度,要求
A(0)=A(1)(0)=A(2)(0)=1 (21)
根据式(21)求解得到如下满足二阶精度要求的时间步进方程的参数关系,即
步骤4.4:在步骤4.1的基础上,以数值耗散可控为目标建立高频段传递因子可调的参数约束关系,如下所示,
对于非耗散格式ρ∞=1,根据式(23)建立如下参数关系,
对于耗散格式0≤ρ∞<1,根据式(23)建立如下参数关系
步骤4.5:在步骤4.1的基础上,以无条件稳定性为目标建立满足BN稳定性条件的参数约束关系;根据BN稳定性定义,建立如下不等式约束条件,即
S1=b1(2c1-b1)≥0 (26)
S2=b2(2c2α-b2)≥0 (27)
S3=(1-b1-b2)[2c3η-(1-b1-b2)]≥0 (28)
S5=S1S3-(1-b1-b2)2[c3(1-β-η)-b1]2≥0 (29)
S6=S2S3-(1-b1-b2)2(c3β-b2)2≥0 (30)
上述不等式约束条件,即公式(26)至公式(31),能够保证在任意时间步长下稳定地计算包含非线性几何项、非线性阻尼项以及它们耦合项的结构动力学系统;根据公式(26)至(31),可建立β和c3之间的关系;对于非耗散格式ρ∞=1,根据公式(26)至公式(31),建立如下参数关系,
对于耗散格式0≤ρ∞<1,根据公式(26)—公式(31),建立如下参数关系,
式中,c3=f(ρ∞)为关于高频段耗散量的已知函数。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202211126469.0A CN115453873B (zh) | 2022-09-16 | 非线性结构动力学系统瞬态响应无条件稳定时间积分方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202211126469.0A CN115453873B (zh) | 2022-09-16 | 非线性结构动力学系统瞬态响应无条件稳定时间积分方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN115453873A true CN115453873A (zh) | 2022-12-09 |
CN115453873B CN115453873B (zh) | 2024-07-09 |
Family
ID=
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116882195A (zh) * | 2023-07-26 | 2023-10-13 | 哈尔滨工业大学 | 一种求解二阶非线性动力学问题的单步显式逐步积分方法、应用及系统 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106844914A (zh) * | 2017-01-09 | 2017-06-13 | 西北工业大学 | 一种空天飞行器机翼振动响应的快速仿真方法 |
CN112147896A (zh) * | 2020-09-28 | 2020-12-29 | 中国科学院数学与系统科学研究院 | 一种非标准型mimo离散非线性系统的自适应控制方法及系统 |
CN112528411A (zh) * | 2020-12-10 | 2021-03-19 | 中国运载火箭技术研究院 | 一种基于模态减缩的几何非线性结构噪声振动响应计算方法 |
US20210207546A1 (en) * | 2018-09-28 | 2021-07-08 | Southeast University | Nonlinear disturbance rejection control apparatus and method for electronic throttle control systems |
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106844914A (zh) * | 2017-01-09 | 2017-06-13 | 西北工业大学 | 一种空天飞行器机翼振动响应的快速仿真方法 |
US20210207546A1 (en) * | 2018-09-28 | 2021-07-08 | Southeast University | Nonlinear disturbance rejection control apparatus and method for electronic throttle control systems |
CN112147896A (zh) * | 2020-09-28 | 2020-12-29 | 中国科学院数学与系统科学研究院 | 一种非标准型mimo离散非线性系统的自适应控制方法及系统 |
CN112528411A (zh) * | 2020-12-10 | 2021-03-19 | 中国运载火箭技术研究院 | 一种基于模态减缩的几何非线性结构噪声振动响应计算方法 |
Non-Patent Citations (1)
Title |
---|
张小青;王晓力;刘韧;: "微气体螺旋槽推力轴承转子系统非线性动力学分析", 北京理工大学学报, no. 11, 15 November 2012 (2012-11-15) * |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116882195A (zh) * | 2023-07-26 | 2023-10-13 | 哈尔滨工业大学 | 一种求解二阶非线性动力学问题的单步显式逐步积分方法、应用及系统 |
CN116882195B (zh) * | 2023-07-26 | 2024-03-22 | 哈尔滨工业大学 | 一种求解二阶非线性动力学问题的显式方法及系统 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Do et al. | Material optimization of tri-directional functionally graded plates by using deep neural network and isogeometric multimesh design approach | |
Lund | Finite element based design sensitivity analysis and optimization | |
CN105975645B (zh) | 一种基于多步的含激波区域飞行器流场快速计算方法 | |
Pen | A high-resolution adaptive moving mesh hydrodynamic algorithm | |
Peng et al. | A novel fast model predictive control with actuator saturation for large-scale structures | |
CN104866692A (zh) | 一种基于自适应代理模型的飞行器多目标优化方法 | |
CN111125963B (zh) | 基于拉格朗日积分点有限元的数值仿真系统及方法 | |
CN113821887A (zh) | 基于无网格efgm和plsm的各向异性结构热力耦合拓扑优化方法 | |
Gang et al. | Mesh deformation on 3D complex configurations using multistep radial basis functions interpolation | |
CN112001004B (zh) | 一种解析中高频振动结构能量密度场的nurbs等几何分析方法 | |
CN114564787A (zh) | 用于目标相关翼型设计的贝叶斯优化方法、装置及存储介质 | |
Park et al. | Multidisciplinary design optimization of a structurally nonlinear aircraft wing via parametric modeling | |
CN115453873A (zh) | 非线性结构动力学系统瞬态响应无条件稳定时间积分方法 | |
CN115453873B (zh) | 非线性结构动力学系统瞬态响应无条件稳定时间积分方法 | |
CN117390778A (zh) | 一种考虑薄板应变强化效应的塑性安定上限载荷计算方法 | |
Ma et al. | An explicit-implicit mixed staggered asynchronous step integration algorithm in structural dynamics | |
Ding et al. | An improved explicit-implicit precise integration method for nonlinear dynamic analysis of structures | |
CN110348158B (zh) | 一种基于分区异步长解的地震波动分析方法 | |
Roy et al. | Advanced heavy water reactor control with the aid of adaptive second-order sliding mode controller | |
CN110460060B (zh) | 一种电力系统中连续潮流计算方法 | |
Shern et al. | The Effects of Weightage Values with Two Objective Functions in iPSO for Electro-Hydraulic Actuator System | |
Huang et al. | Reanalysis-based fast solution algorithm for flexible multi-body system dynamic analysis with floating frame of reference formulation | |
CN113868853A (zh) | 一种梯度增强变保真代理模型建模方法 | |
Liu et al. | A new approach to stabilization of uncertain nonlinear systems | |
CN105893646A (zh) | 一种基于增量Kriging的序列优化实验设计方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant |