CN110008635A - 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法 - Google Patents

利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法 Download PDF

Info

Publication number
CN110008635A
CN110008635A CN201910320149.0A CN201910320149A CN110008635A CN 110008635 A CN110008635 A CN 110008635A CN 201910320149 A CN201910320149 A CN 201910320149A CN 110008635 A CN110008635 A CN 110008635A
Authority
CN
China
Prior art keywords
point
formula
equation
newmark
elasto
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
Application number
CN201910320149.0A
Other languages
English (en)
Other versions
CN110008635B (zh
Inventor
张瑞杰
贾波
孙佳杰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shaanxi Xinxi Commercial Engineering Technology Co Ltd
Original Assignee
Shaanxi Xinxi Commercial Engineering Technology Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Shaanxi Xinxi Commercial Engineering Technology Co Ltd filed Critical Shaanxi Xinxi Commercial Engineering Technology Co Ltd
Priority to CN201910320149.0A priority Critical patent/CN110008635B/zh
Publication of CN110008635A publication Critical patent/CN110008635A/zh
Application granted granted Critical
Publication of CN110008635B publication Critical patent/CN110008635B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • 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)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明公开了一种利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，包括以下步骤：S1：以有限元原则构建弹塑性结构二阶地震方程；S2：引入Newmark法对加速度的基本假定对二阶地震方程进行降阶，转变为一阶状态方程；S3：依据精细积分法的基本原理求解指数矩阵T；S4：采用级数解对非齐次项进行求解；S5：对弹塑性结构的切线刚度变化采用全量迭代算法或增量迭代算法进行积分计算。本发明的方法结合了Newmark法和精细积分法的优势，避免了传统精细积分法“增维降阶”所带来的系统矩阵扩容以及非齐次项解析解所需要的系统矩阵求逆。

Description

