CN112906250A - 一种复杂***模块分类方法 - Google Patents

Abstract

本发明属于复杂***模块分类技术领域，涉及一种复杂***模块分类方法，包括下述步骤：1、采集复杂***现场故障数据；2、将复杂***划分为N个单元，确定各单元的故障分布函数；3、对各单元故障相关性分析及计算；4、计算各单元故障的综合关联强度，建立单元间故障关联矩阵；5、关联矩阵变换，求强连通集合，得到初始聚类模块；6、建立模块分类量化数学模型，构建目标函数；7、使用分组遗传算法对目标函数进行优化，得到最佳模块划分。本发明复杂***模块分类方法使用Copula函数量化故障关联单元之间的关联强弱，避免了人为因素造成的结果不准确，提高了故障关联单元关联强度的精确程度，使用聚合度耦合度作为量化指标，对模块内部以及模块之间的关系进行分析，通过构建目标函数并使用遗传算法优化，得到的最佳划分模块相对符合实际。

Description

一种复杂***模块分类方法

技术领域

本发明属于复杂***模块分类技术领域，涉及一种复杂***模块分类方法，具体涉及单元故障分布函数确定、单元故障相关性分析及计算、单元间故障关联矩阵建立、初始聚类模块划分、模块分类量化数学模型及目标函数建立、使用分组遗传算法对目标函数优化，输出最佳划分模块。

技术背景

复杂***内在的复杂结构和多重功能决定了其各个单元之间故障存在一定的相关性，某一单元的故障，有可能会引起其他单元的故障，因此探究复杂***故障之间的关联性和独立性，并进行模块划分对分析其故障传播和进行故障诊断有着重要的意义。

现有的模块划分方法主要有基于产品功能进行模块划分；基于产品功能和子功能之间的相关程度进行模块划分；基于零部件相关度进行模块划分；基于智能算法进行模块划分；基于ISM(Interpretative Structural Modeling)进行层次划分等。以产品功能和零部件为基础的模块分类，其强调的是产品设计或者产品使用过程中的关联关系，模块的划分倾向于关联功能和零部件的聚类，并没有将故障数据这种反馈考虑在内；ISM方法在构建关系矩阵的过程中，多数依据有向图的连接关系采用0或1进行赋值，仅考虑关联的存在与否，并没有量化关联强度的大小，并且过于主观。

复杂***运行期间，涉及部件众多，运行机理复杂，导致故障传播网络隐藏性较强，传统的方法并没有对故障发生单元的影响区域进行深入分析，不适用于实际复杂***的模块分类。

发明内容

针对现有复杂***模块分类主要以产品功能为基础，并且大都集中在纵向层次划分，忽略了各层级之间各单元故障关联以及各单元故障可能辐射到的区域，本发明提供了一种结合Copula函数、解析结构模型以及相应量化指标的复杂***模块划分方法，利用该方法对复杂***进行模块划分，考虑了各单元的故障相关性，更加符合工程实际。

为解决上述技术问题，本发明采用如下技术方案实现，结合附图说明如下：

一种复杂***模块分类方法，包括下述步骤：

步骤一、采集复杂***现场故障数据；

步骤二、将复杂***划分为N个单元，确定各单元的故障分布函数；

步骤三、各单元故障相关性分析及计算；

步骤四、计算各单元故障的综合关联强度，建立单元间故障关联矩阵；

步骤五、关联矩阵变换，求强连通集合，得到初始聚类模块；

步骤六、建立模块分类量化数学模型，构建目标函数；

步骤七、使用分组遗传算法对目标函数进行优化，得到最佳模块划分。

步骤二中所述确定各单元故障分布函数具体指：

根据复杂***结构和功能，将复杂***划分为N个单元。

假设各单元故障间隔时间样本服从威布尔Weibull分布，各单元概率密度函数为：

可以得到故障分布函数：

对上式进行线性变换得到

设一元线性回归方程为y＝a+bx，按照一元线性回归方程的形式，将威布尔分布进行线性变换为

x＝lnt，t表示按增序排列的故障间隔时间，并且a＝-βlnα，b＝β。由此可以看到，只要得到a，b的值，就可以对形状参数β和尺寸参数α进行估计，很明显a，b是线性回归方程y＝a+bx的截距和斜率，其中x可以直接通过故障间隔时间按x＝lnt得到，y通过式

可以计算，因此需要提前估计F(t_i)的值，本发明中采用中位秩估计法进行估计，即

其中n是样本总数，t_i是样本数据按增序排列的秩次，i＝1～n。按照最小二乘法原理可以得到a，b的值，其中

x_i，y_i表示不同时间下的变量数据，其中

x_i＝(lnt)_i，(lnt)_i是第i个故障间隔时间所求结果。计算检验统计量D_n，最后根据故障时间数据量n及显著性水平α,查表确定D检验临界值D_n(α)，当D_n≤D_n(α)时，则认为lnt与

是线性相关的，故障数据服从假设分布，否则拒绝假设。

步骤三中所说的各单元故障相关性分析及计算具体指：

