matlab杜哈梅尔积分,非平稳地震作用下桥梁系统碰撞间隙宽度的概率评价方法与流程…

基于式(1-6)可得到不同地震强度水平作用下桥梁梁体之间和梁体与桥台 之间的相对位移(△)的响应峰值，从而得到地震作用下梁体间的间隙宽度需求均 值和均方差，为超越碰撞间隙宽度的条件概率需求模型奠定基础。

2)碰撞间隙宽度条件碰撞概率模型

据研究发现相对位移峰值响应Xpeak(t)为随机过程且服从极值I型分布。 Xpeak(t)概率分布函数可表示为

G(Xpeak(t))＝exp{-exp{-αn(Xpeak(t)-un)}} (7)

式中

不同震级下的相对位移响应峰值均值和均方差由式(5)和式(6)得到并表示为 各震级条件下峰值均值和均方差分别为和那么在特定某一震级 下超过间隙宽度并使得桥梁发生碰撞的概率被称为条件碰撞概率。条件碰撞超 越概率表示如下

Pp/m(t)＝1-exp{-exp{-αn(Xpk/m(t)-un)}}(8) 式中

从式(8)中可知，震级条件下峰值均值和均方差已知情况下， Xpk/m(t)的分布可以判断任意间隙宽度时的条件碰撞超越概率。

3)基于矩阵的系统可靠度

假设桥梁结构第i个碰撞点有di个设计间隙宽度，i＝1,…,n，则样本空间被分 为个穷尽互斥基本事件，且用ej表示，j＝1,…,m。c表示任意事件系数 向量且由0和1组成。ej发生的概率可表示为pj＝P(ej)，j＝1,…,m。由于事件ej互斥，那么非平稳地震作用下系统事件的概率Esys(t)可表示为ej事件之和，则系 统条件概率表示为两向量的内积

式中p为概率矩阵的列向量。

桥梁系统每一个碰撞点在多个设计间隙宽度的情况下，系统碰撞事件概率 可表示为Psys(t)＝CTP(t)，Psys(t)表示一个随时间变化的概率矩阵，其中i行j列的 元素表示在第j个条件下第i个事件在时间t的概率。只要确定C和P(t)中各元 素，则碰撞系统事件概率可方便得到。系数矩阵C通过以下迭代方式得到

C[n]的第i列是第i个事件向量，记为有如下的运算关系

同理，p也可由矩阵迭代计算得到

其中，Pi表示第i个事件概率，且

步骤1)中式(2)的虚拟力P由激励功率谱分解得到且表征地震动空 间性。

式中，子矩阵元素Si’j'(iω,t)为3×3维矩阵，参照附图1，相应两个水平向(x，y) 和竖向(z)的地震动分量组成的互功率谱矩阵表示为

三维地面运动的水平向(x,y)分量取相同的功率谱密度函数，而竖向(z)的功率谱 密度取其水平向的0.6倍。故式(15)的非对角线元素为

式中，Si’j’,xx(iω,t)，Si’j’,yy(iω,t)，Si’j’,zz(iω,t)由场地传递函数从基岩功率谱密度矩 阵函数得到。

式中，HU,i(iω)，HV,i(iω)，HW,i(iω)和HU,j(iω)，HV,j(iω)，HW,j(iω)分别为空间 点i-i′、j-j′在U，V，W方向的场地传递函数。i，j点x，y，z方向的自功率 谱密度函数为

式中，Gi,x(ω,t)为第i个空间点x向在时域和频域内的调制函数，Sii,xx(iω)是

基岩x向自功率谱密度函数，γij(iω)为空间点i，j的相干函数，形式如下

γij(iω)＝|γij(iω)|exp[-iωdij/vapp(ω)] (19)

其中，dij为空间点i，j距离在地震波传播方向的投影，vapp(ω)表示视波速。

