对于环境激励下复杂结构的响应分析, 利用限宽白噪声激振手段获取准确的响应幅值估计是较为困难的. 如果使用理想白噪声模型, 利用李雅普诺夫方程求解结构响应的协方差矩阵^[10], 可以很方便地得到白噪声激励下结构的振动信号的方差, 从而估计结构的响应水平.本节以基底输入白噪声激励的三层剪切框架(见图11)为例, 利用理想白噪声模型和李雅普诺夫方程求取框架顶层的相对位移及 绝对加速度响应的理论方差值, 并与数值仿真结果进行对比.

图11

	Figure Option ViewDownloadNew Window
	图11 基底输入白噪声的三层剪切框架

结构的质量矩阵\({\pmb M}\)与刚度矩阵\({\pmb K}\)分别为\[{\pmb M} = \left[ \!\! \begin{array}{ccc} {m_1 } & 0 & 0 0 & {m_2 } & 0 0 & 0 & {m_3 }\end{array} \right] {\pmb K} = \left[\!\! \begin{array}{ccc} {k_1 + k_2 } & { - k_2 } & 0 { - k_2 } & {k_2 + k_3 } & { - k_3 } 0 & { - k_3 } & {k_3 }\end{array} \right] \]结构运动方程为\[{\pmb M}\ddot{\pmb x} + {\pmb C}\dot{\pmb x} + {\pmb K}{\pmb x} = - {\pmb M}\ddot{\pmb u}_{\rm g} \ \ (10)\]式中, \({\pmb x}\)为结构各层相对位移, \(\ddot{\pmb u}_{\rm g}\)为基底激励输入.令状态向量为\({\pmb y} = \left[\!\! \begin{array}{c} {{\pmb x}(t)_{3\times 1} } {\dot{\pmb x}(t)_{3\times 1} }\end{array}\!\! \right]\), 输入向量为\({\pmb w}(t) = \ddot{\pmb u}_{\rm g} (t)\), 构造状态方程\[\dot{\pmb y} ={\pmb A}{\pmb y} + {\pmb B}{\pmb w} \ \ (11)\]式中, 系统矩阵\({\pmb A}\)与控制矩阵\({\pmb B}\)分别为\[{\pmb A }= \left[\!\! \begin{array}{cc} { {\pmb O}_{3\times 3} } & { {\pmb I}_{3\times 3} } { - {\pmb M}^{ - 1}{\pmb K}} & { - {\pmb M}^{ - 1}{\pmb C}}\end{array} \!\! \right]\, {\pmb B} = \left[ \begin{array}{c} { {\pmb O}_{3\times 1} } { - 1} { - 1} { - 1}\end{array} \right]\]式(11)左右同乘以\({\pmb y}^{\rm T}\)有\[\dot{\pmb y}{\pmb y}^{\rm T} = {\pmb A}{\pmb y}{\pmb y}^{\rm T} + {\pmb B}{\pmb w}{\pmb y}^{\rm T} \ \ (12)\]求其期望\[ \begin{array}{l} E[ \dot{\pmb y}{\pmb y}^{\rm T}] = E[ A{\pmb y}{\pmb y}^{\rm T}] + E[ B{\pmb w}{\pmb y}^{\rm T}] = A E[{\pmb y}{\pmb y}^{\rm T}] + B E[{\pmb w}{\pmb y}^{\rm T}] \end{array} \ \ (13)\]令\(\varGamma = E[{\pmb y}{\pmb y}^{\rm T}]\), 有\[\dfrac{d{\pmb\varGamma} }{d t} = d( E[{\pmb y}{\pmb y}^{\rm T}])= E[ \dot {\pmb y}{\pmb y}^{\rm T} + {\pmb y}\dot {\pmb y}^{\rm T}] \ \ (14)\]代入式(13)\[\dfrac{d{\pmb\varGamma} }{d t} = E[ \dot{\pmb y}{\pmb y}^{\rm T} + {\pmb y}\dot {\pmb y}^{\rm T}] = {\pmb A} E[{\pmb y}{\pmb y}^{\rm T}] + {\pmb B} E[{\pmb w}{\pmb y}^{\rm T}] + E[{\pmb y}{\pmb y}^{\rm T}]{\pmb A}^{\rm T} +E[{\pmb y}{\pmb w}^{\rm T}]{\pmb B}^{\rm T} = {\pmb A}{\pmb\varGamma} + {\pmb\varGamma}^{\rm T}{\pmb A}^{\rm T} + {\pmb B} E[{\pmb w}{\pmb y}^{\rm T}] + E[{\pmb y}{\pmb w}^{\rm T}]{\pmb B}^{\rm T} \ \ (15)\]