根据复杂***历史故障数据、***结构及功能等综合考虑分析各单元故障之间的相关性。存在故障相关性的两个单元称为相关单元，考虑到相关单元之间存在的故障形式以及涉及的相关单元的数量，用关联故障连接图可以表示相关单元之间存在的故障形式以及涉及的相关单元的数量关系，其中圆圈表示单元，有向线段表示故障传递方向，常见的关联故障关系图如图2所示。图中1、3、4、5的i单元只影响其他单元，该单元定义为相关故障的起点单元；同理，j单元只被其他单元影响，定义为相关故障的终点；5中k单元既是影响组件又是被影响组件，定义为相关故障的中间点。通过建立关联故障连接图可以很清楚的得到单元故障之间的相关关系，明确具有故障关联的单元。

本发明使用阿基米德Copula(Archimedean Copula)分布函数簇中的GumbelCopula函数结合故障分布函数对具有故障关联单元的相关性进行计算，对函数相关参数θ进行估计，并用1-θ作为关联单元的关联强度。Gumbel Copula函数具体表达式为：

式中：F_i(x)为第i个单元的故障分布函数，i＝1～n；C(F₁(x),…,F_i(x),…,F_n(x))为n个单元故障分布函数的联合分布函数。

当1-θ＝0时，表明这n个单元故障相互独立，当1-θ→1时，意味这n个单元故障趋于完全相关。分别用0、0.2、0.4、0.6、0.8、1表示单元故障相互独立、单元故障基本独立、单元故障间关联微弱、单元故障间存在一定的关联、单元故障间关联较强、单元故障完全相关6个关联强度等级。

步骤四中所述的综合关联强度以及故障关联矩阵是指：

关联强度结合各关联因素权重计算各单元故障的综合关联强度。

用任意两个单元i,j之间故障关联强度的加权值g_ij表示单元故障的综合关联强度：

式中：1-θ_ij为i单元故障与j单元故障的关联强度。L为关联因素的数量。w_l为第l个关联因素的权重，该权重根据具体实际情况定量的给出，0≤w_l≤1并且

(如考虑各关联因素相等的话，w_l可以取1/L)。通过以上建立组件单元间故障关联矩阵：

G_n＝[g_ij]_N×N……………(5)

其中：N为单元数量。

步骤五中所述关联矩阵变换，求强连通集合，得到初始聚类模块按照以下计算步骤：

(1)、故障关联矩阵的布尔化

选取ε作为截取值，根据公式：

式中ε的选取参考相关实际经验和具体情况参照关联强度等级，通过计算对比选取。

将故障关联矩阵布尔化为矩阵A＝[a_ij]，a_ij表示故障i与故障j之间的影响，如果a_ij＝1，则认为故障i与故障j有直接的影响；如果a_ij＝0，则认为故障i与故障j不存在直接影响。

(2)、求可达矩阵

用G_r表示可达矩阵，运用布尔运算法则计算(I∪A)^k，令(I∪A)^k＝A_k，若有A_k＝A_k+1＝…＝A_n，则可达矩阵G_r＝A_k，其中I为单元矩阵，A为布尔化后的关联矩阵，k表示幂次，算子∪为逻辑和，即γ∪λ＝max{γ,λ}，布尔运算法则如下：

0+0＝0,0+1＝1,1+0＝1,1+1＝1

当

时，

A₃＝A₂，所以G_r＝A₂。

在图论中，可达矩阵G_r是有向图中节点之间的关系矩阵，该矩阵表示节点与节点之间的到达关系。在复杂***中用来表示某个单元节点故障与其他单元节点故障之间直接或间接的影响关系。

(3)、求强连通集合

设可达矩阵G_r变形为G_e，

其中b_i为N维行向量(i＝1,2,…,N)，算子∩为逻辑积，即γ∩λ＝min{γ,λ}。这些行向量中所有互不相等的行向量组成的集合为{b′₁,…,b′_i,…,b′_m}(1≤m≤N)，其中b′_i是相同向量的集合，设b′_i(1≤i≤m)中所有值为1的分量是

1≤t≤N，则这些分量组成为一强连通子集。其中λ_t是分量值为1的向量的定位下标。如b_i′表示的向量为[10001001]，则值为1的b′_i分量是

定位下标分别为λ₁,λ₈。通过对G_e作行列交换，使得相同的行聚集，可更直观的看到连通子集。故障之间存在强连通关系，则表明其存在相互影响，通过对强连通集合进行识别,分析单元故障之间存在的相互影响关系。

(4)、进行初始聚类模块划分

把每一个强连通子集看作一个模块，模块内部故障单元编号构成一个模块组合，剩余单独的故障单元编号统一为一个模块，按照这种方式划分初始聚类模块。

步骤六中建立数学模型和目标函数步骤如下：

假定某复杂***存在Q个单元，这Q个单元的故障集用C＝{c₁,…,c_u,…,c_Q}表示，其中c_u,u＝1,2,…,Q表示第u个单元故障；将该***划分为K个模块，模块集合表示为M＝{m₁,…,m_v,…,m_K}，其中m_v,v＝1,2,…,K表示第v个模块，设单元的归属矩阵X＝[x_uv]_Q×K

其中:

用J，E分别表示聚合度、耦合度，采用线性加权的方式得到目标函数

及目标函数的约束条件：

目标函数：

约束：

式中：w_J为J的权重，w_E为E的权重；