利用Newmark精细积分结合法对弹塑性结构地震反应分析的 方法
技术领域
本发明涉及土木工程领域的工程结构的弹塑性地震时程反应分析领域，具体涉及一种利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法。
背景技术
获得工程结构在地震中的反应是工程抗震的重要基础工作。对于弹塑性结构地震反应分析，建立修正刚度的二阶动力微分方程，并采用直接积分法求解是目前最常用的方法。现阶段已发展的直接积分法有中心差分法、Newmark-β法、wilson-θ法、Houbolt法等。上述传统的隐式或显式算法仅具有二阶或三阶精度。钟万勰教授于1994年提出了一种求解结构动力响应的精细积分法。该方法通过引入哈密顿体系对偶变量来实现方程降阶，并且基于时间步内荷载线性化的假定给出了状态方程非齐次项的解析解。相对于传统直接积分法，精细积分法具有精度高和无条件稳定的优点，而且该方法是显式递推方法，积分步长是灵活可变的，结果精度与积分步长无关。传统精细积分方法有两个缺点：其一，在降阶的同时，系统矩阵H的维数增加至刚度矩阵维数的两倍，存储规模为四倍。而且精细积分法的指数矩阵一般为满阵，丧失了刚度矩阵的稀疏性特征。精细积分法在获得高精度、无条件稳定等优势的同时，是以牺牲计算机的存储和运算规模为代价的。因而，传统的精细积分法不利于在大规模多自由度系统中应用。其二，基于时间步内荷载线性化假定的解析解需要对系统矩阵求逆，而对于弹塑性问题，由于状态矩阵是可变的，对状态矩阵求逆不仅工作量大，而且可能会碰到逆矩阵不存在的问题，这样会给求解带来困难。
为克服系统矩阵维数增加的缺点以及避免系统矩阵求逆，在吸收前人研究成果的基础上，本发明提出了Newmark法与精细积分法的结合方法，具体包括全量动力平衡方程和增量动力平衡方程的Newmark-精细积分结合法。其实质是利用Newmark法对加速度的基本假定，对动力方程进行降阶，避免系统矩阵维数增加。针对非齐次项的求解，为避免系统矩阵求逆，提出了非齐次项按线性变化和常量时的级数解，分别适用于全量方程和增量方程。
发明内容
本发明的目的在于提供一种利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，解决传统精细积分法系统矩阵维数增加和求逆的缺点问题。
为解决上述的技术问题，本发明采用以下技术方案：
一种利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，包括以下步骤：
S1：以有限元原则构建弹塑性结构二阶地震方程；
S2：引入Newmark法对加速度的基本假定对二阶地震方程进行降阶，转变为一阶状态方程；
S3：依据精细积分法的基本原理求解指数矩阵T；
S4：采用级数解对非齐次项进行求解；
S5：对弹塑性结构的切线刚度变化采用全量迭代算法或增量迭代算法进行积分计算。
更进一步地，所述S1步骤中以有限元原则构建弹塑性结构二阶地震方程，包括全量形式的地震运动方程和增量形式的地震运动方程：
全量形式地震运动方程：
增量形式地震运动方程：
更进一步地，所述S2步骤中引入Newmark法对加速度的基本假定对二阶地震方程进行降阶的具体步骤如下：
对于全量形式的地震运动方程：
第一步：Newmark法引入γ参数对t+Δt时刻的速度进行修正：
式中，为t+Δt时刻的速度、加速度，为待求的未知量；为t时刻的速度和加速度，为已知量；
第二步：根据(1)式，t+Δt时刻的加速度可以用速度表示为：
第三步：将(2)右端代入全量形式地震运动方程，消去二阶项：
第四步：理(2)式后，得到一阶状态方程：
方程(4)中：
对于增量形式的地震运动方程：
第一步：Newmark法中，[t,t+Δt]范围的速度增量为：
式中，为[t,t+Δt]的加速度增量、速度增量，为待求的未知量；为t时刻的加速度，为已知量。
第二步：利用速度增量式(7)，加速度增量可表达为如下形式：
第三步：将(8)式代入增量形式地震运动方程，消去二阶项：
第四步：(9)式可整理成状态方程形式：
式中：
更进一步地，所述S3步骤中依据精细积分法的基本原理求解指数矩阵T前得到状态方程通解，具体如下：
全量状态方程通解：
增量状态方程通解：
更进一步地，所述S3步骤中依据精细积分法的基本原理求解指数矩阵T的具体方法是：
全量状态方程和增量状态方程通解中都包含指数矩阵T＝eH·Δt，求解指数矩阵采用2N类算法，具体步骤如下：
Ta矩阵采用泰勒级数展开求解：
式中，l为截断阶数；
利用指数函数的加法原理，指数矩阵具有如下公式；
T＝I+Ta (17)；
Ta采用N次循环迭代计算，迭代公式为Ta＝2Ta+Ta×Ta；Ta的最终值代入式(15)，即求出指数矩阵T。
更进一步地，所述S4步骤中采用级数解对非齐次项进行求解的具体方法为：
若为全量形式时，假定非齐次项在[tk tk+1]按照线性变化，级数解公式为：
若为增量形式时，假定非齐次项在[tk tk+1]为常量，级数解公式为：
更进一步地，所述S5步骤中对弹塑性结构的切线刚度变化采用全量迭代算法或增量迭代算法进行计算的方法为：
若采用全量状态方程递推时，采用全量迭代算法，具体步骤如下：
(1)第i积分点的割线刚度和该积分点的结构反应ui、为已知，根据第i积分点的割线刚度，按照全量精细积分递推公式递推第i+1积分点的结构反应ui+1、
(2)根据弹塑性模型滞回规律确定第i+1积分点的割线刚度矩阵Kg(i+1)；
(3)按照下式计算第i+1积分点的不平衡力：
(4)应提前确定不平衡力的容差，当不平衡力R取绝对值后的最大元素超过设定容差时，进行迭代循环；
(5)根据第(2)步骤确定的第i+1积分点的割线刚度Kg(i+1)和第i积分点的结构反应ui、按照全量精细积分递推公式递推第i+1积分点的结构反应ui+1、
(6)根据第(5)步骤计算结果，修正第i+1积分点的割线刚度Kg(i+1)；
(7)执行步骤(3)～(6)，直至不平衡力小于设定容差；
若采用增量状态方程递推时，采用增量迭代算法，具体步骤如下：
(1)第i-1积分点～第i积分点的结构增量Δui、以及第i积分点的结构反应ui、已经在前步计算获得；根据第i积分点的切线刚度Kt，按照增量精细积分递推公式递推第i积分点～第i+1积分点的结构增量反应Δui+1、并计算第i+1积分点的结构反应ui+1、
(2)根据弹塑性模型滞回规律确定第i+1积分点的割线刚度矩阵Kg(i+1)和切线刚度矩阵Kt(i+1)；
(3)按照公式(17)计算第i+1积分点的不平衡力；
(4)当不平衡力R取绝对值后的最大元素超过设定容差时，进行迭代循环。按下述公式对第i+1积分点的结构反应进行修正：
计算等效切线刚度：
计算位移修正量：
计算速度修正量：
计算加速度修正量：
修正第i+1积分点的位移：ui+1＝ui+1+Δu；
修正第i+1积分点的速度：
修正第i+1积分点的加速度：
(5)按修正后的第i+1积分点的结构反应，修正i+1积分点的割线刚度矩阵Kg(i+1)和切线刚度矩阵Kt(i+1)；
(6)执行步骤(3)～(5)，直至不平衡力小于设定容差。
更进一步地，所述S5步骤积分计算完成后，提取结构的位移反应、速度反应或加速度反应并绘制时程曲线。
与现有技术相比，本发明的有益效果是：
本发明的基本原理是利用Newmark法对加速度的基本假定，将计算时刻的未知加速度或加速度增量用已知量表达，代入全量形式或增量形式的地震运动方程，从而将二阶方程降阶为一阶，并转化成状态方程形式。一阶状态方程的求解，关键是指数矩阵和非齐次项的求解。指数矩阵的求解利用精细积分法理论2N类算法计算，非齐次项的求解采用级数解公式，分别导出了非齐次项在时间步内按线性变化和常量时的级数解公式。
与传统的精细积分法相比，首先，本方法的降阶方式没有增维，系统矩阵的维数和刚度矩阵维数相同，避免了运算和存储规模增加，其次，非齐次项的级数解避免了对系统矩阵求逆；与Newmark法相比，一阶状态方程的求解采用精细算法，保留了传统精细积分法计算精度高的优势。
附图说明
图1为本发明Newmark精细积分法步骤流程图。
图2为本发明的弹塑性结构地震响应分析Newmark精细积分结合法增量算法流程图。
图3为本发明的弹塑性结构地震响应分析Newmark精细积分结合法全量算法流程图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
1.按照有限元方法，组集结构的刚度矩阵、质量矩阵，形成阻尼矩阵，全量形式和增量形式的地震运动方程分别如下：
全量形式地震运动方程：
增量形式地震运动方程：
全量形式的地震运动方程采用结构的割线刚度建立，弹性阶段直接递推求解，进入塑性阶段后每步均应迭代求解。增量形式的地震运动方程采用结构的切线刚度建立，无论在弹性阶段还是塑性阶段，只要结构切线刚度不变均可直接递推求解，切线刚度改变时需迭代消除不平衡力。结构进入塑性阶段后，结构的割线刚度时变，采用切线刚度可减少计算量，但应注意切线刚度为0时不适用增量法。
2.根据Newmark法对加速度的基本假定，对上述二阶结构地震运动方程进行降阶，具体步骤如下：
(1)全量动力方程：
第一步：Newmark法引入γ参数对t+Δt时刻的速度进行修正：
式中，为t+Δt时刻的速度、加速度，为待求的未知量；为t时刻的速度和加速度，为已知量；
第二步：根据(1)式，t+Δt时刻的加速度可以用速度表示为：
第三步：将(2)右端代入全量形式地震运动方程，消去二阶项：
第四步：整理(2)式后，得到一阶状态方程：
方程(4)中：
(2)增量动力方程
第一步：Newmark法中，[t,t+Δt]范围的速度增量为：
式中，为[t,t+Δt]的加速度增量、速度增量，为待求的未知量；为t时刻的加速度，为已知量。
第二步：利用速度增量式(7)，加速度增量可表达为如下形式：
第三步：将(8)式代入增量形式地震运动方程，消去二阶项：
第四步：(9)式可整理成状态方程形式：
式中：
3.状态方程的通解
全量状态方程通解：
根据全量状态方程通解，可得到t+Δt时刻的位移，然后根据状态方程(4)可得到t+Δt的速度，再代入全量结构运动方程可求解t+Δt的加速度。
增量状态方程通解：
式中，Δu0为[t-Δt,t]的位移增量，为已知量。根据增量状态方程通解，可得到[t,t+Δt]的位移增量，然后根据状态方程(10)可得到[t,t+Δt]的速度增量，再代入增量动力平衡方程可求解[t,t+Δt]的加速度增量。经过累加可求得t+Δt时刻的加速度、速度、位移。
4.指数矩阵的求解
全量状态方程和增量状态方程通解中都包含指数矩阵T＝eH·Δt，求解指数矩阵采用2N类算法，具体步骤如下：
Ta矩阵采用泰勒级数展开求解：
式中，l为截断阶数，在求解中可取l＝4，N＝15即可得到非常精确的解答。求出Ta矩阵后，由于Ta矩阵中的元素非常小，若直接代入式(13)求矩阵T，则由于计算机的舍入误差会使得指数矩阵T的精度丧失殆尽。利用指数函数的加法原理，指数矩阵T的求解最终采用下式进行。
T＝I+Ta (17)；
Ta采用N次循环迭代计算，迭代公式为Ta＝2Ta+Ta×Ta。Ta的最终值代入式(13)，就可精确的求出指数矩阵T。
5.非齐次项的求解。
全量状态方程和增量状态方程通解中都包含非齐次项的求解。非齐次项的求解采用级数解。对于全量状态方程，可假定非齐次项在[tk tk+1]按照线性变化，级数解公式如下：
对于增量状态方程，可假定非齐次项在[tk tk+1]为常量，级数解公式如下：
6.界点以及时变刚度的迭代计算。
对于弹性结构和弹塑性结构的弹性阶段，由于割线刚度、切线刚度相等且不变，指数矩阵T不变；结构进入塑性状态后，若采用折线型滞回模型，塑性阶段未改变时，切线刚度不变，割线刚度时变。弹塑性结构地震响应分析可采用全量法或增量法或其混合算法进行，弹塑性状态转换后，若时间步内切线刚度或割线刚度产生变化，则需迭代计算，消除不平衡力。
当前阶段采用全量状态方程递推时，全量迭代算法步骤：
(1)第i积分点的割线刚度和该积分点的结构反应ui、为已知，根据第i积分点的割线刚度，按照全量精细积分递推公式递推第i+1积分点的结构反应ui+1、
(2)根据弹塑性模型滞回规律确定第i+1积分点的割线刚度矩阵Kg(i+1)；
(3)按照下式计算第i+1积分点的不平衡力：
(4)应提前确定不平衡力的容差，当不平衡力R取绝对值后的最大元素超过设定容差时，进行迭代循环；
(5)根据第(2)步骤确定的第i+1积分点的割线刚度Kg(i+1)和第i积分点的结构反应ui、按照全量精细积分递推公式递推第i+1积分点的结构反应ui+1、
(6)根据第(5)步骤计算结果，修正第i+1积分点的割线刚度Kg(i+1)；
(7)执行步骤(3)～(6)，直至不平衡力小于设定容差；
当前阶段采用增量状态方程递推时，增量迭代算法步骤：
(1)第i-1积分点～第i积分点的结构增量Δui、以及第i积分点的结构反应ui、已经在前步计算获得；根据第i积分点的切线刚度Kt，按照增量精细积分递推公式递推第i积分点～第i+1积分点的结构增量反应Δui+1、并计算第i+1积分点的结构反应ui+1、
(2)根据弹塑性模型滞回规律确定第i+1积分点的割线刚度矩阵Kg(i+1)和切线刚度矩阵Kt(i+1)；
(3)按照公式(18)计算第i+1积分点的不平衡力；
(4)当不平衡力R取绝对值后的最大元素超过设定容差时，进行迭代循环。按下述公式对第i+1积分点的结构反应进行修正：
计算等效切线刚度：
计算位移修正量：
计算速度修正量：
计算加速度修正量：
修正第i+1积分点的位移：ui+1＝ui+1+Δu；
修正第i+1积分点的速度：
修正第i+1积分点的加速度：
(5)按修正后的第i+1积分点的结构反应，修正i+1积分点的割线刚度矩阵Kg(i+1)和切线刚度矩阵Kt(i+1)；
(6)执行步骤(3)～(5)，直至不平衡力小于设定容差。
本发明的具体实施步骤概述如下：
1.有限元建模
按照有限元原则建立工程结构的有限元模型，具体包括结点信息、单元信息、单元截面属性、单元材料属性、结构边界条件、滞回曲线模型信息；组集工程结构弹性阶段刚度矩阵、质量矩阵，根据预定义的阻尼比等参数信息组集阻尼矩阵。
2.控制参数输入及初始化
输入地震波文件，定义地震波的峰值、地震波步长；积分步长、不平衡力容差等控制信息；输入Newmark法的控制参数γ，一般取0.5；初始化积分步计数变量s，时长累加变量t，弹塑性阶段标识数组等。
3.时程积分
时程积分的方法有全量法和增量法两种方法，切线刚度为0等特殊情况也可以混合使用两种方法。对于折线型模型，采用增量法可减少迭代工作量，且能避免割线刚度为0以及无穷大等特殊情况，因而增量法具有更广泛的用途。
全量法的具体实施方式如下：
1)结构的初始状态为静止，位移、速度、加速度向量均为0，结构处于弹性状态；
2)根据弹性阶段的割线刚度矩阵、质量矩阵、阻尼矩阵等，按照(5)式形成系统矩阵H，并进行指数矩阵T的求解，结构的弹性阶段指数矩阵T不变；
3)以积分时长t是否超出地震波时长为判断条件，积分步长设定为Δt，对后续积分步递推求解，具体实施步骤4)；
4)对当前积分步进行计算，首先形成当前时刻的荷载P0和积分部内的荷载变化率P1，按照(18)式计算非齐次项的级数解；然后，根据全量方程通解递推公式，从前一时刻状态递推下一时刻的状态；
5)每一步计算完成后对结构的弹塑性状态进行判断。若仍为弹性阶段返回步骤4)实施下一步递推；若进入塑性阶段，此后割线刚度时变，则实施步骤6)、7)；
6)以全量迭代算法迭代求解，消除不平衡力。
7)修正计算时刻的割线刚度矩阵，重新形成系统矩阵H，并进行指数矩阵T的求解，下一积分步使用；
8)返回步骤3)实施下一步递推。
增量法的具体实施方式如下：
1)结构的初始状态为静止，位移、速度、加速度向量均为0，结构处于弹性状态；
2)根据弹性阶段的割线或切线刚度矩阵、质量矩阵、阻尼矩阵等，按照(5)式形成系统矩阵H，并进行指数矩阵T的求解，结构的弹性阶段指数矩阵T不变；
3)第一步采用全量法求解，求解出的Δt时刻结构状态与初始结构状态的差值，即为初始增量状态；
4)以积分时长t是否超出地震波时长为判断条件，积分步长设定为Δt，对后续积分步递推求解，具体实施步骤5)；
5)对当前积分步进行计算，首先形成当前时刻的P0，按照(19)式计算非齐次项的级数解；然后，根据增量方程通解递推公式，从前一步增量状态递推下一步的增量状态，并通过累加获得计算时刻的结构状态；
6)每一步计算完成后对结构的弹塑性状态进行判断，并提取计算时刻的切线刚度。若切线刚度阵未变化，返回步骤4)实施下一步递推；若切线刚度产生改变，意味着弹塑性阶段产生改变，则实施步骤7)、8)；
7)以增量迭代算法迭代求解，消除不平衡力。
8)修正积分步的切线刚度矩阵，重新形成系统矩阵H，并进行指数矩阵T的求解，下一积分步使用；
9)返回步骤4)实施下一步递推。
4.结果后处理
积分计算完成后，可提取结构的位移反应、速度反应或加速度反应并绘制时程曲线；根据结构位移，可计算任意截面的内力、应力等信息；积分计算过程中，对结构的恢复力进行存储，可绘制弹塑性结构的滞回曲线。
尽管这里参照本发明的多个解释性实施例对本发明进行了描述，但是，应该理解，本领域技术人员可以设计出很多其他的修改和实施方式，这些修改和实施方式将落在本申请公开的原则范围和精神之内。更具体地说，在本申请公开、附图和权利要求的范围内，可以对主题组合布局的组成部件和/或布局进行多种变型和改进。除了对组成部件和/或布局进行的变形和改进外，对于本领域技术人员来说，其他的用途也将是明显的。

