EI、Scopus 收录
中文核心期刊

含多维随机变量的广义概率密度演化方程解析解: 以Euler-Bernoulli梁为例

周永峰, 李杰

周永峰, 李杰. 含多维随机变量的广义概率密度演化方程解析解: 以Euler-Bernoulli梁为例. 力学学报, 2024, 56(9): 2659-2668. DOI: 10.6052/0459-1879-24-001
引用本文: 周永峰, 李杰. 含多维随机变量的广义概率密度演化方程解析解: 以Euler-Bernoulli梁为例. 力学学报, 2024, 56(9): 2659-2668. DOI: 10.6052/0459-1879-24-001
Zhou Yongfeng, Li Jie. Analytical solutions of the generalized probability density evolution equation with multidimensional random variables: the case of Euler-Bernoulli beam. Chinese Journal of Theoretical and Applied Mechanics, 2024, 56(9): 2659-2668. DOI: 10.6052/0459-1879-24-001
Citation: Zhou Yongfeng, Li Jie. Analytical solutions of the generalized probability density evolution equation with multidimensional random variables: the case of Euler-Bernoulli beam. Chinese Journal of Theoretical and Applied Mechanics, 2024, 56(9): 2659-2668. DOI: 10.6052/0459-1879-24-001
周永峰, 李杰. 含多维随机变量的广义概率密度演化方程解析解: 以Euler-Bernoulli梁为例. 力学学报, 2024, 56(9): 2659-2668. CSTR: 32045.14.0459-1879-24-001
引用本文: 周永峰, 李杰. 含多维随机变量的广义概率密度演化方程解析解: 以Euler-Bernoulli梁为例. 力学学报, 2024, 56(9): 2659-2668. CSTR: 32045.14.0459-1879-24-001
Zhou Yongfeng, Li Jie. Analytical solutions of the generalized probability density evolution equation with multidimensional random variables: the case of Euler-Bernoulli beam. Chinese Journal of Theoretical and Applied Mechanics, 2024, 56(9): 2659-2668. CSTR: 32045.14.0459-1879-24-001
Citation: Zhou Yongfeng, Li Jie. Analytical solutions of the generalized probability density evolution equation with multidimensional random variables: the case of Euler-Bernoulli beam. Chinese Journal of Theoretical and Applied Mechanics, 2024, 56(9): 2659-2668. CSTR: 32045.14.0459-1879-24-001

含多维随机变量的广义概率密度演化方程解析解: 以Euler-Bernoulli梁为例

基金项目: 国家自然科学基金资助项目 (51538010)
详细信息
    通讯作者:

    李杰, 教授, 中国科学院院士, 主要研究方向为随机力学和工程可靠性理论. E-mail: lijie@tongji.edu.cn

  • 中图分类号: O34