K为模块数，MU为最大模块数，即强连通集合矩阵的维数，K以及每个模块内的故障模式数为变量。

模块内聚合度衡量的是模块内故障关系的密切程度，高度互联的故障单元应该划分为一个模块。因此，以故障间的关联度为基础，定义第m_v(v＝1,2,…,K)个模块内的聚合度为：

其中，g_max为关联矩阵G_n中的最大值。

表示第i个单元属于第m_v个模块；

表示第j个单元属于第m_v个模块。

复杂***所有模块的聚合度为：

注：当某个模块内只含有一个单元时，该模块聚合度为0。

定义两个模块

(v₁,v₂＝1,2,…,K,v₁≠v₂)之间的耦合度为：

其中，g_max为关联矩阵G_n中的最大值。

复杂***模块间的耦合度为：

步骤七中使用分组遗传进行目标函数优化过程如下：

(1)、编码方案

分组遗传算法的编码方式概括为：定义目标部分e和组部分h，这两部分组成个体z＝[e|h]。其中e包含n个待分组的基因；假设分成了π组，则h是由π个组的各个顺序标号组成。e中待分组的基因由π个组的标号组成，表示基因属于对应的组。在聚类问题中e表示样本，h表示聚类。

例如，存在某样本e中n＝10，将其划分为4组，即聚类h中π＝4，该样本一个可能的个体编码为1 2 4 3 2 2 1 4 1 6|1 2 3 4，该个体表示将10个样本划分为{1,7,9}{2,5,6}{4,10}{3,8}4组。

(2)、种群初始化

由编码方案可以看出，e的长度固定不变，不同z所得到的聚类数π不相同，因此，随着π的变化，h的长度随之变化。如果π未知，则π取[π_min,π_max]内的整数，其中设定π_min＝2，π_max根据实际情况确定。按照以上变化关系，在进行种群初始化的过程中，每一个个体的生成方法为：π取随机数π∈[π_min,π_max]表示为该个体聚类过程产生的聚类数；接着将随机数序列[1,2,…,π]作为该个体组部分的顺序标号；最后随机用某一组的顺序标号表示e中的一个随机基因，对每个基因确定所属分组。重复ρ次此过程，就生成了规模为ρ的初始化种群。

(3)、交叉

遗传操作中的交叉算子是分组遗传算法中的核心。交叉算子以一定的概率交换两个随机个体的部分基因片段，从而产生对父辈的优良特性进行继承和传递的新的个体基因。

(4)、变异

变异算子对群体中随机个体中的某些基因结构上某点的值进行变动，产生变异基因组合，通过变异加快收敛。

(5)、设计适应度函数

为了更加凸显聚类模块的合理，采用模块内部的聚合度以及模块之间的耦合度通过线性加权的方式建立适应度函数，并求适应度函数最大值，本发明所述适应度函数即式(8)所构建的目标函数：

通过分组遗传算法不断优化，当适应度函数值越大，表明分组越优。

至此得到最佳的模块划分结果，实现复杂***模块划分。

与现有技术相比本发明的有益效果是：

本发明复杂***模块分类方法使用Copula函数量化故障关联单元之间的关联强弱，避免了人为因素造成的结果不准确，提高了故障关联单元关联强度的精确程度，使用聚合度耦合度作为量化指标，对模块内部以及模块之间的关系进行分析，通过构建目标函数并使用遗传算法优化，得到的最佳划分模块相对符合实际。

附图说明

下面结合附图对本发明作进一步的说明：

图1为本发明所述复杂***模块分类方法流程图；

图2为关联故障连接图；

图3为J、E、

函数值迭代图；

图4为可靠度函数对比曲线图；

图5为误差对比曲线图。

具体实施方式

下面结合附图对本发明作详细的描述：

参阅图1所示，本发明的复杂***模块分类方法包括下述步骤：采集故障数据；将复杂***划分为N个单元；确定各单元的故障分布函数；各单元故障相关性分析及计算；建立单元间故障关联矩阵；划分初始聚类模块；建立模块分类量化数学模型，构建目标函数；优化目标函数，输出最佳划分模块。

一、采集复杂***现场故障数据

1、故障数据：针对特定复杂***，通过现场调查，收集复杂***故障数据，包括具体故障发生时间、故障持续时间、维修时间以及故障具体信息等。

2、故障信息：明确各故障的故障模式、故障类型、故障传递关系、故障原因、故障部位、导致后果等数据。

二、将复杂***划分为N个单元，确定各单元的故障分布函数；

根据复杂***结构和功能，将其划分为N个单元。梳理采集到的各单元故障间隔数据，假设其服从二参数的威布尔分布，通过最小二乘法进行参数估计，并用d检验进行假设检验，最终确定各单元的故障分布函数。

各单元故障间隔时间服从二参数Weibull分布，则故障分布函数为

经过线性变换得到

设一元线性回归方程为y＝a+bx，按照该形式，将威布尔分布进行线性变换为

x＝lnt，并且a＝-βlnα，b＝β。由此可以看到，只要得到a，b的值，就可以对形状参数β和尺寸参数α进行估计，很明显a，b是线性回归方程y＝a+bx的截距和斜率，其中x可以直接通过故障间隔时间按x＝lnt得到，y通过式

