基于面向方面灰色GSPN的CPS建模与分析方法

基于面向方面灰色gspn的cps建模与分析方法
技术领域
1.本发明涉及一种cps的建模与分析方法，特别是一种基于面向方面灰色gspn(广义随机petri网)的cps建模与分析方法。


背景技术：

2.信息物理系统(cyber-physical system,cps)是一种集合计算、通信和控制于一体的多维复杂系统，融合了离散的计算过程和连续的物理过程，并基于反馈循环实现物理世界与虚拟世界之间的信息融合以及行为交互。
3.在实际运行过程中，cps利用传感器实时感知物理世界，通过控制器处理获取的信息从而生成指令，使用执行器执行生成的指令，三者之间相互配合保证了其有效性和实时性。cps具有独特的优越性，被广泛应用于智能交通、航空航天、智能医疗等领域，更是在工业4.0中处于核心地位，成为了智能制造领域关注的焦点。
4.工业cps运行在高度动态的复杂环境中，须具备较强的鲁棒性和适应性，这导致了其具有较为复杂的结构，并对其建模与精确分析提出了额外的挑战。一个结构清晰、直观易懂的模型以及对应的性能分析方法不但有助于系统的前期研发、可靠度评估等，而且能保障系统的安全运行，更能为后期维护升级提供理论基础。
5.目前，国内外已经提出了一些方法与分析工具用于cps的建模、仿真以及各种定性、定量分析。例如，提出随机建模框架，基于框架制定并解决联合状态和攻击评估问题；提出cps结构与动态行为协同建模方法等。然而，这些理论模型及分析工具中，一些模型较为抽象、复杂，不易理解，并且无法与实际生产过程进行直观对比；一些方法的建模过程比较繁琐，并且模型的可重用性较差，若后续生产流程调整，模型需大幅更改；还有一些模型很少考虑到工业场景下cps运行的不确定性，例如忽略活动耗时，导致无法准确的进行工业cps性能指标分析。
6.cn111722599a《基于面向对象广义随机petri网的cps建模与分析方法》是发明人早期提出的一种用于cps的建模、仿真以及各种定性、定量分析方法。其采用面向对象的方法，将cps中的各模块抽象成不同的对象进行建模，并引入模糊数学的方法进行性能分析，但是此方法会忽略一些cps流程中的细节，导致准确性不足。
7.基于以上分析，为了解决现有模型较为抽象、复杂度高、描述分析能力不足等问题，针对工业cps的结构流程以及时间属性，提出一种面向方面灰色广义随机petri网，该模型以广义随机petri网为基础，结合了面向方面技术与颜色集概念，再引入灰色系统理论以表征工业cps的不确定性，从而能够对工业cps进行准确地分析。


技术实现要素：

