不完全维修条件下的备件多级库存优化

寇贞贞，李苏剑，顾 涛，吴秀丽

(北京科技大学 机械工程学院，北京 100083)

0 引言

在装备保障研究中，可维修备件库存研究是非常重要的一部分，多等级保障站点和多层级装备结构的复杂情况增加了备件的库存分配难度，而备件维修情况也影响着各级基地的库存。如何在保障装备系统正常运行的情况下合理安排备件库存、提高其可用度并最大幅度地降低费用是研究人员面临的难题。

许多研究人员对可修件多级库存优化问题从不同方面进行了研究[1-2]，其中SHERBROOKE[3]建立了最广泛应用的可维修备件多级管理技术(Multi-Echelon Technology for Recoverable Item Control, METRIC)模型，但模型有很多假设条件过于理想化，如最多考虑两级保障站点结构、要求备件需求稳定、不考虑横向转运、不考虑报废等条件。后来的研究人员在假设条件上有所放宽，所建模型更加贴合实际：模型由单保障等级结构发展为多级库存保障模型[4-5]；王慎等[6]放宽了“完全串件系统”假设，建立了动态管理模型；RUAN等[7]放松了对备件稳定需求的假设，建立了动态配置模型；RICCARDO[8]建立了带有横向转运的多级库存优化模型；在备件的报废研究方面，虽然目前对备件报废量已经有所考量，但都是简单使用报废率和故障率衡量，计算较为粗略[9-11]。另外，METRIC模型中假设备件每次维修后效果“修复如新”，且可以进行无数次维修，然而实际维修过程往往是不完全维修[12]，即“修复如新”和“修复如旧”之间的某一状态，并且在多次不完全维修后备件达到最大维修次数，可能会出现报废的情况，文献[13-17]考虑了设备维修的不完全维修，但上述文献未探讨在不完全维修条件下因有限维修次数限制导致备件报废后，可维修备件多级库存优化模型的建立。

在求解算法方面，边际优化算法是求解该问题最常见的一种方法[18]，该算法虽然简单易懂、应用广泛，但实际上是一种贪婪算法，由于其自身的原因极易陷入局部最优，并不一定能得到最优解。目前，用于求解该问题的算法研究还不是很多，遗传算法、粒子群算法虽有应用[19-20]，但应用时容易陷入局部最优或收敛缓慢，效果不佳。

综上所述，对于可修件备件多级库存优化问题，虽已建立多种模型，但尚未考虑不完全维修条件下有限维修次数约束导致报废以及与备件故障分布结合的情况，求解算法也较为单一。为此，本文提出考虑不完全维修条件下有限维修次数约束的可维修备件库存优化模型，并设计粒子群—边际禁忌混合算法(Particle Swarm Optimization-Marginal Tabu Search hybrid optimization algorithm，PSO-MTS)求解模型。

1 可修件多级库存优化问题

本问题适于各类大型复杂装备维修系统，尤其是对于军用装备来说，可提供合理的库存管理决策依据。可修件多级库存优化问题为：已知装备层级结构和保障站点等级结构，合理安排各个站点j各类备件i的库存Sij，使得系统在保证一定的可用度Am前提下总费用最低。为便于研究，作如下假设：

(1)模型针对价值贵重的装备备件；

(2)装备备件保障独立，不存在任何借用关系；

(3)故障件的维修时间相互独立；

(4)各保障等级的备件库存对策均为(s-1，s)，即one-for-one订购策略；

(5)备件只在最高级别库存保障等级站点处报废。

1.1 总体流程

假设保障结构等级有多级，备件层级有两级，其中外场更换件(Line Replaceable Unit, LRU)是主要的部件即第一层备件，内场更换件(Shop Replaceable Unit, SRU)是LRU的组件，即第二层备件，每个LRU的SRU数量视装备实际情况而定。

以如图1所示的三级库存保障结构为例说明其流程。保障站点有基层级、中继级、基地级共3级，整个维修过程开始于LRU发生故障并送到基层仓库，若基层仓库有一个LRU备件，则将其替换，否则基层级发生一次短缺。LRU故障件在基层级修理只占很小的比例，如果可以在基层级修理，若现有库存有一件备份的SRU，则将其安装到LRU完成修理，经过不完全维修后，进入仓库成为新的备件储备，但之后的工作时间会逐次缩短，故障间隔时间逐渐减少，当备件经过多次维修后，其可靠性降低到一定程度，已无修理必要，则进行报废处理，此时消耗一件备件；如果LRU不能维修，就发往中继级送修，同时安排向中继级申请该LRU一件。中继级和基地级流程同理。

1.2 可修件多级库存优化模型