可以计算，因此需要提前估计F(t_i)的值，本发明中采用中位秩估计法进行估计，即

其中n是样本总数，t_i是样本数据按增序排列的秩次，i＝1～n。按照最小二乘法原理可以得到a，b的值，其中

最后根据故障时间数据量n及显著性水平α,计算D检验临界值D_n(α)，当D_n≤D_n(α)时，则认为lnt与

是线性相关的，故障数据服从假设分布，否则拒绝假设。

三、各单元故障相关性分析及计算

考虑到相关单元之间存在的故障模式以及涉及的相关单元的数量，用关联故障连接图可以表示这种关系，其中圆圈表示单元，有向线段表示故障传递方向，常见的关联故障关系图如图2所示。图中1、3、4、5的i单元只影响其他单元，该单元定义为相关故障的起点单元；同理，j单元只被其他单元影响，定义为相关故障的终点；5中k单元既是影响组件又是被影响组件，定义为相关故障的中间点。通过建立关联故障连接图可以很清楚的得到单元故障之间的相关关系，明确具有故障关联的单元。

使用Gumbel Copula函数对具有故障关联的单元进行相关性计算，对函数相关参数θ进行估计，并用1-θ作为关联单元的关联强度。Gumbel Copula函数具体表达式为：

式中F_i(x)为第i个单元的故障分布函数，i＝1～n；C(F₁(x),F₂(x),…,F_n(x))为n个单元故障分布函数的联合分布函数；当1-θ＝0时，表明这n个单元故障相互独立，当1-θ→1时，意味这n个单元故障趋于完全相关。分别用0、0.2、0.4、0.6、0.8、1表示单元故障无关联、基本独立、关联微弱、存在一定的关联、关联较强、不可分，6个关联强度等级。

四、计算各单元故障的综合关联强度，建立单元间故障关联矩阵

用任意两个单元故障i,j之间关联强度的加权值g_ij表示单元故障的综合关联强度，当i≠j时，

当i＝j时，g_ij＝0。通过以上建立组件单元间故障关联矩阵G_n＝[g_ij]_N×N。

五、关联矩阵变换，求强连通集合，得到初始聚类模块

根据相关实际经验和具体情况参考关联强度等级，选取ε作为截取值，根据

将故障关联矩阵布尔化为矩阵A＝[a_ij]，如果a_ij＝1，则认为故障i与故障j有直接的影响；如果a_ij＝0，则认为故障i与故障j不存在直接影响。

设A_k＝(I∪A)^k，运用布尔运算法则，当k＞n₀时，若有A_k＝A_k+1＝…＝A_n，则可达矩阵G_r＝A_k，其中I为单元矩阵，k表示幂次，n₀指某一正整数，算子∪为逻辑和，布尔运算法则为：0+0＝0,0+1＝1,1+0＝1,1+1＝1。在图论中，可达矩阵G_r是有向图中节点之间的关系矩阵，该矩阵表示节点与节点之间的到达关系。在复杂***中用来表示某个单元节点故障与其他单元节点故障之间直接或间接的影响关系。

设可达矩阵G_r变形为G_e，

其中b_i为N维行向量(i＝1,2,…,N)，算子∩为逻辑积，即a∩b＝min{a,b}。这些行向量中所有互不相等的行向量组成的集合为{b′₁,b′₂,…,b′_m}(1≤m≤N)，其中b′_i是相同向量的集合，设b′_i(1≤i≤m)中所有值为1的分量是

1≤t≤N，则这些分量组成为一强连通子集。其中λ_t是分量值为1的向量的定位下标。如b′_i表示的向量为[10001001]，则值为1的b′_i分量是

定位下标分别为λ₁,λ₈。通过对G_e作行列交换，使得相同的行聚集，可更直观的看到连通子集。故障之间存在强连通关系，则表明其存在相互影响，通过对强连通集合进行识别,分析单元故障之间存在的相互影响关系。

把每一个强连通子集看作一个模块，模块内部故障单元编号构成一个模块组合，剩余单独的故障单元编号统一为一个模块，按照这种方式划分初始聚类模块。

六、建立模块分类量化数学模型，构建目标函数

假定某复杂***存在Q个单元，这Q个单元的故障集用C＝{c₁,c₂,…,c_Q}表示，其中c_u,u＝1,2,…,Q表示第u个单元故障；将该***划分为K个模块，模块集合表示为M＝{m₁,m₂,…,m_K}，其中m_v,v＝1,2,…,K表示第v个模块，设单元的归属矩阵X＝[x_uv]_Q×K，x_uv表示单元的模块归属，当第u个单元故障属于第v个模块时x_uv＝1，当第u个单元故障不属于第v个模块时，x_uv＝0。

用J，E分别表示聚合度、耦合度，采用线性加权的方式得到目标函数

及约束条件

u＝1,2,...,Q，2≤K≤MU。

模块内聚合度衡量的是模块内故障关系的密切程度，高度互联的故障单元应该划分为一个模块。因此，以故障间的关联度为基础，定义第m_v(v＝1,2,…,K)个模块内的聚合度为

复杂***所有模块的聚合度为

定义两个模块

(v₁,v₂＝1,2,…,K,v₁≠v₂)之间耦合度为

复杂***模块间的耦合度为：