8.本发明的目的是针对目前的技术的不足，提出一种基于面向方面灰色广义随机petri网的cps建模与分析方法。
9.本发明方法的大体思想是：
10.在cps的建模阶段，为了说明cps中各基础组件和流程，考虑模型在工业cps中的有效性和可行性，将面向方面技术与灰色广义随机petri网结合，提出了面向方面灰色广义随机petri网(aspect-oriented grey generalized stochastic petri net,aggspn)，作为一种新的cps形式化建模方法。由于引入了面向方面技术，使得一个工业cps被建模成一个基本网和若干个方面网，而非一个整体，使得无法对工业cps的完整流程进行定性、定量的分析。因此，本发明提出一种基于aggspn的编织算法将所有的方面网都正确地织入到基本网中，形成一个完整的编织网，之后即可对其进行性能分析。
11.在cps的性能分析阶段，考虑到cps工作过程中一些生产参数具有“小样本、贫信息”的特征，并且会随着时间的推移而有所变化，无法简单地用一个固定的统计值来表示。为了提高所建立模型的真实性和准确性，在模型中引入了灰色系统理论，利用区间灰数来表示模型中的变迁激发率以及随机开关的触发概率，并在求解模型对应的马尔科夫链时使用数学规划方法，得到最终的灰稳态概率分布。
12.本发明解决其技术问题所采用的方法的步骤如下：
13.基于面向方向灰色广义随机petri网的cps建模与分析方法，包括步骤如下：
14.步骤(1)、根据所建模cps的具体情境，将工业cps抽象为一个基本网和若干个方面网；具体是将工业cps的总运行流程抽象为一个基本网，再将部分子流程细分成若干个方面网；
15.所述工业cps的总运行流程包括多个子流程；
16.步骤(2)、基本网和方面网以gspn为基础，并引入颜色集的概念，即对其中的元素赋色，从而得到aggspn模型(面向方面灰色广义随机petri网，aspect-oriented grey generalized stochastic petri net)；
17.2-1针对基本网，以总运行流程中工业cps生产线上的工业设备状态作为gspn的库所，以总运行流程中工业cps生产线上的资源或者信号作为库所的token，以总运行流程中工业cps生产线上的工业设备动态行为作为gspn的活动变迁；根据总运行流程中库所与活动变迁间的关系作为弧；并对不同的token标注不同的颜色进行区别；(上述网中库所、活动变迁、token、弧均为gspn的元素)
18.2-2针对方面网，以子流程中工业cps生产线上的工业设备具体状态作为gspn的库所，以子流程中工业cps生产线上的资源或者信号作为库所的token，以子流程中工业cps生产线上的工业设备具体动态行为作为gspn的活动变迁；根据子流程中库所与活动变迁间的关系作为弧；并对不同的token标注不同的颜色进行区别；
19.步骤(3)、在保证模型一致性的前提下，将基本网和所有方面网编织成一个编织网，便于后续模型的性能分析；具体是：
20.(3.1)对于任意一个方面网，首先根据该方面网的通知类型，将其与基本网连接到一起，并删除多余辅助元素以及与多余辅助元素相关联的弧；
21.所述方面网的多余辅助元素包括起始元素start、终止元素end和重复的库所；
22.(3.2)若方面网与基本网含有相同库所p，则需将方面网中库所p的输入弧和输出弧转化为基本网中库所p的输入弧和输出弧，并且删除方面网中的库所p。
23.步骤(4)、引入灰色系统理论，利用灰化处理方法求解步骤(3)编织后的aggspn模型中的变迁激发率以及随机开关的触发概率；具体是通过引入一个扰动灰元δ将aggspn模
型中的变迁激发率和触发概率转化为区间灰数，为每个变迁设置一个系数k，即将时间变迁的扰动灰元系数设置为k1，随机开关的扰动灰元系数设置为k2，具体过程如下：
24.(4.1)若aggspn中的活动变迁t1为时间变迁，变迁激发率为λ1，则其区间灰数表示为[λ
1-k1δ,λ1 k1δ]；
[0025]
(4.2)同理，若aggspn中的活动变迁t2为一个随机开关，触发概率为pr(t2)，其区间灰数表示为[pr(t2)-k2δ,pr(t2) k2δ]。
[0026]
步骤(5)、使用马尔科夫方法求解aggspn的灰稳态概率；具体过程如下：
[0027]
(5.1)分析步骤(3)编织网，获取可达标识图以及对应的可达标识集s，s分为实存状态集ts和消失状态集vs，即s＝ts∪vs，其元素个数分别为k
t
和kv，ks＝k
t
kv；
[0028]
(5.2)根据可达标识图得到嵌入式马尔科夫链(embedded markov chain,emc)的转移概率灰矩阵其计算公式如下：
[0029][0030]
由从消失状态向消失状态集转移的转移概率矩阵与从消失状态向实存状态集转移的转移概率矩阵构成；
[0031]
由从实存状态向消失状态集转移的转移概率矩阵与从实存状态向实存状态集转移的转移概率矩阵构成；
[0032]
emc由消失状态集向消失状态集的转移概率矩阵；
[0033]
emc由消失状态集向实存状态集的转移概率矩阵；
[0034]
emc由实存状态集向消失状态集的转移概率矩阵；
[0035]
emc由实存状态集向实存状态集的转移概率矩阵；
[0036]
(5.3)将消失状态从emc中除去，得到remc转移状态矩阵用式(2)表示；
[0037][0038]
移除所有emc中的消失状态，仅剩下实存状态，就可得到一个压缩的emc(reduced emc,remc)，为remc的转移概率矩阵；
[0039]
表示元素的集合，其中表示gspn从给定消失状态r出发，经任意步首次到达实存状态k的概率；
[0040]
(5.4)通过式(3)的线性方程组，求解remc的实存状态灰稳定概率；
[0041]
[0042]
表示一维行向量，元素表示为remc的实存状态的灰稳态概率；
[0043]
(5.5)通过式(4)计算各个实存状态的平均驻留时间
[0044][0045]
(5.6)由各个状态的平均驻留时间以及对应的灰稳定概率可求得系统的循环周期计算方式为式(5)；
[0046][0047]
(5.7)实存状态的循环周期可通过和计算得出；
[0048][0049]
(5.8)模型的灰稳态概率最终可用式(7)表示，为实存状态的驻留时间与其循环周期之比；
[0050][0051]
得到转移状态概率矩阵后，通过式(3)-(7)的gspn稳态概率计算方法可求得的复合灰数表示，之后根据区间灰数运算法则，可得到关于扰动灰元δ的区间灰数表示。由于区间灰数之间的运算会增大结果的不确定度，随着δ的不断增大，的上界或下界必然会落到区间[0,1]之外。
[0052]
步骤(6)、在得到的灰稳态概率的基础上构建数学规划模型，求解最终区间灰数表示；最后利用最终区间灰数表示对该工业cps的性能指标进行分析，具体过程如下：
[0053]
(6.1)在满足实际生产条件(变迁激发速率大于0、随机开关触发概率在0到1之间)以及所有均为可行解的情况下，求出δ的最大值；
[0054]
max δ
ꢀꢀꢀ
(8)
[0055][0056]
0≤δ
[0057][0058][0059]
(6.2)将(8)式所求得的δ代入即可求得的最终区间灰数表示；
[0060]
(6.3)根据具体cps情境，利用求得的最终区间灰数表示对该工业cps的性能指标进行分析。
[0061]
本发明的另一个目的是提供一种电子设备，其特征在于，包括处理器和存储器，所
述存储器存储有能够被所述处理器执行的机器可执行指令，所述处理器执行所述机器可执行指令以实现所述的方法。
[0062]
本发明的又一个目的是提供一种机器可读存储介质，其特征在于，该机器可读存储介质存储有机器可执行指令，该机器可执行指令在被处理器调用和执行时，机器可执行指令促使处理器实现所述的方法。
[0063]
本发明的有益效果：
[0064]
本发明用于诸如工业智能制造、航空航天、智慧医疗等的cps系统的建模与分析。该发明提出了面向方面灰色广义随机petri网模型，并使用其进行cps的形式化建模，能够明确表述cps的逻辑结构及功能流程。本发明采用了灰色理论、马尔科夫链和数学规划的分析方法，有效地处理了cps中的不确定性，从而使结果具有高度的精确性。
附图说明
[0065]
图1 iprpl-cps中的部分设备图，其中(a)为贴膜机，(b)为折弯机以及传感器，(c)为压铆机，(d)为码垛机器人；
[0066]
图2 iprpl-cps的核心生产流程图；
[0067]
图3 iprpl-cps的基本网bn
[0068]
图4 iprpl-cps的两个方面网an1和an2，其中(a)为方面网an1，(b)为方面网an2；
[0069]
图5时间变迁的冲突和解决，其中(a)为时间变迁发生冲突，(b)为瞬时变迁优先发生，(c)为随机开关；
[0070]
图6通知类型为after的方面网的织入过程；
[0071]
图7 iprpl-cps的编织网wn；
[0072]
图8 wn的可达标识图；
[0073]
图9灰稳态概率的图形表示；其中(a)为(b)为(c)为(d)为(e)为(f)为(g)为(h)为图10设备利用率。
具体实施方式
[0074]
下面结合具体实施例子对本发明进行进一步的描述。
[0075]
杭州西奥电梯有限公司(www.xiolift.com，xiolift)的智能板材加固产线(intelligent plate reinforcement production line,iprpl)是一个典型的工业cps，其中包含了贴膜机、折弯机、压铆机、码垛机器人等，部分设备如图1所示。其具体工作过程如下：
[0076]
第一步：贴膜机对板材进行贴膜，同时折弯机对加强筋原料进行折弯；
[0077]
第二步：等待板材贴膜完成和加强筋折弯完成后，压铆机对板材进行压铆。同时，用激光传感器和工业相机对加强筋进行检测以判断是否符合标准，若符合标准，则用压铆机对加强筋进行压铆；否则，折弯机需要依据偏差对其进行校；。
[0078]
第三步：码垛机器人将压铆完成的板材进行码垛，然后该产线将从第一步开始重新迭代。
[0079]
上述步骤可形式化表述为图2所示。
[0080]
本发明所述的基于面向方面灰色广义随机petri网的cps建模与分析方法，在其具体实施之前，相关概念定义及符号说明如下：
[0081]
(1)aggspn网系统定义为一个二元组aggspn＝(bn,ans)。其中：bn表示基本网；ans＝{an1,an2,
…
,ann}则表示一个有限的方面网集合。
[0082]
(2)由于需要描述工业cps中各活动的耗时以及资源、信号的流动，因此，aggspn网系统中bn和ans各子网均以gspn为基础，并结合了颜色集的概念，其中，基本网bn对应着工业cps中的核心生产流程，其定义如下：
[0083]
2-1aggspn网系统中基本网bn定义为一个六元组bn＝(p,t,f,c,m0,λ)，其中：
[0084]
p＝{p0,p1,
…
,pn}表示库所的有限集合；
[0085]
t＝{t0,t1,
…
,tn}表示活动变迁的有限集合，活动变迁t被划分为时间变迁集t
t
＝{t0,t1,
…
,tk}和瞬时变迁集ti＝{t
k 1
,t
k 2
,
…
,tn}两个子集，即t＝t
t
∪ti，时间变迁的发生时延服从于指数的随机分布，瞬时变迁的发生时延为零；
[0086]
f＝{f0,f1,
…
,fn}表示连接库所与活动变迁之间弧的有限集合；
[0087]
c＝{c0,c1,
…
,cn}表示颜色的有限集合；
[0088]
m0：p
→
{0,1,2,
…
}表示petri网的初始标识，是从p到非负整数的映射，对于m0(pi)表示m0标识下库所pi中token的数量；
[0089]
λ＝{λ0,λ1,
…
,λn}表示与时间变迁集t
t
关联的变迁激发率集合。
[0090]
虽然基本网bn描述了工业cps的总体运行流程，但其中某些步骤需要进一步细分(例如码垛过程需要进行细分)或者在其前后进行额外处理(例如在压铆前需对折弯后的加强筋进行检测)，因此，引入方面网ani代表工业cps核心功能流程中的某个步骤的具体细节或者额外处理，其定义如下：
[0091]
2-2方面网集合ans中的方面网ani定义为一个三元组ani＝(pc,ad,vn)，其中：
[0092]
pc表示ani对应于bn的切入点，可为库所或者活动变迁；
[0093]
ad＝{before,after,around}表示ani对应于bn的通知类型，before表示ani需在切入点前执行，after表示ani需在切入点后执行，around则表示用ani来替代切入点；
[0094]
vn表示一个六元组vn＝(p,t,f,c,m0,λ)，其定义与bn相同。
[0095]
2-3每一个方面网ani都对应着基本网中某一步骤，为了方便后续将其织入基本网，在这方面网中额外增加两个辅助元素：起始元素start与终止元素end，表示该方面网织入基本网时的起止点。
[0096]
2-4工业cps中涉及到资源(例如iprpl-cps中的加强筋等)以及信号(例如iprpl-cps中加强筋折弯的检测结果)的流动，因此在aggspn中，对token进行赋色以区分其类型，并借由活动变迁的发生规则来描述工业cps中的动态行为。在定义变迁的发生规则之前，元素的前集与后集定义如下：
[0097]
设pn＝(p,t,f)为一个petri网，对于任意的元素x∈p∪t，若存在元素y∈p∪t，使得弧(y,x)∈f，则将所有满足条件的元素y称之为x的前集，记为
·
x＝{y|(y,x)∈f}。同样的，使得满足条件(x,y)∈f的所有元素y称之为x的后集，记为x
·
＝{y|(x,y)∈f}。
[0098]
2-5为每条弧都关联一个表达式w：f
→
c，表示该弧所对应的颜色，弧的颜色与相连库所的颜色相同，即w(pi,ti)/w(ti,pi)＝c(pi)，则变迁ti的发生规则如下：
[0099]
对于若库所pi含有token且颜色为w(pi,ti)，则满足ti的发生条件，可将
token从pi中移除，并使ti使能。等待ti完成后，对于将token的颜色改为w(ti,pi)并将其放入库所pi中。
[0100]
(3)为了提高所建立模型的真实性和准确性，在模型中引入了灰色系统理论，利用区间灰数来表示模型中的变迁激发率以及随机开关的触发概率。区间灰数是最常见的灰数，并且其用法比较简单，可以表明实际生产中某个变量的波动范围，其定义如下：
[0101]
区间灰数是指既有上界又有下界的灰数，被记为a≤b。其中a是区间灰数的下界，b是区间灰数的上界，若a＝b，则该区间灰数就是一个白数，即白数是一种特殊的区间灰数。
[0102]
设有两个区间灰数与它们之间相互运算后，结果仍为区间灰数。区间灰数的运算法则如下所示，其中运算时，c≠0,d≠0；
[0103][0104][0105][0106][0107][0108][0109][0110]
δ：灰化处理算法中的扰动灰元；
[0111]
kj：活动变迁系数，时间变迁的扰动灰元系数设置为k1，随机开关的扰动灰元系数设置为k2；
[0112]
pr＝{pr
k 1
,pr
k 2
,
…
,prn}表示瞬时变迁ti为随机开关时的触发概率集合；
[0113]
τ＝{τ0,τ1,
…
,τk}表示时间变迁对应的平均实施时间集合；
[0114]
(7)马尔科夫链相关概念：
[0115]
mi：为可达标识集s中某一个状态，若此状态下有瞬时变迁处于使能状态，则mi为消失状态，否则为实存状态；
[0116]
s：gspn的状态空间集，元素个数ks，s可分为实存状态集ts和消失状态集vs，即s＝ts∪vs，其元素个数分别为k
t
和kv，ks＝k
t
kv；
[0117]
gspn的同构马尔科夫链(emc)的转移概率灰矩阵；
[0118]
由从消失状态向消失状态集转移的转移概率矩阵与从消失状态向实存状态集转移的转移概率矩阵构成；
[0119]
由从实存状态向消失状态集转移的转移概率矩阵与从实存状态向实存状态集转移的转移概率矩阵构成；
[0120]
emc由消失状态集向消失状态集的转移概率矩阵；
[0121]
emc由消失状态集向实存状态集的转移概率矩阵；
[0122]
emc由实存状态集向消失状态集的转移概率矩阵；
[0123]
emc由实存状态集向实存状态集的转移概率矩阵；
[0124]
将消失状态从emc中移除后，仅含有实存状态的emc(reduced emc,remc)，为remc的转移概率矩阵；
[0125]
其中的元素表示gspn从某个消失状态r开始，经任意步首次到达实存状态k的概率；
[0126]
表示一个一维行向量，元素表示remc中实存状态i的灰稳态概率分布；
[0127]
实存状态的平均驻留时间；
[0128]hi
：状态i下可实施的变迁集；
[0129]
实存状态的循环周期；
[0130]
系统的循环周期；
[0131]
模型的灰稳态概率。
[0132]
本发明具体实施步骤如下：
[0133]
步骤(1)根据所建模cps的具体情境，将工业cps抽象为一个基本网和若干个方面网：
[0134]
(1.1)根据工业cps的生产流程，定义一个基本网和若干方面网,基本网对应核心生产运行流程；
[0135]
(1.2)由于其中某些步骤需要进一步细分或者在其前后进行额外处理，定义若干方面网代表工业cps核心功能流程中的某个步骤的具体细节或者额外处理；
[0136]
步骤(2)构建aggspn模型，具体过程如下：
[0137]
(2.1)根据iprpl-cps的结构流程，将其划分为一个基本网bn与两个方面网an1和an2；针对基本网，以总运行流程中工业cps生产线上的工业设备状态作为gspn的库所，以总运行流程中工业cps生产线上的资源或者信号作为库所的token，以总运行流程中工业cps生产线上的工业设备动态行为作为gspn的活动变迁；根据总运行流程中库所与活动变迁间的关系作为弧；并对不同的token标注不同的颜色进行区别；(上述网中库所、活动变迁、token、弧均为gspn的元素)；针对方面网，以子流程中工业cps生产线上的工业设备具体状态作为gspn的库所，以子流程中工业cps生产线上的资源或者信号作为库所的token，以子流程中工业cps生产线上的工业设备具体动态行为作为gspn的活动变迁；根据子流程中库所与活动变迁间的关系作为弧；并对不同的token标注不同的颜色进行区别；
[0138]
(2.2)aggspn模型中，每个子网都以gspn为基础，基本网bn如图3所示，其库所与变迁的具体含义如表1所示，两个方面网an1和an2如图4所示，其库所与变迁的具体含义如表2所示；
[0139]
表1 bn中的库所与变迁及其含义
[0140][0141][0142]
表2 an1和an2中的库所与变迁及其含义
[0143][0144]
步骤(3)提出一种编织算法，以便于性能分析；
[0145]
对于任意一个方面网，其织入过程主要分为两步：首先将该方面网与基本网连接到一起，并删除多余的辅助元素；之后合并该方面网与基本网中的相同元素。所述方面网的多余辅助元素包括起始元素start、终止元素end和重复的库所。具体步骤如下：
[0146]
(3.1)若该方面网的通知类型为before，则将该方面网插入到基本网中切入点与其前集之间，替代其原先的输入弧集；若该方面网的通知类型为after，则将方面网插入到基本网中切入点与其后集之间，替代其原先的输出弧集；若方面网的通知类型为around，则将该方面网替代基本网中的切入点，并与其前集和后集相连。之后删除辅助元素以及相关
的弧。例如图6所示的是一个通知类型为ofter，切入点为t1的方面网的织入过程。
[0147]
(3.2)初步连接方面网与基本网后，若方面网与基本网含有相同库所p，则需将方面网中库所p的输入弧和输出弧转化为基本网中库所p的输入弧和输出弧，并且删除方面网中的库所p。例如图3的bn与图4(a)的an1均含有库所p1，则删除an1.p1以及弧(an1.t
10
,an1.p1)与(an1.p1,an1.t9)，并增加新弧(an1.t
10
,bn.p1)和(bn.p1,an1.t9)。
[0148]
按照该编织算法将iprpl-cps的基本网与所有方面网编织到一起后，得到编织网如图7所示，其中an1的通知类型为before，所以织入至切入点p6之前，an2的通知类型为around，则作为切入点t6的替代。因此最终的编织网wn不含t6，bn与an1均含有p1，所以进行了合并。
[0149]
图8即为iprpl-cps编织网对应的可达标识图，其中深色代表实存状态，浅色代表消失状态，椭圆表示该状态为一个中间状态，矩形则表示该状态为初始状态。
[0150]
表4 wn的可达标识集s
[0151][0152]
在aggspn中，如果多个时间变迁同时使能，则会发生冲突，如图5(a)所示，t1、t2和t3三个时间变迁发生冲突。但若有瞬时变迁存在并使能时，瞬时变迁会优先发生，如图5(b)所示，t2优先于t1和t3发生。因此，可以在各个时间变迁前增加一个瞬时变迁以消除冲突，这种结构被定义为随机开关，如图5(c)所示，t1、t2和t3构成了一个随机开关，触发概率分别为pr1、pr2和1-pr
1-pr2，解决了t4、t5和t6产生的冲突。这种结构可在cps建模中用于描述一些事件的逻辑选择概率。例如在iprpl-cps中，需要检测加强筋折弯和板材压铆是否合格，合格率则需要根据实际的生产数据统计确定。因此，iprpl-cps中共定义了两组随机开关，分别为an1中的t8与t9以及an2中的t
12
与t
13
。
[0153]
步骤(4)引入灰色系统理论，利用灰化处理方法求解步骤(3)编织后的aggspn模型中的变迁激发率以及随机开关的触发概率，具体做法是通过引入一个扰动灰元δ将模型的变迁激发率和触发概率转化为区间函数，为每个变迁设置一个系数k。具体做法是将时间变迁的扰动灰元系数设置为k1，随机开关的扰动灰元系数设置为k2，具体过程如下：
[0154]
(4.1)若aggspn中的活动变迁t1为时间变迁，变迁激发率为λ1，则其区间灰数表示[λ
1-k1δ,λ1 k1δ]；
[0155]
(4.2)同理，若aggspn中的活动变迁t2为一个随机开关，触发概率为pr(t
12
)，其区间灰数表示[pr(t
12
)-k2δ,pr(t
12
) k2δ]；
[0156]
通过调阅iprpl-cps的生产数据以及对某些生产步骤耗时进行实地测量，可得到模型中时间变迁的变迁激发率以及随机开关的触发概率。其中时间变迁{t2,t3,t5,t
10
,t
14
,t
15
}的平均实施延时为：τ2＝1，τ3＝1/3，τ5＝1/3，τ
10
＝1/4，τ
14
＝1/6，τ
15
＝1/6(单位为min)，即对应的变迁激发率为：λ2＝1，λ3＝3，λ5＝3，λ
10
＝4，λ
14
＝6，λ
15
＝6；随机开关t8和t9的触发概率分别为pr1＝0.95，pr2＝0.05，随机开关t
12
和t
13
的触发概率分别为pr3＝0.98，pr4＝0.02。结合实际生产情况，将时间变迁激发率的扰动灰元系数k1设置为其本身的10％，随机开关触发概率的扰动灰元系数k2设置为0.01，之后使用该灰化处理算法对这些值进行灰化处理，得到其区间灰数表示，如表3所示。
[0157]
表3时间变迁激发率与随机开关概率的区间灰数表示
[0158][0159][0160]
步骤(5)使用马尔科夫方法求解模糊gspn的稳态概率。具体过程如下：
[0161]
使用马尔科夫方法求解aggspn的灰稳态概率；具体过程如下：
[0162]
(5.1)分析编制网得到可达标识图以及对应的可达标识集；
[0163]
(5.2)根据可达标识图可得到emc的转移概率灰矩阵其计算公式如下：
[0164][0165]
(5.3)计算remc的转移概率矩阵u'，用式(2)表示；
[0166][0167]
以iprpl-cps为例，可得到该gspn同构的remc的转移概率矩阵
[0168][0169]
(5.4)通过式(3)的线性方程组，求解remc的实存状态灰稳定概率；
[0170][0171]
表示一维行向量，元素表示remc的实存状态的灰稳态概率；
[0172]
(5.5)通过式(4)计算各个实存状态的平均驻留时间
[0173][0174]
(5.6)由各个状态的平均驻留时间以及对应的灰稳定概率可求得系统的循环周期计算方式为式(5)；
[0175][0176]
(5.7)实存状态的循环周期可通过和计算得出；
[0177][0178]
(5.8)模型的灰稳态概率最终可用式(7)表示，为实存状态的驻留时间与其循环周期之比；
[0179][0180]
得到转移状态概率矩阵后，通过式(3)-(7)的gspn稳态概率计算方法可求得的复合灰数表示，之后根据区间灰数运算法则，可得到关于扰动灰元δ的区间灰数表示。由于区间灰数之间的运算会增大结果的不确定度，随着δ的不断增大，
的上界或下界必然会落到区间[0,1]之外。
[0181]
步骤(6)在得到的稳态概率的基础上构建数学规划模型，并求解最终准确的数值结果，具体过程如下：
[0182]
(6.1)在满足实际生产条件(变迁激发速率大于0、随机开关触发概率在0到1之间)以及所有均为可行解的情况下，求出δ的最大值；
[0183]
max δ
[0184][0185]
0≤δ
[0186][0187][0188]
(6.2)将(8)式所求得的δ代入即可求得的最终区间灰数表示；
[0189]
通过以上步骤可求得其关于扰动灰元δ的区间灰数表示，其几何表示如图9所示。
[0190]
通过图10可发现，只需的上界则其余均为可行解。则根据步骤4中的数学规划模型，可求得δ的最大值为1.473，将其代入各个中得到其区间灰数表示。将δ＝0代入各个中即得到其白化值，最终得到的灰稳态概率区间灰鼠数表示及其白化值如表5所示。
[0191]
表5灰稳态概率的区间灰数表示及其白化值
[0192][0193]
(6.3)根据具体cps情境，利用求得的灰稳态概率对该工业cps的性能指标进行分析。
[0194]
得到aggspn的灰稳态概率后，即可进行所建模系统的多种指标的分析，例如设备利用率、变迁吞吐率、生产率等。
[0195]
以设备利用率为例。在aggspn中，某些库所pi指代具体的设备是否可用，当pi中含有token时，则表示该设备正在空闲中，若pi中不含token时，则表示该设备正在工作中。因此，该设备的利用率可表示所有m(pi)＝0的稳态概率之和。
[0196]
在iprpl-cps中，折弯机的利用率为cps中，折弯机的利用率为贴膜机的利用率为压铆机的利用率为码垛机器人的利用率为码垛机器人的利用率为各设备的利用率如图10所示。由图可知，码垛机器人的利用率最低，折弯机的利用率最高，贴膜机与压铆机的利用率相差无几，但均远低于折弯机。折弯机的利用率最高是因为其对加强筋进行折弯以及校准所需时间最多，若长期维持高利用率可能会因负荷过高导致机器故障。