1.2.1 备件报废量的确定

(1)不完全维修过程

在实际应用中，一个备件的生命周期如图2所示，一个新品备件经历第一次投入使用，工作一段时间后，出现了第一次故障，随后进行故障件的更换拆卸,故障件进行第一次维修，并入库作为下一次更换的备件存储，储存一段时间后再次投入使用，历经相同的维修、存储过程，直至备件达到最大维修次数Nr，在第Nr+1次故障后不再维修，进行报废处理。

役龄回退因子η表示对某一个可修复系统进行维修后，实际服役年龄回退的程度，0≤η≤1。一次修复后，回退役龄为装备备件通过维修可以让服役年龄向前回退维修间隔期的η倍，Te为平均维修时间间隔，回退役龄te为：

te=η·Te。

(1)

(2)备件报废量的计算

对一个装备单元来说，第一次修复之前平均使用时间θ1为:

(2)

式中f(t)为故障概率密度函数，常见的故障分布有泊松分布、威布尔分布、指数分布、正态分布等等[21-22]。由于后面每一次的故障间隔时间均与前一次有关，通过z维空间取值范围内的z次概率密度函数进行多重积分，可以得到该备件在第z-1次修复和第z次修复之间的平均使用时间θz[23]为：

(3)

(4)

其中：Nr为最大维修次数，该备件经历了Nr+1次故障后报废，至此一个备件经历了一个完整的生命周期，备件消耗量加1。整个生命周期内的总平均使用时间为θ0，Nm为装备数量，Z0为装备单元单机安装数量，当装备工作了T时间后，该备件报废量md为：

(5)

1.2.2 优化模型

用i表示SRU的项目编号i=1,2,3,…,I，i=0表示LRU；用k表示保障站点编号，h表示该保障站点的上级站点，j表示该保障站点的下级站点。其他参数定义如表1所示。

表1 参数定义表

由于系统的特殊性，停机损失巨大，将可用度作为约束条件，成本费用为优化目标(由于备件价格较高，运输、库存成本在本模型中可忽略)，则所建立的备件库存优化模型为:

(6)

(7)

设装备由N个LRU部件组成，任何LRU失效都会导致装备出现故障，因此装备可用度为:

(8)

式(8)反映了系统可用度与备件短缺数期望值的关系，求最大可用度等同于求备件短缺数期望值最小。式中BO(s0j|E[X0j],Var[X0j])/NjZ0表示供应渠道均值为E[X0j]、方差为Var[X0j]、基层级库存为S0j时的备件短缺数期望值。由基层级开始对各个站点的各类备件依次计算其需求，顺序递推到中继级、基地级等，再计算最高级别的供应渠道短缺数均值和方差，反向递推到基层级即可得到EBO(s0j|E[X0j],Var[X0j])/NjZ0，具体如下：

(1)各级备件需求的确定

各站点的需求由两部分组成：该保障站点负责的所有下级站点不能维修的故障件之和；该站点维修LRU时所需的SRU需求。

(9)

(10)

(11)

(12)

式(9)为站点k处LRU的需求率，仅包含下级站点不能维修的故障件数量之和；式(10)～式(12)为站点k处SRU的需求数量，除不能维修的故障件数量之和还要加上该站点维修LRU时所需的SRU数量。

(2)各级供应渠道备件数的期望和方差的确定

供应渠道备件数主要由4部分构成：①该保障站点的在修件数量；②正在对该站点进行供应和补给的数量；③故障件LRU在该站点维修时，因等待SRU维修而造成延误的数量；④因报废需要补充的数量。具体如下：

fik=mik(1-rik)/mih;

(13)

ρik=m0k·qik·r0k/mik;

(14)

f0k=m0k(1-r0k)/m0h;

(15)

E[Xik]=mik·[(1-rik)·tik+rik·Tik]+

fik·EBO(sih|mih·Tih);

(16)

Var[Xik]=mik·[(1-rik)·tik+rik·Tik]+

fik·(1-fik)·EBO(sih|mih·Tih)+

(17)

E[X0k]=m0k·[(1-r0k)·t0k+r0k·T0k]+

(18)

Var[X0k]=m0k·[(1-r0k)·t0k+r0k·T0k]+

ρik·(1-ρik)·EBO(sik|E[Xik],Var[Xik])+

(19)

其中：式(13)～式(15)用于计算因故障件维修造成的需求占总需求的比例;式(16)～式(17)用于计算站点k处SRU的期望供应渠道短缺数和方差;式(18)～式(19)为站点k处LRU的供应渠道期望短缺数和方差。

2 粒子群—边际禁忌混合优化算法

