CN110909500A

CN110909500A - 一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法

Info

Publication number: CN110909500A
Application number: CN201911139241.3A
Authority: CN
Inventors: 冯德成; 程浩然; 王偲; 陈科玮; 泰子煜
Original assignee: Southeast University
Current assignee: Southeast University
Priority date: 2019-11-19
Filing date: 2019-11-19
Publication date: 2020-03-24

Abstract

本发明公开了一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法。建立结构的有限元模型并定义纤维截面，然后将结构做均布荷载转化成集中质量处理，运用KR算法得到地震响应图像。本发明中所运用的无条件稳定显式积分算法能够在很大程度上弥补上述如隐式积分法由于迭代导致不收敛、加载路径或显示积分算法中步长过小、效率精度不高等问题；本发明的分析方法无时间步长要求，即可得到稳定解，在进行计算大跨结构和复杂结构时，效率较高。

Description

一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法

技术领域

本发明涉及一种大跨度桥梁的地震响应分析方法，特别是一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法。

背景技术

通过研究桥梁在不同等级地震下的地震响应分析，可以直观地通过数据了解到桥梁在三维立体空间不同方向的位移、加速度以及形变程度等，从而可以针对性提出相应解决办法，以便加强各种级别的桥梁在不同程度地震下的抗震能力。

目前仿真模拟采用较多的是New-markβ算法和中心差分法，但这两种算法需要对时间步长的大小Δt进行限制。New-markβ中，只有加速度影响速度权重β≥0.5和加速度影响位移权重γ≥0.25(0.5+β)²时，时间步长才不受限制。因此这两种算法在解决复杂结构时，会出现解的不稳定性。且计算效率不高。

发明内容

发明目的：针对现有技术中的不足，本发明的目的是提供一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法，解决了现有技术中地震响应分析方法中的计算收敛性问题和计算稳定性问题。

技术方案：本发明的无条件稳定的超大跨桥梁多点激励倒塌显式分析方法，包括如下步骤：

步骤1、选取地震波；结合实际情况与所选取的模型，选取在各种情境下所可能遇到的地震波等类似扰动，通过其地震响应分析来确定本发明的分析方法的高效性。

步骤2、根据目标结构的结构类型和结构参数，进行有限元建模。

定义节点，每个节点有自己的质量、刚度、阻尼，再将每个节点拼装起来得到结构的质量矩阵M，刚度矩阵K和阻尼矩阵C，其表达式如下所示：

上式中，i为定义节点个数，m_i为第i号节点的质量；K₁₁表示在1处产生单位位移在1处产生的力，K₁₂表示在2处产生单位位移在1处产生的力，以此类推，K_ij表示在j处产生单位位移在i处产生的力；其中i、j均为正整数，j≤i。

其中C＝α[M]+β[K]，α、β由任意二阶振型自拟比，由给定的任意二阶振型的阻尼比ξ_i、ξ_j反算求得。

步骤3、定义纤维截面：将混凝土划分成n个区域；将混凝土划分成钢纤维和混凝土纤维，每个纤维赋予响应的本构关系。

混凝土纤维可采用Mander单轴混凝土本构模型，钢材纤维和钢筋纤维采用基于Giuffr-Menegotto-Pinto的钢材本构关系。

构件截面纤维网格的划分方式、疏密程度和高斯点积分的个数均会影响程序求解的精度和效率。纤维点个数越多，求解越精确，但效率越慢。

并定义材料的参数、定义梁柱；

步骤4、均布荷载转化成集中质量模型；将杆件均布荷载等效为杆件两端节点质量，每个转化的节点质量为杆件均布荷载的1/2。

利用节点与杆件的连接情况，将荷载进行处理，使其能够适应于KR算法。

例如：第j号节点与k个节点相连，j号节点与这k个节点分别构成k个单元，第i个单元的均布荷载为q_i，每个单元的长度为l_i，则每个节点承担的质量

为：

步骤5、导入地震波，利用显式方法求解如下结构动力方程：

其中，M为质量矩阵，C为阻尼矩阵，

为t时刻的加速度，

为t时刻的速度，

为t时刻的结构抗力，F(t)为t时刻外力。

在求解过程中，构建速度方程，位移更新方程以及动力方程离散形式，分别如下式(2)、(3)、(4)所示：

式(2)、式(3)中，X_i、

以及

分别表示第i时间步的位移、速度以及加速度向量；α₁和α₂为积分参数矩阵；Δt为积分时间步长；同样地，X_i+1、

以及

分别表示第i+1时间步的位移、速度以及加速度向量。

式(4)中的