七、使用分组遗传算法对目标函数进行优化，得到最佳模块划分

(1)、编码方案

分组遗传算法的编码概括方式为：定义目标部分e和组部分h，这两部分组成个体z＝[e|h]。其中e包含n个待分组的基因；假设分成了π组，则h是由π个组的各个顺序标号组成。e中待分组的基因由π个组的标号组成，表示基因属于对应的组。在聚类问题中e表示样本，h表示聚类。

(2)、种群初始化

由编码方案可以看出，e的长度固定不变，不同z所得到的聚类数π不相同，因此，随着π的变化，h的长度随之变化。如果π未知，则π取[π_min,π_max]内的整数，其中设定π_min＝2，π_max根据实际情况确定。按照以上变化关系，在进行种群初始化的过程中，每一个个体的生成方法为：π取随机数π∈[π_min,π_max]表示为该个体聚类过程产生的聚类数；接着将其序列[1,2,…,π]作为该个体组部分的顺序标号；最后随机用某一组的顺序标号表示e中的一个随机基因，对每个基因确定所属分组。重复ρ次此过程，就生成了规模为ρ的初始化种群。

(3)、交叉

遗传操作中的交叉算子是分组遗传算法中的核心。交叉算子以一定的概率交换两个随机个体的部分基因片段，从而产生对父辈的优良特性进行继承和传递的新的个体基因。

(4)、变异

变异算子对群体中随机个体中的某些基因结构上某点的值进行变动，产生变异基因组合，通过变异加快收敛。

(5)、设计适应度函数

为了更加凸显聚类模块的合理，采用模块内部的聚合度以及模块之间的耦合度通过线性加权的方式建立适应度函数，并求其最大值，因此建立如下适应度函数：

通过分组遗传算法不断优化，当适应度函数值越大，表明分组越优。

实施例

加工中心***模块分类

采集40台某型号加工中心的实际运行历史故障数据，依据其结构与功能将其划分为Q＝11个子***单元，通过分析确定其子***单元的故障规律及其线性相关系数，如表1所示。

表1子***单元故障分布函数及其线性相关系数

从表1可知，几乎所有的子***单元故障分布函数的线性相关系数都偏小，可认为这些子***单元不具备故障独立性。

在进行关联因素划分时需要考虑到各子***单元内在属性，确定各关联因素的权重，再逐一确定每个子***单元故障的关联强度。将加工中心划分为机械、电子、辅助三个关联因素，统计加工中心各关联因素发生的故障频率，将其作为个关联因素的权重。加工中心各个关联因素及其权重如下表2所示：

表2加工中心故障关联因素及其权重

经过整理得到加工中心子***单元部分关联故障如下表3所示：

表3部分加工中心子***单元关联故障统计表

在得到各子***单元故障分布函数的基础上，依照表3确定故障关联子***单元，并根据Copula函数求其关联参数，利用式(3)进行参数估计，再根据式(4)得到综合关联强度，形成加工中心故障关联矩阵G_n。

参照关联强度等级，当故障关联等级为0时，两个单元之间的关联关系可以忽略，此处取截取值ε＝0为例展开过程论述。

如得到故障关联矩阵G_n为，

按照等式(6)将关联矩阵布尔化为如下矩阵A：

通过布尔运算，求得矩阵A的可达矩阵G_r：

将可达矩阵G_r变形为矩阵G_e，通过对G_e作行列交换，使得相同的行聚集，得到强连通子集，G_e变形为矩阵G_e′，方括号内表示强连通子集。

按照强连通子集的划分，将11个加工中心子***单元初始聚类为{S1,S10}、{S7,S9}、{S3,S11}、{S2,S4,S5,S6,S8}等4模块；据此根据式(7)对各故障子***单元进行分类，得到故障子***单元归属矩阵X为：

按照式(9～12)计算得出加工中心各个模块内部子***单元的聚合度、模块整体聚合度以及模块之间的耦合度、模块整体耦合度，最后通过式(8)经过线性加权得到初始聚类模块的目标函数值(设置聚合度和耦合度的权重都为0.5)。得到聚合度为0.1201，耦合度为0.0132，函数值为0.5534。

同理当截取值ε＝0.2时，按照以上步骤将加工中心分为{S3,S11}、{S7,S9}、{S1,S2,S4,S5,S6,S8,S10}等3个模块，聚合度为0.0709，耦合度为0.0131，函数值为0.5289；当截取值ε＝0时，加工中心分为

和{S1,S2,S3,S4,S5,S6,S8,S10,S11}两个模块，聚合度为0.0413，耦合度为0.0227，函数值为0.5093；当截取值ε＝0.6时，结果与ε＝0.4时一致；根据实际，ε不能取0.8及以上。故聚类结果如表4所示：

表4初始聚类模块量化指标以及目标函数值

由表4得出，ε＝0时聚合度、函数值及耦合度都优于其他取值，因此按ε＝0的取值，***初始聚类模块为{S1,S10}、{S7,S9}、{S3,S11}、{S2,S4,S5,S6,S8}4个模块，并以此作为分组遗传算法的输入，设聚合度J和耦合度E的权重均为0.5，设聚类模块数K最小值为2，最大值为8，即2≤K≤8，根据式(8)得到目标函数及约束条件为：