2.1 算法描述

可修件多级库存优化问题是NP难题，传统方法采用边际优化算法进行求解[24]。但边际优化算法实际上是一种贪婪算法，考虑函数在边际点上的极值取舍最佳点，很容易陷入局部最优，且循环迭代次数多、计算量大，因此设计了PSO-MTS算法解决此问题。首先用粒子群算法快速得到一个全局较优的解，一旦判断粒子群连续迭代过程中结果无较大改进后，在算法后期引入禁忌搜索开发局部搜索能力；同时考虑到边际优化算法在本问题中的应用广泛、求解方便，在禁忌搜索的邻域构建部分，也引入了边际优化算法以保证邻域的合理性。本算法结合粒子群算法快速的收敛速度和禁忌搜索算法的局部开发能力，同时避免了粒子群后期收敛速度变慢的问题并提供给禁忌搜索一个较好的初始解。

算法流程如图3所示，具体步骤如下：

步骤1初始化。包括粒子群规模N1；最大迭代次数maxgen；维度D；初始惯性权重w1、迭代到最大次数时的惯性权重w2；粒子位置取值范围[Xmin，Xmax]；粒子最大速度Vmax、最小速度Vmin；学习因子c1、c2；禁忌表长度list；聚集距离阈值Borderdist；连续不变化阈值maxstep。随机生成粒子并初始化位置和速度。

步骤2根据适应度函数评估各粒子的初始适应度即备件成本费用，并将其赋值给pBest以记录每个粒子的历史最优值，同时将最佳适应度的粒子赋值给gBest即全局最优值。此时迭代次数gen=0。

步骤3若满足迭代条件gen

步骤4评估新粒子群的适应度值，依次比较各个粒子与其历史最优和全局最优的适应度，若更优则替换pBest、gBest，否则保持不变。

步骤5为防止算法陷入局部最优，对种群最佳位置的粒子重新初始化赋值。

步骤6计算粒子的平均聚集距离和连续不变化步数，以此判断是否可以进入禁忌搜索，若无，依据公式计算各个粒子的速度和位置并更新，转步骤3；否则，转步骤7。

步骤7进入禁忌搜索，给定前述步骤的结果gBest作为禁忌搜索部分的初始解，置禁忌表为空。

步骤8判断算法是否满足终止条件，若是则结束算法并输出优化结果gBest；否则，转步骤9。

步骤9利用当前解的邻域函数结合边际优化算法产生满足条件的邻域解，并从中确定若干候选解。

步骤10判断候选解是否满足藐视准则，若成立，则用满足藐视准则的最佳状态替代gBest成为新的当前解，并替换禁忌对象，然后转步骤12；否则，转步骤11。

步骤11判断候选解对应的各对象的禁忌属性，选择候选解集中非禁忌对象对应的最佳状态为新的当前解，同时用与之对应的禁忌对象替换最早进入禁忌表的禁忌对象元素。

步骤12若未达到结束条件，转步骤8，否则输出gBest并结束。

2.2 算法详细步骤

2.2.1 速度和位置更新

粒子的速度和位置更新公式如下：

(20)

(21)

动态调整的惯性权重和学习因子更利于算法的收敛和寻优，惯性权重将由式(22)计算得出，其中w1=0.9，w2=0.4。

(22)

粒子在搜索过程中，初期希望学习因子速度大一些，尽量搜索到整个空间，后期希望粒子速度变小，从而到达精确位置，对学习因子的改进公式如下[25]：

(23)

2.2.2 平均聚集距离的判定

利用平均聚集距离来判断群内粒子聚集程度，平均聚集距离为群内所有粒子到历史最优位置的欧几里得空间距离：

(24)

聚集距离阈值BorderDist为判断群内粒子聚集程度的距离阈值，计算如式(25)所示：

(25)

式中XMaxd是粒子位置X每一维度对应的最大值，Dist为一常数，可取1 000、100、10等。

在算法运算过程中，记录粒子历史最优位置连续不变化或者变化极小的迭代次数logjamstep，同时设置连续不变化次数的阈值maxstep。在粒子群的历史最优粒子位置连续无变化或变化极小时，若粒子的聚集情况较为严重、迭代过程中解连续没有改进，即logjamstep>maxstep同时Meandist

2.2.3 禁忌搜索邻域设计

粒子群算法在后期极易陷入局部最优，若迭代过程中解连续没有改进，则进入禁忌搜索算法。对于邻域的设计，传统的禁忌搜索算法通常用两点交换等方法产生邻域，因其操作简单易于实现，但在本算法中，需要保证系统可用度大于目标可用度，而可用度计算并非简单的线性计算，因此在边际优化的思想基础上进行邻域的设计，这样在保证系统可用度的情况下既能增加邻域设置的合理性，又能增强局部搜索能力。下面设计5种邻域的产生办法：