已知结构状态向量\({\pmb y}\)在零初值下的积分解\[{\pmb y} = \int_0^t {\rm e}^{{A}(t - \tau )} {\pmb B}{\pmb w}(\tau )d\tau \ \ (16)\]因此有\[ E[ {\pmb w}(t) {\pmb y}^{\rm T}(t)] = E[{\pmb w} (t)\int _0^t {{\pmb w}^{\rm T}(\tau )} {\pmb B}^{\rm T}{\rm e}^{{A}^{\rm T}(t - \tau )}d \tau ] =\ \ E[ \int _0^t {{\pmb w} (t) {\pmb w}^{\rm T}(\tau )} {\pmb B}^{\rm T}{\rm e}^{{A}^{\rm T}(t - \tau )}d \tau ] =\ \ \int _0^t {E[{\pmb w} (t) {\pmb w}^{\rm T}(\tau )} ] {\pmb B}^{\rm T}{\rm e}^{{A}^{\rm T}(t - \tau )}d \tau \ \ (17)\]由于\({\pmb w} (t)\)为理想白噪声, 参考式(1), 有\[E[{\pmb w} (t) {\pmb w}^{\rm T}(\tau )] = 2\pi S_0 \delta (\tau - t) \ \ (18)\]故\[\begin{array}{c} E[{\pmb w} (t) {\pmb y}^{\rm T}(t)] =\ \ E \Big[ \int _0^t {2\pi S_0 } \delta (\tau - t) {\pmb B}^{\rm T}{\rm e}^{{A}^{\rm T}(t - \tau )} d \tau \Big] =\ \ E[ \pi S_0 {\pmb B}^{\rm T}] = \pi S_0 {\pmb B}^{\rm T} \end{array} \ \ (19)\]

注意, 上式积分过程中, \(\int _0^t {f(\tau )\delta (\tau - t)} d(\tau ) =\dfrac{1}{2}f(t)\), 其原因在于\(\delta (\tau - t)\)处于积分右边界上, \(\delta(\tau - t)\)的积分只能计入一半的范围, 见图12所示.

图12

	Figure Option ViewDownloadNew Window
	图12 函数积分边界上的狄拉克\(\delta\)函数

代入式(15), 有\[ \dfrac{d{\pmb\varGamma} }{d t} = {\pmb A}{\pmb\varGamma} + {\pmb\varGamma}^{\rm T}{\pmb A}^{\rm T} + \pi S_0{\pmb B}{\pmb B}^{\rm T} + \pi S_0 {\pmb B}{\pmb B}^{\rm T}= {\pmb A}{\pmb\varGamma} +{\pmb\varGamma}^{\rm T}{\pmb A}^{\rm T} + 2\pi S_0 {\pmb B}{\pmb B}^{\rm T} \ \ (20)\]

对于时不变系统, 有\[\dfrac{d{\pmb\varGamma}}{d t} = {\bf 0} \ \ (21)\]

故\[{\pmb A}{\pmb\varGamma}+{\pmb\varGamma}^{\rm T}{\pmb A}^{\rm T} + 2\pi S_0 {\pmb B}{\pmb B}^{\rm T} = {\bf 0} \ \ (22)\]