ANALYTICAL SOLUTIONS OF THE GENERALIZED PROBABILITY DENSITY EVOLUTION EQUATION WITH MULTIDIMENSIONAL RANDOM VARIABLES: THE CASE OF EULER-BERNOULLI BEAM

  • 摘要: 广义概率密度演化方程的解析解, 不仅具有重要理论价值, 而且具有校验数值解、进而标定数值算法误差的作用. 以Euler-Bernoulli简支梁为例, 推导给出了梁受迫振动时跨中位移响应所对应的广义概率密度演化方程解析解. 包括非平稳非高斯随机载荷作用下的解(包含2维随机变量)以及同时考虑载荷随机性和结构参数随机性时的解(分别包含2维、4维和5维随机变量). 分析结果表明, 真实的概率密度演化是一个十分复杂的过程, 远不能用简单的概率分布函数加以描述. 这一进展, 可为概率密度演化理论的进一步深入研究提供一个方面的基础.
    Abstract: The analytical solution of the generalized probability density evolution equation not only holds significant theoretical value, but also serves the purpose of validating numerical solutions and subsequently calibrating the errors in numerical algorithms. However, the analytical solution for this equation is limited to a small number of simple systems under a single random variable. Therefore, taking the Euler-Bernoulli simply supported beam as an example, the analytical solution of the generalized probability density evolution equation corresponding to the mid-span displacement response of the beam under forced vibration is derived. The solutions include those under non-stationary and non-Gaussian random excitations (involving 2-dimensional random variables) as well as those considering both the randomness of excitations and structural parameters (involving 2-dimensional, 4-dimensional and 5-dimensional random variables, respectively). The analysis results indicate that the real evolution of probability density is a highly intricate process that cannot be described through simple probability distribution functions. This advancement can provide a foundational aspect for further in-depth research into the theory of probability density evolution.
  • 随机系统分析是兼具理论价值与工程应用前景的重要研究领域. 从随机性在物理系统中传播这一全新视角[1]研究随机物理系统, 为随机系统研究确立了基本路线. 矩传递和概率密度演化, 被认为是两种主要的随机特征传播形式. 其中, 基于矩传递的传统分析方法, 例如经典随机振动理论、随机摄动理论和正交多项式展开理论等已发展成熟[2-4], 但难以应用于非线性系统. 概率密度演化理论的建立, 为求解复杂随机系统, 尤其是非线性系统奠定了基础[5-8]. 相比经典的概率密度方程, 如Liouville方程[9]、FPK方程[10-12]和Dostupov-Pugachev方程[13], 广义概率密度演化方程的维数可以降低至一维, 从而使得概率密度演化方程的求解大为简化.

    在一般意义上, 广义概率密度演化方程具有解析解[7,14]. 但受其形式——包含Dirac函数的多维积分的限制, 已获得的解析解尚局限于仅含有单一随机变量的情形[15-17]. 对于含有多维随机变量系统的广义概率密度演化方程解析解, 尚无系统研究. 这不仅限制了理论的发展, 也给在多维概率空间内研究广义概率密度演化方程的数值解并定量分析其误差带来了困难. 为此, 本文以Euler-Bernoulli梁受迫振动为例, 根据将多维积分转化为一维积分的基本思想, 研究了多维随机变量系统的广义概率密度演化方程解析解.

    考虑一般的随机动力系统, 其状态方程可写为

    $$ {\boldsymbol{\dot X}} = {\boldsymbol{G}}\left( {{\boldsymbol{X}},{\boldsymbol{\varTheta }},t} \right),{\text{ }}{\boldsymbol{X}}\left( 0 \right) = {{\boldsymbol{x}}_0} $$ (1)

    式中, ${\boldsymbol{X}} = {\left( {{X_1},{X_2}, \cdots ,{X_n}} \right)^{\text{T}}}$为$n$维状态向量; ${\boldsymbol{\varTheta }} = \left( {{\varTheta }_1}, {{\varTheta }_2}, \cdots , {{\varTheta }_s} \right)^{\text{T}}$为$s$维源随机向量; ${\boldsymbol{G}}$为系统函数向量; ${{\boldsymbol{x}}_0} = {\left( {{x_{1,0}},{x_{2,0}}, \cdots ,{x_{n,0}}} \right)^{\text{T}}}$为确定性初始条件.

    对于适定的动力系统, 在给定初始条件后方程(1)存在唯一的连续解答

    $$ {\boldsymbol{X}} = {\boldsymbol{H}}\left( {{\boldsymbol{\varTheta }},t} \right) $$ (2)

    根据概率守恒原理的随机事件描述, 随机系统演化过程中的某一随机事件$\left\{ {\left[ {{\boldsymbol{X}}\left( t \right),{\boldsymbol{\varTheta }}} \right] \in {{\varOmega }_t} \times {{\varOmega }_{\boldsymbol{\varTheta }}}} \right\}$, 在经历微小时间${\text{d}}t$之后, 应满足如下关系[1,7]

    $$ \begin{split} & \int_{{{\varOmega }_{t + {\text{d}}t}} \times {{\varOmega }_{\boldsymbol{\varTheta }}}} {{p_{{\boldsymbol{X\varTheta }}}}\left( {{\boldsymbol{x}},{\boldsymbol{\theta }},t + {\text{d}}t} \right){\text{d}}{\boldsymbol{x}}{\text{d}}{\boldsymbol{\theta }}} = \\ &\qquad \int_{{{\varOmega }_t} \times {{\varOmega }_{\boldsymbol{\varTheta }}}} {{p_{{\boldsymbol{X\varTheta }}}}\left( {{\boldsymbol{x}},{\boldsymbol{\theta }},t} \right){\text{d}}{\boldsymbol{x}}{\text{d}}{\boldsymbol{\theta }}} \end{split} $$ (3)

    式中, ${p_{{\boldsymbol{X\varTheta }}}}\left( {{\boldsymbol{x}},{\boldsymbol{\theta }},t} \right)$为$t$时刻$\left[ {{\boldsymbol{X}}\left( t \right),{\boldsymbol{\varTheta }}} \right]$的联合概率密度函数; ${{\varOmega }_t} \times {{\varOmega }_{\boldsymbol{\varTheta }}}$为$t$时刻${\boldsymbol{X}}\left( t \right)$与${\boldsymbol{\varTheta }}$所在的值域; $ \times $表示${{\varOmega }_t}$和${{\varOmega }_{\boldsymbol{\varTheta }}}$共同张成概率空间.

    由式(3), 可以导出广义概率密度演化方程[1,7]

    $$ \frac{{\partial {p_{{\boldsymbol{X\varTheta }}}}\left( {{\boldsymbol{x}},{\boldsymbol{\theta }},t} \right)}}{{\partial t}} + \sum\limits_{i = 1}^m {{{\dot X}_i}\left( {{\boldsymbol{\theta }},t} \right)} \frac{{\partial {p_{{\boldsymbol{X\varTheta }}}}\left( {{\boldsymbol{x}},{\boldsymbol{\theta }},t} \right)}}{{\partial {x_i}}} = 0 $$ (4)

    通常, 感兴趣的系统物理量可以表示为系统状态量的函数

    $$ {\boldsymbol{Z}} = h\left( {{\boldsymbol{\dot X}},{\boldsymbol{X}},{\boldsymbol{\varTheta }},t} \right) $$ (5)

    因此, 对于一般状态向量, 存在如下广义概率密度演化方程

    $$ \frac{{\partial {p_{{\boldsymbol{Z\varTheta }}}}\left( {{\boldsymbol{z}},{\boldsymbol{\theta }},t} \right)}}{{\partial t}} + \sum\limits_{j = 1}^m {{{\dot Z}_j}\left( {{\boldsymbol{\theta }},t} \right)} \frac{{\partial {p_{{\boldsymbol{Z\varTheta }}}}\left( {{\boldsymbol{z}},{\boldsymbol{\theta }},t} \right)}}{{\partial {z_j}}} = 0 $$ (6)

    当仅考察某一具体物理量$Z$时, 上述偏微分方程退化为一维形式

    $$ \frac{{\partial {p_{Z{\boldsymbol{\varTheta }}}}\left( {z,{\boldsymbol{\theta }},t} \right)}}{{\partial t}} + \dot Z\left( {{\boldsymbol{\theta }},t} \right)\frac{{\partial {p_{Z{\boldsymbol{\varTheta }}}}\left( {z,{\boldsymbol{\theta }},t} \right)}}{{\partial z}} = 0 $$ (7)

    并有初始条件

    $$ {p_0} = {p_{Z{\boldsymbol{\varTheta }}}}\left( {z,{\boldsymbol{\theta }},{t_0}} \right) = \delta \left( {z - {z_0}} \right){p_{\boldsymbol{\varTheta }}}\left( {\boldsymbol{\theta }} \right) $$ (8)

    采用特征线法, 可以给出广义概率密度演化方程(7)的解析解

    $$ \begin{split} & {p_{Z{\boldsymbol{\varTheta }}}}\left( {z,{\boldsymbol{\theta }},t} \right) = {p_0}\left[ {z - \int_0^t {\dot Z\left( {{\boldsymbol{\theta }},\tau } \right){\mathrm{d}}\tau } } \right] = \\ &\qquad \delta \left[ {z - {h_Z}\left( {{\boldsymbol{\theta }},t} \right)} \right]{p_{\boldsymbol{\varTheta }}}\left( {\boldsymbol{\theta }} \right) \end{split} $$ (9)

    关于源随机向量${\boldsymbol{\varTheta }}$积分, 即可得到$Z$随$t$变化的概率密度函数

    $$ {p_Z}\left( {z,t} \right) = \int_{{{\varOmega }_{\boldsymbol{\varTheta }}}} {\delta \left[ {z - {h_Z}\left( {{\boldsymbol{\theta }},t} \right)} \right]{p_{\boldsymbol{\varTheta }}}\left( {\boldsymbol{\theta }} \right){\text{d}}{\boldsymbol{\theta }}} $$ (10)

    对于绝大多数系统, 式(10)难以直接求解. 因此, 发展出了一套完整的数值求解流程, 关键步骤包括概率空间剖分和偏微分方程的数值求解. 前者形成了以数论选点[18-19]或Sobol点集[20]为基础的初始选点策略及基于GF偏差最小的点集调整策略[21]. 后者则包含有限差分[7,22]、有限元[23]、再生核函数配点[24]及群演化[25]在内的多种算法. 这些数值方法的适用性均已得到验证, 但其对比依据, 仍然是计算效率低且具有随机误差的蒙特卡洛模拟[26-27]. 本文后续给出的多维随机变量下的广义概率密度演化方程解析解, 则有望成为新的数值算法参考基准.

    图1所示等截面Euler-Bernoulli简支梁, 其受迫振动的运动方程和边界条件为[28]

    图  1  Euler-Bernoulli简支梁跨中集中载荷模型
    Figure  1.  Model of Euler-Bernoulli simply supported beam subjected to midspan concentrated load
    $$\left. \begin{split} & \frac{{{\partial ^4}y\left( {x,t} \right)}}{{\partial {x^4}}} + \frac{{\bar m}}{{EI}}\frac{{{\partial ^2}y\left( {x,t} \right)}}{{\partial {t^2}}} = P\left( t \right) \\ & y\left( {0,t} \right) = 0,{\text{ }}{y{''}}\left( {0,t} \right) = 0 \\ & y\left( {L,t} \right) = 0,{\text{ }}{y{''}}\left( {L,t} \right) = 0 \end{split}\right\} $$ (11)

    式中, $L$为梁的长度; $\bar m$为单位长度梁的质量; $E$为弹性模量; $I$为截面惯性矩; $y\left( {x,t} \right)$为$x$处$t$时刻的横向位移; $P\left( t \right)$为随时间$t$变化的横向载荷.

    依据结构动力学知识, 当作用图2(a)所示的阶跃载荷$P\left( t \right) = {A_0}\delta \left( {t - {t_0}} \right)$时, 梁有位移响应

    图  2  不同的载荷形式
    Figure  2.  Different forms of load
    $$ y\left( {x,t} \right) = \frac{{2{A_0}{L^3}}}{{{{\text{π}}^4}EI}}\sum\limits_{n = 1}^\infty {\frac{{{\alpha _n}}}{{{n^4}}}\left( {1 - \cos {\omega _n}t} \right)\sin \frac{{n\text{π} x}}{L}} $$ (12)

    当作用图2(b)所示的简谐载荷$P\left( t \right) = {A_0}\cos \left( {\bar \omega t + \phi } \right)$时, 梁有位移响应

    $$ \begin{split} & y\left( {x,t} \right) = \frac{{2{A_0}{L^3}}}{{{{\text{π}}^4}EI}}\sum\limits_{n = 1}^\infty {\frac{{{\alpha _n}{\omega _n}}}{{{n^4}}}} \sin \frac{{n{\text{π}}x}}{L} \cdot \\ &\qquad \frac{{{\omega _n}\left[ {\cos {\omega _n}t\cos \phi - \cos \left( {\bar \omega t + \phi } \right)} \right] - \bar \omega \sin {\omega _n}t\sin \phi }}{{{{\bar \omega }^2} - {\omega _n}^2}}\end{split} $$ (13)

    式中, ${\omega _n}$为梁的自振频率

    $$ {\omega _n} = {n^2}{{\text{π}}^2}\sqrt {\frac{{EI}}{{\bar m{L^4}}}} ,{\text{ }}n = 1,2, \cdots $$ (14)

    ${\alpha _n}$为与$n$相关的系数

    $$ {\alpha }_{n} = \left\{\begin{split} &1,\quad n = 1,5,9,\cdots \\ &-1,\quad n = 3,7,11,\cdots \\ &0,\quad n = {\mathrm{even}}\;{\mathrm{number}}\end{split}\right. $$ (15)

    对于图3所示等截面Euler-Bernoulli简支梁的受迫振动, 其运动方程和边界条件同式(11).

    图  3  Euler-Bernoulli简支梁3分点集中载荷模型
    Figure  3.  Model of Euler-Bernoulli simply supported beam subjected to third-point concentrated load

    当同时作用阶跃载荷${P_1}\left( t \right) = {A_1}\delta \left( {t - {t_0}} \right)$和简谐载荷${P_2}\left( t \right) = {A_2}\cos \left( {\bar \omega t + \phi } \right)$时, 梁的位移响应为

    $$ \begin{split} & y\left( {x,t} \right) = \frac{{2{A_1}{L^3}}}{{{{\text{π}}^4}EI}}\sum\limits_{n = 1}^\infty {\frac{{{\alpha _{1n}}}}{{{n^4}}}\left( {1 - \cos {\omega _n}t} \right)\sin \frac{{n\text{π} x}}{L}} + \\ &\qquad \frac{{2{A_2}{L^3}}}{{{{\text{π}}^4}EI}}\sum\limits_{n = 1}^\infty {\frac{{{\alpha _{2n}}{\omega _n}}}{{{n^4}}}\sin \frac{{n{\text{π}}x}}{L}} \cdot \\ &\qquad \frac{{{\omega _n}\left[ {\cos {\omega _n}t\cos \phi - \cos \left( {\bar \omega t + \phi } \right)} \right] - \bar \omega \sin {\omega _n}t\sin \phi }}{{{{\bar \omega }^2} - {\omega _n}^2}} \end{split} $$ (16)

    式中, ${\alpha _{in}}\left( {i = 1,2} \right)$为与$n$相关的系数

    $$ {\alpha _{1n}} = \left\{\begin{aligned} & {\sqrt 3 /2,{\text{ }}n = 6k - 5 \cup 6k - 4{\text{ }}\left( {k = 1,2, \cdots } \right)} \\[-2pt] & { - \sqrt 3 /2{\text{, }}n = 6k - 2 \cup 6k - 1{\text{ }}\left( {k = 1,2, \cdots } \right)} \\[-2pt] & {0,{\text{ }}n = 6k - 3 \cup 6k{\text{ }}\left( {k = 1,2, \cdots } \right)} \end{aligned}\right. $$ (17)
    $$ {\alpha _{2n}} = \left\{\begin{aligned} & \sqrt 3 /2,{\text{ }}n = 3k - 2{\text{ }}\left( {k = 1,2, \cdots } \right) \\[-2pt] & - \sqrt 3 /2,{\text{ }}n = 3k - 1{\text{ }}\left( {k = 1,2, \cdots } \right) \\[-2pt] & 0,{\text{ }}n = 3k{\text{ }}\left( {k = 1,2, \cdots } \right) \end{aligned} \right. $$ (18)

    (1) 随机激励

    考虑跨中简谐载荷$P\left( t \right) = {A_0}\cos \left( {\bar \omega t + \phi } \right)$, 取载荷幅值${A_0}$和相位$\phi $为相互独立的随机变量, 且分别服从均匀分布${U_{{A_0}}}\left( {{A_{0,1}},{A_{0,2}}} \right)$和${U_\phi }\left( {{\phi _1},{\phi _2}} \right)$(记为${\theta _1}$和${\theta _2}$).

    关于跨中位移响应$Y$的广义概率密度演化方程为

    $$ \frac{{\partial {p_{Y{\boldsymbol{\varTheta }}}}\left( {y,{\boldsymbol{\theta }},t} \right)}}{{\partial t}} + \dot Y\left( {{\boldsymbol{\theta }},t} \right)\frac{{\partial {p_{Y{\boldsymbol{\varTheta }}}}\left( {y,{\boldsymbol{\theta }},t} \right)}}{{\partial y}} = 0 $$ (19)

    根据式(10)和式(13)并利用Dirac函数的积分法则, 可得$Y$的概率密度函数

    $$\begin{split} & {p_Y}\left( {y,t} \right) = \int_{{{\varOmega }_{\boldsymbol{\varTheta }}}} {\delta \left[ {y - {H_Y}\left( {{\boldsymbol{\theta }},t} \right)} \right]} {p_{\boldsymbol{\varTheta }}}\left( {\boldsymbol{\theta }} \right){\text{d}}{\boldsymbol{\theta }} = \\ &\qquad \int_{{{\varOmega }_{{{{{\varTheta}} }_2}}}} {\left| {{J}} \right|} {p_{\boldsymbol{\varTheta}} }\left( {{\theta _1},{\theta _2}} \right)\left| {_{{\theta _1} = {H_Y}^{ - 1}}} \right.{\text{d}}{\theta _2}\end{split} $$ (20)

    式中

    $$ {\theta _1} = {H_Y}^{ - 1} = \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{y}{{a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2}}} $$ (21)
    $$ {{J}} = \frac{{\partial {\theta _1}}}{{\partial y}} = \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{1}{{a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2}}} $$ (22)

    其中

    $$ a\left( t \right) = \sum\limits_{n = 1}^\infty {\frac{{{\alpha _n}{\omega _n}\left( {{\omega _n}\cos {\omega _n}t - {\omega _n}\cos \bar \omega t} \right)}}{{{n^4}\left( {{{\bar \omega }^2} - {\omega _n}^2} \right)}}} \sin \frac{{n{\text{π}}}}{2} $$ (23)
    $$ b\left( t \right) = \sum\limits_{n = 1}^\infty {\frac{{{\alpha _n}{\omega _n}\left( {{\omega _n}\sin \bar \omega t - \bar \omega \sin {\omega _n}t} \right)}}{{{n^4}\left( {{{\bar \omega }^2} - {\omega _n}^2} \right)}}} \sin \frac{{n{\text{π}}}}{2} $$ (24)

    将式(21) ~ 式(24)代入式(20)可以得到

    $$ \begin{split} & {p_Y}\left( {y,t} \right) =\\ &\int_{{{\varOmega }_{{{\varTheta }_2}}}} {\left| {\frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{1}{{a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2}}}} \right|} \frac{1}{{{A_{0,2}} - {A_{0,1}}}}\frac{1}{{{\phi _2} - {\phi _1}}}{\text{d}}{\theta _2} =\\ & \frac{{{{\text{π}}^4}EI}}{{2{L^3}\left( {{A_{0,2}} - {A_{0,1}}} \right)\left( {{\phi _2} - {\phi _1}} \right)}}\int_{{{\varOmega }_{{{\varTheta }_2}}}} {\frac{1}{{\left| {a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2}} \right|}}} {\text{d}}{\theta _2} \end{split} $$ (25)

    由于

    $$ \begin{split} & \int_{{{\varOmega }_X}} {\frac{1}{{a\cos x + b\sin x}}} {\text{d}}x = \\ &\qquad \frac{2}{{\sqrt {{a^2} + {b^2}} }}{\text{arctanh}}\left( {\frac{{ - b + a\tan \dfrac{x}{2}}}{{\sqrt {{a^2} + {b^2}} }}} \right)\Biggr| {_{{{\varOmega }_X}}} \end{split} $$ (26)

    $$ c = \frac{{{{\text{π}}^4}EI}}{{2{L^3}\left( {{A_{0,2}} - {A_{0,1}}} \right)\left( {{\phi _2} - {\phi _1}} \right)}} $$ (27)

    则式(25)可以改写为

    $$\begin{split} & {p_Y}\left( {y,t} \right)= \\ &\quad \frac{{2c}}{{\sqrt {a{{\left( t \right)}^2} + b{{\left( t \right)}^2}} }}{\text{arctanh}}\left[ {\frac{{ - b\left( t \right) + a\left( t \right)\tan \dfrac{{{\theta _2}}}{2}}}{{\sqrt {a{{\left( t \right)}^2} + b{{\left( t \right)}^2}} }}} \right]\Biggr| {_{{{\varOmega }_{{{\varTheta }_2} + }}}} - \\ &\quad \frac{{2c}}{{\sqrt {a{{\left( t \right)}^2} + b{{\left( t \right)}^2}} }}{\text{arctanh}}\left[ {\frac{{ - b\left( t \right) + a\left( t \right)\tan \dfrac{{{\theta _2}}}{2}}}{{\sqrt {a{{\left( t \right)}^2} + b{{\left( t \right)}^2}} }}} \right]\Biggr| {_{{{\varOmega }_{{{\varTheta }_2} - }}}} \end{split} $$ (28)

    式中的积分限可表达为

    $$ \begin{split} & {{\varOmega }_{{{\varTheta }_2} + }} = a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2} > 0 \cap {\phi _1} \leqslant {\theta _2} \leqslant {\phi _2} \cap \\ &\qquad {A_{0,1}} \leqslant \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{y}{{a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2}}} \leqslant {A_{0,2}} \end{split} $$ (29)
    $$ \begin{split} & {{\varOmega }_{{{\varTheta }_2} - }} = a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2} < 0 \cap {\phi _1} \leqslant {\theta _2} \leqslant {\phi _2}\cap \\ &\qquad {A_{0,1}} \leqslant \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{y}{{a\left( t \right)\cos {\theta _2} + b\left( t \right)\sin {\theta _2}}} \leqslant {A_{0,2}} \end{split} $$ (30)

    式中, $ \cap $表示集合的交集, 其定义为$ A\cap B = \left\{x|x\in A 且x\in B\right\} $.

    式(28) ~ 式(30)给出了仅考虑随机激励时具有2个随机变量的一个典型随机系统的概率密度解析解. 不妨对式(23)和式(24)取一阶截断, 由此可得到梁位移响应的概率密度一阶近似解. 值得指出的是, 在概率密度演化理论中, 对物理方程求解所作的近似并不影响对广义概率密度演化方程的求解. 事实上, 只要在广义概率密度演化方程求解中取相同的物理近似, 不同截断阶数所得到的概率密度解答均能作为数值解的标定基准. 对梁位移响应采用一阶截断, 函数$a\left( t \right)$和$b\left( t \right)$简化为

    $$\qquad a\left( t \right) = \frac{{{\omega _1}\left( {{\omega _1}\cos {\omega _1}t - {\omega _1}\cos \bar \omega t} \right)}}{{{{\bar \omega }^2} - {\omega _1}^2}} $$ (31)
    $$\qquad b\left( t \right) = \frac{{{\omega _1}\left( {{\omega _1}\sin \bar \omega t - \bar \omega \sin {\omega _1}t} \right)}}{{{{\bar \omega }^2} - {\omega _1}^2}} $$ (32)

    若令随机变量的分布参数为$\left( {{A_{0,1}},{A_{0,2}}} \right) = \left( 4 \times {{10}^4}{\text{ N}},6 \times {{10}^4}{\text{ N}} \right)$和$\left( {{\phi _1},{\phi _2}} \right) = \left( {{\text{π}}/4,3{\text{π}}/4} \right)$, 而频率为确定性变量$\bar \omega = 0.6{\omega _1}$, 则$P\left( t \right)$有均值函数$\mu \left( t \right)$和自相关函数$R\left( {t,t + \tau } \right)$, 分别为

    $$ \mu \left( t \right) = E\left[ {{A_0}} \right]E\left[ {\cos \left( {\bar \omega t + \phi } \right)} \right] = - \frac{{\sqrt 2 \times {{10}^5}}}{{\text{π}}}\sin \bar \omega t $$ (33)
    $$ \begin{split} & R\left( {t,t + \tau } \right) = E\left[ {{A_0}^2} \right]E\left[ {\cos \left( {\bar \omega t + \phi } \right)\cos \left( {\bar \omega t + \bar \omega \tau + \phi } \right)} \right]= \\ &\qquad \frac{{38 \times {{10}^8}}}{3}\left[ { - \frac{2}{{\text{π}}}\cos \left( {2\bar \omega t + \bar \omega \tau } \right) + \cos \bar \omega \tau } \right]\end{split} $$ (34)

    同时, 对任意实数$\left( {{l_1},{l_2}, \cdots ,{l_n}} \right)$有[29]

    $$ \begin{split} & \sum\limits_{k = 1}^n {{l_k}P\left( {{t_k}} \right)} = \sum\limits_{k = 1}^n {{l_k}{A_0}\left( {\cos \bar \omega {t_k}\cos \phi - \sin \bar \omega {t_k}\sin \phi } \right)} = \\ &\qquad \left( {\sum\limits_{k = 1}^n {{l_k}\cos \bar \omega {t_k}} } \right){A_0}\cos \phi - \left( {\sum\limits_{k = 1}^n {{l_k}\sin \bar \omega {t_k}} } \right){A_0}\sin \phi \end{split} $$ (35)

    显然, ${A_0}\cos \phi $和${A_0}\sin \phi $均不服从正态分布. 因此, 由式(33) ~ 式(35)可知, $P\left( t \right)$为非平稳非高斯随机过程.

    对梁位移响应取一阶截断时, 式(29)和式(30)所表达的积分限可进一步简化为

    $$ {{\varOmega }_{{{\varTheta }_2} + }} = \left\{\begin{aligned} &\emptyset\;\; \left( {y > 0} \right) \\ & \left\{ {\max \left[ {{\phi _1},{d_1}\left( {y,{A_{0,2}},t} \right)} \right],\min \left[ {{\phi _2},{d_1}\left( {y,{A_{0,1}},t} \right)} \right]} \right\} \cup \\ &\quad \left\{ {\max \left[ {{\phi _1},{d_2}\left( {y,{A_{0,1}},t} \right)} \right],\min \left[ {{\phi _2},{d_2}\left( {y,{A_{0,2}},t} \right)} \right]} \right\}\\ &\quad \left( {y < 0,{\text{ }}b\left( t \right) > 0} \right) \\ & \left\{ {\max \left[ {{\phi _1},{d_1}\left( {y, - {A_{0,1}},t} \right)} \right],\min \left[ {{\phi _2},{d_1}\left( {y, - {A_{0,2}},t} \right)} \right]} \right\} \cup \\ &\quad \left\{ {\max \left[ {{\phi _1},{d_2}\left( {y, - {A_{0,2}},t} \right)} \right],\min \left[ {{\phi _2},{d_2}\left( {y, - {A_{0,1}},t} \right)} \right]} \right\}\\ &\quad \left( {y < 0,{\text{ }}b\left( t \right) < 0} \right) \end{aligned} \right. $$ (36)
    $$ {{\varOmega }_{{{\varTheta }_2} - }} = \left\{\begin{aligned} & \left\{ {\max \left[ {{\phi _1},{d_1}\left( {y,{A_{0,1}},t} \right)} \right],\min \left[ {{\phi _2},{d_1}\left( {y,{A_{0,2}},t} \right)} \right]} \right\} \cup \\ &\quad \left\{ {\max \left[ {{\phi _1},{d_2}\left( {y,{A_{0,2}},t} \right)} \right],\min \left[ {{\phi _2},{d_2}\left( {y,{A_{0,1}},t} \right)} \right]} \right\}\\ &\quad \left( {y > 0,{\text{ }}b\left( t \right) > 0} \right) \\ & \left\{ {\max \left[ {{\phi _1},{d_1}\left( {y, - {A_{0,2}},t} \right)} \right],\min \left[ {{\phi _2},{d_1}\left( {y, - {A_{0,1}},t} \right)} \right]} \right\} \cup \\ &\quad \left\{ {\max \left[ {{\phi _1},{d_2}\left( {y, - {A_{0,1}},t} \right)} \right],\min \left[ {{\phi _2},{d_2}\left( {y, - {A_{0,2}},t} \right)} \right]} \right\}\\ &\quad \left( {y > 0,{\text{ }}b(t) < 0} \right) \\ & \emptyset\;\; \left( {y < 0} \right) \end{aligned} \right. $$ (37)

    式中, $ \cup $表示集合的并集, 其定义为$ A\cup B = \Bigr\{x|x\in A\text{ }或x\in B\Bigr\} $, 并有

    $$ \begin{split} & {d_1}\left( {y,A,t} \right) = \\ &\qquad \arcsin \frac{{{\text{π} ^4}EIy}}{{2{L^3}A\sqrt {a{{\left( t \right)}^2} + b{{\left( t \right)}^2}} }} - \arctan \frac{{a\left( t \right)}}{{b\left( t \right)}} \end{split} $$ (38)
    $$\begin{split} & {d_2}\left( {y,A,t} \right) = \\ &\qquad {\text{π}} - \arcsin \frac{{{\text{π} ^4}EIy}}{{2{L^3}A\sqrt {a{{\left( t \right)}^2} + b{{\left( t \right)}^2}} }} - \arctan \frac{{a\left( t \right)}}{{b\left( t \right)}} \end{split} $$ (39)

    据此, 可给出图1所示梁跨中位移响应$Y$的概率密度解答, 如图4所示. 图中, 同时给出了联合概率密度演化曲面及典型时刻的概率分布密度. 可见, 即使对这样一个简单问题, 真实的概率密度演化过程也是十分复杂的: 对于具有均匀分布的源随机变量, 系统响应的瞬时概率密度既非均匀、也非正态, 而是具有不断变化的、有时双峰、有时单峰的复杂演化过程. 联系到文献[7]中所述单变量系统的响应概率密度特征, 不难推测随机系统响应的复杂性.

    图  4  两变量系统的概率密度解答(随机激励)
    Figure  4.  Analytical solution of probability density for the two-variable system (random excitations)

    (2)随机激励和随机结构参数

    考虑跨中阶跃载荷$P\left( t \right) = {A_0}\delta \left( {t - {t_0}} \right)$, 取载荷幅值${A_0}$和结构单位质量$\bar m$为相互独立的随机变量. 为简单计, 取结构响应的一阶截断, 则一阶频率${\omega _1}$为与$\bar m$完全相关的随机变量, 令${A_0}$和${\omega _1}$服从均匀分布${U_{{A_0}}}\left( {{A_{0,1}},{A_{0,2}}} \right)$和${U_{{\omega _1}}}\left( {{\omega _{1,1}},{\omega _{1,2}}} \right)$(记为${\theta _1}$和${\theta _2}$).

    关于跨中位移响应$Y$的广义概率密度演化方程同式(19). 根据式(10)和式(12)并利用Dirac函数的积分法则, 可得$Y$的概率密度函数

    $$ \begin{split} & {p_Y}\left( {y,t} \right) = \int_{{{\varOmega }_{\boldsymbol{\varTheta }}}} {\delta \left[ {y - {H_Y}\left( {{\boldsymbol{\theta }},t} \right)} \right]} {p_{\boldsymbol{\varTheta }}}\left( {\boldsymbol{\theta }} \right){\text{d}}{\boldsymbol{\theta }} =\\ &\qquad \int_{{{\varOmega }_{{{\varTheta }_2}}}} {\left| {{J}} \right|} {p_{\boldsymbol{\varTheta}} }\left( {{\theta _1},{\theta _2}} \right)\left| {_{{\theta _1} = {H_Y}^{ - 1}}} \right.{\text{d}}{\theta _2} \end{split} $$ (40)

    式中

    $$ \qquad\qquad {\theta _1} = {H_Y}^{ - 1} = \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{y}{{1 - \cos {\theta _2}t}} $$ (41)
    $$\qquad\qquad {{J}} = \frac{{\partial {\theta _1}}}{{\partial y}} = \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{1}{{1 - \cos {\theta _2}t}} $$ (42)

    将式(41) ~ 式(42)代入式(40)可以得到

    $$ \begin{split} & {p_Y}\left( {y,t} \right) = \\ &\qquad \int_{{{\varOmega }_{{{\varTheta }_2}}}} {\left| {\frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{1}{{1 - \cos {\theta _2}t}}} \right|} \frac{1}{{{A_{0,2}} - {A_{0,1}}}}\frac{1}{{{\omega _{1,2}} - {\omega _1}}}{\text{d}}{\theta _2} =\\ &\qquad \frac{{{{\text{π}}^4}EI}}{{2{L^3}\left( {{A_{0,2}} - {A_{0,1}}} \right)\left( {{\omega _{1,2}} - {\omega _{1,1}}} \right)}}\int_{{{\varOmega }_{{{\varTheta }_2}}}} {\frac{1}{{1 - \cos {\theta _2}t}}} {\text{d}}{\theta _2} \end{split} $$ (43)

    由于

    $$ \int_{{\varOmega _X}} {\frac{1}{{1 - \cos xt}}} {\mathrm{d}}x = - \frac{{\cot xt/2}}{t}\left| {_{{\varOmega _X}}} \right. $$ (44)

    $$ c = - \frac{{{{\text{π}}^4}EI}}{{2{L^3}\left( {{A_{0,2}} - {A_{0,1}}} \right)\left( {{\omega _{1,2}} - {\omega _{1,1}}} \right)}} $$ (45)

    则式(43)可以改写为

    $$ {p_Y}\left( {y,t} \right) = c\frac{{\cot ({{\theta _2}t/2}) }}{t}\left| {_{{\varOmega _{{\varTheta _2}}}}} \right. $$ (46)

    式中的积分限可表达为

    $$ {\varOmega _{{\varTheta _2}}} = {\omega _{1,2}} \leqslant {\theta _2} \leqslant {\omega _{1,1}} \cap {A_{0,1}} \leqslant \frac{{{{\text{π}}^4}EI}}{{2{L^3}}}\frac{y}{{1 - \cos {\theta _2}t}} \leqslant {A_{0,2}} $$ (47)

    式(46)和式(47)给出了同时考虑随机激励和随机结构参数时的典型随机系统的概率密度解析解. 为给出具体分析结果, 取随机变量的分布参数为$\left( {{A_{0,1}},{A_{0,2}}} \right) = \left( {4 \times {{10}^4}{\text{ N}},6 \times {{10}^4}{\text{ N}}} \right)$, $\left( {{\omega _{1,1}},{\omega _{1,2}}} \right) = \left( {9.75{\text{ rad/s, }} 19.5{\text{ rad/s}}} \right)$, 则式(47)所表达的积分限可进一步简化为

    $$ \begin{split} & {{\varOmega }_{{{\varTheta }_2}}} = \left\{ {\max \left[ {{\omega _{1,1}},{d_1}\left( {y,{A_{0,2}},t} \right)} \right],\min \left[ {{\omega _{1,2}},{d_1}\left( {y,{A_{0,1}},t} \right)} \right]} \right\} \cup \\ &\quad \left\{ {\max \left[ {{\omega _{1,1}},{d_2}\left( {y,{A_{0,1}},t} \right)} \right],\min \left[ {{\omega _{1,2}},{d_2}\left( {y,{A_{0,2}},t} \right)} \right]} \right\} \end{split} $$ (48)

    式中

    $$ {d_1}\left( {y,A,t} \right) = \frac{{\arccos \left( {1 - \dfrac{{{{\text{π}}^4}EIy}}{{2{L^3}A}}} \right)}}{t} $$ (49)
    $$ {d_2}\left( {y,A,t} \right) = \frac{{{{2\text{π} }} - \arccos \left( {1 - \dfrac{{{{\text{π}}^4}EIy}}{{2{L^3}A}}} \right)}}{t} $$ (50)

    由此得到图5所示梁跨中位移响应$Y$的概率密度解答. 可见, 不同时刻的概率密度函数变化复杂, 难以用经典的概率分布形式加以描述.

    图  5  两变量系统的概率密度解答(随机激励和随机结构参数)
    Figure  5.  Analytical solution of probability density for the two-variable system (random excitations and random structural parameters)

    考虑三分点作用阶跃载荷${P_1}\left( t \right) = {A_1}\delta \left( {t - {t_0}} \right)$和简谐载荷${P_2}\left( t \right) = {A_2}\cos \left( {\bar \omega t + \phi } \right)$. 取载荷幅值${A_1}$, ${A_2}$和相位$\phi $, 以及结构一阶频率${\omega _1}$为相互独立的随机变量, 且分别服从均匀分布${U_{{A_1}}}\left( {{A_{1,1}},{A_{1,2}}} \right)$, $ {U_{{A_2}}}\left( {{A_{2,1}},{A_{2,2}}} \right) $, ${U_\phi }\left( {{\phi _1},{\phi _1}} \right)$和${U_{{\omega _1}}}\left( {{\omega _{1,1}},{\omega _{1,2}}} \right)$(记为${\theta _1}$, ${\theta _2}$, ${\theta _3}$和${\theta _4}$).

    关于跨中位移响应$Y$的广义概率密度演化方程仍同式(19). 根据式(10)和式(16)并利用Dirac函数的积分法则, 可得$Y$的概率密度函数

    $$ \begin{split} & {p_Y}\left( {y,t} \right) = \int_{{{\varOmega }_{\boldsymbol{\varTheta }}}} {\delta \left[ {y - {H_Y}\left( {{\boldsymbol{\theta }},t} \right)} \right]} {p_{\boldsymbol{\varTheta }}}\left( {\boldsymbol{\theta }} \right){\text{d}}{\boldsymbol{\theta }} = \\ &\qquad \int_{{{\varOmega }_{{{\varTheta }_2}}}} {\int_{{{\varOmega }_{{{\varTheta }_3}}}} {\int_{{{\varOmega }_{{{\varTheta }_4}}}} {\left| {{J}} \right|{p_\varTheta }\left( {{\theta _1},{\theta _2},{\theta _3},{\theta _4}} \right)\left| {_{{\theta _1} = {H_Y}^{ - 1}}} \right.} } } {\text{d}}{\theta _2}{\text{d}}{\theta _3}{\text{d}}{\theta _4} \end{split} $$ (51)

    式中

    $$ \begin{split} & {\theta _1} = {H_Y}^{ - 1} = \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}}}\frac{1}{{1 - \cos {\theta _4}t}} \cdot \\ &\qquad \left\{ {y - \frac{{\sqrt 3 {\theta _2}{L^3}}}{{{{\text{π}}^4}EI}}\left[ {a\left( {{\theta _4},t} \right)\cos {\theta _3} + b\left( {{\theta _4},t} \right)\sin {\theta _3}} \right]} \right\} \end{split} $$ (52)
    $$ {{J}} = \frac{{\partial {\theta _1}}}{{\partial y}} = \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}}}\frac{1}{{1 - \cos {\theta _4}t}} $$ (53)

    其中

    $$ a\left( {{\theta _4},t} \right) = \frac{{{\theta _4}\left( {{\theta _4}\cos {\theta _4}t - {\theta _4}\cos \bar \omega t} \right)}}{{{{\bar \omega }^2} - {\theta _4}^2}} $$ (54)
    $$ b\left( {{\theta _4},t} \right) = \frac{{{\theta _4}\left( {{\theta _4}\sin \bar \omega t - \bar \omega \sin {\theta _4}t} \right)}}{{{{\bar \omega }^2} - {\theta _4}^2}} $$ (55)

    将式(52) ~ 式(55)代入式(51)可以得到

    $$ \begin{split} & {p_Y}\left( {y,t} \right) = \int_{{{\varOmega }_{{{\varTheta }_2}}}} {\int_{{{\varOmega }_{{{\varTheta }_3}}}} {\int_{{{\varOmega }_{{{\varTheta }_4}}}} {\left| {\frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}}}\frac{1}{{1 - \cos {\theta _4}t}}} \right|} } } \frac{1}{{{A_{1,2}} - {A_{1,1}}}}\cdot \\ &\quad \frac{1}{{{A_{2,2}} - {A_{2,1}}}}\frac{1}{{{\phi _2} - {\phi _1}}}\frac{1}{{{\omega _{1,2}} - {\omega _{1,1}}}}{\text{d}}{\theta _2}{\text{d}}{\theta _3}{\text{d}}{\theta _4} = \\ &\quad \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}\left( {{A_{1,2}} - {A_{1,1}}} \right)\left( {{A_{2,2}} - {A_{2,1}}} \right)\left( {{\phi _2} - {\phi _1}} \right)\left( {{\omega _{1,2}} - {\omega _{1,1}}} \right)}}\cdot \\ &\quad \int_{{{\varOmega }_{{{\varTheta }_2}}}} {\int_{{{\varOmega }_{{{\varTheta }_3}}}} {\int_{{{\varOmega }_{{{\varTheta }_4}}}} {\frac{1}{{1 - \cos {\theta _4}t}}} } } {\text{d}}{\theta _2}{\text{d}}{\theta _3}{\text{d}}{\theta _4} \end{split} $$ (56)

    由于

    $$ \int_{{{\varOmega }_X}} {\int_{{{\varOmega }_Y}} {\int_{{{\varOmega }_Z}} {\frac{1}{{1 - \cos xt}}} } } {\text{d}}x{\text{d}}y{\text{d}}z = - \frac{{\cot xt/2}}{t}yz\left| {_{{{\varOmega }_X} \times {{\varOmega }_Y} \times {{\varOmega }_Z}}} \right. $$ (57)

    $$ c = - \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}\left( {{A_{1,2}} - {A_{1,1}}} \right)\left( {{A_{2,2}} - {A_{2,1}}} \right)\left( {{\phi _2} - {\phi _1}} \right)\left( {{\omega _{1,2}} - {\omega _{1,1}}} \right)}} $$ (58)

    则式(55)可改写为

    $$ {P_Y}\left( {y,t} \right) = c\frac{{\cot \left( {{\theta _4}t/2} \right)}}{t}{\theta _2}{\theta _3}\left| {_{{{\varOmega }_{{\varTheta _2}}} \times {{\varOmega }_{{\varTheta _3}}} \times {{\varOmega }_{{\varTheta _4}}}}} \right. $$ (59)

    式中的积分限可表达为

    $$ \begin{split} & {{\varOmega }_{{{\varTheta }_2}}} \times {{\varOmega }_{{{\varTheta }_3}}} \times {{\varOmega }_{{{\varTheta }_4}}} = {A_{2,1}} \leqslant {\theta _2} \leqslant {A_{2,2}} \cap {\phi _1} \leqslant {\theta _3} \leqslant {\phi _2} \cap \\ &\qquad {\omega _{1,1}} \leqslant {\theta _4} \leqslant {\omega _{1,2}} \cap {A_{1,1}} \leqslant \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}}}\frac{1}{{1 - \cos {\theta _4}t}} \cdot \\ &\qquad \left\{ {y - \frac{{\sqrt 3 {\theta _2}{L^3}}}{{{{\text{π}}^4}EI}}\left[ {a\left( {{\theta _4},t} \right)\cos {\theta _3} + b\left( {{\theta _4},t} \right)\sin {\theta _3}} \right]} \right\} \leqslant {A_{1,2}}\end{split} $$ (60)

    式(59)和式(60)给出了同时考虑激励随机性和结构参数随机性的、具有4个随机变量的典型随机系统的响应概率密度解答. 为给出具体分析结果, 取随机变量的分布参数分别为$\left( {{A_{1,1}},{A_{1,2}}} \right) = \left( 4 \times {{10}^4}{\text{ N}}, 6 \times {{10}^4}{\text{ N}} \right)$, $\left( {{A_{2,1}},{A_{2,2}}} \right) = \left( {3 \times {{10}^4}{\text{ N}},7 \times {{10}^4}{\text{ N}}} \right)$, $\left( {\phi _1}, {\phi _2} \right) = \left( {{\text{π}}/2,3{\text{π}}/2} \right)$和$\left( {{\omega _{1,1}},{\omega _{1,2}}} \right) = \left( {9.75{\text{ rad/s, }} 19.5{\text{ rad/s}}} \right)$, 同时, 取载荷频率为确定性变量$ {\bar \omega _1} = 29.24{\text{ rad/s}} $. 由此, 可得图6所示的概率密度解答. 可见, 在不同时刻, 概率密度分布可由单峰分布逐渐演化为双峰分布.

    图  6  四变量系统的概率密度解析解
    Figure  6.  Analytical solution of probability density for the four-variable system

    类似地, 取三分点载荷的幅值${A_1}$和${A_2}$、相位$\phi $、频率$\bar \omega $以及结构一阶频率${\omega _1}$为相互独立且服从均匀分布的随机变量(记为${\theta _1}$, ${\theta _2}$, ${\theta _3}$, ${\theta _4}$和${\theta _5}$), 可给出五变量系统的解答为

    $$ {p_Y}\left( {y,t} \right) = c\frac{{\cot \left( {{\theta _4}t/2} \right)}}{t}{\theta _2}{\theta _3}{\theta _5}\left| {_{{{\varOmega }_{{\varTheta _2}}} \times {{\varOmega }_{{\varTheta _3}}} \times {{\varOmega }_{{\varTheta _4}}} \times {{\varOmega }_{{\varTheta _5}}}}} \right. $$ (61)

    式中

    $$ \begin{split} & c = - \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}}}\frac{1}{{{A_{1,2}} - {A_{1,1}}}}\frac{1}{{{A_{2,2}} - {A_{2,1}}}} \cdot \\ &\qquad \frac{1}{{{\phi _2} - {\phi _1}}}\frac{1}{{{\omega _{1,2}} - {\omega _{1,1}}}}\frac{1}{{{{\bar \omega }_2} - {{\bar \omega }_1}}}\end{split} $$ (62)

    相应的积分限为

    $$ \begin{split} & {{\varOmega }_{{{\varTheta }_2}}} \times {{\varOmega }_{{{\varTheta }_3}}} \times {{\varOmega }_{{{\varTheta }_4}}} \times {{\varOmega }_{{{\varTheta }_5}}} = {A_{2,1}} \leqslant {\theta _2} \leqslant {A_{2,2}} \cap {\phi _1} \leqslant {\theta _3} \leqslant {\phi _2} \cap \\ & {\omega _{1,1}} \leqslant {\theta _4} \leqslant {\omega _{1,2}} \cap {{\bar \omega }_1} \leqslant {\theta _5} \leqslant {{\bar \omega }_2} \cap {A_{1,1}} \leqslant \frac{{{{\text{π}}^4}EI}}{{\sqrt 3 {L^3}}}\frac{1}{{1 - \cos {\theta _4}t}} \cdot \\ & \left\{ {y - \frac{{\sqrt 3 {\theta _2}{L^3}}}{{{{\text{π}}^4}EI}}\left[ {a\left( {{\theta _4},{\theta _5},t} \right)\cos {\theta _3} + b\left( {{\theta _4},{\theta _5},t} \right)\sin {\theta _3}} \right]} \right\} \leqslant {A_{1,2}} \end{split} $$ (63)

    其中

    $$ a\left( {{\theta _4},{\theta _5},t} \right) = \frac{{{\theta _4}\left( {{\theta _4}\cos {\theta _4}t - {\theta _4}\cos {\theta _5}t} \right)}}{{{\theta _5}^2 - {\theta _4}^2}} $$ (64)
    $$ b\left( {{\theta _4},{\theta _5},t} \right) = \frac{{{\theta _4}\left( {{\theta _4}\sin {\theta _5}t - {\theta _5}\sin {\theta _4}t} \right)}}{{{\theta _5}^2 - {\theta _4}^2}} $$ (65)

    由此得到的典型分析结果如图7所示. 在这一分析中, 与图7相对应的分布参数列于表1.

    图  7  五变量系统的概率密度解析解
    Figure  7.  Analytical solution of probability density for the five-variable system
    表  1  五变量系统中随机变量的分布参数
    Table  1.  Distribution parameters of random variables in the five-variable system
    P1/N P2/N φ ω1/s−1 $\bar \omega $/s−1
    upper boundary 4 × 104 3 × 104 π / 2 9.75 9.75
    lower boundary 6 × 104 7 × 104 3π / 2 19.5 43.86
    下载: 导出CSV 
    | 显示表格

    对四变量和五变量系统的概率密度解答进行积分[30], 可得到均值函数和标准差函数, 如图8所示. 可以看出, 虽然四变量系统中的确定性激励频率$\bar \omega $恰好为五变量系统中随机频率的均值, 但二者响应的均值曲线与标准差曲线均不相同、且经过一定时间后, 四变量系统标准差逐渐趋于平稳, 而五变量系统的标准差曲线并不趋于平稳, 且远大于四变量系统.

    图  8  四变量和五变量系统的矩函数对比
    Figure  8.  Comparison of moment functions between four-variable and five-variable systems

    本文以Euler-Bernoulli梁为例, 研究了含多维随机变量的广义概率密度演化方程解析解. 给出了含2维、4维和5维随机变量系统(含随机激励与随机结构参数)的解析解及具体分析结果. 这些分析结果, 为在多维概率空间内研究广义概率密度演化方程的各类数值求解方法奠定了基础.

    值得特别提出的是, 上述解析解的最终形式均为在多维函数空间内寻找与分布参数相对应的积分限, 这样的求解方式给出了在更高维系统中求取广义概率密度演化方程解析解的可能性.

  • 图  1   Euler-Bernoulli简支梁跨中集中载荷模型

    Figure  1.   Model of Euler-Bernoulli simply supported beam subjected to midspan concentrated load

    图  2   不同的载荷形式

    Figure  2.   Different forms of load

    图  3   Euler-Bernoulli简支梁3分点集中载荷模型

    Figure  3.   Model of Euler-Bernoulli simply supported beam subjected to third-point concentrated load

    图  4   两变量系统的概率密度解答(随机激励)

    Figure  4.   Analytical solution of probability density for the two-variable system (random excitations)

    图  5   两变量系统的概率密度解答(随机激励和随机结构参数)

    Figure  5.   Analytical solution of probability density for the two-variable system (random excitations and random structural parameters)

    图  6   四变量系统的概率密度解析解

    Figure  6.   Analytical solution of probability density for the four-variable system

    图  7   五变量系统的概率密度解析解

    Figure  7.   Analytical solution of probability density for the five-variable system

    图  8   四变量和五变量系统的矩函数对比

    Figure  8.   Comparison of moment functions between four-variable and five-variable systems

    表  1   五变量系统中随机变量的分布参数

    Table  1   Distribution parameters of random variables in the five-variable system

    P1/N P2/N φ ω1/s−1 $\bar \omega $/s−1
    upper boundary 4 × 104 3 × 104 π / 2 9.75 9.75
    lower boundary 6 × 104 7 × 104 3π / 2 19.5 43.86
    下载: 导出CSV
  • [1] 李杰. 工程结构可靠性分析原理. 北京: 科学出版社, 2021 (Li Jie. Fundamental of Structural Reliability Analysis. Beijing: Science Press, 2021 (in Chinese)

    Li Jie. Fundamental of Structural Reliability Analysis. Beijing: Science Press, 2021 (in Chinese)

    [2]

    Lin YK. Probabilistic Theory of Structural Dynamics. New York: McGraw-Hill Book Company, 1967

    [3] 李杰. 随机结构系统——分析与建模. 北京: 科学出版社, 1996 (Li Jie. Stochastic Structural System—Analysis and Modeling. Beijing: Science Press, 1996 (in Chinese)

    Li Jie. Stochastic Structural System—Analysis and Modeling. Beijing: Science Press, 1996 (in Chinese)

    [4]

    Spanos PD, Ghanem R. Stochastic finite element expansion for random media. Journal of Engineering Mechanics, 1989, 115(5): 1035-1053 doi: 10.1061/(ASCE)0733-9399(1989)115:5(1035)

    [5]

    Li J, Chen JB. The probability density evolution method for dynamic response analysis of non-linear stochastic structures. International Journal for Numerical Methods in Engineering, 2006, 65(6): 882-903 doi: 10.1002/nme.1479

    [6]

    Li J, Chen JB. The principle of preservation of probability and the generalized density evolution equation. Structural Safety, 2008, 30(1): 65-77 doi: 10.1016/j.strusafe.2006.08.001

    [7]

    Li J, Chen JB. Stochastic Dynamics of Structures. Hoboken: John Wiley & Sons, 2009

    [8] 李杰, 陈建兵. 随机动力系统中的概率密度演化方程及其研究进展. 力学进展, 2010, 40(2): 170-188 (Li Jie, Chen Jianbing. Advances in the research on probability density evolution equations of stochastic dynamical systems. Advances in Mechanics, 2010, 40(2): 170-188 (in Chinese)

    Li Jie, Chen Jianbing. Advances in the research on probability density evolution equations of stochastic dynamical systems. Advances in Mechanics, 2010, 40(2): 170-188 (in Chinese)

    [9]

    Kozin F. On the probability densities of the output of some random systems. Journal of Applied Mechanics, 1961, 28(2): 161-164 doi: 10.1115/1.3641646

    [10]

    Fokker AD. Die mittlere energie rotierender elektrischer dipole im strahlungsfeld. Annalen der Physik, 1914, 348(5): 810-820 doi: 10.1002/andp.19143480507

    [11]

    Planck M. Ueber einen satz der statistichen dynamik und eine erweiterung in der quantumtheorie//Preussischen Akad Wiss, 1917: 324-341

    [12]

    Kolmogoroff A. Über die analytischen methoden in der wahrscheinlichkeitsrechnung. Mathematische Annalen, 1931, 104(1): 415-458 doi: 10.1007/BF01457949

    [13]

    Dostupov B, Pugachev V. The equation to define a probability distribution of the integral of a system of ordinary differential equations with random parameters. Avtomat i Telemekh, 1957, 18(7): 620-630

    [14]

    Lyu MZ, Feng DC, Chen JB, et al. A decoupled approach for determination of the joint probability density function of a high-dimensional nonlinear stochastic dynamical system via the probability density evolution method. Computer Methods in Applied Mechanics and Engineering, 2024, 418: 116443 doi: 10.1016/j.cma.2023.116443

    [15]

    Li J, Chen JB. Probability density evolution method for dynamic response analysis of structures with uncertain parameters. Computational Mechanics, 2004, 34(5): 400-409 doi: 10.1007/s00466-004-0583-8

    [16] 蒋仲铭, 李杰. 三类随机系统广义概率密度演化方程的解析解. 力学学报, 2016, 48(2): 413-421 (Jiang Zhongming, Li Jie. Analytical solutions of the generalized probability density evolution equation of three classes stochastic systems. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(2): 413-421 (in Chinese)

    Jiang Zhongming, Li Jie. Analytical solutions of the generalized probability density evolution equation of three classes stochastic systems. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(2): 413-421 (in Chinese)

    [17]

    Pourtakdoust SH, Khodabakhsh AH. A deep learning approach for the solution of probability density evolution of stochastic systems. Structural Safety, 2022, 99: 16

    [18] 华罗庚, 王元. 数论在近似分析中的应用. 北京: 科学出版社, 1978 (Hua Luogeng, Wang Yuan. Applications of Number Theory in Approximation Analysis. Beijing: Science Press, 1978 (in Chinese)

    Hua Luogeng, Wang Yuan. Applications of Number Theory in Approximation Analysis. Beijing: Science Press, 1978 (in Chinese)

    [19] 陈建兵, 李杰. 结构随机响应概率密度演化分析的数论选点法. 力学学报, 2006, 38(1): 134-140 (Chen Jianbing, Li Jie. Strategy of selecting points via number theoretical method in probability density evolution analysis of stochastic response of structures. Chinese Journal of Theoretical and Applied Mechanics, 2006, 38(1): 134-140 (in Chinese)

    Chen Jianbing, Li Jie. Strategy of selecting points via number theoretical method in probability density evolution analysis of stochastic response of structures. Chinese Journal of Theoretical and Applied Mechanics, 2006, 38(1): 134-140 (in Chinese)

    [20]

    Dick J, Pillichshammer F. Digital Nets and Sequences: Discrepancy Thoory and Quasi-Monte Carlo Integeration. Cambridge: Cambridge University Press, 2010

    [21]

    Chen JB, Yang JY, Li J. A GF-discrepancy for point selection in stochastic seismic response analysis of structures with uncertain parameters. Structural Safety, 2016, 59: 20-31 doi: 10.1016/j.strusafe.2015.11.001

    [22] 石晟, 杜东升, 王曙光等. 概率密度演化方程TVD格式的自适应时间步长技术及其初值条件改进. 力学学报, 2019, 51(4): 1223-1234 (Shi Sheng, Du Dongsheng, Wang Shuguang, et al. Non-uniform time step TVD scheme for probability density evolution function with improvement of initial condition. Chinese Journal of Theoretical and Applied Mechanics, 2019, 51(4): 1223-1234 (in Chinese)

    Shi Sheng, Du Dongsheng, Wang Shuguang, et al. Non-uniform time step TVD scheme for probability density evolution function with improvement of initial condition. Chinese Journal of Theoretical and Applied Mechanics, 2019, 51(4): 1223-1234 (in Chinese)

    [23]

    Papadopoulos V, Kalogeris I. A Galerkin-based formulation of the probability density evolution method for general stochastic finite element systems. Computational Mechanics, 2016, 57: 701-716 doi: 10.1007/s00466-015-1256-9

    [24]

    Wang D, Sun WL, Li J. An RKPM-based formulation of the generalized probability density evolution equation for stochastic dynamic systems. Probabilistic Engineering Mechanics, 2021, 66: 103152 doi: 10.1016/j.probengmech.2021.103152

    [25]

    Tao WF, Li J. An ensemble evolution numerical method for solving generalized density evolution equation. Probabilistic Engineering Mechanics, 2017, 48: 1-11 doi: 10.1016/j.probengmech.2017.03.001

    [26]

    Shinozuka M. Monte Carlo solution of structural dynamics. Computers & Structures, 1972, 2(5-6): 855-874

    [27]

    Li J, Wang D. Comparison of PDEM and MCS: Accuracy and efficiency. Probabilistic Engineering Mechanics, 2023, 71: 103382 doi: 10.1016/j.probengmech.2022.103382

    [28]

    Clough RW, Penzien J. Dynamics of Structures. New York: McGraw-Hill, 1975

    [29]

    Tong YL. Fundamental properties and sampling distributions of the multivariate normal distribution//The Multivariate Normal Distribution. New York: Springer, 1990

    [30]

    Ross SM. Introduction to Probability Models. New York: Academic Press, 2014

图(8)  /  表(1)
计量
  • 文章访问数:  152
  • HTML全文浏览量:  32
  • PDF下载量:  57
  • 被引次数: 0
出版历程
  • 收稿日期:  2024-01-01
  • 录用日期:  2024-03-17
  • 网络出版日期:  2024-03-17
  • 发布日期:  2024-03-18
  • 刊出日期:  2024-09-17

目录

/

返回文章
返回