表示等效加速度向量，

表示等效速度向量，

表示等效抗力向量，

表示等效外力向量；其具体形式分别如下：

上式中，I为单位矩阵；α₃为第三个积分参数矩阵；α_f为积分参数标量；R_i、F_i分别为第i步的抗力和外力向量；R_i+1、F_i+1分别为第i+1步的抗力和外力向量。

步骤6、求解α₁、α₂和α₃；通过设定上述KR方法的特征值及放大矩阵与广义-α相同，可得前述积分参数矩阵，如下式所示：

α₁＝[M+γΔtC+β(Δt)²K]^-1M (9)

α₃＝[M+γΔtC+β(Δt)²K]^-1[α_mM+α_fγΔtC+α_fβ(Δt)²K] (11)

上式中，K为初始刚度矩阵；参数γ，β，α_m，α_f为高频谱半径ρ_m的函数，其中ρ_m为人为定义的参数，其计算公式如下：

由上式(9)～式(14)可见，KR算法仅有一个用户定义参数，即高频谱半径ρ_m。该参数的取值范围为0～1。ρ_m＝1意味着不引入数值阻尼(即无数值能量耗散)，随着ρ_m的减小，该方法的数值阻尼增大，直至ρ_m＝0时，数值阻尼达到最大。ρ_m的取值与实际问题以及时间步长相关，较大的ρ_m可能引起结构内力的高频振荡，这种振荡可以通过减小时间步长或者减小ρ_m来消除。但减小时间步长又会丧失该显式方法的计算效率。通常，ρ_m值可以从1开始减小直至得到稳定的计算结果。需要指出的是，一旦得到稳定的结果后，ρ_m便不可减小，因为太多的数值阻尼会影响低阶模态的结果。

步骤7，高频谱半径ρ_m从1开始往0衰减，得到稳定解。从而可以得到位移-时间响应结果，力-位移响应结果。

本发明采用动力积分算法，该算法是无条件稳定的，避免了收敛问题，且可以采用较大的步长进行计算，大大提高了计算的效率。因此是一种无条件稳定的超大跨桥梁多点激励倒塌显式分析方法。

本发明的方法的主要步骤为：首先利用OpenSees软件建立大跨度桥梁的结构模型，并完成材料和截面的定义，并转换成集中质量模型，输入相应的地震波，利用KR算法，得到相应的地震响应结果。

本发明采用无条件稳定显式积分算法，首先针对大跨桥梁在地震下的响应分析用到的方法是直接积分法，将空间变量进行有限元离散，而对时间变量有限差分离散，则得到“直接积分法”，又分为隐式动力分析法与显式动力分析法，这是目前广泛采用分析结构在荷载作用下力学表现的方法之一。在这些方法里，隐式积分算法一般情况下无条件稳定，比如Newmarkβ法、Wilson法，但是其计算过程需要迭代，对于强非线性的试验体牵涉到加载路径的问题，与迭代过程对应的加-卸载是不允许的。而显式动力分析算法则是利用已知信息求解位移，例如本研究项目中将会对比的中心差分法，虽然没有隐式算法迭代问题，但大部分算法条件稳定，受制于积分步长的限制，步长又和模型单元的尺寸相关，在应用于高阶振型显著的结构体系时，积分步长很小，模拟效率不高，所以在大跨桥梁中一般不会使用。

因此，本发明中所运用的无条件稳定显式积分算法能够在很大程度上弥补上述如隐式积分法由于迭代导致不收敛、加载路径或显示积分算法中步长过小、效率精度不高等问题。无时间步长要求即可得到稳定解，在进行计算大跨结构和复杂结构时，效率较高。

有益效果：与现有技术相比，

(1)本发明采用仿真技术，可以避免物理试验的过程，可以在非常短的时间内完成分析和预测过程，大幅度提高了计算的效率和精度；

(2)本发明的无条件稳定积分建立在算法的选择上，可以很好地避免不收敛问题；

(3)本发明中的时间步长Δt不影响解的稳定性，仅只有高频谱半径ρ_m影响计算结果的稳定性；对大型复杂结构，计算效率更高，计算稳定性更高。

(4)本发明考虑到大跨度桥梁偏于柔性，且尺度大的特点，避免了常见积分算法的计算收敛性、稳定性问题，便于大跨桥梁地震响应分析；结合相应的工程结构分析软件并运用此方法能够进行有效的各中小跨桥梁乃至大跨桥梁的高效率分析。

附图说明

图1是本发明的流程示意图；