Claims (8)

1.一种利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：包括以下步骤：
S1：以有限元原则构建弹塑性结构二阶地震方程；
S2：引入Newmark法对加速度的基本假定对二阶地震方程进行降阶，转变为一阶状态方程；
S3：依据精细积分法的基本原理求解指数矩阵T；
S4：采用级数解对非齐次项进行求解；
S5：对弹塑性结构的切线刚度变化采用全量迭代算法或增量迭代算法进行积分计算。
2.根据权利要求1所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S1步骤中以有限元原则构建弹塑性结构二阶地震方程，包括全量形式的地震运动方程和增量形式的地震运动方程。
3.根据权利要求2所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S2步骤中引入Newmark法对加速度的基本假定对二阶地震方程进行降阶的具体步骤如下：
对于全量形式的地震运动方程：
第一步：Newmark法引入γ参数对t+Δt时刻的速度进行修正：
式中，为t+Δt时刻的速度、加速度，为待求的未知量；为t时刻的速度和加速度，为已知量；
第二步：根据(1)式，t+Δt时刻的加速度可以用速度表示为：
第三步：将(2)右端代入全量形式地震运动方程，消去二阶项：
第四步：理(2)式后，得到一阶状态方程：
方程(4)中：
对于增量形式的地震运动方程：
第一步：Newmark法中，[t,t+Δt]范围的速度增量为：
式中，为[t,t+Δt]的加速度增量、速度增量，为待求的未知量；为t时刻的加速度，为已知量；
第二步：利用速度增量式(7)，加速度增量可表达为如下形式：
第三步：将(8)式代入增量形式地震运动方程，消去二阶项：
第四步：(9)式可整理成状态方程形式：
式中：
4.根据权利要求1所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S3步骤中依据精细积分法的基本原理求解指数矩阵T前得到状态方程通解，具体如下：
全量状态方程通解：
增量状态方程通解：
5.根据权利要求1所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S3步骤中依据精细积分法的基本原理求解指数矩阵T的具体方法是：
全量状态方程和增量状态方程通解中都包含指数矩阵T＝eH·Δt，求解指数矩阵采用2N类算法，具体步骤如下：
Ta矩阵采用泰勒级数展开求解：
式中，l为截断阶数；
利用指数函数的加法原理，指数矩阵具有如下公式：
T＝I+Ta (17)；
Ta采用N次循环迭代计算，迭代公式为Ta＝2Ta+Ta×Ta；Ta的最终值代入式(15)，即求出指数矩阵T。
6.根据权利要求1所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S4步骤中采用级数解对非齐次项进行求解的具体方法为：
若为全量形式时，假定非齐次项在[tk tk+1]按照线性变化，级数解公式为：
若为增量形式时，假定非齐次项在[tk tk+1]为常量，级数解公式为：
7.根据权利要求1所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S5步骤中对弹塑性结构的切线刚度变化采用全量迭代算法或增量迭代算法进行计算的方法为：
若采用全量状态方程递推时，采用全量迭代算法，具体步骤如下：
(1)第i积分点的割线刚度和该积分点的结构反应ui、为已知，根据第i积分点的割线刚度，按照全量精细积分递推公式递推第i+1积分点的结构反应ui+1、
(2)根据弹塑性模型滞回规律确定第i+1积分点的割线刚度矩阵Kg(i+1)；
(3)按照下式计算第i+1积分点的不平衡力：
(4)应提前确定不平衡力的容差，当不平衡力R取绝对值后的最大元素超过设定容差时，进行迭代循环；
(5)根据第(2)步骤确定的第i+1积分点的割线刚度Kg(i+1)和第i积分点的结构反应ui、按照全量精细积分递推公式递推第i+1积分点的结构反应ui+1、
(6)根据第(5)步骤计算结果，修正第i+1积分点的割线刚度Kg(i+1)；
(7)执行步骤(3)～(6)，直至不平衡力小于设定容差；
若采用增量状态方程递推时，采用增量迭代算法，具体步骤如下：
(1)第i-1积分点～第i积分点的结构增量Δui、以及第i积分点的结构反应ui、已经在前步计算获得；根据第i积分点的切线刚度Kt，按照增量精细积分递推公式递推第i积分点～第i+1积分点的结构增量反应Δui+1、并计算第i+1积分点的结构反应ui+1、
(2)根据弹塑性模型滞回规律确定第i+1积分点的割线刚度矩阵Kg(i+1)和切线刚度矩阵Kt(i+1)；
(3)按照公式(17)计算第i+1积分点的不平衡力；
(4)当不平衡力R取绝对值后的最大元素超过设定容差时，进行迭代循环。按下述公式对第i+1积分点的结构反应进行修正：
计算等效切线刚度：
计算位移修正量：
计算速度修正量：
计算加速度修正量：
修正第i+1积分点的位移：ui+1＝ui+1+Δu；
修正第i+1积分点的速度：
修正第i+1积分点的加速度：
(5)按修正后的第i+1积分点的结构反应，修正i+1积分点的割线刚度矩阵Kg(i+1)和切线刚度矩阵Kt(i+1)；
(6)执行步骤(3)～(5)，直至不平衡力小于设定容差。
8.根据权利要求1所述的利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法，其特征在于：所述S5步骤积分计算完成后，提取结构的位移反应、速度反应或加速度反应并绘制时程曲线。
CN201910320149.0A 2019-04-19 2019-04-19 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法 Active CN110008635B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910320149.0A CN110008635B (zh) 2019-04-19 2019-04-19 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910320149.0A CN110008635B (zh) 2019-04-19 2019-04-19 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法