使用Matlab软件对目标函数进行优化，设置种群规模为50，交叉概率为0.8，变异概率为0.2，迭代代数为200，优化次数为10次，统计结果并进行分析，结果如表5所示：

表5 Matlab优化10次的结果

从表5数据可知，当聚合度增加时，耦合度呈现减少趋势，目标函数值的变化趋势与聚合度一致。聚合度达到最大值时，耦合度达到了最小值，此时目标函数值也达到了最大值。经过模块优化，耦合度稳定到0.01范围内，整体聚合度从0.1201提高到了0.3589，目标函数值从0.5534提高到0.6729，说明经过优化后的模块聚类具有优越性。

迭代图如图3所示。

由图3可知，随着迭代次数的增加，聚合度和目标函数值逐步增加，最后趋于平稳；耦合度随着聚合度的增加呈现减少趋势，最终也趋于平稳。分组遗传算法在迭代30次左右的时候已经接近最优解，收敛速度较快。

由表5和图3可知，最优目标函数值为0.6729，此时将聚类模块为{S1,S10}、{S2,S6}、{S3,S11}、{S4,S5}、{S7,S9}、{S8}共6组。与优化前相比，目标函数值明显提高。

通过基于故障相关的聚类分析，将加工中心划分为6个模块，各模块的联合故障分布函数根据式(3)可以得到，分别为：

C(F₈(x))＝F₈(x)

据此得到各模块可靠度函数分别为：

R_1,10＝1-C(F₁(x),F₁₀(x))

R_2,6＝1-C(F₂(x),F₆(x))

R_3,11＝1-C(F₃(x),F₁₁(x))

R_4,5＝1-C(F₄(x),F₅(x))

R_7,9＝1-C(F₇(x),F₉(x))

R₈＝1-F₈(x)

加工中心在运行过程中任意子***单元出现故障都会影响设备运行，所以，加工中心是典型的串联***。其***可靠度为各子***单元可靠度的乘积。

当考虑故障相关时，加工中心***可靠度为：

R_{opt_group}＝R_1,10·R_2,6·R_3,11·R_4,5·R_7,9·R₈

当假设加工中心各子***单元故障独立时，***可靠度为：

R_indep＝R₁·R₂·R₃·R₄·R₅·R₆·R₇·R₈·R₉·R₁₀·R₁₁

假设加工中心***故障服从二参数Weibull分布，根据故障数据，使用最小二乘法进行参数估计，经线性相关性检验和d检验，得到***可靠度模型为：

结合加工中心故障数据，计算不同情况下***可靠度，并建立可靠度函数曲线。此处选择最优模块划分结果以及其他3个非最优模块划分结果建立加工中心***可靠度对比曲线，如图4所示。

表6不同情况下的可靠度

表6中：R_{opt_group}为最优模块划分***可靠度；R_{ord_group}为非最优模块划分***可靠度，对应分组为：

opt_group：{1,11}{3,10}{2,6}{4,5}{7,9}{8}

indep：{1}{2}{3}{4}{5}{6}{7}{8}{9}{10}{11}

ord_group1：{1,2}{4,7}{3,11}{5,6}{9,10}{8}

ord_group2：{1,10}{2,3}{5,8}{4,7}{9,11}{6}

ord_group3：{1,7}{2,3}{5,8}{6,9}{10,11}{4}

由图4和表6可知，考虑子***单元故障相关的***可靠度曲线较子***单元故障独立的***可靠度曲线下降明显缓慢；最优模块划分结果下***可靠度曲线较非最优结果可靠度曲线更接近***的可靠度曲线。由对比得到的***可靠度曲线可知，考虑故障相关的***可靠度曲线更加合理。

为进一步说明不同聚类模块下的***可靠度和独立串联***可靠度建模误差，建立各可靠度误差对比曲线，如图5所示。

由图5可以发现，在加工中心运行过程中，最优模块划分下的***可靠度误差最小；随着加工中心运行时间的推进，不同情况下的可靠度误差都呈现先递增再递减的变化趋势，并且在150h左右的时候误差达到最大，并在之后不断递减，在1000h左右可靠度误差曲线几乎重合。

最后，通过AIC准则(赤池信息量准则)衡量该优化分组的优越性。AIC＝2k+nln(RSS/n)，其中k是参数的数量，n为观察数，RSS为残差平方和，AIC越小说明模型越优。因对于该加工中心，k与n相同，故AIC最小等同于RSS最小，故比较不同情况下运行时间内***可靠度的RSS即可，具体如表7所示。

表7不同模块下的加工中心***可靠度的AIC与RSS

模块划分	opt_group	indep	ord_group1	ord_group2	ord_group3
						AIC	-455.4984	-212.9278	-392.4765	-398.0267	-413.0047
RSS	1.3889	13.6938	2.5171	2.3887	2.0739

由表7数据可以看出，最优模块划分的AIC、RSS值最小，进一步说明该分组最符合实际情况。

(1)本发明对复杂***各单元之间故障关联进行分析，通过Copula函数结合传统量化解析方法，避免了传统方法人为因素对关联强度造成主观影响。