为了得到式(16)-(18)中的频率传递函数，基于一维波动理论，假设基岩中地 震波由面外SH波或者面内P-SV波组成，场地频率传递函数HU,i(iω)，HV,i(iω)， HW,i(iω)表示为

KSHXSH＝PSH，KP-SVXP-SV＝PP-SV (20)

其中，XSH，PSH和XP-SV，PP-SV是关于SH波和P-SV波的位移和荷载向量；刚度 矩阵XSH和KP-SV主要由土体性质、入射波类型、入射角和圆频率ω确定；动力荷 载向量PSH和PP-SV主要依赖于基岩的特性、类型、频率和入射波的幅值。通过式(20) 可以得到在每一个离散频率点面内或面外的频率传递函数。

为了构造虚拟力，三维非平稳地震激励功率谱矩阵可分解为

式中，P是3m×r维矩阵，r是矩阵的秩，上标*和T分别表示复数共 轭和转置。

式(20)中的W(iω)为3m×3m矩阵，表示了行波效应和场地土条件对相位的 变化，表达式如下

其中，Ti,x表示地震波沿x方向到第i个支撑处所需时间，θi,x，θm,y，θm,z分别表 示在第i个支撑处x，y，z方向的相位角，将由下面的公式得到

式中，θij,xx为空间点i-i′和j-j′在x方向的相位差。

式(11)中的Sg(iω)表达式如下

式中，Sgi,x代表空间点i的x方向上的功率谱密度函数，水平向(x,y)与竖向(z) 的功率谱比值为3:2。

式(11)中的G(ω,t)为三维非平稳地面运动非均匀调制函数，形式如下

一般认为各个空间每一个方向上调制函数一致，故可取

Gi,x(ω,t)＝Gi,y(ω,t)＝Gi,z(ω,t)＝G(ω,t) (27)

式(11)中，矩阵R表示多维多点地面运动的相干矩阵，一般是正定或半正定 矩阵，可通过LDLT分解如下

