CN104504189A - 随机激励下大规模结构设计方法 - Google Patents

Abstract

本发明公开了一种随机激励下大规模结构设计方法，用于解决现有的结构设计方法实用性差的技术问题。技术方案是采用虚拟激励法结合模态加速度法计算随机激励下的位移响应均方根，然后以结构指定位置的位移响应均方根最小为目标，以结构质量为约束对结构进行设计。相比背景技术的设计方法，本发明方法明显提高了随机位移响应均方根的分析精度，最终设计得到清晰有效的结构构型，从而满足工程实际中大规模结构的设计需求。本发明方法实施例中结构的自由度数量为14762，是背景技术中结构自由度数量3782的3.9倍，设计得到的结构构型清晰有效，易于在实际工程中使用。

Description

随机激励下大规模结构设计方法

技术领域

本发明涉及一种大规模结构设计方法，特别是涉及一种随机激励下大规模结构设计方法。

背景技术

工程实际中的结构常常会承受各种随机振动激励，例如自然界中存在的风激励，航空航天飞行器在服役时受到的气动激励。结构在随机激励作用下破坏的事件时有发生，因此在结构构型设计时，考虑结构在随机激励下的性能表现非常重要。

文献“Zhang Q,Zhang WH,Zhu JH,Gao T.Layout optimization of multi-componentstructures under static loads and random excitations.Engineering Structures.2012.43:120-128”公开了一种随机激励下以结构指定位置的位移响应均方根最小化为目标，结构体积为约束的结构构型优化设计方法，该方法采用传统的Complete QuadraticCombination(CQC)方法计算随机激励下的位移响应均方根。CQC方法计算随机响应时计算量巨大，并且在实际应用中存在由截断模态引起的分析误差，更值得注意的是分析误差会随着结构规模的增加而变大，因此文献中的方法只适用于设计结构规模较小的问题。文献所设计的结构中，最大的自由度数量为3782，规模很小，可以得到清晰的结构构型。但是当使用文献中的方法设计大规模结构时，设计得到的结构构型边界不清晰，构型复杂、结构上空洞较多，难以应用于工程实际(见图3)。

发明内容

为了克服现有的结构设计方法实用性差的不足，本发明提供一种随机激励下大规模结构设计方法。该方法采用虚拟激励法结合模态加速度法计算随机激励下的位移响应均方根，然后以结构指定位置的位移响应均方根最小为目标，以结构质量为约束对结构进行设计。相比背景技术的设计方法，本发明方法可以明显提高随机位移响应均方根的分析精度，最终能够设计得到清晰有效的结构构型，从而能满足工程实际中大规模结构的设计需求。

本发明解决其技术问题所采用的技术方案是：一种随机激励下大规模结构设计方法，其特点是采用以下步骤：

(a)建立设计空间有限元模型，设置拓扑设计变量η_h初始值，h是正整数表示单元编号，1≤h≤N_h，N_h表示结构单元总数量。给定材料密度ρ和杨氏模量E。给定质量约束上限

(b)设置激励载荷，给出随机激励f(t)的功率谱密度矩阵S_f(ω)，f(t)为p维列向量，p为载荷中力的个数，t表示时间，S_f(ω)为p维方阵，其下标f表示其为激励f(t)的功率谱矩阵。ω为激励角频率，载荷的激励频段为[ω，]。ω表示激励角频率的下限，表示激励角频率的上限。根据矩阵LDLT分解，存在下式成立

$S_f\left(\mathrm{\ω}\right)=\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}{\left({\mathrm{\γ}}_q\right)}^{*}{\left({\mathrm{\γ}}_q\right)}^T$

其中Q为矩阵S_f(ω)的秩，γ_q为p维列向量表示第q个虚拟简谐激励，1≤q≤Q，上标T表示向量或矩阵的转置。

(c)根据当前设计变量值，采用以下公式分别计算每一个有限元单元的材料密度ρ_h和杨氏模量E_h

ρ_h＝η_hρ

$E_h=\frac{{15\mathrm{\η}}_h^5+{\mathrm{\η}}_h}{16}E$

更新结构有限元模型中的相应材料属性并进行结构有限元分析。