(2)本发明通过求解故障关联单元之间的相关参数，建立故障关联矩阵，划分初始模块；用聚合度和耦合度作为量化指标建立单元故障模式聚类数学模型目标函数，使用分组遗传发算法对目标函数进行优化，确定最佳聚类模块。

(3)以某加工中心这类复杂***为研究对象对本发明方法进行了计算验证，通过对其子***单元故障关联关系进行分析、量化、模块分类、分类优化以及优化验证，证明了所提方法的有效性，该方法在实际生产中可以为设备维修和故障诊断提供一定参考。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，都应涵盖在本发明的保护范围之内。同时本说明书中未作详细描述的内容均属于本领域技术人员公知的现有技术。

Claims (7)

1.一种复杂***模块分类方法，包括下述步骤：

步骤一：采集复杂***现场故障数据；

步骤二：将复杂***划分为N个单元，确定各单元的故障分布函数；

步骤三：各单元故障相关性分析及计算；

步骤四：计算各单元故障的综合关联强度，建立单元间故障关联矩阵；

步骤五：关联矩阵变换，求强连通集合，得到初始聚类模块；

步骤六：建立模块分类量化数学模型，构建目标函数；

步骤七：使用分组遗传算法对目标函数进行优化，得到最佳模块划分。

2.根据权利要求1所述的复杂***模块分类方法，其特征在于：

步骤二中所述确定各单元故障分布函数具体指：

根据复杂***结构和功能，将复杂***划分为N个单元；

假设各单元故障间隔时间样本服从威布尔Weibull分布，各单元概率密度函数为：

可以得到故障分布函数：

对上式进行线性变换得到

设一元线性回归方程为y＝a+bx，按照一元线性回归方程的形式，将威布尔分布进行线性变换为

x＝lnt，t表示按增序排列的故障间隔时间，并且a＝-βlnα，b＝β；

得到a，b的值，对形状参数β和尺寸参数α进行估计，a，b是线性回归方程y＝a+bx的截距和斜率，其中x直接通过故障间隔时间按x＝lnt得到，y通过式

计算；

估计F(t_i)的值，采用中位秩估计法进行估计，即

其中n是样本总数，t_i是样本数据按增序排列的秩次，i＝1～n；

按照最小二乘法原理得到a，b的值，其中

x_i，y_i表示不同时间下的变量数据，其中

x_i＝(lnt)_i，(lnt)_i是第i个故障间隔时间所求结果；计算检验统计量D_n，并根据故障时间数据量n及显著性水平α,查表确定D检验临界值D_n(α)，当D_n≤D_n(α)时，则认为lnt与

是线性相关的，故障数据服从假设分布，否则拒绝假设。

3.根据权利要求2所述的复杂***模块分类方法，其特征在于：

步骤三中所述的各单元故障相关性分析及计算具体指：

根据复杂***历史故障数据、***结构及功能综合考虑分析各单元故障之间的相关性；存在故障相关性的两个单元称为相关单元，用关联故障连接图表示相关单元之间存在的故障形式以及涉及的相关单元的数量关系，其中圆圈表示单元，有向线段表示故障传递方向；

i单元只影响其他单元，该单元定义为相关故障的起点单元；j单元只被其他单元影响，定义为相关故障的终点；k单元既是影响组件又是被影响组件，定义为相关故障的中间点；通过建立关联故障连接图得到单元故障之间的相关关系，明确具有故障关联的单元；

使用阿基米德Copula分布函数簇中的Gumbel Copula函数结合故障分布函数对具有故障关联单元的相关性进行计算，对函数相关参数θ进行估计，并用1-θ作为关联单元的关联强度；Gumbel Copula函数具体表达式为：

式中：F_i(x)为第i个单元的故障分布函数，i＝1～n；C(F₁(x),…,F_i(x),…,F_n(x))为n个单元故障分布函数的联合分布函数；

当1-θ＝0时，表示这n个单元故障相互独立，当1-θ→1时，表示这n个单元故障趋于完全相关。

4.根据权利要求3所述的复杂***模块分类方法，其特征在于：

步骤四中所述的综合关联强度是指：

关联强度结合各关联因素权重计算各单元故障的综合关联强度；

用任意两个单元i,j之间故障关联强度的加权值g_ij表示单元故障的综合关联强度：

式中：1-θ_ij为i单元故障与j单元故障的关联强度；L为关联因素的数量；w_l为第l个关联因素的权重，0≤w_l≤1并且

步骤四中所述的故障关联矩阵是指：

依据单元故障的综合关联强度g_ij，建立组件单元间故障关联矩阵：

G_n＝[g_ij]_N×N……………(5)

其中：N为单元数量。

5.根据权利要求4所述的复杂***模块分类方法，其特征在于：

步骤五中所述关联矩阵变换，求强连通集合，得到初始聚类模块按照以下计算步骤：

(1)、故障关联矩阵的布尔化

选取ε作为截取值，根据公式：

式中ε的选取参考相关实际经验和具体情况参照关联强度等级，通过计算对比选取。

将故障关联矩阵布尔化为矩阵A＝[a_ij]，a_ij表示故障i与故障j之间的影响，如果a_ij＝1，则认为故障i与故障j有直接的影响；如果a_ij＝0，则认为故障i与故障j不存在直接影响。

(2)、求可达矩阵