(1)边际取值法

1)随机选取两点，将其数值置为0；

2)边际优化计算。

对于备件多级库存这类问题，边际优化算法通过对边际单元的效益和费用进行权衡分析，达到对有效资源的合理利用。费效比ε定义为相邻期望短缺数的减少值与相应备件费用之比，

(26)

边际分析的主要思想是在算法中每次只增加一个备件，判断增加哪个备件对边际效益影响最大，即费效比最大，最大的备件数量相应增加一个，其他保持不变。依次类推，直到费用不够购买下一个备件为止。此时装备备件的总短缺数即为各项备件对应的期望短缺数之和。边际优化算法流程如图4所示。

对于处理后的粒子，按照边际优化算法计算其约束函数即可用度，对每一位置依次加1并计算可用度，直到系统可用度达到目标可用度时，可停止计算，此时得到的粒子视为新邻域。

(2)邻域准则取值法

1)设定邻域准则[1 2 3 -1 -2 -3]，选择邻域准则中的一个随机数对随机产生的位置点库存进行加减；

2)边际优化计算以保证邻域合理。

(3)两点交换取值法

1)随机产生两个位置，两点上的数值交换，其他部分不变，得到新粒子；

2)边际优化计算以保证邻域合理。

(4)插入取值法

1)随机产生两个位置，将靠后的位置点插入靠前位置点的前面，其余位置向后顺延；

2)边际优化计算以保证邻域合理。

(5)片段倒序取值法

1)随机产生两个位置，将两点之间的片段倒序；

2)边际优化计算以保证邻域合理。

在每次算法迭代中，均使用了以上5种邻域产生办法，由于邻域的产生存在一定的随机性，每次最优值对应的方法并不固定，需依次利用5种方法产生邻域，比较并保留最优值。

3 实验结果与分析

3.1 实验设计

实验主要验证所建模型的正确性和求解算法的效果，因此分为两个实验进行。PSO-MTS算法在Intel Core i5 1.6 GHz CPU、8 GB RAM、Windows 10操作系统和MATLAB编程环境下编译。

实验1用于验证模型的准确性：收集某个备件故障时间间隔Tei数据：x1,x2,x3,…,xn,将其分为等间隔的kn组，进行初步数据处理，计算各组中值mi及频数fi，可以得到各组平均故障概率密度f(mi)、可靠度R(mi)和累计故障概率F(mi)，对k组mi和f(mi)进行数据拟合，得到对应任一时间t的最佳故障概率密度分布公式f(t)以及其他可靠性指标，随后可计算出备件报废量，带入模型可得到库存矩阵及系统总费用，与完全维修、无维修次数限制的计算结果相比，可验证模型的正确性。

实验数据如下：在一个如图5所示的三级维修体系中，由一个基地级站点H0、两个中继级站点R1、R2、三个基层级站点J1、J2、J3组成，装备系统的层级结构如图6所示。

该装备在各个使用现场数量分别为14、18、12；装备平均每周工作时间H为40 h；站点向上申请备件及其运输时间Oj=0.01 年,Ok=0.02 年，报废后采购所需时间TL=0.05 年；备件单价C=[6 7 1 1 2 1 1]，单位：万元；单机安装数Z0=[2 2]，Zj=[1 2 2 1 2]；LRU1和LRU2最大维修次数分别为3、2次；役龄回退因子为0.8、0.7；其余参数见表2～表4。

表2 实验1：备件年需求数据

表3 实验1：平均维修时间

表4 实验1：维修率

续表4

实验2用于观察算法的效果，实验2-1：利用PSO-MTS算法求解文献[26]并对比结果；实验2-2：将实验1的数据分别利用边际优化算法、粒子群算法、PSO-MTS求解20次，对比优化结果。

3.2 实验1结果

3.2.1 完全维修、无维修次数约束的情况

当不考虑不完全维修和维修次数约束的影响时，利用PSO-MTS算法测试20次。算法参数如下：初始化最大迭代次数maxgen=2 000；种群规模N=20；维度D=6×7=42；种群生成方式为randi函数随机生成整数矩阵，经过多次实验得到库存值均小于10，因此初始化x=randi([0,10],N,D)；禁忌表长度list=⌊42⌋+1=7；聚集距离阈值BorderDist的计算中，经过测试，Dist取18为最佳；连续不变化阈值maxstep=10。