Publications (2)

Publication Number Publication Date
CN110008635A true CN110008635A (zh) 2019-07-12
CN110008635B CN110008635B (zh) 2023-06-02

Family

ID=67173257

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910320149.0A Active CN110008635B (zh) 2019-04-19 2019-04-19 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法

Country Status (1)

Country Link
CN (1) CN110008635B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110378060A (zh) * 2019-07-26 2019-10-25 中国海洋大学 一种顶张式立管随机耦合振动的计算方法
CN111404825A (zh) * 2020-03-13 2020-07-10 重庆第二师范学院 数据传输方法、系统、计算机设备及存储介质
CN115438513A (zh) * 2022-11-07 2022-12-06 人工智能与数字经济广东省实验室(广州) 分数阶阻尼减震结构抗震设计的分析方法、系统、设备和介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109145386A (zh) * 2018-07-24 2019-01-04 张瑞杰 一种用于工程结构的地震碰撞反应分析方法
CN109446715A (zh) * 2018-11-15 2019-03-08 李鲁 一种大跨度桥梁地震响应时程分析方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109145386A (zh) * 2018-07-24 2019-01-04 张瑞杰 一种用于工程结构的地震碰撞反应分析方法
CN109446715A (zh) * 2018-11-15 2019-03-08 李鲁 一种大跨度桥梁地震响应时程分析方法

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110378060A (zh) * 2019-07-26 2019-10-25 中国海洋大学 一种顶张式立管随机耦合振动的计算方法
CN110378060B (zh) * 2019-07-26 2021-02-09 中国海洋大学 一种顶张式立管随机耦合振动的计算方法
CN111404825A (zh) * 2020-03-13 2020-07-10 重庆第二师范学院 数据传输方法、系统、计算机设备及存储介质
CN115438513A (zh) * 2022-11-07 2022-12-06 人工智能与数字经济广东省实验室(广州) 分数阶阻尼减震结构抗震设计的分析方法、系统、设备和介质