用G_r表示可达矩阵，运用布尔运算法则计算(I∪A)^k，令(I∪A)^k＝A_k，若有A_k＝A_k+1＝…＝A_n，则可达矩阵G_r＝A_k，其中I为单元矩阵，A为布尔化后的关联矩阵，k表示幂次，算子∪为逻辑和，即γ∪λ＝max{γ,λ}，布尔运算法则如下：

0+0＝0,0+1＝1,1+0＝1,1+1＝1

当

时，

A₃＝A₂，所以G_r＝A₂；

在图论中，可达矩阵G_r是有向图中节点之间的关系矩阵，该矩阵表示节点与节点之间的到达关系；在复杂***中用来表示某个单元节点故障与其他单元节点故障之间直接或间接的影响关系。

(3)、求强连通集合

设可达矩阵G_r变形为G_e，

其中b_i为N维行向量(i＝1,2,…,N)，算子∩为逻辑积，即γ∩λ＝min{γ,λ}；这些行向量中所有互不相等的行向量组成的集合为{b′₁,…,b′_i,…,b′_m}(1≤m≤N)，其中b′_i是相同向量的集合，设b′_i(1≤i≤m)中所有值为1的分量是

则这些分量组成为一强连通子集；其中λ_t是分量值为1的向量的定位下标；如b′_i表示的向量为[10001001]，则值为1的b′_i分量是

定位下标分别为λ₁,λ₈；通过对G_e作行列交换，使得相同的行聚集，更直观的看到连通子集；故障之间存在强连通关系，则表明其存在相互影响，通过对强连通集合进行识别,分析单元故障之间存在的相互影响关系；

(4)、进行初始聚类模块划分

把每一个强连通子集看作一个模块，模块内部故障单元编号构成一个模块组合，剩余单独的故障单元编号统一为一个模块，按照这种方式划分初始聚类模块。

6.根据权利要求5所述的复杂***模块分类方法，其特征在于：

步骤六中建立数学模型和目标函数步骤如下：

假定某复杂***存在Q个单元，这Q个单元的故障集用C＝{c₁,…,c_u,…,c_Q}表示，其中c_u,u＝1,2,…,Q表示第u个单元故障将该***划分为K个模块，模块集合表示为M＝{m₁,…,m_v,…,m_K}，其中m_v,v＝1,2,…,K表示第v个模块，设单元的归属矩阵X＝[x_uv]_Q×K

其中:

用J，E分别表示聚合度、耦合度，采用线性加权的方式得到目标函数

及目标函数的约束条件：

目标函数：

约束：

式中：w_J为J的权重，w_E为E的权重；

K为模块数，MU为最大模块数，即强连通集合矩阵的维数，K以及每个模块内的故障模式数为变量；

模块内聚合度衡量的是模块内故障关系的密切程度，高度互联的故障单元应该划分为一个模块；以故障间的关联度为基础，定义第m_v(v＝1,2,…,K)个模块内的聚合度为：

其中，g_max为关联矩阵G_n中的最大值；

表示第i个单元属于第m_v个模块；

表示第j个单元属于第m_v个模块；

复杂***所有模块的聚合度为：

当某个模块内只含有一个单元时，该模块聚合度为0；

定义两个模块

之间的耦合度为：

其中，g_max为关联矩阵G_n中的最大值；

复杂***模块间的耦合度为：

7.根据权利要求6所述的复杂***模块分类方法，其特征在于：

步骤七中使用分组遗传算法对目标函数进行优化过程如下：

(1)、编码方案

分组遗传算法的编码方式概括为：定义目标部分e和组部分h，这两部分组成个体z＝[e|h]；其中e包含n个待分组的基因；假设分成了π组，则h是由π个组的各个顺序标号组成；e中待分组的基因由π个组的标号组成，表示基因属于对应的组；在聚类问题中e表示样本，h表示聚类；

(2)、种群初始化

由编码方案可以看出，e的长度固定不变，不同z所得到的聚类数π不相同，随着π的变化，h的长度随之变化；如果π未知，则π取[π_min,π_max]内的整数，其中设定π_min＝2，π_max根据实际情况确定；按照以上变化关系，在进行种群初始化的过程中，每一个个体的生成方法为：π取随机数π∈[π_min,π_max]表示为该个体聚类过程产生的聚类数；接着将随机数序列[1,2,…,π]作为该个体组部分的顺序标号；最后随机用某一组的顺序标号表示e中的一个随机基因，对每个基因确定所属分组；重复ρ次此过程，就生成了规模为ρ的初始化种群；

(3)、交叉

遗传操作中的交叉算子是分组遗传算法中的核心，交叉算子以一定的概率交换两个随机个体的部分基因片段，产生对父辈的优良特性进行继承和传递的新的个体基因；

(4)、变异

变异算子对群体中随机个体中的某些基因结构上某点的值进行变动，产生变异基因组合，通过变异加快收敛；

(5)、设计适应度函数

采用模块内部的聚合度以及模块之间的耦合度通过线性加权的方式建立适应度函数，并求适应度函数最大值，本发明所述适应度函数即式(8)所构建的目标函数：

通过分组遗传算法不断优化，当适应度函数值越大，表明分组越优；

至此得到最佳的模块划分结果，实现复杂***模块划分。

Images