规定整个装备系统所有装备可用度不低于0.95，得到20次的平均费用为344.2 万元，取最优结果得到系统可用度为0.951 5，费用341 万元，库存配置情况如表5所示。

表5 实验1：完全维修时的库存配置表

3.2.2 考虑不完全维修和维修次数约束的情况

对采集到的故障时间间隔数据进行初步处理如表6和表7所示，对数据进行拟合，分别得到两备件最佳分布，LRU1、LRU2故障概率密度分别为：f(x)=0.001 713×e(-0.001 564×x)，f(x)=0.002 859×e(-0.002 765×x)。

表6 实验1：LRU1数据处理

表7 实验1：LRU2数据处理

由式(3)得LRU1和LRU2平均工作时间如表8所示。

表8 实验1：平均工作时间 h

根据式(5)得备件报废量：

(27)

(28)

将报废量带入模型，规定整个装备系统所有装备可用度不低于0.95，用粒子群—边际禁忌混合算法PSO-MTS求解20次，算法参数同3.2.1节，平均费用为515.1 万元，最优解系统可用度为0.950 3，费用508 万元，库存配置情况如表9所示。

表9 实验1：考虑不完全维修的库存配置

库存配置柱状图如图7所示，可以看出，对比不考虑不完全维修和维修次数限制的情况，在同样的目标可用度下，由于后续报废备件的产生，LRU的储备均有所增加，且基地级库存增加幅度显著，总费用由341万元提高到508万元，符合实际情况。

3.3 实验2结果

3.3.1 实验2-1

为了验证PSO-MTS算法的准确性，本节选用文献[26]的算例进行对比验证，具体数据参见文献。算法参数如下：初始化最大迭代次数maxgen=2 000；种群规模N=20；维度D=6×9=54；种群生成方式为randi函数随机生成整数矩阵，经过多次实验得到库存值均小于10，因此初始化x=randi([0,10],N,D)；禁忌表长度list=⌊54⌋+1=8；聚集距离阈值BorderDist的计算中，经过测试，Dist取18为最佳；连续不变化阈值maxstep=10。目标可用度为0.95，用PSO-MTS算法求解20次，取最优结果得到系统可用度为0.950 3，费用315.2 万元，算法收敛曲线如图8所示，库存配置如表10所示。

表10 实验2-1：PSO-MTS计算库存配置

该算法测试20次得到的结果同文献[26]得到的结果对比如表11所示。可以看出，粒子群—边际禁忌算法结果均更优，平均值318.07万元，优化幅度为5.9%，且最优解315.2 万元时改进幅度达到6.74%，表现较好。

表11 实验2-1：PSO-MTS算法及文献算例结果对比

3.3.2 实验2-2

分别利用粒子群(Particle Swarm Optimization, PSO)算法、边际优化算法和PSO-MTS算法来求解实验1中不完全维修条件下有维修次数限制的模型，测试算法20次，结果对比如表12所示。

可以看出，虽然粒子群算法部分解相比传统的边际优化算法有了一定的改进，但结果不稳定浮动较大，粒子群—边际禁忌算法的结果改进幅度更大且更稳定。相比粒子群算法，20次结果均为更优；相比边际优化算法，平均值515.1万，优化幅度5.66%，最优解508万对应改进幅度达到了6.96%。从各自最优解的收敛曲线图9也可以看出，后半段禁忌算法的加入有效地改进了优化结果，该算法在解决可修件多级库存优化问题上表现更加优秀。

表12 实验2-2：边际优化、PSO算法、PSO-MTS算法求解结果对比

4 结束语

本文针对装备可修件多级库存优化问题建立了考虑不完全维修条件下有限维修次数的库存优化模型，并设计了粒子群—边际禁忌混合优化算法PSO-MTS进行求解，得出如下结论：

(1)针对可修件不完全维修条件下、有限维修次数约束的情况，以系统可用度为约束条件，以费用为优化目标建立了多保障等级、两装备层级的库存优化模型，与不考虑不完全维修及维修次数限制的情况相比较，所建模型正确且所得结果符合实际情况。

(2)提出了粒子群—边际禁忌混合算法，该算法前期以粒子群算法为主求得优质解，后期则以禁忌搜索优化，其中利用边际优化算法进行邻域构造是一大亮点，以此扩大搜索范围、保证邻域的合理性。采用该算法求解文献算例并与文献求解结果进行了对比，发现该算法优化效果明显，最优解改进幅度达到6.74%；采用该算法求解实际算例并与粒子群、边际优化算法比较，相比边际优化算法最优解优化幅度达到6.96%，平均改进5.66%，表现优异。

未来将针对备件多级库存优化问题作进一步研究，如考虑非稳定需求条件、多目标优化算法等。