满足李雅普诺夫方程. 注意到李雅普诺夫方程中未考虑观测矩阵, 求出的结果直接为状态向量的统计值. 为得到期望输出的方差估计, 需要在结果中增加相应的观测矩阵\({V}\)

\[{\pmb\varSigma} ={\pmb V}{\pmb\varGamma }{\pmb V}^{\rm T} \ \ (23)\]

对于本算例, 希望估计顶层的相对位移与绝对加速度的方差值, 观测矩阵、转换矩阵(用于在3个自由度的位移、加速度的估计值中提取第3层的数据)、方差估计分别为\[{\pmb V} = \left[ \!\!\begin{array}{cc} { {\pmb I}_{3\times 3} } & { {\pmb O}_{3\times 3} } { - {\pmb M}^{ - 1}{\pmb K}} & { - {\pmb M}^{ - 1}{\pmb C}}\end{array}\!\! \right]_{6\times 6} {\pmb T} = \left[\!\!\begin{array}{cccccc} 0 & 0 & 1 & 0 & 0 & 0 0 & 0 & 0 & 0 & 0 & 1 \end{array}\!\!\right] {\pmb\varSigma} = ({\pmb T}{\pmb V}) {\pmb \varGamma} ({\pmb T}{\pmb V})^{\rm T} \]

令算例中的参数取值为

\[ m_1 = m_2 = m_3 = 1\, 000\, {\rm kg}\\k_1 = k_2 = k_3 = 200\, {\rm MN} \\ {\pmb C} = 0.2 {\pmb M}+ 0.002 {\pmb K} \]

限宽白噪声的功率谱密度\(S_0 = 1\); 代入式(22)和式(23), 得\[{\pmb\varSigma} = \left[ \!\!\begin{array}{cc} {0.42} & {-17.06} \ { - 17.06} & {761.20}\end{array} \right]\]

可知顶层相对位移与绝对加速度的方差值分别为0.42与761.20, 注意该结果是在理想白噪声模型下得到的, 是理论解. 利用MatlabSimulink工具箱, 以及BLWN限宽白噪声工具, 进行信号模拟. 信号采样频率为1 000Hz, 信号时长为10 000s, 可得到结构顶层位移与加速度响应信号. 之后求取信号的协方差矩阵, 即得到顶层相对位移与绝对加速度的信号方差值\[{\pmb\varSigma} _{\rm Simulink} = \left[\!\! \begin{array}{cc} { 0.42} & { - 16.99} { - 16.99} & { 757.75}\end{array} \right]\]

由于结构为线性系统, 顶层响应同样也为白噪声. 对于顶层绝对加速度响应, 利用\(S_0 = \dfrac{\sigma ^2_{\rm s}}{2\pi f}\), 可求得信号的功率谱密度为0.12.也就是说, 对算例中的三层剪切框架, 功率谱密度为1的白噪声基底加速度, 可激励起结构顶层0.12的加速度响应.数值仿真的结果 与利用理想白噪声模型和李雅普诺夫方程计算的理论解十分接近.当降低信号采样频率, 或减少信号时长时, 数值计算的结果与理论值相比会有一定的误差. 注意在Matlab Simulink BLWN的设置中, 噪声为\(2\pi S_0 \).