图2是Mander单轴混凝土本构模型示意图和Giuffr-Menegotto-Pinto的钢材本构关系示意图；其中图(a)为Mander单轴混凝土本构模型示意图，图(b)为Giuffr-Menegotto-Pinto的钢材本构关系示意图；

图3是采用纤维截面方法和opensees软件得到的苏通大桥结构模型；

图4是桥墩处XYZ方向上的力与时间的关系图；

图5是桥墩处XYZ方向上的力矩与时间的关系图；

图6是塔顶处XYZ方向上的位移与时间的关系图。

具体实施方式

下面结合实施例对本发明作进一步地详细描述。

本实施例以苏通大桥为例进行地震响应分析。采用Opensees软件进行模拟，如图1所示为流程示意图。

步骤1、选取el-centrol东西向地震波；

步骤2、用node命令定义600个节点，用element命令定义854个单元杆件，在杆件部位定义恰当的受力节点；计算得到结构质量矩阵M、刚度矩阵K和阻尼矩阵C。

质量矩阵M如下式所示：

m₁为1号节点的质量，m₂为2号节点质量，以此类推，总共600个节点。

刚度矩阵K如下式所示：

其中，K₁₁为1号节点发生单位位移在1处产生的力；K₁₂在2号节点发生单位位移在1处产生的力，以此类推。

阻尼矩阵C如下式所示：

C＝α[M]+β[K]；

其中，α、β由任意二阶振型自拟比，由给定的任意二阶振型的阻尼比ξ_i、ξ_j反算求得。根据振型正交条件，特定常数α和β与振型阻尼比之间的关系应满足：

任意给定两个振型阻尼比ξ_i和ξ_j后，可按下式确定比例常数：

其中，ω_i、ω_j分别为第i、j振型的原频率，取前两阶。

步骤3、定义纤维截面，将纤维点划分成面积为20mm*20mm纤维。

采用section Fiber$secTag<-GJ$GJ>{

patch rect$matTag$numSubdivY$numSubdivZ$yI$zI$yJ$zJ

fiber$yLoc$zLoc$A$matTag

}命令。

定义局部坐标系下y方向的纤维数和z方向的纤维数及顶点坐标。

定义局部坐标系中纤维的坐标与面积及与该纤维相关联的材料标签。

单元设置3～5个积分点就可以取得足够的计算精度。

步骤4、利用Opensees软件采用命令定义选取材料的参数，如uniaxialMaterialSteel01$matTag$Fy$E0$b，定义选用材料的屈服强度，初始弹性切线，应变硬化比(屈服后正切与初始弹性正切之比)，在程序中定义完合适的材料参数有助于后续的数据分析，并能体现一定的真实度。

步骤5、定义梁柱，如采用element elasticBeamColumn$eleTag$iNode$jNode$A$E$G$J$Iy$Iz$transfTa命令。定义梁柱的材料类型，再定义截面面积，杨氏模量，剪切模量，设定杆件惯性矩，坐标转换对象的标识符，以便后续程序的合理分析与计算，能够正确导出在各地震波下的不同地震分析。

步骤6、均布荷载转换成集中质量。删去原来的均布荷载命令。第i号节点采用如下的程序命令：

从而，将结构转化成集中质量模型。

步骤7、导入地震波，利用KR算法命令integrator KRAlphaExplicit。

提前预设ρ_m等于1，并设定分析子结构的初始相对位移向量x₀、初始速度向量

从而获取分析子结构的初始加速度向量

并采用双显式积分法，得到分析子结构第i+1时刻的速度向量

和相对位移向量x_i+1。

其中α₁、α₂、α₃的表达式如下：

α₁＝[M+γΔtC+β(Δt)²K]^-1M

α₃＝[M+γΔtC+β(Δt)²K]^-1[α_mM+α_fγΔtC+α_fβ(Δt)²K]

通过高频谱半径ρ_m从1往0衰减，得到第一个稳定值，便得到地震时刻地震反力与时间的关系图、位移响应-时间关系图。

采用上述算法的分析方法，本桥梁模型所计算出的地震响应分析数据如图4～6所示。图4中，图(a)、(b)、(c)分别为桥墩处X、Y、Z方向上的力与时间的关系图；图5中，图(a)、(b)、(c)分别是桥墩处X、Y、Z方向上的力矩与时间的关系图；图6中，图(a)、(b)、(c)分别是桥墩处X、Y、Z方向上的力矩与时间的关系图。

因此本实施例通过以苏通大桥为例，采用集中质量模型，并运用KR算法得出苏通大桥在东西向的el-centrol波下的塔顶、桥墩响应图像；同时可以采用上述方法进行各种不同类型桥梁的地震响应分析。