Also Published As

Publication number Publication date
CN110008635B (zh) 2023-06-02

Similar Documents

Publication Publication Date Title
CN110008635A (zh) 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法
Rein et al. WHFAST: a fast and unbiased implementation of a symplectic Wisdom–Holman integrator for long-term gravitational simulations
Dokainish A new approach for plate vibrations: combination of transfer matrix and finite-element technique
Monterrubio et al. Proof of convergence for a set of admissible functions for the Rayleigh–Ritz analysis of beams and plates and shells of rectangular planform
CN109255188B (zh) 基于有限元的六轴工业机器人动态性能优化方法
CN106778012A (zh) 一种小天体附着探测下降轨迹优化方法
CN103886160A (zh) 一种基于基础激励响应数据的考虑阻尼的模型修正方法
CN109145386B (zh) 一种用于工程结构的地震碰撞反应分析方法
CN109902404A (zh) 不同阻尼形式的结构时程响应积分的统一递推计算方法
WO2021082959A1 (zh) 一种复杂建筑结构非线性动力时程分析方法
Li et al. Woodbury approximation method for structural nonlinear analysis
CN108108559A (zh) 一种基于子结构的结构响应获取方法及灵敏度获取方法
Chen et al. A universal method for structural static reanalysis of topological modifications
CN106354954A (zh) 一种基于叠层基函数的三维力学模态仿真模拟方法
CN109726506B (zh) 一种基于等效静态载荷法的汽车保险杠面罩尺寸优化方法
CN105260499A (zh) 一种三维柔性梁系统力学的运动仿真方法及系统
CN108376192B (zh) 一种确定模态叠加法计算加速度反应所需振型数目的方法
Haug An index 0 Differential-Algebraic equation formulation for multibody dynamics: Holonomic constraints
CN111554355B (zh) 一种基于非冯诺依曼架构的分子动力学计算方法
CN108427833A (zh) 基于三维弹塑性地震反应谱的性能化抗震设计方法
CN104778309A (zh) 飞机结构强度校核方法和装置
CN112084592A (zh) 折叠式桁架动力学分析系统、方法、装置和存储介质
CN110909500A (zh) 一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法
Kujawski et al. Generalized time finite element algorithm for non‐linear dynamic problems
CN112100863B (zh) 用于桩基工程的桩基负摩擦力解析解分析方法

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
GR01 Patent grant