(d)从有限元分析结果中提取每个单元的刚度矩阵k_h和质量矩阵m_h，结构的前l阶模态频率值ω_i，1≤i≤l，模态振型为n行l列矩阵，n为结构总自由度数目。设置结构前l阶模态的阻尼比ζ_i。采用虚拟激励法结合模态加速度法计算结构自由度r的随机位移响应均方根的公式为

${\mathrm{\σ}}_{u_r}={\left({\∫}_{\underset{\‾}{\mathrm{\ω}}}^{\stackrel{\‾}{\mathrm{\ω}}}\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}{\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^2\mathrm{d\ω}\right)}^{1/2}$

式中u表示位移，||(g_q(t))_r||表示复数(g_q(t))_r的模，g_q(t)为n维列向量表示结构在第q个虚拟简谐激励γ_q下的位移响应，其第r项的计算公式为

式中a为n维列向量，只有第r项为1，其它项均为0。为的第i列。b为n行p列由0、1组成的载荷分布矩阵，假如f(t)中第d个力施加在第z个自由度上，则b的第d列中只有第z个元素值是1，d列中其它元素值均为0。e^jωt表示以自然常数e为底数的指数函数，j²＝-1。上式中，

$H_i={({\mathrm{\ω}}_i^2-{\mathrm{\ω}}^2+{2\mathrm{j\ζ}}_i{\mathrm{\ω}}_i\mathrm{\ω})}^{-1}$

x_q＝k^-1(bγ_q)

式中k为结构有限元整体刚度矩阵，x_q是第q个静力载荷bγ_q下的位移向量。

结构自由度r的随机位移响应均方根对每一设计变量的灵敏度的计算公式为

