CN111783236B - 基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法 - Google Patents

基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法 Download PDF

Info

Publication number
CN111783236B
CN111783236B CN202010407601.XA CN202010407601A CN111783236B CN 111783236 B CN111783236 B CN 111783236B CN 202010407601 A CN202010407601 A CN 202010407601A CN 111783236 B CN111783236 B CN 111783236B
Authority
CN
China
Prior art keywords
sample
failure
function
value
samples
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010407601.XA
Other languages
English (en)
Other versions
CN111783236A (zh
Inventor
王璐
吕震宙
王嘉琪
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Northwestern Polytechnical University
Original Assignee
Northwestern Polytechnical University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Northwestern Polytechnical University filed Critical Northwestern Polytechnical University
Priority to CN202010407601.XA priority Critical patent/CN111783236B/zh
Publication of CN111783236A publication Critical patent/CN111783236A/zh
Application granted granted Critical
Publication of CN111783236B publication Critical patent/CN111783236B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/17—Mechanical parametric or variational design
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Geometry (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • General Engineering & Computer Science (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Complex Calculations (AREA)
  • Testing Of Devices, Machine Parts, Or Other Structures Thereof (AREA)

Abstract

本发明提供一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，包括：训练Kriging模型；将得到的每一层失效事件对应的条件失效概率的估算值做乘法运算，得到最终的涡轮机匣失效概率的估计值
Figure DDA0002491911290000011
得到可靠性局部灵敏度指标和可靠性全局灵敏度指标。本发明提供的基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法具有以下优点：本发明提供一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，是一种有效可行的同时求解可靠性局部灵敏度指标和可靠性全局灵敏度指标的方法。可以克服现有方法所需样本量较大及耗时较长的问题，本发明提高涡轮机匣灵敏度分析的效率。

Description

基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法
技术领域
本发明属于灵敏度分析技术领域，具体涉及一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法。
背景技术
涡轮机匣是航空发动机的重要零件之一，是整个发动机的基座，其外形结构复杂，不同的发动机、发动机不同部位，其机匣形状各不相同，机匣零件的功能决定了机匣的形状，但其基本特征是圆筒形或圆锥形的壳体和支板组成的构件。由于机匣零件设计难度大、周期长，在整个发动机的设计中，机匣的设计占相当大的比重，提高机匣的设计效率对压缩发动机整机的设计周期具有重要意义。
为了查出对机匣失效概率影响较大的设计参数，需要对机匣进行灵敏度分析。传统的灵敏度分析方法需要多次有限元模型，在实际应用中需要过多的分析时间，且对于小失效概率问题需要更多的样本量，具有灵敏度分析效率低的问题，从而制约了其在工程实践中的应用。
发明内容
针对现有技术存在的缺陷，本发明提供一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，可有效解决上述问题。
本发明采用的技术方案如下：
本发明提供一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，包括以下步骤：
步骤1，建立涡轮机匣模型；对所述涡轮机匣模型进行分析，确定影响涡轮机匣灵敏度的n维随机向量x＝{x1,x2,…,xn}T；
步骤2，令k＝1；
步骤3，根据n维随机向量x的联合概率密度函数fX(x)，产生Nk个样本构建样本集合
Figure GDA0003982326970000021
其中，
Figure GDA0003982326970000022
代表样本集合Sk中的第i个样本，为n维随机向量x＝{x1,x2,…,xn}T中各维自变量取具体数值后的n维数组；
步骤4，从样本集合Sk中随机选择T1个样本，形成训练集合
Figure GDA0003982326970000023
调用涡轮机匣失效的实际功能函数y＝g(x)，分别计算训练集合T中每个样本的实际功能函数值，得到实际功能函数值集合
Figure GDA0003982326970000024
将实际功能函数值集合gt中的各个元素从大到小排序，排序结果为
Figure GDA0003982326970000025
取第
Figure GDA0003982326970000026
个实际功能函数值作为第k层失效事件Fk＝{x:g(x)≤bk}的临界值bk，即
Figure GDA0003982326970000027
其中，
Figure GDA0003982326970000028
表示向上取整符号，p0为预设已知值；
步骤5，由训练集合
Figure GDA0003982326970000029
以及相应的实际功能函数值集合
Figure GDA00039823269700000210
构建得到Kriging模型gK(x)如下：
Figure GDA00039823269700000211
其中：
f(x)＝[f1(x),f2(x),...,fp(x)]T是n维随机向量x的p个基函数形成的集合，用于提供设计空间内的全局近似模型；p是基函数的个数；
β＝[β1,β2,...,βp]为回归函数的p个系数值形成的集合；每个回归函数的系数值通过训练集合
Figure GDA00039823269700000212
以及相应的实际功能函数值集合
Figure GDA00039823269700000213
得到；
z(x)代表随机过程；
步骤6，对于样本集合Sk中任意的样本
Figure GDA00039823269700000214
代入Kriging模型gK(x)，得到对应的一个代理功能函数值
Figure GDA00039823269700000215
每个代理功能函数值
Figure GDA00039823269700000216
均服从代理功能函数均值为
Figure GDA00039823269700000217
标准差为
Figure GDA00039823269700000218
的高斯分布：
Figure GDA00039823269700000219
因此，对于任意的样本
Figure GDA00039823269700000220
可得到代理功能函数均值为
Figure GDA00039823269700000221
代理功能函数标准差
Figure GDA0003982326970000031
步骤7，根据下式，计算样本集合Sk中每个样本
Figure GDA0003982326970000032
对应的U学习函数值
Figure GDA0003982326970000033
Figure GDA0003982326970000034
步骤8，因此，由于样本集合Sk中共有Nk个样本，共计算得到Nk个U学习函数值，判断每个U学习函数值是否均大于等于2，即：判断下式是否成立：
Figure GDA0003982326970000035
如果成立，则代表当前Kriging模型已收敛，停止更新当前Kriging模型，得到训练好的Kriging模型，执行步骤9；
否则，根据下式从样本集合Sk选取新的样本，也就是选取具有最小U学习函数值的样本，表示为
Figure GDA0003982326970000036
Figure GDA0003982326970000037
计算样本
Figure GDA0003982326970000038
的实际功能函数值
Figure GDA0003982326970000039
将
Figure GDA00039823269700000310
和
Figure GDA00039823269700000311
分别加入训练集合T和实际功能函数值集合gt中，由此得到更新后的训练集合T和更新后的实际功能函数值集合gt，再返回步骤5；
步骤9，利用训练好的Kriging模型，计算样本集合Sk中每个样本
Figure GDA00039823269700000312
对应的代理功能函数均值
Figure GDA00039823269700000313
根据下述指示函数
Figure GDA00039823269700000314
从样本集合Sk中挑选可使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本：
Figure GDA00039823269700000315
若
Figure GDA00039823269700000316
则样本
Figure GDA00039823269700000317
可使第k层失效事件Fk＝{x:g(x)≤bk}发生；
样本集合Sk中，所有使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本组成的集合记为失效样本集合
Figure GDA0003982326970000041
其中，Mk为使第k层失效事件Fk＝{x:g(x)≤bk}发生的失效样本总数；
根据下式计算条件失效概率Pk的估算值
Figure GDA0003982326970000042
Figure GDA0003982326970000043
其中：
Figure GDA0003982326970000044
表示第k层失效事件Fk发生概率的估算值；
Figure GDA0003982326970000045
表示第k-1层失效事件Fk-1发生的条件下，第k层失效事件Fk发生概率的估算值；
P{Fk-1}表示第k-1层失效事件Fk-1发生概率；
Figure GDA0003982326970000046
表示重要抽样密度函数；
根据下式计算条件失效概率Pk的估算值
Figure GDA0003982326970000047
对
Figure GDA0003982326970000048
的偏导数
Figure GDA0003982326970000049
Figure GDA00039823269700000410
其中：
Figure GDA00039823269700000411
代表n维随机向量x＝{x1,x2,…,xn}T中的第r维自变量xr的第s个分布参数；
Pq代表第q层失效事件Fq发生的概率，q＝1,2,...,k-1；
步骤10，判断bk是否等于0，如果bk＝0，则循环过程结束，执行步骤15；否则，执行步骤11；
步骤11，令k＝k+1，在失效样本集合RSk-1中，选出联合概率密度函数值最大的样本，作为重要抽样密度函数
Figure GDA00039823269700000412
的抽样中心；
然后，根据n维随机向量x的重要抽样密度函数
Figure GDA0003982326970000051
以联合概率密度函数值最大的样本为抽样中心，产生Nk个样本构建样本集合
Figure GDA0003982326970000052
步骤12，利用训练好的Kriging模型，计算样本集合Sk中每个样本
Figure GDA0003982326970000053
对应的代理功能函数值
Figure GDA0003982326970000054
和代理功能函数均值
Figure GDA0003982326970000055
根据下述指示函数
Figure GDA0003982326970000056
从样本集合Sk中挑选可使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本：
Figure GDA0003982326970000057
样本集合Sk中，所有使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本组成的集合记为失效样本集合
Figure GDA0003982326970000058
其中，Mk为使第k层失效事件Fk＝{x:g(x)≤bk}发生的失效样本总数；
步骤13，利用训练好的Kriging模型，计算失效样本集合RSk中每个样本对应的代理功能函数均值；然后，将失效样本集合RSk中每个样本对应的代理功能函数均值从大到小排序，
Figure GDA0003982326970000059
取第
Figure GDA00039823269700000510
个代理功能函数均值作为第k层失效事件Fk＝{x:g(x)≤bk}的临界值bk，即
Figure GDA00039823269700000511
步骤14，判断bk是否小于0；如果小于，则令bk＝0，然后返回步骤7；如果不小于，直接返回步骤7；
步骤15，将至此得到的每一层失效事件对应的条件失效概率的估算值做乘法运算，得到最终的涡轮机匣失效概率的估计值
Figure GDA00039823269700000512
具体的，令m＝k，根据下式计算涡轮机匣失效概率的估计值
Figure GDA00039823269700000513
Figure GDA00039823269700000514
通过下式计算得到最终的n维随机向量x＝{x1,x2,…,xn}T中的第r维自变量xr的第s个分布参数的可靠性局部灵敏度指标
Figure GDA0003982326970000061
Figure GDA0003982326970000062
优选的，影响涡轮机匣灵敏度的n维随机向量中，各分量包括高温区域温度、流动冷却区域温度、冲击冷却区域温度、前段机匣长度及机匣内半径。
优选的，还包括：
步骤16，最后一层失效事件即为第m层失效事件Fm，使第m层失效事件Fm发生的样本组成的失效样本集合为RSm，设失效样本集合RSm中失效样本总数为NC；
将失效样本集合RSm中的所有失效样本放入矩阵C，则矩阵C中的样本表示为
Figure GDA0003982326970000063
其中，b＝1,2,...,NC；
步骤17，初始化
Figure GDA0003982326970000064
并令β＝1；
步骤18，构建并计算根据重要抽样密度函数产生的样本
Figure GDA0003982326970000065
向目标分布中第β+1状态zβ+1转化的接受概率
Figure GDA0003982326970000066
为：
Figure GDA0003982326970000067
其中：
Figure GDA0003982326970000068
为样本
Figure GDA0003982326970000069
的联合概率密度函数；
Figure GDA00039823269700000610
为第β状态zβ的重要抽样密度函数；
fX(zβ)为第β状态zβ的联合概率密度函数；
Figure GDA00039823269700000611
为样本
Figure GDA00039823269700000612
的重要抽样密度函数；
步骤19，根据以下判断式，判断zβ+1是否接受
Figure GDA00039823269700000613
Figure GDA00039823269700000614
其中，random[0,1]表示0与1之间的随机数；
即：如果
Figure GDA00039823269700000615
则将
Figure GDA00039823269700000616
赋值给zβ+1；否则，将zβ赋值给zβ+1；
步骤20，令β＝β+1，判断β是否大于NC，如果大于，则执行步骤21；否则，执行步骤18；
步骤21，由此得到NC个zβ，分别为
Figure GDA0003982326970000071
将
Figure GDA0003982326970000072
放入到矩阵B，矩阵B中样本服从输入变量的条件概率密度函数fX(x|F)；
步骤22，根据矩阵B中的样本估计失效域内xr条件概率密度函数
Figure GDA0003982326970000073
再结合步骤15得到的涡轮机匣失效概率的估计值
Figure GDA0003982326970000074
以及第r维自变量xr的联合概率密度函数
Figure GDA0003982326970000075
代入下式得到输入变量xr的可靠性全局灵敏度指标Sr：
Figure GDA0003982326970000076
其中：
E[·]表示期望运算。
本发明提供的基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法具有以下优点：
本发明提供一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，是一种有效可行的同时求解可靠性局部灵敏度指标和可靠性全局灵敏度指标的方法。可以克服现有方法所需样本量较大及耗时较长的问题，本发明提高涡轮机匣灵敏度分析的效率。
附图说明
图1为本发明提供的基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法的流程示意图。
图2是本发明应用子集模拟引入中间事件的示意图；
图3是本发明涡轮机匣的平面简图；
图4是本发明涡轮机匣的受力示意图。
具体实施方式
为了使本发明所解决的技术问题、技术方案及有益效果更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。
本发明提供一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，是一种有效可行的同时求解可靠性局部灵敏度指标和可靠性全局灵敏度指标的方法。该方法旨在克服现有方法所需样本量较大及耗时较长的问题，本发明提高灵敏度分析的效率。
本发明采用的技术方案主要包含两个部分，分别为：(一)利用子集模拟法产生中间失效事件和分层抽样的样本、利用Kriging模型自适应的选点和代理功能函数；(二)基于贝叶斯公式求解失效概率的可靠性全局灵敏度指标。
参考图1，本发明提供的一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，包括以下步骤：
第一部分：利用子集模拟法和自适应Kriging构造功能函数的代理模型及求解可靠性局部灵敏度指标，其中，各子集的样本由重要抽样法得到的抽样密度函数产生，Kriging模型由U学习函数进行加点并判断是否满足收敛条件。具体步骤包括步骤1-步骤14：
步骤1，建立涡轮机匣模型；对所述涡轮机匣模型进行分析，确定影响涡轮机匣灵敏度的n维随机向量x＝{x1,x2,…,xn}T；
步骤2，令k＝1；
步骤3，根据n维随机向量x的联合概率密度函数fX(x)，产生Nk个样本构建样本集合
Figure GDA0003982326970000081
其中，
Figure GDA0003982326970000082
代表样本集合Sk中的第i个样本，为n维随机向量x＝{x1,x2,…,xn}T中各维自变量取具体数值后的n维数组；
步骤4，从样本集合Sk中随机选择T1个样本，形成训练集合
Figure GDA0003982326970000083
调用涡轮机匣失效的实际功能函数y＝g(x)，分别计算训练集合T中每个样本的实际功能函数值，得到实际功能函数值集合
Figure GDA0003982326970000084
将实际功能函数值集合gt中的各个元素从大到小排序，排序结果为
Figure GDA0003982326970000085
取第
Figure GDA0003982326970000086
个实际功能函数值作为第k层失效事件Fk＝{x:g(x)≤bk}的临界值bk，即
Figure GDA0003982326970000091
其中，
Figure GDA0003982326970000092
表示向上取整符号，p0为预设已知值；
步骤5，由训练集合
Figure GDA0003982326970000093
以及相应的实际功能函数值集合
Figure GDA0003982326970000094
构建得到Kriging模型gK(x)如下：
Figure GDA0003982326970000095
其中：
f(x)＝[f1(x),f2(x),...,fp(x)]T是n维随机向量x的p个基函数形成的集合，用于提供设计空间内的全局近似模型；p是基函数的个数；
β＝[β1,β2,...,βp]为回归函数的p个系数值形成的集合；每个回归函数的系数值通过训练集合
Figure GDA0003982326970000096
以及相应的实际功能函数值集合
Figure GDA0003982326970000097
得到；
z(x)代表随机过程；
步骤6，对于样本集合Sk中任意的样本
Figure GDA0003982326970000098
代入Kriging模型gK(x)，得到对应的一个代理功能函数值
Figure GDA0003982326970000099
每个代理功能函数值
Figure GDA00039823269700000910
均服从代理功能函数均值为
Figure GDA00039823269700000911
标准差为
Figure GDA00039823269700000912
的高斯分布：
Figure GDA00039823269700000913
因此，对于任意的样本
Figure GDA00039823269700000914
可得到代理功能函数均值为
Figure GDA00039823269700000915
代理功能函数标准差
Figure GDA00039823269700000916
步骤7，根据下式，计算样本集合Sk中每个样本
Figure GDA00039823269700000917
对应的U学习函数值
Figure GDA00039823269700000918
Figure GDA00039823269700000919
步骤8，因此，由于样本集合Sk中共有Nk个样本，共计算得到Nk个U学习函数值，判断每个U学习函数值是否均大于等于2，即：判断下式是否成立：
Figure GDA00039823269700000920
如果成立，则代表当前Kriging模型已收敛，停止更新当前Kriging模型，得到训练好的Kriging模型，执行步骤9；
否则，根据下式从样本集合Sk选取新的样本，也就是选取具有最小U学习函数值的样本，表示为
Figure GDA0003982326970000101
Figure GDA0003982326970000102
计算样本
Figure GDA0003982326970000103
的实际功能函数值
Figure GDA0003982326970000104
将
Figure GDA0003982326970000105
和
Figure GDA0003982326970000106
分别加入训练集合T和实际功能函数值集合gt中，由此得到更新后的训练集合T和更新后的实际功能函数值集合gt，再返回步骤5；
步骤9，利用训练好的Kriging模型，计算样本集合Sk中每个样本
Figure GDA0003982326970000107
对应的代理功能函数均值
Figure GDA0003982326970000108
根据下述指示函数
Figure GDA0003982326970000109
从样本集合Sk中挑选可使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本：
Figure GDA00039823269700001010
若
Figure GDA00039823269700001011
则样本
Figure GDA00039823269700001012
可使第k层失效事件Fk＝{x:g(x)≤bk}发生；
样本集合Sk中，所有使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本组成的集合记为失效样本集合
Figure GDA00039823269700001013
其中，Mk为使第k层失效事件Fk＝{x:g(x)≤bk}发生的失效样本总数；
根据下式计算条件失效概率Pk的估算值
Figure GDA00039823269700001014
Figure GDA00039823269700001015
其中：
Figure GDA00039823269700001016
表示第k层失效事件Fk发生概率的估算值；
Figure GDA0003982326970000111
表示第k-1层失效事件Fk-1发生的条件下，第k层失效事件Fk发生概率的估算值；
P{Fk-1}表示第k-1层失效事件Fk-1发生概率；
Figure GDA0003982326970000112
表示重要抽样密度函数；
根据下式计算条件失效概率Pk的估算值
Figure GDA0003982326970000113
对
Figure GDA0003982326970000114
的偏导数
Figure GDA0003982326970000115
Figure GDA0003982326970000116
其中：
Figure GDA0003982326970000117
代表n维随机向量x＝{x1,x2,…,xn}T中的第r维自变量xr的第s个分布参数；
Pq代表第q层失效事件Fq发生的概率，q＝1,2,...,k-1；
步骤10，判断bk是否等于0，如果bk＝0，则循环过程结束，执行步骤15；否则，执行步骤11；
步骤11，令k＝k+1，在失效样本集合RSk-1中，选出联合概率密度函数值最大的样本，作为重要抽样密度函数
Figure GDA0003982326970000118
的抽样中心；
然后，根据n维随机向量x的重要抽样密度函数
Figure GDA0003982326970000119
以联合概率密度函数值最大的样本为抽样中心，产生Nk个样本构建样本集合
Figure GDA00039823269700001110
步骤12，利用训练好的Kriging模型，计算样本集合Sk中每个样本
Figure GDA00039823269700001111
对应的代理功能函数值
Figure GDA00039823269700001112
和代理功能函数均值
Figure GDA00039823269700001113
根据下述指示函数
Figure GDA00039823269700001114
从样本集合Sk中挑选可使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本：
Figure GDA0003982326970000121
样本集合Sk中，所有使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本组成的集合记为失效样本集合
Figure GDA0003982326970000122
其中，Mk为使第k层失效事件Fk＝{x:g(x)≤bk}发生的失效样本总数；
步骤13，利用训练好的Kriging模型，计算失效样本集合RSk中每个样本对应的代理功能函数均值；然后，将失效样本集合RSk中每个样本对应的代理功能函数均值从大到小排序，
Figure GDA0003982326970000123
取第
Figure GDA0003982326970000124
个代理功能函数均值作为第k层失效事件Fk＝{x:g(x)≤bk}的临界值bk，即
Figure GDA0003982326970000125
步骤14，判断bk是否小于0；如果小于，则令bk＝0，然后返回步骤7；如果不小于，直接返回步骤7；
步骤15，将至此得到的每一层失效事件对应的条件失效概率的估算值做乘法运算，得到最终的涡轮机匣失效概率的估计值
Figure GDA0003982326970000126
具体的，令m＝k，根据下式计算涡轮机匣失效概率的估计值
Figure GDA0003982326970000127
Figure GDA0003982326970000128
通过下式计算得到最终的n维随机向量x＝{x1,x2,…,xn}T中的第r维自变量xr的第s个分布参数的可靠性局部灵敏度指标
Figure GDA0003982326970000129
Figure GDA00039823269700001210
第二部分：基于Bayes公式
Figure GDA00039823269700001211
能建立无条件失效概率P{F}与条件失效概率P{F|Xi}之间的关系，由此可以进一步得到基于失效概率的可靠性全局灵敏度指标Si，但在步骤一中得到的失效样本均服从条件重要抽样概率密度函数
Figure GDA0003982326970000131
故为了求解Si需先利用Metropolis-Hastings准则将其转换为服从fX(x|F)的样本。具体包括以下步骤16-步骤22。
步骤16，最后一层失效事件即为第m层失效事件Fm，使第m层失效事件Fm发生的样本组成的失效样本集合为RSm，设失效样本集合RSm中失效样本总数为NC；
由于使最后一层失效事件Fm发生的失效样本较多，故利用最后一层的失效样本进行转换。
将失效样本集合RSm中的所有失效样本放入矩阵C，则矩阵C中的样本表示为
Figure GDA0003982326970000132
其中，b＝1,2,...,NC；
步骤17，初始化
Figure GDA0003982326970000133
并令β＝1；
步骤18，构建并计算根据重要抽样密度函数产生的样本
Figure GDA0003982326970000134
向目标分布中第β+1状态zβ+1转化的接受概率
Figure GDA0003982326970000135
为：
Figure GDA0003982326970000136
其中：
Figure GDA0003982326970000137
为样本
Figure GDA0003982326970000138
的联合概率密度函数；
Figure GDA0003982326970000139
为第β状态zβ的重要抽样密度函数；
fX(zβ)为第β状态zβ的联合概率密度函数；
Figure GDA00039823269700001310
为样本
Figure GDA00039823269700001311
的重要抽样密度函数；
步骤19，根据以下判断式，判断zβ+1是否接受
Figure GDA00039823269700001312
Figure GDA00039823269700001313
其中，random[0,1]表示0与1之间的随机数；
即：如果
Figure GDA00039823269700001314
则将
Figure GDA00039823269700001315
赋值给zβ+1；否则，将zβ赋值给zβ+1；
步骤20，令β＝β+1，判断β是否大于NC，如果大于，则执行步骤21；否则，执行步骤18；
步骤21，由此得到NC个zβ，分别为
Figure GDA0003982326970000141
将
Figure GDA0003982326970000142
放入到矩阵B，矩阵B中样本服从输入变量的条件概率密度函数fX(x|F)；
步骤22，根据矩阵B中的样本估计失效域内xr条件概率密度函数
Figure GDA0003982326970000143
再结合步骤15得到的涡轮机匣失效概率的估计值
Figure GDA0003982326970000144
以及第r维自变量xr的联合概率密度函数
Figure GDA0003982326970000145
代入下式得到输入变量xr的可靠性全局灵敏度指标Sr：
Figure GDA0003982326970000146
其中：
E[·]表示期望运算。
下面介绍一种具体实施例：
针对某型发动机涡轮机匣进行灵敏度分析，首先需要建立参数化的涡轮机匣部件模型。为了简化原始机匣模型，本发明忽略安装边上的孔及支座的影响，将该机匣简化为轴对称的模型，机匣外表面的加强筋不予考虑，机匣的安装边的安装孔也简化。简化后的机匣由外机匣组成，确定性设计参数见表1。假定6个随机输入变量分别为高温区域温度、流动冷却区域温度、冲击冷却区域温度、前段机匣长度及机匣内半径，随机输入变量的分布参数表见表2。之后对其利用ANSYS软件进行有限元建模，并进行参数化处理可供MATLAB调用。
表1涡轮机匣参数表
Figure GDA0003982326970000147
表2随机输入变量的分布参数
Figure GDA0003982326970000151
然后，按本发明提供的基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，对涡轮机匣灵敏度进行分析，结果见表3：
表3可靠性灵敏度指标的计算结果
Figure GDA0003982326970000152
因此，对于6维随机向量x中的每一维分量，分别为高温区域温度T1(℃)、流动冷却区域温度T2(℃)、冲击冷却区域温度T3(℃)、压强p(Mpa)，前段机匣长度l1和机匣内半径r1(mm)，对于每一维分量，基于步骤15的可靠性局部灵敏度指标计算公式，分别取分布参数为均值和标准差，可分别计算得到对均值的可靠性局部灵敏度指标，以及对标准差的可靠性局部灵敏度指标；根据步骤22的公式，可得到每一维分量对应的可靠性全局灵敏度指标。
通过对表3进行分析，可以看出各维分量对均值的可靠性局部灵敏度指标、对标准差的可靠性局部灵敏度指标以及可靠性全局灵敏度指标的排序是一致的，即：对均值的可靠性局部灵敏度指标，高温区域温度T1(℃)、流动冷却区域温度T2(℃)、冲击冷却区域温度T3(℃)、压强p(Mpa)，前段机匣长度l1和机匣内半径r1(mm)的值逐渐增大；同样，对标准差的可靠性局部灵敏度指标，这6维分量的值逐渐增大；同样，对可靠性全局灵敏度指标，这6维分量的值逐渐增大。并且，机匣内半径r1(mm)的各个指标值均为最大，说明其对涡轮机匣部失效有最重要的影响，因此，机匣内半径r1(mm)为对涡轮机匣可靠性最为敏感的设计参数，在对涡轮机匣的可靠性优化设计过程中，应着重考虑机匣内半径的尺寸的优化。以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视本发明的保护范围。

Claims (1)

1.一种基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法，其特征在于，包括以下步骤：
步骤1，建立涡轮机匣模型；对所述涡轮机匣模型进行分析，确定影响涡轮机匣灵敏度的n维随机向量x＝{x1,x2,…,xn}T；
步骤2，令k＝1；
步骤3，根据n维随机向量x的联合概率密度函数fX(x)，产生Nk个样本构建样本集合
Figure FDA0003982326960000011
其中，
Figure FDA0003982326960000012
代表样本集合Sk中的第i个样本，为n维随机向量x＝{x1,x2,…,xn}T中各维自变量取具体数值后的n维数组；
步骤4，从样本集合Sk中随机选择T1个样本，形成训练集合
Figure FDA0003982326960000013
调用涡轮机匣失效的实际功能函数y＝g(x)，分别计算训练集合T中每个样本的实际功能函数值，得到实际功能函数值集合
Figure FDA0003982326960000014
将实际功能函数值集合gt中的各个元素从大到小排序，排序结果为
Figure FDA0003982326960000015
取第
Figure FDA0003982326960000016
个实际功能函数值作为第k层失效事件Fk＝{x:g(x)≤bk}的临界值bk，即
Figure FDA0003982326960000017
其中，
Figure FDA0003982326960000018
表示向上取整符号，p0为预设已知值；
步骤5，由训练集合
Figure FDA0003982326960000019
以及相应的实际功能函数值集合
Figure FDA00039823269600000110
构建得到Kriging模型gK(x)如下：
Figure FDA00039823269600000111
其中：
f(x)＝[f1(x),f2(x),...,fp(x)]T是n维随机向量x的p个基函数形成的集合，用于提供设计空间内的全局近似模型；p是基函数的个数；
β＝[β1,β2,...,βp]为回归函数的p个系数值形成的集合；每个回归函数的系数值通过训练集合
Figure FDA00039823269600000112
以及相应的实际功能函数值集合
Figure FDA00039823269600000113
得到；
z(x)代表随机过程；
步骤6，对于样本集合Sk中任意的样本
Figure FDA0003982326960000021
代入Kriging模型gK(x)，得到对应的一个代理功能函数值
Figure FDA00039823269600000221
每个代理功能函数值
Figure FDA0003982326960000022
均服从代理功能函数均值为
Figure FDA0003982326960000023
标准差为
Figure FDA0003982326960000024
的高斯分布：
Figure FDA0003982326960000025
因此，对于任意的样本
Figure FDA0003982326960000026
得到代理功能函数均值为
Figure FDA0003982326960000027
代理功能函数标准差
Figure FDA0003982326960000028
步骤7，根据下式，计算样本集合Sk中每个样本
Figure FDA0003982326960000029
对应的U学习函数值
Figure FDA00039823269600000210
Figure FDA00039823269600000211
步骤8，因此，由于样本集合Sk中共有Nk个样本，共计算得到Nk个U学习函数值，判断每个U学习函数值是否均大于等于2，即：判断下式是否成立：
Figure FDA00039823269600000212
如果成立，则代表当前Kriging模型已收敛，停止更新当前Kriging模型，得到训练好的Kriging模型，执行步骤9；
否则，根据下式从样本集合Sk选取新的样本，也就是选取具有最小U学习函数值的样本，表示为
Figure FDA00039823269600000213
Figure FDA00039823269600000214
计算样本
Figure FDA00039823269600000215
的实际功能函数值
Figure FDA00039823269600000216
将
Figure FDA00039823269600000217
和
Figure FDA00039823269600000218
分别加入训练集合T和实际功能函数值集合gt中，由此得到更新后的训练集合T和更新后的实际功能函数值集合gt，再返回步骤5；
步骤9，利用训练好的Kriging模型，计算样本集合Sk中每个样本
Figure FDA00039823269600000219
对应的代理功能函数均值
Figure FDA00039823269600000220
根据下述指示函数
Figure FDA0003982326960000031
从样本集合Sk中挑选使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本：
Figure FDA0003982326960000032
若
Figure FDA0003982326960000033
则样本
Figure FDA0003982326960000034
使第k层失效事件Fk＝{x:g(x)≤bk}发生；
样本集合Sk中，所有使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本组成的集合记为失效样本集合
Figure FDA0003982326960000035
其中，Mk为使第k层失效事件Fk＝{x:g(x)≤bk}发生的失效样本总数；
根据下式计算条件失效概率Pk的估算值
Figure FDA0003982326960000036
Figure FDA0003982326960000037
其中：
Figure FDA0003982326960000038
表示第k层失效事件Fk发生概率的估算值；
Figure FDA0003982326960000039
表示第k-1层失效事件Fk-1发生的条件下，第k层失效事件Fk发生概率的估算值；
P{Fk-1}表示第k-1层失效事件Fk-1发生概率；
Figure FDA00039823269600000310
表示重要抽样密度函数；
根据下式计算条件失效概率Pk的估算值
Figure FDA00039823269600000311
对
Figure FDA00039823269600000312
的偏导数
Figure FDA00039823269600000313
Figure FDA00039823269600000314
其中：
Figure FDA00039823269600000315
代表n维随机向量x＝{x1,x2,…,xn}T中的第r维自变量xr的第s个分布参数；
Pq代表第q层失效事件Fq发生的概率，q＝1,2,...,k-1；
步骤10，判断bk是否等于0，如果bk＝0，则循环过程结束，执行步骤15；否则，执行步骤11；
步骤11，令k＝k+1，在失效样本集合RSk-1中，选出联合概率密度函数值最大的样本，作为重要抽样密度函数
Figure FDA0003982326960000041
的抽样中心；
然后，根据n维随机向量x的重要抽样密度函数
Figure FDA0003982326960000042
以联合概率密度函数值最大的样本为抽样中心，产生Nk个样本构建样本集合
Figure FDA0003982326960000043
步骤12，利用训练好的Kriging模型，计算样本集合Sk中每个样本
Figure FDA0003982326960000044
对应的代理功能函数值
Figure FDA0003982326960000045
和代理功能函数均值
Figure FDA0003982326960000046
根据下述指示函数
Figure FDA0003982326960000047
从样本集合Sk中挑选使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本：
Figure FDA0003982326960000048
样本集合Sk中，所有使第k层失效事件Fk＝{x:g(x)≤bk}发生的样本组成的集合记为失效样本集合
Figure FDA0003982326960000049
其中，Mk为使第k层失效事件Fk＝{x:g(x)≤bk}发生的失效样本总数；
步骤13，利用训练好的Kriging模型，计算失效样本集合RSk中每个样本对应的代理功能函数均值；然后，将失效样本集合RSk中每个样本对应的代理功能函数均值从大到小排序，
Figure FDA00039823269600000410
取第
Figure FDA00039823269600000411
个代理功能函数均值作为第k层失效事件Fk＝{x:g(x)≤bk}的临界值bk，即
Figure FDA00039823269600000412
步骤14，判断bk是否小于0；如果小于，则令bk＝0，然后返回步骤7；如果不小于，直接返回步骤7；
步骤15，将至此得到的每一层失效事件对应的条件失效概率的估算值做乘法运算，得到最终的涡轮机匣失效概率的估计值
Figure FDA0003982326960000051
具体的，令m＝k，根据下式计算涡轮机匣失效概率的估计值
Figure FDA0003982326960000052
Figure FDA0003982326960000053
通过下式计算得到最终的n维随机向量x＝{x1,x2,…,xn}T中的第r维自变量xr的第s个分布参数的可靠性局部灵敏度指标
Figure FDA0003982326960000054
Figure FDA0003982326960000055
其中，影响涡轮机匣灵敏度的n维随机向量中，各分量包括高温区域温度、流动冷却区域温度、冲击冷却区域温度、前段机匣长度及机匣内半径；
其中，还包括：
步骤16，最后一层失效事件即为第m层失效事件Fm，使第m层失效事件Fm发生的样本组成的失效样本集合为RSm，设失效样本集合RSm中失效样本总数为NC；
将失效样本集合RSm中的所有失效样本放入矩阵C，则矩阵C中的样本表示为
Figure FDA0003982326960000056
其中，b＝1,2,...,NC；
步骤17，初始化
Figure FDA0003982326960000057
并令β＝1；
步骤18，构建并计算根据重要抽样密度函数产生的样本
Figure FDA0003982326960000058
向目标分布中第β+1状态zβ+1转化的接受概率
Figure FDA0003982326960000059
为：
Figure FDA00039823269600000510
其中：
Figure FDA00039823269600000511
为样本
Figure FDA00039823269600000512
的联合概率密度函数；
Figure FDA00039823269600000513
为第β状态zβ的重要抽样密度函数；
fX(zβ)为第β状态zβ的联合概率密度函数；
Figure FDA0003982326960000061
为样本
Figure FDA0003982326960000062
的重要抽样密度函数；
步骤19，根据以下判断式，判断zβ+1是否接受
Figure FDA0003982326960000063
Figure FDA0003982326960000064
其中，random[0,1]表示0与1之间的随机数；
即：如果
Figure FDA0003982326960000065
则将
Figure FDA0003982326960000066
赋值给zβ+1；否则，将zβ赋值给zβ+1；
步骤20，令β＝β+1，判断β是否大于NC，如果大于，则执行步骤21；否则，执行步骤18；
步骤21，由此得到NC个zβ，分别为
Figure FDA0003982326960000067
将
Figure FDA0003982326960000068
放入到矩阵B，矩阵B中样本服从输入变量的条件概率密度函数fX(x|F)；
步骤22，根据矩阵B中的样本估计失效域内xr条件概率密度函数
Figure FDA0003982326960000069
再结合步骤15得到的涡轮机匣失效概率的估计值
Figure FDA00039823269600000610
以及第r维自变量xr的联合概率密度函数
Figure FDA00039823269600000611
代入下式得到输入变量xr的可靠性全局灵敏度指标Sr：
Figure FDA00039823269600000612
其中：
E[·]表示期望运算。
CN202010407601.XA 2020-05-14 2020-05-14 基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法 Active CN111783236B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010407601.XA CN111783236B (zh) 2020-05-14 2020-05-14 基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010407601.XA CN111783236B (zh) 2020-05-14 2020-05-14 基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法

Publications (2)

Publication Number Publication Date
CN111783236A CN111783236A (zh) 2020-10-16
CN111783236B true CN111783236B (zh) 2023-03-14

Family

ID=72754091

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010407601.XA Active CN111783236B (zh) 2020-05-14 2020-05-14 基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法

Country Status (1)

Country Link
CN (1) CN111783236B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112528511B (zh) * 2020-12-17 2023-12-12 中国矿业大学 一种刮板输送机链轮磨损可靠性评估方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB0916094D0 (en) * 2007-03-14 2009-10-28 Halliburton Energy Serv Inc Neural-network based surrogate model construction methods and applications thereof
CN108287808A (zh) * 2018-01-29 2018-07-17 华中科技大学 一种面向结构可靠性分析的通用动态追踪序列采样方法
CN108491284A (zh) * 2018-02-13 2018-09-04 西北工业大学 多失效模式复杂机构可靠性及全局灵敏度分析方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9740917B2 (en) * 2012-09-07 2017-08-22 Stone Lock Global, Inc. Biometric identification systems and methods
CN107563067A (zh) * 2017-09-06 2018-01-09 电子科技大学 基于自适应代理模型的结构可靠性分析方法
CN108304679A (zh) * 2018-03-07 2018-07-20 西北工业大学 一种自适应可靠性分析方法
CN109257749B (zh) * 2018-09-05 2021-05-18 中国人民解放军陆军工程大学 面向动态拓扑的无线传感器网络自适应分层入侵检测方法
CN110532723B (zh) * 2019-09-06 2021-06-22 北京航空航天大学 一种基于egra的涡轮盘多失效模式可靠性优化方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB0916094D0 (en) * 2007-03-14 2009-10-28 Halliburton Energy Serv Inc Neural-network based surrogate model construction methods and applications thereof
GB2462380A (en) * 2007-03-14 2010-02-10 Halliburton Energy Serv Inc Neural-network based surrogate model construction methods and applications thereof
CN108287808A (zh) * 2018-01-29 2018-07-17 华中科技大学 一种面向结构可靠性分析的通用动态追踪序列采样方法
CN108491284A (zh) * 2018-02-13 2018-09-04 西北工业大学 多失效模式复杂机构可靠性及全局灵敏度分析方法

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
An efficient Kriging-based subset simulation method for hybrid reliability analysis under random and interval variables with small failure probability;Mi Xiao 等;《Structural and Multidisciplinary Optimization》;20190118;第59卷;第2077-2092页 *
可靠性全局灵敏度求解的Meta-IS-AK算法;周苏婷 等;《航空学报http://kns.cnki.net/kcms/detail/11.1929.V.20191010.0931.002.html》;20191010;第1-11页 *
基于Kriging模型的子集模拟方法;闫阳 等;《计算机仿真》;20130215;第30卷(第02期);第217-220页 *
基于主动学习Kriging模型和子集模拟的可靠度分析;黄晓旭等;《应用力学学报》;20161031;第33卷(第05期);第866-871页 *
基于子集模拟和重要抽样的可靠性灵敏度分析方法;宋述芳 等;《力学学报》;20080930;第40卷(第5期);第654-662页 *

Also Published As

Publication number Publication date
CN111783236A (zh) 2020-10-16

Similar Documents

Publication Publication Date Title
CN111814956B (zh) 一种基于多维度二次特征提取的多任务学习的空气质量预测方法
CN112580228B (zh) 一种涡扇发动机混合结构风扇叶片结构设计优化方法
CN111597631B (zh) 基于自适应代理模型的汽车风阻系数优化方法
CN105893669A (zh) 一种基于数据挖掘的全局仿真性能预测方法
CN114741946B (zh) 一种航空发动机典型机构多构件失效相关的可靠性仿真方法
CN111797535A (zh) 一种面向多种代理模型的结构可靠性分析自适应加点方法
CN113504768B (zh) 面向装配质量预测的高精密产品数字孪生可计算性方法
CN117494584B (zh) 基于神经网络对抗迁移学习的高维可靠性设计优化方法
CN111783236B (zh) 基于自适应模型和子集模拟的涡轮机匣灵敏度分析方法
CN114254915A (zh) 轴类零件全流程加工质量合格状态决策及优化方法
CN107622144A (zh) 基于序贯方法的区间不确定性条件下多学科可靠性优化设计方法
Gal et al. Using deep neural networks and derivative free optimization to accelerate coverage closure
CN118094977A (zh) 一种新能源汽车的翼子板结构强度评价方法及系统
CN106874561B (zh) 基于牛顿迭代的多学科不确定性传播分析方法
CN116296239A (zh) 基于卷积神经网络的风洞试验方法、装置、设备和介质
CN114417653B (zh) 基于不确定性的卫星组件电路板热可靠性分析方法
CN110826182A (zh) 一种基于顶点法和序贯优化策略的飞行器结构气动弹性鲁棒设计方法
CN114329805A (zh) 基于自适应代理模型的连接件多学科协同设计优化方法
CN114492112A (zh) 基于采样分位数回归的卫星框架结构可靠性分析方法
CN114880792A (zh) 一种基于形变预测的全方位多角度优化方法
CN114996863A (zh) 基于特征提取的涡扇发动机t-s模糊建模方法
CN114330104A (zh) 一种航天器在轨空间环境异常风险预测方法
Antonova et al. Multicriteria Selection of Algorithms in Grid Methods for the Study of Large-Scale Systems
CN116644603B (zh) 基于区域划分自适应采样的复杂装备可靠性优化设计方法
CN116595883B (zh) 数值反应堆实时在线系统状态修正方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant