基于响应面法的机械可靠性设计与仿真方法研究

刘征宇

（海装驻上海地区第七军事代表室，上海 200000）

1 机械可靠性特点

随着可靠性理论的日益发展，研究对象逐渐由相关理论较为成熟的电子产品转向机械产品，机械产品可靠性与产品的结构设计、制造工艺、工作条件以及运行环境等因素密切相关，导致机械产品故障模式具有随机性与不确定性[1]。机械可靠性具有如下特点：失效机理多样化；失效模式多样化，且存在相关性；寿命分布复杂；失效及故障积累数据不足；可靠性试验费用大；模拟试验困难等[2]。

目前已经颁布并投入使用的GJB 899A-2009《可靠性鉴定与验收试验》适用对象仅限于寿命服从指数分布的电子产品，故机械产品暂时没有相关标准进行可靠性试验方面的研究。由于机械可靠性设计理论方法的不断发展，技术人员更多的是在机械结构设计阶段，提高机械结构的可靠性，设计出可靠性更高的产品，满足日益激烈的产品质量竞争环境。

2 机械可靠性设计方法

2.1 不确定性理论

安全系数法是指根据经验选取合适的安全系数，用计算应力除以材料的许用应力，当比值大于安全系数时即认为结构强度满足设计要求能够安全工作[3]。安全系数法视结构分析中各参数为定值,忽视了不确定因素产生的变异性,因而具有一定的局限性，容易造成过设计或欠设计。因此在机械结构可靠性设计中，引入概率统计学，考虑结构分析中各参数的不确定性与随机性，这种不确定性主要包括机械结构材料属性的不确定性、承受载荷的不确定性以及结构尺寸的不确定性等。

2.2 应力-强度干涉模型

A.M.Freudenthal提出的“应力—强度干涉模型”是从可靠性角度考虑机械结构的设计，影响机械结构可靠性的各种因素可以概括为“应力”和“强度”[4]。广义的应力是指包含各种影响因素产生的力，例如结构尺寸、材料属性、载荷等对结构有影响的变量；材料强度在材料力学各向同性的假设下是确定值，但实际上强度也服从特定的随机分布规律，当应力大于强度的时候，结构失效；应力小于强度的时候，结构可靠。

“应力—强度干涉模型”具体如图1所示，fs(S)表示应力的概率密度分布曲线，fr(R)表示强度的概率密度分布曲线，随时间的增加，由于损耗会导致强度下降，故二者曲线产生相交区域，相交区域即结构发生失效的概率[5]。

设XRi是与结构强度有关的变量；XSi是与结构应力有关的变量，由R、S两个基本随机变量表示结构的极限状态函数，大大简化了极限状态函数的形式，结构的极限状态函数如式（1）所示。

图1 应力-强度干涉模型

当Z＞0时，结构可靠，Ps表示结构的可靠度；当Z＜0时，结构失效，Pf表示失效概率。若R和S相互独立，设fr(r)、fS(s)为随机变量R、S的概率密度函数，则结构的可靠度Ps如式（2）所示[5]。

2.3 响应面法

由上节介绍的结构可靠度计算公式可知，当结构的极限状态函数确定的时候，可以通过计算得出可靠度及失效概率等参数，但并不是所有的结构极限状态函数都可以通过显式函数表示，当结构极限状态函数为隐式函数时，可通过响应面法计算结构可靠度。

响应面法通过在样本空间采集一定数量的样本点，利用多项式函数近似结构隐式功能函数，拟合出极限状态方程，利用拟合出的极限状态方程，替代有限元计算。这样一来，原本耗时很久的有限元模型计算就化简为对拟合极限状态方程的计算，大大减少了计算时间，并且没有计算次数的限制，通常仅需几秒钟就可以对拟合方程进行上万次的计算[6]。

响应面法通常用二次多项式函数来拟合结构的极限状态函数。假设输入变量与输出变量的函数关系符合二次多项式函数，一般选取下式进行数据拟合：

上式中，a是常数项，b（i=1，2，3…k）是线性项系数，c（i=1，2，3…k）是二次项系数，利用最小二乘回归分析的方法得到这些系数。

具体过程如下：

1）样本空间采样

类似蒙特卡洛法，在样本空间抽取一定数量的数据点，作为拟合函数的样本值。

2）拟合响应面

根据选定的样本点数据，选取二次多项式函数模型拟合响应面可以提高计算效率。

3）计算可靠度

通过拟合的响应面作为模型进行计算，利用蒙特卡洛法的计算方法，求解结构的可靠度。

3 机械结构可靠性仿真

3.1 ANSYS-PDS

ANSYS中自带的PDS（概率分析模块）可分析一个或多个变量对分析结果的影响，通过定义各参数的概率分布类型和特性参数，解决结构分析中各参数不确定性对计算结果的影响[7]。应用PDS模块进行结构可靠性分析的步骤如下：

1）创建分析文件

应用APDL语言参数化建模的方法建立计算有限元模型，包括前处理中的定义参数、选择单元类型、建模、划分网格、设置边界条件、加载、计算求解等常规步骤，在文件的结尾添加参数提取命令，建立极限状态方程[7]。

2）定义输入变量。

定义输入参数的分布类型，包括高斯分布、威布尔分布等；选择计算方法，PDS中主要的计算方法有蒙特卡洛法和响应面法；设置抽样次数，开始计算。

3）后处理

后处理环节可以计算得出可靠度和实效概率的数值，分析各输入参数对结果的影响及输出参数灵敏度等参数。

3.2 法兰算例

以法兰为算例，介绍运用PDS模块进行结构可靠性分析的过程。首先利用APDL语言参数化建立有限元模型，全过程采用APDL命令流语言编程生成，材料属性见表1，边界条件设置为螺栓孔位置全约束，网格在ANSYS中自动生成，在法兰端面加载压力P，几何模型见图2，有限元模型如图3所示。

图4、图5分别是求解结果的等效应力云图和位移云图，由云图可以看出，法兰最大应力为127.2 MPa，最大应力位置集中在螺栓孔处；最大位移为0.01 mm，最大变形位置集中在法兰外侧边缘位置。

在编制好的完整的APDL命令流结尾处添加如下参数提取命令，建立极限状态方程Z。

接下来进入PDS模块，将已生成的文件名为falan的分析文件导入PDS模块，定义随机输入参数与随机输出参数，见表2。

表1 法兰材料参数

图2 几何模型

图3 有限元模型

定义好输入输出参数后，选择响应面法进行模拟抽样，如图6所示，这里采用CCD实验设计法，ANSYS系统自动确定需要拟合的样本数为283个。

通过输入变量与输出变量的拟合函数关系拟合生成的响应面如图7、图8所示，图中给出了283个样本点分布情况，这里仅给出R3、R8、S作为输入变量，Z作为输出变量的响应面。对响应面即拟合函数的循环模拟，这里进行10 000次的模拟循环，得到计算结果。

图4 等效应力云图

图5 位移云图

表2 随机输入变量与随机输出变量

图6 拟合样本点

图7 Z-S-R3响应面

图8 Z-S-R8响应面

后处理阶段首先查看随机输入变量的历史抽样曲线，这里给出弹性模量E的抽样曲线，见图9、10，从图中可以看出，样本抽样点曲线密集，样本抽样统计分布直方图光滑，说明抽样次数足够多。

经过10 000次模拟循环，极限状态函数Z的均值为39.2 MPa，最大值137.2 MPa，最小值-75.9 MPa，在95 %的置信区间下，法兰的失效概率（见图11，即Z＜0）为7.1 %，可靠度为92.9 %。

在后处理阶段还可以查看输出变量的灵敏度，分析各参数对结构可靠性的影响，结构可靠性的灵敏度见图12，可以看出对法兰可靠性影响最大的变量是屈服极限S与压力P，S的灵敏度为正，与法兰的结构强度可靠性为正相关，P的灵敏度为负，与法兰结构强度可靠性为负相关，图13给出了Z-S散点图，可以看出二者符合正相关的函数关系。

图9 样本点抽样曲线(E)

图10 样本点统计分布直方图(E)

图11 法兰结构失效概率

图12 灵敏度(Z)

图13 Z-S散点图

4 结论

本文介绍了机械结构可靠性设计的基本理论，主要包括应力-强度干涉模型及响应面法。介绍了一种仿真方法的分析流程，基于ANSYS-PDS模块对法兰结构强度可靠性进行分析，得出结构可靠度为92.9 %，通过灵敏度分析，屈服极限S和压力P对法兰结构强度的影响最大，为结构优化提供一种思路和参考，具有一定的工程应用价值。