$\frac{{\∂\mathrm{\σ}}_{u_r}}{{\∂\mathrm{\η}}_h}=\frac{1}{{2\mathrm{\σ}}_{u_r}}{\∫}_{\underset{\‾}{\mathrm{\ω}}}^{\stackrel{\‾}{\mathrm{\ω}}}\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}\frac{{\∂\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^2}{{\∂\mathrm{\η}}_h}\mathrm{d\ω}$

式中为偏微分符号，上式中

$\frac{\∂{\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^2}{{\∂\mathrm{\η}}_h}=2\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|\frac{\∂{\left(g_q\left(t\right)\right)}_r\left|\right|}{{\∂\mathrm{\η}}_h}$

式中

$\frac{\∂\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}{{\∂\mathrm{\η}}_h}=(\mathrm{Re}\left({\left(g_q\left(t\right)\right)}_r\right)\mathrm{Re}\left(\frac{\∂{\left(g_q\left(t\right)\right)}_r}{{\∂\mathrm{\η}}_h}\right)+\mathrm{Im}\left({\left(g_q\left(t\right)\right)}_r\right)\mathrm{Im}\left(\frac{\∂{\left(g_q\left(t\right)\right)}_r}{{\∂\mathrm{\η}}_h}\right)){\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^{-1}$

式中Re表示复数的实部，Im表示复数的虚部，而

式中

$a^T\frac{{\∂x}_q}{{\∂\mathrm{\η}}_h}=-{\left({\mathrm{\Λ}}^T\right)}_h\frac{{\∂k}_h}{{\∂\mathrm{\η}}_h}{\left(x_q\right)}_h$

$\frac{{\∂H}_i}{{\∂\mathrm{\η}}_h}=-H_i^2({2\mathrm{\ω}}_i\frac{{\∂\mathrm{\ω}}_i}{{\mathrm{\η}}_h}+{2\mathrm{j\ζ}}_i\mathrm{\ω}\frac{{\∂\mathrm{\ω}}_i}{{\mathrm{\η}}_h})$

式中(x_q)_h表示结构在第q个静力载荷bγ_q下单元h的位移向量，(Λ^T)_h表示结构在静力载荷向量a下单元h的位移向量的转置。g表示模态阶数，表示单元h的第i阶模态振型向量。上式中

Λ＝k^-1a

$\frac{{\∂k}_h}{{\∂\mathrm{\η}}_h}=\frac{{75\mathrm{\η}}_h^4+1}{{15\mathrm{\η}}_h^5+{\mathrm{\η}}_h}{k}_h$

$\frac{{\∂m}_h}{{\∂\mathrm{\η}}_h}=\frac{1}{{\mathrm{\η}}_h}{m}_h$

读取每个单元的体积V_h，计算结构整体质量M及其对每一设计变量的灵敏度计算式分别为

$M=\underset{h=1}{\overset{N_h}{\mathrm{\Σ}}}{\mathrm{\ρ}}_h{V}_h$

$\frac{\∂M}{{\∂\mathrm{\η}}_h}={\mathrm{\ρV}}_h$

(e)根据当前设计变量值和灵敏度值，以结构自由度r的随机位移响应均方根为目标函数，结构整体质量M为设计约束，采用数学规划法对优化问题进行求解得到新的设计变量值；

(f)重复步骤(c)至步骤(e)，直至最近两次迭代计算得到目标函数相对误差小于1％或达到预设的最大迭代次数，得到设计结果。

本发明的有益效果是：该方法采用虚拟激励法结合模态加速度法计算随机激励下的位移响应均方根，然后以结构指定位置的位移响应均方根最小为目标，以结构质量为约束对结构进行设计。相比背景技术的设计方法，本发明方法明显提高了随机位移响应均方根的分析精度，最终设计得到清晰有效的结构构型，从而满足工程实际中大规模结构的设计需求。

本发明方法经过实施例160步迭代后得到设计结果。初始结构指定自由度r＝304的随机位移响应均方根为0.06m，设计得到结构的指定自由度r＝304的随机位移响应均方根为0.0015m，随机位移响应均方根降幅达到99.75％。由图2可见设计得到的结构构型清晰有效，易于在实际工程中使用。本发明方法实施例中结构的自由度数量为14762，是背景技术中结构自由度数量3782的3.9倍。

下面结合附图和具体实施方式对本发明作详细说明。

附图说明

图1是本发明方法实施例中几何模型与边界条件示意图。

图2是本发明方法实施例的设计结果图。

图3是背景技术方法的设计结果图。

具体实施方式

参照图1-2。本发明随机激励下大规模结构设计方法具体步骤如下：

(a)建立设计空间有限元模型：将长宽厚分别为0.8m，0.4m和0.001m的矩形平面结构划分为120×60的正方形网格，结构左边界固定。设置拓扑设计变量η_h初始值均为0.5，结构单元总数量N_h＝7200。给定材料密度ρ＝7800kg/m³和杨氏模量E＝200GPa。给定质量约束上限

(b)设置激励载荷为作用在结构右边界中点竖直向下的随机载荷f(t)，结构只有p＝1个力载荷，因此随机载荷的功率谱密度矩阵S_f(ω)为1维矩阵，其值设为2500N²/(rad/s)，载荷的激励频段为[0，1884]rad/s。根据矩阵LDLT分解，存在下式成立。

$S_f\left(\mathrm{\ω}\right)=\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}{\left({\mathrm{\γ}}_q\right)}^{*}{\left({\mathrm{\γ}}_q\right)}^T$

此处Q＝1为矩阵S_f(ω)的秩，γ₁＝50表示第1个虚拟简谐激励。

(c)根据当前设计变量值，采用以下公式分别计算每一个有限元单元的材料密度ρ_h和杨氏模量E_h

ρ_h＝η_hρ

$E_h=\frac{{15\mathrm{\η}}_h^5+{\mathrm{\η}}_h}{16}E$

更新结构有限元模型中的相应材料属性并进行结构有限元分析。

(d)从有限元分析结果中提取每个单元的刚度矩阵k_h和质量矩阵m_h，前l＝30阶模态频率值ω_i，模态振型为n＝14762行l＝30列矩阵。设置结构前l＝30阶模态的阻尼比ζ_i均为0.03。采用虚拟激励法结合模态加速度法计算结构自由度r＝304的随机位移响应均方根的公式为

${\mathrm{\σ}}_{u_{304}}={\left({\∫}_0^{1884}\underset{q=1}{\overset{1}{\mathrm{\Σ}}}{\left|\right|{\left(g_q\left(t\right)\right)}_{304}\left|\right|}^2\mathrm{d\ω}\right)}^{1/2}$

式中u表示位移，||(g_q(t))₃₀₄||表示复数(g_q(t))₃₀₄的模，g_q(t)为14762维列向量表示结构在第q个虚拟简谐激励γ_q下的位移响应，其第r＝304项的计算公式为

式中a为14762维列向量，只有第r＝304项元素值为1，其它项均为0。b为14762行1列由0，1组成的载荷分布矩阵，f(t)中只有1个力施加在第z＝304个自由度上，因此b的第1列中只有第304个元素值是1，其它元素值均为0。e^jωt表示以自然常数e为底数的指数函数，j²＝-1。上式中

$H_i={({\mathrm{\ω}}_i^2-{\mathrm{\ω}}^2+{2\mathrm{j\ζ}}_i{\mathrm{\ω}}_i\mathrm{\ω})}^{-1}$

x_q＝k^-1(bγ_q)

式中k为结构有限元整体刚度矩阵。x_q是第q个静力载荷bγ_q下的位移向量。

结构自由度r＝304的随机位移响应均方根对每一设计变量的灵敏度的计算公式为

$\frac{{\∂\mathrm{\σ}}_{u_{304}}}{{\∂\mathrm{\η}}_h}=\frac{1}{{2\mathrm{\σ}}_{u_{304}}}{\∫}_0^{1884}\underset{q=1}{\overset{1}{\mathrm{\Σ}}}\frac{{\∂\left|\right|{\left(g_q\left(t\right)\right)}_{304}\left|\right|}^2}{{\∂\mathrm{\η}}_h}\mathrm{d\ω}$

式中

$\frac{\∂{\left|\right|{\left(g_q\left(t\right)\right)}_{304}\left|\right|}^2}{{\∂\mathrm{\η}}_h}=2\left|\right|{\left(g_q\left(t\right)\right)}_{304}\left|\right|\frac{\∂{\left(g_q\left(t\right)\right)}_{304}\left|\right|}{{\∂\mathrm{\η}}_h}$

式中

$\frac{\∂\left|\right|{\left(g_q\left(t\right)\right)}_{304}\left|\right|}{{\∂\mathrm{\η}}_h}=(\mathrm{Re}\left({\left(g_q\left(t\right)\right)}_{304}\right)\mathrm{Re}\left(\frac{\∂{\left(g_q\left(t\right)\right)}_{304}}{{\∂\mathrm{\η}}_h}\right)+\mathrm{Im}\left({\left(g_q\left(t\right)\right)}_{304}\right)\mathrm{Im}\left(\frac{\∂{\left(g_q\left(t\right)\right)}_{304}}{{\∂\mathrm{\η}}_h}\right)){\left|\right|{\left(g_q\left(t\right)\right)}_{304}\left|\right|}^{-1}$

式中Re表示复数的实部，Im表示复数的虚部，而

式中

$a^T\frac{{\∂x}_q}{{\∂\mathrm{\η}}_h}=-{\left({\mathrm{\Λ}}^T\right)}_h\frac{{\∂k}_h}{{\∂\mathrm{\η}}_h}{\left(x_q\right)}_h$

$\frac{{\∂H}_i}{{\∂\mathrm{\η}}_h}=-H_i^2({2\mathrm{\ω}}_i\frac{{\∂\mathrm{\ω}}_i}{{\mathrm{\η}}_h}+{2\mathrm{j\ζ}}_i\mathrm{\ω}\frac{{\∂\mathrm{\ω}}_i}{{\mathrm{\η}}_h})$

式中(x_q)_h表示结构在第q个静力载荷bγ_q下单元h的位移向量，(Λ^T)_h表示结构在静力载荷向量a下单元h的位移向量的转置。g表示模态阶数，表示单元h的第i阶模态振型向量。上式中

Λ＝k^-1a

$\frac{{\∂k}_h}{{\∂\mathrm{\η}}_h}=\frac{{75\mathrm{\η}}_h^4+1}{{15\mathrm{\η}}_h^5+{\mathrm{\η}}_h}{k}_h$

$\frac{{\∂m}_h}{{\∂\mathrm{\η}}_h}=\frac{1}{{\mathrm{\η}}_h}{m}_h$

读取每个单元的体积V_h，计算结构整体质量M及其对每一设计变量的灵敏度计算式分别为

$M=\underset{h=1}{\overset{N_h}{\mathrm{\Σ}}}{\mathrm{\ρ}}_h{V}_h$

$\frac{\∂M}{{\∂\mathrm{\η}}_h}={\mathrm{\ρV}}_h$

(e)根据当前设计变量值和灵敏度值，以结构自由度r＝304的随机位移响应均方根为目标函数，结构整体质量M为设计约束，采用数学规划法对优化问题进行求解得到新的设计变量值；

(f)重复步骤(c)至步骤(e)，直至最近两次迭代计算得到目标函数相对误差小于1％或达到预设的最大迭代次数，得到设计结果。

本实施例经过160步迭代后得到设计结果。初始结构指定自由度r＝304的随机位移响应均方根为0.06m，设计得到结构的指定自由度r＝304的随机位移响应均方根为0.0015m，随机位移响应均方根降幅达到99.75％。由图2可见设计得到的结构构型清晰有效，易于在实际工程中使用。本实施例中结构的自由度数量为14762，是背景文献中结构自由度数量3782的3.9倍，对此大规模结构使用文献中的技术方法得到的设计结构如图3所示，其边界模糊、构型复杂、结构上空洞较多，无法实际应用。本实施例表明了本发明方法在处理随机激励下大规模结构设计问题上的有效性。

Claims (1)

1.一种随机激励下大规模结构设计方法，其特征在于包括以下步骤：

(a)建立设计空间有限元模型，设置拓扑设计变量η_h初始值，h是正整数表示单元编号，1≤h≤N_h，N_h表示结构单元总数量；给定材料密度ρ和杨氏模量E；给定质量约束上限

(b)设置激励载荷，给出随机激励f(t)的功率谱密度矩阵S_f(ω)，f(t)为p维列向量，p为载荷中力的个数，t表示时间，S_f(ω)为p维方阵，其下标f表示其为激励f(t)的功率谱矩阵；ω为激励角频率，载荷的激励频段为 表示激励角频率的下限，表示激励角频率的上限；根据矩阵LDLT分解，存在下式成立

$S_f\left(\mathrm{\ω}\right)=\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}{\left({\mathrm{\γ}}_q\right)}^{*}{\left({\mathrm{\γ}}_q\right)}^T$

其中Q为矩阵S_f(ω)的秩，γ_q为p维列向量表示第q个虚拟简谐激励，1≤q≤Q，上标T表示向量或矩阵的转置；

(c)根据当前设计变量值，采用以下公式分别计算每一个有限元单元的材料密度ρ_h和杨氏模量E_h

ρ_h＝η_hρ

$E_h=\frac{15{\mathrm{\η}}_h^5+{\mathrm{\η}}_h}{16}E$

更新结构有限元模型中的相应材料属性并进行结构有限元分析；

(d)从有限元分析结果中提取每个单元的刚度矩阵k_h和质量矩阵m_h，结构的前l阶模态频率值ω_i，1≤i≤l，模态振型 为n行l列矩阵，n为结构总自由度数目；设置结构前l阶模态的阻尼比ζ_i；采用虚拟激励法结合模态加速度法计算结构自由度r的随机位移响应均方根的公式为

${\mathrm{\σ}}_{u_r}={\left({\∫}_{\underset{\‾}{\mathrm{\ω}}}^{\stackrel{\‾}{\mathrm{\ω}}}\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}{\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^2\mathrm{d\ω}\right)}^{1/2}$

式中u表示位移，||(g_q(t))_r||表示复数(g_q(t))_r的模，g_q(t)为n维列向量表示结构在第q个虚拟简谐激励γ_q下的位移响应，其第r项的计算公式为

式中a为n维列向量，只有第r项为1，其它项均为0；为的第i列；b为n行p列由0、1组成的载荷分布矩阵，假如f(t)中第d个力施加在第z个自由度上，则b的第d列中只有第z个元素值是1，d列中其它元素值均为0；e^jωt表示以自然常数e为底数的指数函数，j²＝-1；上式中，

$H_i={({\mathrm{\ω}}_i^2-{\mathrm{\ω}}^2+2j{\mathrm{\ζ}}_i{\mathrm{\ω}}_i\mathrm{\ω})}^{-1}$

x_q＝k^-1(bγ_q)

式中k为结构有限元整体刚度矩阵，x_q是第q个静力载荷bγ_q下的位移向量；

结构自由度r的随机位移响应均方根对每一设计变量的灵敏度的计算公式为

$\frac{{\∂\mathrm{\σ}}_{u_r}}{{\∂\mathrm{\η}}_h}=\frac{1}{{2\mathrm{\σ}}_{u_r}}{\∫}_{\underset{\‾}{\mathrm{\ω}}}^{\stackrel{\‾}{\mathrm{\ω}}}\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}\frac{{\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^2}{{\∂\mathrm{\η}}_h}\mathrm{d\ω}$

式中为偏微分符号，上式中

$\frac{\∂{\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^2}{{\∂\mathrm{\η}}_h}=2\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|\frac{\∂\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}{{\∂\mathrm{\η}}_h}$

式中

$\frac{\∂\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}{{\∂\mathrm{\η}}_h}=(\mathrm{Re}\left({\left(g_q\left(t\right)\right)}_r\right)\mathrm{Re}\left(\frac{\∂{\left(g_q\left(t\right)\right)}_r}{{\∂\mathrm{\η}}_h}\right)+\mathrm{Im}\left({\left(g_q\left(t\right)\right)}_r\right)\mathrm{Im}\left(\frac{\∂{\left(g_q\left(t\right)\right)}_r}{{\∂\mathrm{\η}}_h}\right)){\left|\right|{\left(g_q\left(t\right)\right)}_r\left|\right|}^{-1}$

式中Re表示复数的实部，Im表示复数的虚部，而

式中

$a^T\frac{{\∂x}_q}{{\∂\mathrm{\η}}_h}=-{\left({\mathrm{\Λ}}^T\right)}_h\frac{{\∂k}_h}{{\∂\mathrm{\η}}_h}{\left(x_q\right)}_h$

$\frac{{\∂H}_i}{{\∂\mathrm{\η}}_h}=-H_i^2({2\mathrm{\ω}}_i\frac{{\∂\mathrm{\ω}}_i}{{\mathrm{\η}}_h}+{2\mathrm{j\ζ}}_i\mathrm{\ω}\frac{{\∂\mathrm{\ω}}_i}{{\mathrm{\η}}_h})$

式中(x_q)_h表示结构在第q个静力载荷bγ_q下单元h的位移向量，(Λ^T)_h表示结构在静力载荷向量a下单元h的位移向量的转置；g表示模态阶数，1≤g≤l，表示单元h的第i阶模态振型向量；上式中

Λ＝k^-1a

$\frac{{\∂k}_h}{{\∂\mathrm{\η}}_h}=\frac{75{\mathrm{\η}}_h^4+1}{15{\mathrm{\η}}_h^5+{\mathrm{\η}}_h}{k}_h$

$\frac{{\∂m}_h}{{\∂\mathrm{\η}}_h}=\frac{1}{{\mathrm{\η}}_h}{m}_h$

读取每个单元的体积V_h，计算结构整体质量M及其对每一设计变量的灵敏度计算式分别为

$M=\underset{h=1}{\overset{N_h}{\mathrm{\Σ}}}{\mathrm{\ρ}}_h{V}_h$

$\frac{\∂M}{{\∂\mathrm{\η}}_h}={\mathrm{\ρV}}_h$

(e)根据当前设计变量值和灵敏度值，以结构自由度r的随机位移响应均方根为目标函数，结构整体质量M为设计约束，采用数学规划法对优化问题进行求解得到新的设计变量值；

(f)重复步骤(c)至步骤(e)，直至最近两次迭代计算得到目标函数相对误差小于1％或达到预设的最大迭代次数，得到设计结果。