R＝[Q]3m×r[QT]r×3m＝[|γij|]3m×3m(当i＝j时，|γij|＝1；当i≠j时，0≤|γij|

式中，R被分解为非零特征值αj与其相应特征向量(j＝1,2,…,r,r≤3m)的和，其 形式如下

部分相干情况下，第阶虚拟激励表示为

从式(30)可知，相干效应、行波效应、局部场地效应和非平稳性分别通过虚拟激 励项中的W(iω)、Sg(iω)和G(ω,t)来体现。

完全相干和完全不相干虚拟力可表示为

根据式(20)和式(21)可得第j个特征值下的绝对位移响应为

式中

式中，Ij(iω,t)是杜哈梅尔积分，h(t-τ)单位脉冲函数，在部分相干情况的绝对 位移响应功率谱密度函数表示为

同样，可以得到完全相干和完全不相干情况下的响应为

完全相干

式中

完全不相干

式中

由公式(33)、(36)、(38)可以看出，需要在每一个频率点下做一次瞬态分 析计算，虽然相对于传统的随机振动理论来说，直接求解的虚拟激励法已经非 常节约时间了，但是如果离散的频率点太细，那么计算的瞬态分析也就越多， 这样也就非常耗时，故下面引入精细积分方法并在通用有限元软件实现，这样 只需要两个瞬态分析就能确定结构任意一个关注节点的响应，这样大大缩短了 计算时间。

在有限元计算中，式(2)可另写成如下形式

式中

为确定频率点下的非平稳随机激励时程，Mb为桥梁基础处附加的大质量。

将每一个确定频点下的动力方程式(39)写成状态方程形式

其中

状态方程(41)的一般解为

式中eHt为指数矩阵，对式(43)进行离散积分，设时间步长为Δt，由递推法 可知，在ti时刻的响应V(ti)＝Vi可用ti-1时刻的响应V(ti-1)＝Vi-1来表示如下

式中，T＝eH△t为指数矩阵，改写成如下形式

T＝eH△t＝[eH△t/m]m (45)

令τ＝△t/m，当取m＝2N很大的，△τ非常小，由泰勒级数展开式有

当N＝20时，与一般计算机的舍入误差相比较，泰勒级数的截断误差要小的多， 不会由于截断带来数值误差，故T(τ)的计算精度实际上已经给出了计算机上的精 确解。

[I+Tai]≡[I+Ta,i-1]2＝[I+2Ta,i-1+Ta,i-1×Ta,i-1](i＝1,2,…,N) (47)

因此依次类推

[I+TaN]＝[I+Ta,N-1]2＝[I+Ta,N-2]4＝…＝[I+Ta0]m＝T(τ) (48)

每次运算单元矩阵[I]不参与计算，因为Tai很小，当它与单位矩阵[I]相加时 就成为其尾数，在计算机的舍入操作中就会自动舍去。为了避免这种情况的发 生，采用下列递推方式求解

[Tai]＝2[Ta，i-1]+[Ta，i-1][Ta，i-1](i＝1,2,…,N) (49)

假定在每一个时间积分步(ti-1,ti)中，荷载为线性变换，则非平稳虚拟力 F(ω,t)可离散为t0,t1,t2,…,tk处的随机变量F0,F1,F2,…,Fk,非平稳虚拟力可表示为

F(ω,τ)＝r0+r1(τ-ti-1) (50)

式中，r0，r1是一个是不变向量。

将式(50)代入式(44)可得

式中H-1的计算可根据式(42)可得

将式(51)进一步写成下式

令S1＝(I-T)(H-2/△t)+TH-1,S2＝(T-I)(H-2/△t)-H-1 (54)

式中，H-2＝H-1·H-1，则式(53)简化为

Vi＝T△Vi-1+S1Fi-1+S2Fi(i＝1,2,…,k) (55)

每一固定频点时，V0＝0，则通过式(55)可得到固定频率点各时刻的响应表达式为

式中

S3＝TS2+S1 (57)

若用Ai,0，Ai,1，…，Ai,i分别表示F0，F1，…，Fi前面的系数，Ai,i中第一 个i表示ti时刻，第二个表示第i个荷载离散点位置，则式(56)可表示为

Vi＝Ai,0F0+Ai,1F1+…+Ai,i-1Fi-1+Ai,iFi(i＝1,2,…,k) (58)

式中，Ai,0，Ai,1，…，Ai,i只与结构本身有关，而与外荷载无关，体现了结构 本身在外荷载激励下所表现出来的固有属性。Ai,i的计算过程是一个递推过程， 用Vi-1的系数Ai-1,0，Ai-1,1，…，Ai-1,i-1来推导下一步的系数为

将式(59)代入式(58)写成矩阵的形式，则

由式(60)可知，系数矩阵中所有其他所有的元素都是第一列和第二列元素的 重复，也就是说只要能得到第一列和第二列的所有元素，那么整个系数矩阵就 能确定。令F0＝1，Fk＝0(k>0)，进行一个瞬态分析就可以得到系数矩阵的第一列 元素。同理，令F0＝0，F1＝1，Fk＝0(k>2)，进行一个瞬态分析就可以得到系数矩 阵的第二列元素，求解加载方式如附图2。

只要得到系数矩阵，那么任意频率点处的非平稳时程荷载作用下的响应就 能根据式(60)得到。当运用绝对位移直接求解时，对于多维多点激励的非平稳地 震动荷载，需要对每一个激励点三个方向依次进行激励，分别得到所需的结构 关键点响应的系数矩阵A，然后将每一个方向的激励时程与对应的系数矩阵A 相乘后叠加就得到多维多点在固定频率点的结构响应。最后，将所有的频率点 下的激励时程都对应代入后就得到整个频率范围的结构响应。

本发明的技术方案不限于上述具体实施例的限制，凡是根据本发明的技术方案 做出的技术变形，均落入本发明的保护范围之内。

相关资源：凯利公司 KDZ 系列有刷串励、永磁、他励电机控制器设置软件.rar