CN105512404B - 基于混沌多项式展开的时变可靠性全局灵敏度分析方法 - Google Patents
基于混沌多项式展开的时变可靠性全局灵敏度分析方法 Download PDFInfo
- Publication number
- CN105512404B CN105512404B CN201510923656.5A CN201510923656A CN105512404B CN 105512404 B CN105512404 B CN 105512404B CN 201510923656 A CN201510923656 A CN 201510923656A CN 105512404 B CN105512404 B CN 105512404B
- Authority
- CN
- China
- Prior art keywords
- pce
- reliability
- time
- distribution
- performance
- 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.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 34
- 238000010206 sensitivity analysis Methods 0.000 title claims abstract description 20
- 238000006731 degradation reaction Methods 0.000 claims abstract description 20
- 230000035945 sensitivity Effects 0.000 claims abstract description 10
- 230000007246 mechanism Effects 0.000 claims abstract description 7
- 230000008569 process Effects 0.000 claims abstract description 6
- 238000013461 design Methods 0.000 claims description 66
- 238000009826 distribution Methods 0.000 claims description 64
- 238000012360 testing method Methods 0.000 claims description 36
- 238000004364 calculation method Methods 0.000 claims description 21
- 230000000739 chaotic effect Effects 0.000 claims description 13
- 238000004088 simulation Methods 0.000 claims description 13
- 239000000758 substrate Substances 0.000 claims description 12
- 241001463014 Chazara briseis Species 0.000 claims description 11
- 238000012614 Monte-Carlo sampling Methods 0.000 claims description 9
- 238000005070 sampling Methods 0.000 claims description 5
- 230000008859 change Effects 0.000 claims description 4
- 230000003993 interaction Effects 0.000 claims description 4
- 238000012417 linear regression Methods 0.000 claims description 4
- 238000005215 recombination Methods 0.000 claims description 4
- 238000009827 uniform distribution Methods 0.000 claims description 4
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims description 3
- 230000004044 response Effects 0.000 claims description 3
- 102100029469 WD repeat and HMG-box DNA-binding protein 1 Human genes 0.000 claims description 2
- 101710097421 WD repeat and HMG-box DNA-binding protein 1 Proteins 0.000 claims description 2
- 238000000354 decomposition reaction Methods 0.000 claims description 2
- 239000011159 matrix material Substances 0.000 claims description 2
- 239000011541 reaction mixture Substances 0.000 claims description 2
- 230000006798 recombination Effects 0.000 claims description 2
- 230000008521 reorganization Effects 0.000 claims description 2
- 230000009466 transformation Effects 0.000 claims 1
- 230000003412 degenerative effect Effects 0.000 abstract 1
- 230000000694 effects Effects 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 4
- 230000015556 catabolic process Effects 0.000 description 3
- 238000006243 chemical reaction Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 230000010354 integration Effects 0.000 description 2
- 238000005457 optimization Methods 0.000 description 2
- IQVNEKKDSLOHHK-FNCQTZNRSA-N (E,E)-hydramethylnon Chemical compound N1CC(C)(C)CNC1=NN=C(/C=C/C=1C=CC(=CC=1)C(F)(F)F)\C=C\C1=CC=C(C(F)(F)F)C=C1 IQVNEKKDSLOHHK-FNCQTZNRSA-N 0.000 description 1
- 238000012351 Integrated analysis Methods 0.000 description 1
- 238000000342 Monte Carlo simulation Methods 0.000 description 1
- 241000286209 Phasianidae Species 0.000 description 1
- 102100035115 Testin Human genes 0.000 description 1
- 101710070533 Testin Proteins 0.000 description 1
- 241000082085 Verticillium <Phyllachorales> Species 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 238000012938 design process Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000013076 uncertainty analysis Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/30—Circuit design
- G06F30/36—Circuit design at the analogue level
- G06F30/367—Design verification, e.g. using simulation, simulation program with integrated circuit emphasis [SPICE], direct methods or relaxation methods
Landscapes
- Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Microelectronics & Electronic Packaging (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Complex Calculations (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,包括以下步骤:一,在基于故障机理的可靠性与性能一体化模型的基础上,计算产品退化过程离散时刻的性能混沌多项式展开;二,基于正交实验设计方案,计算产品退化过程离散时刻的可靠性混沌多项式展开;三,根据移动最小二乘原理,计算产品退化过程时变可靠性混沌多项式展开;四,将上述时变可靠性混沌多项式展开的系数按照索博尔(Sobol’)分解重组,计算全局灵敏度Sobol’指标。该方法能够高效地计算复杂工程模型在退化过程中可靠性全局灵敏度分析问题,具有精度较高、适用范围广等特点。
Description
(一)技术领域:
本发明提供一种基于混沌多项式展开的可靠性全局灵敏度分析方法,注重于解决复杂工程模型的可靠全局灵敏度分析问题,属于***优化设计与可靠性设计的交叉技术领域。
(二)背景技术:
***可靠性与性能一体化技术是一项在产品设计过程中考虑故障和环境扰动,采用可靠性优化和不确定性分析等方法来实现可靠性与性能综合设计分析的新技术。由于产品单元退化机理在时间轴上的发展及其在不同单元间的传播,产品性能输出随着时间的增长表现为一个逐渐退化的随机过程,时变可靠性作为可实时反映***可靠性及质量特性的指标,逐渐得到一体化分析设计人员的重视。曾声奎、陈云霞等人对时变可靠性建模仿真技术开展了大量的工作。
灵敏度分析是在一项产品设计阶段研究***设计变量对产品特性影响的方法。它包括局部灵敏度分析和全局灵敏度分析。局部灵敏度分析主要考虑由单个设计变量在标称值处的线性梯度所引起产品性能变化的比率;而全局灵敏度分析则能够衡量所有设计变量的不确定性对产品性能的综合作用,为设计人员优化设计方案提供一种有效手段。1990年,Sobol’提出一种基于方差的全局灵敏度分析指标(Sobol’指标),它能够快速简便地计算出设计变量高阶交叉影响项,在工程界得到广泛应用。但是该方法需要利用蒙特卡洛积分,面对具有复杂耦合关系和随机过程的复杂工程模型,计算负担很大。后来Saltelli、吕震宙等人将Sobol’指标向可靠性领域进行扩展,但是为获得某些小失效概率的结果,往往也需要大量蒙特卡洛仿真作为支撑,进一步增加了计算成本。混沌多项式展开(Polynomial chaosexpanse,PCE)采用相互正交的多项式作为基底,它将***性能输出投影到概率空间,用标准随机变量来表示设计变量不确定性,从而建立起产品性能与设计变量之间关系的代理模型,具有较高的精度且较小的计算量。2007年,Bruno Sudret首次提出将PCE按照Sobol’分解、重组,然后利用重组后的PCE系数直接求解Sobol’指标,引起学术界和工程界的广泛关注。但是这些研究的对象没有涉及产品的退化过程,不能反映产品在整个寿命周期内的真实情况。
如果能提供一种针对时变可靠性灵敏度分析方法,即使对于复杂工程模型,也能就其退化过程进行高效、精确的分析,将显著地提高灵敏度分析的工程实用性,扩大其应用范围。
(三)发明内容:
(1)目的:
针对上述问题,本发明提出了一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,为可靠性与性能一体化设计提供一种客观、合理的分析技术手段。
(2)技术方案:
本发明是一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,该方法包括如下四个步骤:
步骤一:根据PCE与基于故障机理的可靠性与性能一体化仿真模型之间的关系,计算产品退化过程中各离散时刻的性能PCE,用来表示性能输出与关键设计变量的关系;
步骤二:根据蒙特卡洛抽样原理,计算各离散时刻的可靠性PCE,用来描述可靠性与关键设计变量的关系;
步骤三:根据移动最小二乘原理,计算退化过程的时变可靠性PCE,其中,所述的时变是指PCE系数随时间变化;
步骤四:在上述时变可靠性PCE的基础上,计算时变全局灵敏度的Sobol’s指标。
其中,在步骤一中所述的“计算产品退化过程中各离散时刻的性能PCE”,其计算步骤可细分为两步。
步骤1)确定性能PCE的基底和阶数
首先根据关键设计变量的分布类型选择相应的标准正交基底,如下列表1所示,其中,标准正交基底是相应标准随机变量的多项式函数。注意,当模型设计变量存在多种分布类型时,只能根据实际情况选择一种主要的正交基底形式。然后根据工程经验确定阶数p,并分别计算p=k和p=k+1两种情况下PCE及其相应的误差估计值,如两者误差差别不大,则可将阶数最终确定为k+1,否则再计算k+2,直至相邻两阶PCE的误差估计值基本一致,取最高阶PCE作为最终结果。这里k为大于1的正整数。
表1混沌多项式类型及对应的随机变量
表中,Hermite表示厄米多项式正交基底,N(0,1)表示均值为0、方差为1的正态分布;Legendre表示勒让德多项式正交基底,U[-1,1]表示上下界分别为1和-1的均匀分布;Laguerre表示拉盖尔多项式正交基底;Generalized Laguerre表示广义拉盖尔多项式正交基底,Γ(α+1,1)表示分布参数为α+1和1的伽马分布。上述正交基底的具体形式可以通过公开文献获得,这里不再赘述。
由此可以得到时刻t性能输出y(ξ;t)的PCE的一般形式为:
式(1)中n是关键设计变量个数;p是PCE展开阶数;N为PCE所包含系数的总个数,它由n和p来确定,即:
是PCE在各离散时刻t的系数,其中是标准正交多项式基底,ξi~N(0,1),i=1,…,n。
步骤2)计算性能PCE的系数
根据待定系数个数N-1,在比PCE阶数度高一阶的标准正交多项式基底的根中选择合适数量的配点其中,ξk是PCE基底对应的标准随机变量的一组样本点;然后将配点转化变为仿真模型的变量输入,并求解***响应y(ξk;t)。以Hermit正交基底为例,由于Hermit正交基底对已的标准随机变量ξh服从标准正态分布N(0,1),则配点与仿真模型设计变量输入之间的转换关系如下列表2所示。
表2常见分布与标准正态分布关系
其中,是高斯误差函数,y=exp(x)=ex是指数函数。
再将配点代入到式(1)中多项式部分,就能得到多组计算样本
利用多元线性回归求取系数
即:
式中M为配点数目,为保证系数矩阵的条件数,样本数量要求不小于未知系数个数的两倍(M≥2N),同时还要增加一种配点为零的选择方案,并且在布置配点时,应尽量关于原点对称。
按上述步骤不断重复,获得各离散时刻的性能PCE。
式中,Mt为离散时刻数目。
其中,在步骤二中所述的“计算各离散时刻的可靠性PCE”,其目的是获得关键设计变量的分布参数与产品可靠度之间的关系。其中,可以认为关键设计变量的分布参数在取值范围内的任何取值都是等可能的,即关键设计变量分布参数在取值范围内服从均匀分布。根据表1选择可靠性PCE的基底为Legendre基底;而可靠性PCE的阶数确定方法与性能PCE的阶数确定方法相同,但根据实际经验,一般和性能PCE同阶数。各离散时刻可靠性PCE的形式为:
式中, 是Legendre正交多项式,是可靠性PCE系数,Nr是可靠性PCE系数的总个数。其具体计算步骤如下:
步骤1)离散时刻可靠性PCE配点设计
根据工程经验,确定关键设计变量的分布参数的取值范围。然后将各关键设计变量的分布参数水平均划分为m等,进行随机搭配,得到可靠性PCE的试验方案。其中,总试验方案数目M′=mn,且M′不小于2N。然后将配点映射到区间[-1,1]上,得到均匀分布标准随机变量与试验方案的关系。一般地,令关键设计变量的分布参数μi服从区间[a,b]的均匀分布,即,μi~U[a,b];均匀分布标准随机变量ζi~U[-1,1]。则有
步骤2)各配点对应可靠度的计算
依据上述试验方案进行蒙特卡洛抽样,以获得相应的可靠度。该步骤需要利用到各离散时刻性能PCE,即y(ξ;t)=PCEt。该计算过程包括以下两个步骤:
I)建立试验方案中各个配点的关键设计变量的分布参数与步骤一中配点的关键设计变量的初始分布参数以及均匀分布标准随机变量之间的关系,得到反映试验方案的性能PCE修正随机变量和均匀分布标准随机变量为后续计算提供抽样分布。
例如,在以Hermit正交多项式为基底的PCE中,设某设计变量的初始分布参数为即均值为μi、方差为的正态分布,对应离散时刻性能PCE的标准正态变量ξi=(X-μ)/σ;在第k个试验方案中则该试验方案对应的性能PCE修正随机变量为:
又如,在以Hermit正交多项式为基底的PCE中,设某设计变量的初始分布参数为Xi~U[ai,bi],即上下限为ai和bi的均匀分布;在第k个试验方案中则该试验方案对应的性能PCE修正随机变量为:
式中,
此时,根据式(6),第k个试验方案中对应的均匀分布标准随机变量值为:
II)利用各离散时刻的性能PCE和各试验方案的修正随机变量,进行蒙特卡洛抽样。修正随机变量反映了设计变量分布参数变化对失效概率的影响。例如,根据式(7)中对配点进行抽样时,PCE中随机变量的分布不再服从标准正态分布,而是服从参数为的正态分布。
然后根据抽样的结果计算可靠度:
式(10)中,y(ξ′(k);t|MC)表示在t时刻对第k个试验方案的性能PCE修正随机变量进行蒙特卡洛抽样,Nsim是抽样次数,num(y(ξ′;t)∈Ωf)表示仿真结果落入失效域的数目,Rk(t)表示对应于第k个试验方案的可靠度。注意,与均匀分布标准随机变量对应。
步骤3)计算离散时刻可靠性PCE的系数
在各离散时刻,根据上述步骤1)和步骤2),对各试验方案对应的修正随机变量进行蒙特卡罗抽样计算可靠度的同时,将该试验方案对应的均匀分布标准随机变量代入(5)中的多项式部分,可以得到多组样本点然后类似步骤二中求解各离散时刻性能PCE系数的方法,利用多元线性回归求取可靠性PCE的系数。
式(11)中,为各试验方案对应的均匀分布标准随机变量值。
按上述步骤不断重复,获得各离散时刻的可靠性
式中,Mt为离散时间点数目。
其中,在步骤三中所述的“退化过程的时变可靠性PCE”,其计算过程如下。
这里的“时变”是指可靠性PCE的所有系数分别随时间变化,即展开系数分别是时间的函数。该函数利用已知的各离散时刻的值,通过移动最小二乘拟合,即
式(13)中,p(t)是移动最小二乘拟合基函数,一般可取p(t)=(1,t,t2)T,bj是移动最小二乘拟合系数。其中移动最小二乘的计算原理可从公开资料中获得,其计算结果为:
式(14)中,Mt为离散时间点数目;w(t)=w(||t-ti||)是权函数,比较常用的有样条函数,径向基函数,高斯函数等。代入式(13)得到时变可靠性PCE:
式中,为通过MLS估计得到的时变可靠性PCE系数。
其中,在步骤四中所述的“计算时变全局灵敏度的Sobol’指标”,其具体计算步骤如下:
首先根据时变可靠性PCE计算时变可靠性方差D(R(ζ,t)):
式中,
然后将式(15)的系数按照Sobol’分解的形式进行重组,通过重组后的PCE的系数直接计算,具体步骤如下:
步骤1)PCE系数重组
系数重组按照如下规则进行:
(1)被加项不展开,将单独的各个变量的一次项、二次项……放在一堆。
(2)被加项不展开,将具有两个变量(交互作用)的被加项的一次项、二次项……放在一堆。
(3)被加项不展开,将具有三个变量(交互作用)的被加项的一次项、二次项……放在一堆。
重组后数学表达式为:
式(17)中,α=(α1,…,αn)是一个整数序列序列,是的具体值,满足
步骤2)Sobol’s指标
然后利用PCE基底的正交性等优良性质,可以直接得到Sobol’s指标:
(3)本发明的优点及功效:
该方法能够高效地计算复杂工程模型在退化过程中可靠性全局灵敏度分析问题,具有精度较高、适用范围广等特点。
(四)附图说明:
图1是连杆机构结构示意图。
图2本发明所述方法流程图。
图3时变可靠性PCE系数曲线。
图4各关键设计变量一阶时变可靠性Sobol’指标曲线。
图5各关键设计变量总的时变可靠性Sobol’指标曲线。
图中符号说明如下:
(Xi,Yi) 第i点的横纵坐标
PCEt t时刻性能混沌多项式展开
RPCEt t时刻可靠性混沌多项式展开
Rtk t时刻第k个试验方案对应的可靠度值
时变可靠性混沌多项式展开系数
Si 第i个关键设计变量一阶的可靠性Sobol’指标
STi 第i个关键设计变量总的可靠性Sobol’指标
(五)具体实施方式:
见图1-5,下面将结合某连杆装置对本发明作进一步的详细说明。该收放装置包括滑动杆、后主推臂、后支撑臂、铰链、前主推臂、前支撑臂、液压装置及负载构成,如下列表3所示,主要设计变量数据如表3所示。故障判据为最大摩擦阻力Max-Force不超过了液压装置所能输出力的最大值。研究目标是分析2000次收放过程中,连杆装置的时变可靠性对设计变量均值的灵敏度。
表3主要设计变量参数值
符号 | 含义 | 名义值 | 符号 | 含义 | 名义值 |
XA | A点的横坐标 | -910.49mm | Gd | 负载重量 | 2600N |
XB | B点的横坐标 | -150mm | W1 | 铰链上连杆宽 | 50 |
YC/YD | C点的纵坐标 | -264.67mm | D1 | 铰链上连杆深 | 20 |
Rlug | 铰链耳片半径 | 15mm | R1 | 作动筒半径 | 50 |
L1 | 铰链上连杆长 | 70mm | R2 | 作动筒滑动杆半径 | 25 |
L2 | 作动筒长 | 500mm | W2/3 | 前/后支撑臂宽度 | 50 |
L3 | 作动筒滑动杆 | 1580.16mm | D2 | 前/后支撑臂深度 | 20 |
L4/5 | 前/后支撑臂 | 50mm | W3/4 | 前/后主推臂宽度 | 50 |
L6/7 | 前/后主推臂长 | 260mm | D3/4 | 前/后主推臂深度 | 20 |
本发明一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,如下列表4所示,其中已经建立基于故障机理的可靠性与性能一体化仿真模型。它包括如下步骤:
步骤一,计算产品退化过程中各离散时刻的性能PCE
根据工程经验,在本案例中关键设计变量为A点横坐标XA、C点纵坐标YC、负载重量Gd以及铰链半径Rlug四项,它们的不确定性输入参数的概率分布如表4所示:
表4各输入变量的概率分布
符号 | 描述 | 设计均值/ai | 标准差/bi | 分布类型 |
XA | A点横坐标 | -910.49 | 3 | 正态分布 |
YC | C点纵坐标 | -264.67 | 1 | 正态分布 |
Gd | 负载重量 | 2600 | 10 | 正态分布 |
Rlug | 铰链半径 | 14.9 | 15.1 | 均匀分布 |
根据表4,选择Hermit多项式作为PCE的基底,其第零、一、二、三、四阶Hermit多项式的形式分别为1,ξ,ξ2-1,ξ3-3ξ和ξ4-6ξ2+3。相应的标准正交随机变量为标准正态随机变量。然后选择阶数p=2和p=3时分别构造PCE,其对应的系数个数N分别为15和35。
然后设计配点以2阶PCE为例,其配点的坐标从三阶Hermit多项式ξ3-3ξ的根中选,每个配点是这三个根的组合,配点方案的总数目M至少为30。由于根中已经包含0点,无需再额外增加配点为零的选择方案。考虑到配点尽量关于原点对称,其设计方案如下列表5所示。
表5二阶性能PCE配点方案
根据设计变量与标准正态随机变量之间的转化关系如下列表6所示,其中ξi~N(0,1),i=1,2,3,4,将然后将配点转化变为仿真模型的变量输入。以配点方案10为例,即仿真模型输入为(905.2938,-262.9379,2600,14.9)T。
表6设计变量与标准变量之间转换关系
根据各配点所对应的仿真模型输入,可以计算得到各离散时刻***的响应输出其中Mt为离散时刻的数目,对本案例而言,将一次收放结束定义为离散时间间隔,即Mt=2000。再将配点代入到式(1)中多项式部分得到回归样本其中,以二阶PCE的配点方案10在初始时刻的回归样本为为例,具体形式为
即
利用式(3),可以求出各离散时刻性能PCE。以在初始时刻三阶PCE为例,其回归方程组为
计算结果为:
同时,由2阶PCE、3阶PCE以及MCS得到结果的误差在0.1%左右,如下列表7所示,满足精度要求,故最终选择PCE的阶数为p=3,基底个数即PCE系数总个数为35。
表7各阶PCE对最大力估计值的对比
方法 | 调用模型次数 | 均值 | 方差(×105) |
MCS | 10000 | 8074.602 | 6.3366 |
2阶PCE | 30 | 8038.331 | 5.9674 |
相对误差 | --- | 0.0045 | 0.0583 |
3阶PCE | 70 | 8056.497 | 6.2845 |
相对误差 | --- | 0.0022 | 0.0082 |
步骤二,计算产品退化过程中各离散时刻的可靠性PCE。
可靠性PCE的基底为Legendre正交多项式,其第零、一、二、三、四阶Hermit多项式的形式分别为1,ζ,和相应的标准正交随机变量为均匀分布标准随机变量ζ~U[-1,1],可靠性PCE的阶数取3阶。
首先构造试验方案。根据工程经验,设计变量的设计均值变化范围为设计均值的±1%。将各实验因素(关键设计变量均值)在其变化范围内等分为3个水平,进行随机搭配,得到共计34=81个试验方案,达到了PCE基底个数的2倍,满足要求;其次,将各试验方案映射到[-1,1],获得其对应的均匀分布标准随机变量值;同时计算各试验方案对应性能PCE修正随机变量值,如下列表8所示。
表8试验方案及对应均匀分布标准随机变量值
然后,计算试验方案对应的性能PCE修正随机变量,以试验方案1为例,
关键设计变量均值的取值为(μ1,μ2,μ3,μ4)T=(905.94,-263.35,2587.00,14.83)T,均值增量分别为以及则性能PCE修正随机变量分别为:
通过对第k个试验方案对应的性能PCE修正随机变量进行蒙特卡洛抽样,结合故障判据,按照式(10)计算一系列的可靠度数值其中Nsim取10000。再将对应的均匀分布标准随机变量值代入到式(5)中多项式部分得到回归样本
其中,以试验方案1在初始时刻的回归样本为为例,回归样本的具体形式为:
即
最后,根据回归方程组计算可靠性PCE。以初始时刻为例,其余不再赘述。
得可靠性PCE:
步骤三,计算产品退化过程时变可靠性PCE。
选择三次样条函数作为移动最小二乘的权重函数,即:
式(26)中,s=||t-ti||,ti是离散时刻。
然后结合上述离散时刻可靠性PCE的系数,根据式(14),计算时变系数表达式。以式(15)中系数为例,其计算结果如图3所示。图中,离散的点是经过各离散时刻PCE的系数,曲线是利用移动最小二乘原理估计得到的曲线。
步骤四,计算时变全局灵敏度的Sobol’指标。
首先,根据式(16)求出时变可靠性PCE计算时变可靠性方差D(R(ζ,t))。
然后,将上述时变可靠性PCE按照式(17)重新展开。以初始时刻为例时变可靠性PCE为例,其系数重组后形式如式(27)所示。
最后,根据式(19)和式(20)可直接求出Sobol’指标,如图4和图5所示。
Claims (4)
1.一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,其特征在于:该方法包括如下四个步骤:
步骤一:根据PCE与基于故障机理的可靠性与性能一体化仿真模型之间的关系,计算产品退化过程中各离散时刻的性能PCE,用来表示性能输出与关键设计变量的关系;
步骤二:根据蒙特卡洛抽样原理,计算各离散时刻的可靠性PCE,用来描述可靠性与关键设计变量的关系;
步骤三:根据移动最小二乘原理,计算退化过程的时变可靠性PCE,其中,所述的时变是指PCE系数随时间变化;
步骤四:在上述时变可靠性PCE的基础上,计算时变全局灵敏度的Sobol’s指标;
其中,在步骤一中所述的“计算产品退化过程中各离散时刻的性能PCE”,计算步骤细分为两步:
步骤1.1确定性能PCE的基底和阶数
首先,根据关键设计变量的分布类型选择相应的标准正交基底,如下所列;
表1所示,标准正交基底是相应标准随机变量的多项式函数,然后确定阶数p,并分别计算p=k和p=k+1两种情况下PCE及其相应的误差估计值,如两者误差差别不大,则将阶数最终确定为k+1,否则再计算k+2,直至相邻两阶PCE的误差估计值一致,取最高阶PCE作为最终结果,这里k为大于1的正整数;
表1 混沌多项式类型及对应的随机变量
表中,Hermite表示厄米多项式正交基底,N(0,1)表示均值为0、方差为1的正态分布;Legendre表示勒让德多项式正交基底,U[-1,1]表示上下界分别为1和-1的均匀分布;Laguerre表示拉盖尔多项式正交基底;Generalized Laguerre表示广义拉盖尔多项式正交基底,Γ(α+1,1)表示分布参数为α+1和1的伽马分布;上述正交基底的具体形式能通过公开文献获得;
由此能得到时刻t性能输出y(ξ;t)的PCE的形式为:
式(2)中n是关键设计变量个数;p是PCE展开阶数;N为PCE所包含系数的总个数,它由n和p来确定,即:
是PCE在各离散时刻t的系数,其中,是标准正交多项式基底,ξi~N(0,1),i=1,…,n;
步骤1.2计算性能PCE的系数
根据待定系数个数N-1,在比PCE阶数度高一阶的标准正交多项式基底的根中选择合适数量的配点其中,ξk是PCE基底对应的标准随机变量的一组样本点;然后将配点转化变为仿真模型的变量输入,并求解***响应y(ξk;t);由于Hermit正交基底对已的标准随机变量ξh服从标准正态分布N(0,1),则配点与仿真模型设计变量输入之间的转换关系如下列表2所示
表2 常见分布与标准正态分布关系
其中,是高斯误差函数,y=exp(x)=ex是指数函数;
再将配点代入到式(1)中多项式部分,
就能得到多组计算样本{ψ0(ξk),ψ1(ξk),…,ψN-1(ξk);利用多元线性回归求取系数即:
式中M为配点数目,为保证系数矩阵的条件数,样本数量要求不小于未知系数个数的两倍,即M≥2N,同时还要增加一种配点为零的选择方案,并且在布置配点时,应关于原点对称;
按上述步骤不断重复,获得各离散时刻的性能PCE;
式中,Mt为离散时刻数目。
2.根据权利要求1所述的一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,其特征在于:在步骤二中所述的“计算各离散时刻的可靠性PCE”,其目的是获得关键设计变量的分布参数与产品可靠度之间的关系;其中,认为关键设计变量的分布参数在取值范围内的任何取值都是等可能的,即关键设计变量分布参数在取值范围内服从均匀分布;
根据表1选择可靠性PCE的基底为Legendre基底;而可靠性PCE的阶数确定方法与性能PCE的阶数确定方法相同,各离散时刻可靠性PCE的形式为:
式中, 是Legendre正交多项式,是可靠性PCE系数,Nr是可靠性PCE系数的总个数;具体计算步骤如下:
步骤2.1离散时刻可靠性PCE配点设计
确定关键设计变量的分布参数的取值范围,然后将各关键设计变量的分布参数水平均划分为m等,进行随机搭配,得到可靠性PCE的试验方案;其中,总试验方案数目M′=mn,且M′不小于2N;然后将配点映射到区间[-1,1]上,得到均匀分布标准随机变量与试验方案的关系;令关键设计变量的分布参数μi服从区间[a,b]的均匀分布,即,μi~U[a,b];均匀分布标准随机变量ζi~U[-1,1],则有
步骤2.2各配点对应可靠度的计算
该步骤需要利用到各离散时刻性能PCE,即y(ξ;t)=PCEt;该计算过程包括以下两个步骤:
步骤2.21建立试验方案中各个配点的关键设计变量的分布参数与步骤一中配点的关键设计变量的初始分布参数以及均匀分布标准随机变量之间的关系,得到反映试验方案的性能PCE修正随机变量和均匀分布标准随机变量为后续计算提供抽样分布;
在以Hermit正交多项式为基底的PCE中,设某设计变量的初始分布参数为即均值为μi、方差为的正态分布,对应离散时刻性能PCE的标准正态变量ξi=(X-μ)/σ;在第k个试验方案中则该试验方案对应的性能PCE修正随机变量为:
在以Hermit正交多项式为基底的PCE中,设某设计变量的初始分布参数为Xi~U[ai,bi],即上下限为ai和bi的均匀分布;在第k个试验方案中则该试验方案对应的性能PCE修正随机变量为:
式中,
根据式(6),第k个试验方案中对应的均匀分布标准随机变量值为:
步骤2.22利用各离散时刻的性能PCE和各试验方案的修正随机变量,进行蒙特卡洛抽样;修正随机变量反映了设计变量分布参数变化对失效概率的影响,根据式(7)中对配点进行抽样时,PCE中随机变量的分布不再服从标准正态分布,而是服从参数为的正态分布;
然后计算可靠度:
式(10)中,y(ξ′(k);t|MC)表示在t时刻对第k个试验方案的性能PCE修正随机变量进行蒙特卡洛抽样,Nsim是抽样次数,num(y(ξ′;t)∈Ωf)表示仿真结果落入失效域的数目,Rk(t)表示对应于第k个试验方案的可靠度;注意,与均匀分布标准随机变量对应;
步骤2.3计算离散时刻可靠性PCE的系数
在各离散时刻,根据上述步骤2.1和步骤2.2,对各试验方案对应的修正随机变量进行蒙特卡罗抽样计算可靠度的同时,将该试验方案对应的均匀分布标准随机变量代入(5)中的多项式部分,得到多组样本点{ψL0(ζk),ψL1(ζk),…,ψL,N-1(ζk);然后根据步骤二中求解各离散时刻性能PCE系数的方法,利用多元线性回归求取可靠性PCE的系数;
式(11)中,为各试验方案对应的均匀分布标准随机变量值;
获得各离散时刻的可靠性
式中,Mt为离散时间点数目。
3.根据权利要求1所述的一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,其特征在于:在步骤三中所述的“退化过程的时变可靠性PCE”,计算过程如下:
这里的“时变”是指可靠性PCE的所有系数分别随时间变化,即展开系数分别是时间的函数,该函数利用已知的各离散时刻的值,通过移动最小二乘拟合,即
式(13)中,p(t)是移动最小二乘拟合基函数,取p(t)=(1,t,t2)T,bj是移动最小二乘拟合系数;其中,移动最小二乘的计算原理能从公开资料中获得,计算结果为:
式(14)中,Mt为离散时间点数目;w(t)=w(||t-ti||)是权函数,代入式(13)得到时变可靠性PCE:
式中,为通过MLS估计得到的时变可靠性PCE系数。
4.根据权利要求1所述的一种基于混沌多项式展开的时变可靠性全局灵敏度分析方法,其特征在于:在步骤四中所述的“计算时变全局灵敏度的Sobol’指标”,具体计算步骤如下:
首先,根据时变可靠性PCE计算时变可靠性方差D(R(ζ,t)):
式中,
然后将式(15)的系数按照Sobol’分解的形式进行重组,通过重组后的PCE的系数直接计算,具体步骤如下:
步骤4.1PCE系数重组
系数重组按照如下规则进行:
(1)被加项不展开,将单独的各个变量的一次项、二次项……放在一堆;
(2)被加项不展开,将具有两个变量交互作用的被加项的一次项、二次项……放在一堆;
(3)被加项不展开,将具有三个变量交互作用的被加项的一次项、二次项……放在一堆;
重组后数学表达式为:
式(17)中,α=(α1,…,αn)是一个整数序列序列,是的具体值,满足
步骤4.2Sobol’s指标
然后利用PCE基底的正交性,直接得到Sobol’s指标:
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510923656.5A CN105512404B (zh) | 2015-12-14 | 2015-12-14 | 基于混沌多项式展开的时变可靠性全局灵敏度分析方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510923656.5A CN105512404B (zh) | 2015-12-14 | 2015-12-14 | 基于混沌多项式展开的时变可靠性全局灵敏度分析方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105512404A CN105512404A (zh) | 2016-04-20 |
CN105512404B true CN105512404B (zh) | 2018-08-28 |
Family
ID=55720384
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201510923656.5A Expired - Fee Related CN105512404B (zh) | 2015-12-14 | 2015-12-14 | 基于混沌多项式展开的时变可靠性全局灵敏度分析方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105512404B (zh) |
Families Citing this family (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106021772B (zh) * | 2016-05-31 | 2018-10-02 | 北京航空航天大学 | 一种基于区间集合包络函数的服役结构静态等效可靠性评价方法 |
CN113673120B (zh) * | 2018-02-05 | 2024-01-30 | 西北工业大学 | 一种包括退化模型的机构可靠性建模方法 |
CN108491284A (zh) * | 2018-02-13 | 2018-09-04 | 西北工业大学 | 多失效模式复杂机构可靠性及全局灵敏度分析方法 |
CN109472048B (zh) * | 2018-09-28 | 2022-09-27 | 广西电网有限责任公司电力科学研究院 | 基于稀疏多项式混沌扩展评估智能电表结构可靠度的方法 |
CN109583111A (zh) * | 2018-12-05 | 2019-04-05 | 中国航空工业集团公司西安飞行自动控制研究所 | 一种基于任意型混沌多项式的加速度计不确定性分析方法 |
CN110941881A (zh) * | 2019-10-16 | 2020-03-31 | 北京航空航天大学 | 一种基于混沌多项式展开的混合不确定性结构疲劳寿命分析方法 |
CN111141974A (zh) * | 2019-12-26 | 2020-05-12 | 吉林大学 | 一种多导体传输线辐射敏感度分析方法 |
CN111125928B (zh) * | 2019-12-31 | 2021-02-12 | 湖南大学 | 一种概率-区间混合的加速度计偏值不确定性分析方法 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP1739582A1 (en) * | 2005-06-29 | 2007-01-03 | Siemens Aktiengesellschaft | Probabilistic design tool for optimizing a technical system |
CN102651054A (zh) * | 2012-04-10 | 2012-08-29 | 北京航空航天大学 | 一种基于贝叶斯理论的电子产品寿命模型概率化方法 |
CN103559398A (zh) * | 2013-11-04 | 2014-02-05 | 湖南大学 | 汽车盘式制动器***振动稳定性的不确定性分析方法 |
-
2015
- 2015-12-14 CN CN201510923656.5A patent/CN105512404B/zh not_active Expired - Fee Related
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP1739582A1 (en) * | 2005-06-29 | 2007-01-03 | Siemens Aktiengesellschaft | Probabilistic design tool for optimizing a technical system |
CN102651054A (zh) * | 2012-04-10 | 2012-08-29 | 北京航空航天大学 | 一种基于贝叶斯理论的电子产品寿命模型概率化方法 |
CN103559398A (zh) * | 2013-11-04 | 2014-02-05 | 湖南大学 | 汽车盘式制动器***振动稳定性的不确定性分析方法 |
Non-Patent Citations (4)
Title |
---|
《Non-intrusive Polynomial Chaos Based Mechanism》;Shaohua Du;《2014 Prognostics and System Health Management Conference (PHM-2014 Hunan)》;20140427;全文 * |
《Time-Dependent Global Sensitivity Analysis for Long-Term》;Jianbin Guo et;《2014 Prognostics and System Health Management Conference》;20140827;正文5.1,8.2,表5,表6 * |
《Time-Variant Reliability Analysis using Polynomial Chaos Expansion》;Lara Hawchar et;《12th International Conference on Applications of Statistics and Probability in Civil Engineering》;20150715;全文 * |
《动态***故障的混杂传播特征及建模方法》;郭健彬等;《***工程与电子技术》;20150115;第37卷(第1期);224-228 * |
Also Published As
Publication number | Publication date |
---|---|
CN105512404A (zh) | 2016-04-20 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105512404B (zh) | 基于混沌多项式展开的时变可靠性全局灵敏度分析方法 | |
Gu et al. | State space model identification of multirate processes with time-delay using the expectation maximization | |
Eldred et al. | Comparison of non-intrusive polynomial chaos and stochastic collocation methods for uncertainty quantification | |
Huang et al. | Two finite difference schemes for time fractional diffusion-wave equation | |
Alsuyuti et al. | Modified Galerkin algorithm for solving multitype fractional differential equations | |
Vreman et al. | COMPARISION of numerical schemes in large‐eddy simulation of the temporal mixing layer | |
Ghasemi et al. | Numerical solution of nonlinear delay differential equations of fractional order in reproducing kernel Hilbert space | |
Pulch et al. | Sensitivity analysis and model order reduction for random linear dynamical systems | |
Wang et al. | Uncertainty propagation of heat conduction problem with multiple random inputs | |
Sezer et al. | Rational Chebyshev collocation method for solving higher‐order linear ordinary differential equations | |
Bougoffa et al. | Solving Cauchy integral equations of the first kind by the Adomian decomposition method | |
Diao et al. | Mixed, componentwise condition numbers and small sample statistical condition estimation of Sylvester equations | |
Licea et al. | Analytic and numerical solutions of a Riccati differential equation with random coefficients | |
Abdullayev | Identification of the functions of response to loading for stationary systems | |
Cai et al. | Probabilistic load flow algorithms considering correlation between input random variables: A review | |
Costin et al. | Numerical study of radial basis function interpolation for data transfer across discontinuous mesh interfaces | |
Kumar et al. | Optimal input signal design for identification of interactive and ill-conditioned systems | |
Ali et al. | An inverse problem for a family of two parameters time fractional diffusion equations with nonlocal boundary conditions | |
Pekalp et al. | Power series expansions for the probability distribution, mean value and variance functions of a geometric process with gamma interarrival times | |
Yang et al. | Differential quadrature method for fractional logistic differential equation | |
Maleknejad et al. | Numerical solution of volterra functional integral equation by using cubic B‐spline scaling functions | |
Aiello et al. | Scattering asymptotic conditions in Euclidean relativistic quantum theory | |
CN107133418B (zh) | 基于交替tvd差分算法的地球流体物质平流输运模拟方法 | |
Wu et al. | Probabilistic load flow by generalized polynomial chaos method | |
Deka et al. | Tractable structure learning in radial physical flow networks |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20180828 |
|
CF01 | Termination of patent right due to non-payment of annual fee |