调频质量阻尼器(tuned mass damper, TMD)通常简化称为TMD系统, 其由质量、弹簧、阻尼器构成, 通过改善结构共振特性来减振.自1909年Frahm首次提出后已在减振结构中广为应用\(^{[11\hbox{-}12]}\).本节以广州新电视塔TMD系统为研究对象, 利用随机地震金井清模型\(^{[13\hbox{-}14]}\)对结构及TMD系统的响应统计进行分析, 并对TMD系统进行参数优化. 结构模型选已标定的38节点模型^[15](见图13), 其中1号节点嵌固于地面, 1\(\sim\)27单元为塔体, \n 28\(\sim \)37单元为天线桅杆.每个节点有5个自由度(依次为2个平动、3个转动), 结构共185个自由度, 其质量、阻尼、刚度阵分别为\({\pmb M}\), \({\pmb C}\), \({\pmb K}\). TMD系统安装于塔体顶层, 即第28个自由度. TMD系统为2个600\, t (\(m_{\rm TMD}\))的消防水箱, 设计TMD弹簧刚度\(k_{\rm TMD} = 541\)\, kN/m, 阻尼系数\(c_{\rm TMD} = 2.04\times 10^5\). 本节仅讨论结构\({\pmb x}\)向平动.

图13

	Figure Option ViewDownloadNew Window
	图13 广州新电视塔37结点模型

构造结构与TMD系统整体的质量阵\(\bar {\pmb M}\)、阻尼阵\(\bar {\pmb C}\)、刚度阵\(\bar {\pmb K}\), 有\[ \bar {\pmb M} = \left[ \!\! \begin{array}{cc} {\pmb M} & {\bf 0} \ {\bf 0 } & {\pmb m}_{\rm TMD }\end{array} \!\! \right]\]

$$ \bar {\pmb C} = \left[\!\! \begin{array}{cc} \left[\!\! \begin{array}{ccc} \cdots & \cdots & \cdots \\ \cdots & { \cdots + c_{\rm TMD} } & \cdots \\ \cdots & \cdots & \cdots\end{array} \!\! \right]_C & \begin{array}{c} \vdots \\ { - c_{\rm TMD } } \\ \vdots \end{array} \\ \begin{array}{ccc} \cdots & - c_{\rm TMD } & \cdots\end{array} & c_{\rm TMD } \\\end{array} \right]$$

\[ \bar {\pmb K} = \left[\!\!\begin{array}{cc} \left[\!\!\begin{array}{ccc} \cdots & \cdots & \cdots \\ \cdots & \cdots + k_{\rm TMD } & \cdots \\ \cdots & \cdots & \cdots\end{array} \!\! \right]_K & \begin{array}{c} \vdots \\{ - k_{\rm TMD } } \\ \vdots\end{array} \\ \begin{array}{ccc} \cdots & { - k_{\rm TMD} } & \cdots\end{array} & {k_{\rm TMD } } \end{array} \right]\]

\(\bar{\pmb C}\)和\(\bar{\pmb K}\)中省略号表示不含TMD系统的原结构阻尼和刚度阵的参数, 在此省略.金井清模型基岩加速度为\(\ddot {\pmb U}\). 地表土层相对加速度为\(\ddot {\pmb x}_{\rm g} \), 绝对加速度 为\(\ddot{\pmb u}_{\rm g} \)

\[\ddot {\pmb u}_{\rm g} = \ddot {\pmb x}_{\rm g} + \ddot {\pmb U} = - (2\zeta _{\rm g} \omega _{\rm g}\dot{\pmb x}_{\rm g} - \omega _{\rm g}^2 {\pmb x}_{\rm g} )\]

式中取\(\omega _{\rm g} = 4.5\)\, rad/s, \(\zeta _{\rm g}= 0.1\), 为软土地基.结构整体的运动方程\[\bar{\pmb M}\ddot{\pmb x} + \bar{\pmb C}\dot{\pmb x} + \bar{\pmb K}{\pmb x} = - \bar {\pmb M}{\pmb L}\ddot{\pmb u}_{\rm g} \ \ (24)\]

式中\({\pmb L}\)为单位阵, 提取结构\({\pmb x}\)向的平动自由度. 构造状态向量与状态方程

\[{\pmb y} = \left[ {\pmb x} \ \ {{\pmb x}_{\rm g} } \ \ \dot {\pmb x} \ \ {\dot{\pmb x}_{\rm g} } \right]^{\rm T}\ \ (25)\]

\[\dot {\pmb y} = {\pmb A}{\pmb y} + {\pmb B}{\pmb w}\ \ (26)\]

式中

\[{\pmb A} = \left[ \begin{array}{cccc}{\bf 0} & {\bf 0} & {\bf 1} & {\bf 0} \\ 0 & 0 & 0 & 1 \\ { - \bar {\pmb M}^{ - 1}\bar{\pmb K}} & {\omega _{\rm g}^2 \cdot {\pmb I}} & { -\bar {\pmb M}^{ - 1}\bar {\pmb C}} & { - 2\zeta _{\rm g} \omega _{\rm g} \cdot {\pmb I}} \\ 0 & { - \omega _{\rm g}^2 } & 0 & {2\zeta _{\rm g} \omega _{\rm g} } \end{array} \right]\]

\[{\pmb B} = \left[ {\bf 0} \ \ 0 \ \ {\pmb L} \ \ { - 1} \right]^{\rm T} \]

\({\pmb w}\)为白噪声, 功率谱密度\(S_0 = 1\).同上节, 令\({\pmb\varGamma} = E[ {\pmb y}{\pmb y}^{\rm T}]\), 则有\[{\pmb A}{\pmb \varGamma} + {\pmb \varGamma}^{\rm T}{\pmb A}^{\rm T} + 2\pi S_0 {\pmb B}{\pmb B}^{\rm T} = {\bf 0} \ \ (27)\]

提取结构及TMD系统的相对位移(相对地面)方差、绝对加速度方差进行响应分析. 以塔体(不包括天线桅杆)顶层相对位移标准差为控制目标, 对TMD参数进行敏感性分析, 见图14所示.图中, 横坐标为TMD阻尼比, 纵坐标为TMD自振频率与塔体1阶自振频 率之比. 图中横向虚线自下向上分别表示结构1\(\sim\)10阶自振频率, 横向点划线表示塔体设计TMD频率所在位置. 可以看出, 设计TMD的频率基本处于对塔体顶层位移的最优控制范围.

仍然以塔体顶层相对位移为控制目标, 利用约束置信域优化方法对TMD参数进行优化计算, 得到最优的参 数为\(k_{\rm TMD} = 554\)\, kN/m, 阻尼比\(\zeta _{\rm TMD}= 0.088\) (或阻尼系数\(c_{\rm TMD} = 1.43\times 10^5\)). 此时塔体及TMD系统的整体位移状态见图15所示.与TMD系统的设计参数(主要针对风载荷激励进行优化)相比, 刚度参数较为近似, 阻尼有一定差别.与文献[16]相比, 本文优化刚度更接近实际取值(弱轴方向), 且二者均存在对天线振动控制不佳的现象; 与文献[17]相比, TMD在地震作用下优化计算得到的阻尼比接近.另一方面, 相对时程分析方法, 随机振动方法计算速度更快, 在复杂结构的质量调谐阻尼器的参数优化中可凸显其计算效率高的优点.

图14

	Figure Option ViewDownloadNew Window
	图14 TMD参数敏感性分析

图15

	Figure Option ViewDownloadNew Window
	图15 结构TMD系统优化后的结点位移图

本文对理想白噪声模型与限宽白噪声模型的模型概念、时频特征等进行了详细分析讨论, 其中包括信号的自相关函数与功率谱密度函数分布、信号的水平和带宽控制等. 对于连续系统, 推导了连续限宽白噪声的自相关函数分布, 并利用数值仿真对信号的频域特征进行了近似实现. 详细推导了利用李雅普诺夫方程与理想白噪声模型预测结构振动信号的方差特性的方法, 并结合三层框架进行了算例展示. 结合广州新电视塔及其TMD系统, 对利用随机振动方法进行结构控制优化设计的问题进行了研究. 本文为结构振动分析中的白噪声模型的使用提供了较为完整的理论思路.