在制造业生产系统中20%～50%的加工费用于物料运输，而合理的设施布局可使这一费用减少10%～30%[1]．设计一种优良的设施布置方案不仅可以减少生产系统的总体运行费用，还能加快物料处理效率，减少在制品的停留时间，提高企业生产效率，增强企业的市场竞争力，因此开展生产企业设施布局优化设计的研究具有重要的实际意义．在全球经济一体化的进程中，企业竞争愈发激烈，市场需求变化也越来越快．动态设施布局将整个生产周期划分为多个生产阶段，不同生产阶段设施之间物流量随着市场的变化而发生变化，综合考虑各个生产阶段的设施布局，因其对市场变化的快速响应能力而引起研究者的广泛关注[2]．近年来，研究者对不等面积动态设施布局问题(UA-DFLP)进行了研究．文献[3]将UA-DFLP转化为重布局过程和多个子生产阶段的静态设施布局问题，并采用带精英策略的非支配遗传算法(NSGA-II)进行求解．文献[4]采用多边形内近似方法将设施面积非线性约束线性化，提出了一种融合线性规划方法的问题演化算法．文献[5]将一些启发式的格局更新策略与基于直方图函数的Wang-Landau抽样算法相结合，为UA-DFLP提出了一种启发式设施布置方法．文献[6]基于柔性区带结构(FBS)对UA-DFLP进行了研究，提出了一种遗传算法(GA)．文献[7]提出了一种改进的粒子群(PSO)算法，通过两种局部搜索方法和设施交换方法进一步改进布局的质量．文献[8]建立了基于水平和垂直区域的动态设施布置模型，提出了一种结合模拟退火(SA)算法和变邻域搜索的混合算法．UA-DFLP传统的干涉性处理(即确保设施之间不重叠)方法主要包括设施分行放置策略[3]、启发式方法[4-5]、罚函数方法[7]和FBS策略[6,8]．但目前各种干涉性处理方法均存在明显缺陷：分行放置容易导致设施超出纵向的布局区域，启发式方法得到的布局结构并不紧凑，罚函数方法难以设置合理的惩罚因子故效果并不理想，FBS策略由于其特殊的区带结构而限制了布局的多样性．本研究针对UA-DFLP传统的干涉性处理方法的不足提出了一种拟物方法，即将设施与车间外部区域均想象为具有弹性的光滑实体，通过模拟弹性物体在挤压弹性力作用下不断运动而解决设施间的干涉性约束问题；同时，为提高布局优化效果，对传统禁忌搜索算法中禁忌对象与解的接收准则进行改进，提出了一种改进的禁忌搜索算法，实验证明该算法有效．1 数学模型动态设施布局问题是在对未来T个生产阶段产品需求可预测的前提下进行布局设计的，主要根据各阶段设施间物流量的变化调整，使得整体布局的成本最小．对于具有T个生产阶段的动态设施布局，每一阶段可看作静态设施布局，在布局时假设：a. 所要布局的设施为矩形，忽略其细节形状；b. 布局所在车间F为长和宽已知的矩形，并足以布置下所有设施；c. 重置设施所造成的重布局费用为固定值；d. 设施只能沿水平方向或垂直方向在车间内任意放置；e. 物料搬运的位置为设施的中心，设施间的物料搬运距离为曼哈顿距离．假设笛卡尔坐标系的原点位于矩形车间的左下角，给出描述模型的相关参数：lti和wti分别为设施i在t阶段的实际长和宽，其中，t∈{1，2，…，T}，i，j∈{1，2，…，N}，T为阶段总数，N为设施总数；dti为二元变量，dti=1表示设施i在t阶段水平放置，dti=0表示设施i在t阶段垂直放置；L和W分别为车间的长和宽．由以上问题假设，建立如下数学模型．a. 目标函数minE(X)=∑t=1T∑i=1N-1∑j=1,jiNctijftij(|xti-xtj| +|yti-ytj|)+∑t=2T∑i=1Nktirti.目标函数E(X)是最小化整个布局方案的物料搬运费用与设施重布局费用之和．X=(x11，y11，d11，x12，y12，d12，…，x1N，y1N，d1N；x21，y21，d21，x22，y22，d22，…，x2N，y2N，d2N；…；xt1，yt1，dt1，xt2，yt2，dt2，…，xtN，ytN，dtN ；…；xT1，yT1，dT1，xT2，yT2，dT2，…，xTN，yTN，dTN)为构型(解)，即布局方案；ctij和ftij分别为在阶段t中设施i和j之间的单位物料搬运费用和物流量；xti和yti为在t阶段中设施i的中心坐标；kti为t阶段设施i的重布局成本；rti为二元变量，rti=1表示设施i在t阶段初始时被重新布局(xt-1,i≠xti或yt-1,i≠yti或dt-1,i≠dti)，rti=0表示设施i的位置没有变化．b. 约束条件0.5(lti+ltj)-(xti-xtj)≤L(Stij+Vtij); (1)0.5(lti+ltj)-(xtj-xti)≤L(1-Stij+Vtij); (2)0.5(wti+wtj)-(yti-ytj)≤W(1+Stij-Vtij); (3)0.5(wti+wtj)-(ytj-yti)≤W(2-Stij-Vtij)；(4)xti-0.5lti≥0；(5)xti+0.5lti≤L；(6)yti-0.5wti≥0；(7)yti+0.5wti≤W；(8)max{lti,wti}/min{lti,wti}≤d；(9)minLti≤lti≤maxLti；(10)minWti≤wti≤maxWti；(11)xti,yti≥0；(12)Stij∈{0,1}；(13)Vtij∈{0,1}．(14)式(1)～(4)表示设施i和j之间互不干涉(嵌入)；式(5)～(8)确保设施不会放置在车间之外；式(9)是设施最大长宽纵横比约束，其中d=2；式(10)和(11)表示设施的长宽范围，其中min Lti，min Wti，max Lti，max Wti分别为设施i在阶段t的最小长宽和最大长宽；式(12)是设施中心坐标的范围约束；式(13)和(14)中二元变量Stij和Vtij保证式(1)～(4)在任何时候只有一个约束起作用．事实上，当Stij=Vtij=0时，设施i位于设施j的右边，式(1)起作用；当Stij=1且Vtij=0时，设施i位于设施j的左边，式(2)起作用；当Stij=0且Vtij=1时，设施i位于设施j的上面，式(3)起作用；当Stij=Vtij=1时，设施i位于设施j的下面，式(4)起作用．2 干涉性处理如何处理不等面积设施之间的干涉性约束一直是解决UA-DFLP的难点问题．本研究采用拟物方法[9-10]来解决设施之间干涉性约束问题，即首先将所有设施和矩形车间F的外部区域O(含边界)视为光滑的弹性实体，在设施布置过程中，一旦设施之间或设施与车间外部区域之间发生干涉，它们之间必定存在挤压弹性势能，该设施必将在挤压弹性力作用下不断运动，最终有望达到不干涉的目的，具体步骤如下．在车间布置系统中引入挤压弹性势能U(X)=∑t=1T∑i=1N-1∑j=i+1Nutij+∑t=1T∑i=1NutiO=∑t=1T∑i=1N-1∑j=i+1Nkmtij2+∑t=1T∑i=1NkmtiO2,   式中：utij为在t阶段设施i与j之间的挤压弹性势能；utiO为在t阶段车间外部区域O与设施i之间的挤压弹性势能；k为弹性系数，本研究设置k=1；mtij为在t阶段设施i与j之间的干涉量，计算方法参见文献[11]；mtiO为在t阶段设施i与车间外部区域O之间的干涉量，计算方法参见文献[12]．为了使挤压弹性势能U(X)近似为0，将运用梯度法(GM)[9]求解优化问题min F(X)=E(X)+U(X)．梯度法又称最速下降法，在梯度法中，若已得到迭代点X1且X1不是极小点，则从X1出发，沿着目标函数F(X1)的负梯度方向找到新的迭代点X2=X1-▽(F(X1))h，其中：h为迭代步长；-▽(F(X1))为迭代方向，即函数F(X1)在X1处的负梯度方向，负梯度方向是目标函数值下降得最快的方向．对于步长h，本研究采用自适应策略，即在当前构型X1附近，若新产生的构型X2的目标函数值F(X2)高于构型X1的目标函数值F(X1)，则缩小步长h，令h=0.8h，否则继续按照步长h迭代．在GM算法中越接近极小点目标函数值下降得越慢，为了缩短迭代时间，提出一个“提前结束”策略：若目标函数值F(X2)非常接近目标函数值F(X1)，即当|F(X2)-F(X1)|1×10-3时，则认为迭代可能到达了局部极小点，此时若目标函数值F(X2)δ(δ的值根据具体动态设施布局问题多次优化的实验结果而定)，则继续执行梯度法，直至可行解出现或步长h小于给定的值(1×10-4)，否则提前结束梯度法的搜索．将融合自适应步长和提前结束策略的梯度法称为自适应梯度算法．自适应梯度算法的执行过程实际上是设施的运动过程，运动结果为所有设施之间及设施与车间外部区域之间不发生干涉或整个布局系统处于各方向作用力平衡的状态．若后者情形出现，则执行基于邻域构型的启发式布局更新策略，再从新构型出发继续执行自适应梯度算法，交替执行布局更新策略和梯度法直至得到合法构型为止．3 禁忌搜索算法本研究对禁忌搜索算法进行改进，结合基于邻域构型的启发式布局更新策略和自适应梯度算法，提出求解UA-DFLP的启发式禁忌搜索(HTS)算法．3.1　初始构型的产生禁忌搜索算法是一种基于局部搜索的全局寻优方法，通常从一个合法的好的初始解出发，以利于提高算法搜索性能，但即使从所有设施均相互干涉(重叠)的初始解出发，在交替执行基于局部搜索的自适应梯度算法和基于邻域构型的启发式布局更新策略后，仍可使设施之间的干涉性在弹性力作用下逐渐消除并搜索到优良的合法构型，因此本研究为随机产生初始解(构型)．3.2　邻域构型的产生在启发式禁忌搜索算法中，每次迭代都须更新构型，根据车间布局的特点提出一种基于邻域构型的启发式布局更新策略．定义1 动作．从当前构型X的任意一个布局阶段t中挑选最大单位物料搬运费用的设施，并将该设施的中心放置在该阶段车间内任意一个空位点(即空白区域内的点)，放置方向可以采用水平放置或垂直放置，称这一挑选设施并将其放置在空位点位置的一系列操作为做一个动作．布局阶段t中设施i的单位物料搬运费用为Mti=∑j=1Nctijftij(|xti-xtj|+|yti-ytj|)/∑j=1Nftij．定义2 邻域构型与邻域．对于当前构型X=(x11，y11，d11，x12，y12，d12，…，x1i，y1i，d1i，…，x1N，y1N，d1N；x21，y21，d21，x22，y22，d22，…，x2l，y2l，d2l，…，x2N，y2N，d2N；…；xt1，yt1，dt1，xt2，yt2，dt2，…，xtj，ytj，dtj，…，xtN，ytN，dtN ；…；xT1，yT1，dT1，xT2，yT2，dT2，…，xTk，yTk，dTk，…，xTN，yTN，dTN)，X′=(x11，y11，d11，x12，y12，d12，⋯，x′1i，y′1i，d′1i，⋯，x1N，y1N，d1N；x21，y21，d21，x22，y22，d22，⋯，x′2l，y′2l，d′2l，⋯，x2N，y2N，d2N；⋯；xt1，yt1，dt1，xt2，yt2，dt2，⋯，x′tj，y′tj，d′tj，⋯，xtN，ytN，dtN；⋯；xT1，yT1，dT1，xT2，yT2，dT2，⋯，x′Tk，y′Tk，d′Tk，⋯，xTN，yTN，dTN)为X的一个邻域构型．假设第一阶段、第t阶段和第T阶段中最大单位物料搬运费用的设施分别为i，l，j，k，(x′1i，y′1i，d′1i)，(x′2l，y′2l，d′2l)，(x′tj，y′tj，d′tj)，(x′Tk，y′Tk，d′Tk)分别为第一阶段中设施i、第二阶段中设施l、第t阶段中设施j及第T阶段中设施k在做一个动作后设施的中心坐标与方向．X的所有邻域构型的集合称为X的邻域．定义3 候选邻域．对于当前构型X=(x11，y11，d11，x12，y12，d12，…，x1i，y1i，d1i，…，x1N，y1N，d1N；x21，y21，d21，x22，y22，d22，…，x2l，y2l，d2l，…，x2N，y2N，d2N；…；xt1，yt1，dt1，xt2，yt2，dt2，…，xtj，ytj，dtj，…，xtN，ytN，dtN ；…；xT1，yT1，dT1，xT2，yT2，dT2，…，xTk，yTk，dTk，…，xTN，yTN，dTN)，将第1，2，…，t，…，T个阶段中设施i，l，…，j，…，k随机执行10个动作后得到的10个邻域构型为X的候选邻域．X的候选邻域是X的邻域的子集．基于邻域构型的启发式布局更新策略如下：首先在每个阶段布局区域内随机产生100个空位点，对每个阶段均采用空位点策略做100个动作；然后计算每个阶段挑选的设施在各空位点的单位物料搬运费，确定具有最小费用的空位点并将该设施正式放置于此，从而得到一个新的构型．对于当前构型X，假设第1，2，…，t，…，T阶段执行动作的设施分别为i，l，…，j，…，k，将这些设施序列记为R，为避免算法重复搜索，将R设置为禁忌对象．3.3　改进的禁忌搜索算法禁忌搜索算法的基本过程为：首先产生一个初始解(即为当前解)；然后在当前解的邻域中确定若干候选邻域解，若最佳候选邻域解的目标值优于当前最优解的目标值，则采用藐视准则(即当禁忌表中某一禁忌对象对应的最佳候选邻域解比当前最优解目标值更优时，忽略其禁忌属性的准则)，忽视其禁忌属性，用其替换当前解和当前最优解，并将相应的对象加入禁忌表，同时修改禁忌表中各对象的任期，释放任期为0的对象，若不存在上述候选邻域解，则在候选邻域中选择非禁忌的最佳解为新的当前解(忽视其与当前解的优劣比较)，同时将相应的对象加入禁忌表，并修改禁忌表中各对象的任期，释放任期为0的对象．如此重复上述迭代搜索过程，直至满足停止条件．改进的禁忌搜索算法基本过程为：首先产生一个初始解(即为当前解)；然后在当前解的邻域中确定若干候选邻域解，若最佳候选邻域解对应的对象是禁忌对象且其目标值优于当前最优解的目标值，则采用藐视准则，忽视其禁忌属性，将该禁忌对象从禁忌表中删除，用其替换当前解和当前最优解，同时修改禁忌表中各对象的任期，释放任期为0的对象，否则当前解和当前最优解保持不变；若最佳候选邻域解对应的对象不是禁忌对象且其目标值优于前一解的目标值，则仍将其替换当前解，同时将禁忌表中各对象的任期减1，释放任期为0的对象，否则不接受该最佳候选邻域解；若多次产生候选邻域解均不被接受，则恢复前一解为当前解，同时将相应的对象加入禁忌表(任期不减)．如此重复上述迭代搜索过程，直至满足停止条件．改进的禁忌搜索算法与传统的禁忌搜索算法有以下区别．a. 禁忌对象不同．前者是若某对象在若干次更新布局后得到的候选邻域解均不被接受，则将之设置为禁忌对象；后者是一旦某一对象对应的候选邻域解被接受，就将之设置为禁忌对象．这种改进策略能够给当前解多次寻优的机会，即从其出发经过基于邻域构型的启发式布局更新策略更新布局后，再进一步通过局部搜索策略(自适应梯度算法)对全局最优解进行有效搜索．b. 接受准则不同．对于非禁忌对象，前者是只有当其对应候选邻域解的目标值比前一解好时才被接受，否则就不被接受；后者是对于非禁忌对象，无论其相应候选邻域解目标值是否比前一解目标值好，均挑选最好的一个候选邻域解接受．这种改进策略可以给当前有前景(即目标值较好)的构型足够的机会，从其出发进一步搜索有望搜索到全局最优解，可以克服传统禁忌搜索算法从目标值较差的解出发而不利于全局最优解的搜索问题．3.4　启发式禁忌搜索算法本研究通过将改进的禁忌搜索算法与基于邻域构型的启发式布局更新策略及自适应梯度算法相结合，提出了一种求解UA-DFLP的HTS算法．在HTS算法中，首先随机产生一个初始构型初始化禁忌表．然后判断算法终止条件是否满足，若满足则结束算法，输出最优构型及其费用，否则继续执行算法．随后从当前构型出发，基于邻域构型的启发式布局更新策略得到最佳候选邻域构型，对该最佳候选邻域构型通过执行自适应梯度算法进行合法化操作，得到一个新的构型．判断新构型是否满足藐视准则，若满足则将该构型替换当前构型和当前最优构型，同时修改禁忌表中各对象的任期，释放任期为0的对象，并将新构型对应的禁忌对象从禁忌表中删除；若不满足则不接受新构型且当前构型保持不变；若新构型相应的对象不是禁忌对象，且其目标值优于前一构型的目标值，则仍接受新构型，同时修改禁忌表中各对象的任期，释放任期为0的对象，否则不接受新构型且当前构型保持不变；若多次产生的新构型均不被接受，则恢复前一构型为当前构型，并将相应对象加入禁忌表．重复上述步骤继续迭代，直至满足结束条件．HTS算法的具体步骤如下．步骤1 随机产生初始构型X1，计算F(X1)=E(X1)+U(X1)．令l=1，M=1×105，禁忌表H={}，最优构型Xopt=X1，最优目标值Fopt=F(X1)，Fmin=6 000(根据不同算例而定)．步骤2 判断算法是否满足终止准则(lM)，若满足则输出X1；否则令l=l+1，转步骤3．步骤3 找出各个布局阶段中最大单位物料搬运费用的设施，将这些设施序列记为R．备份X1′=X1，令i=1．步骤4 基于X1，执行基于邻域构型的启发式布局更新策略，得到新构型X2'．步骤5 从X2'出发，执行自适应梯度算法进行干涉性处理，得到合法构型X2．步骤6 若F(X2)Fmin，则输出X2，迭代结束；否则分两种情况考虑．a. 当R是禁忌对象时，若F(X2)Fopt，则将R从禁忌表H中删除，更新Fopt=F(X2)，转步骤7；否则不接受X2，恢复X1为当前构型，令i=i+1，转步骤8．b. 当R不是禁忌对象时，若F(X2)Fopt，则更新Xopt=X2，X1=X2，Fopt=F(X2)．若F(X2)F(X1)，则转步骤7；否则不接受X2，恢复X1为当前构型，令i=i+1，转步骤8．步骤7 接受X2为当前构型，令X1=X2，F(X1)=F(X2)，将H中各对象任期减1，释放任期为0的对象，转步骤2．步骤8 若i6，则将R中各设施单位物料搬运费设置为0，并将R加入H中，转步骤3；否则转步骤4．4 实验结果和分析为了检验启发式禁忌搜索算法的效率，引用文献[5]中的3组算例对HTS算法的性能进行测试：第一组包含两个经典动态设施布局测试算例P6-6和P12-4；第二组包含一个实际生产应用中的测试算例；第三组包含两个已知最优布局的测试算例．为了对比实验结果，所有算法均采用Java语言编写，并在配置为Intel Core 2 Duo，2.94 GHz CPU，2.0 GiB内存的电脑上运行．4.1　经典动态设施布局算例的测试对于UA-DFLP，由于文献中带设施长宽纵横比约束的算例非常少，因此为了测试本研究模型和算法性能，将经典动态设施布局问题测试算例中的约束条件进行修改：将设施长宽的约束条件由原来的固定长宽修改为长宽可在一定范围内变化，并且满足最大长宽比在2以内，每个设施添加最小面积约束．有关约束修改如表1和表2所示．算例P6-6和P12-4的设施重置费用分别设置为19元和50元，分别运行启发式禁忌搜索算法、传统的禁忌搜索算法(TS)及文献[5]中启发式Wang-Landau抽样算法(HWL)各10次，将获得的最优结果列于表3中．10.13245/j.hust.210206.T001表1算例P6-6长宽约束值设施编号原长原宽长范围宽范围最小面积1542～81～7162986～125～11573653～92～8244643～91～7195441～71～7126532～81～61210.13245/j.hust.210206.T002表2算例P12-4长宽约束值设施编号原长原宽长范围宽范围最小面积1542～81～7162754～102～8283653～92～8244441～71～7125663～93～9286542～81～71671077～134～10568754～102～8289653～92～82410552～82～82011552～82～82012643～91～71910.13245/j.hust.210206.T003表3计算结果比较算例HWLTSHTSE(X)/元s/minE(X)/元s/minE(X)/元s/minP6-66 857.258.617 086.2412.426 787.3310.23P12-428 751.5362.8330 542.4872.4127 614.6979.78从表3可以看出：对于这两个算例，HTS算法得到的目标值E(X)(物料搬运费和设施重置费之和)均优于HWL算法和TS算法得到的目标值．对于算例P6-6，HTS算法所得目标值比HWL算法所得目标值减少了1.0%，比TS算法所得目标值减少了4.21%．对于算例P12-4，HTS算法所得目标值比HWL算法所得目标值减少了3.95%，比TS算法所得目标值减少了9.58%．算例P6-6和P12-4的计算结果验证了HTS算法的改进策略比TS算法有效．另外，由于HTS算法采用自适应梯度算法进行干涉性处理，HWL算法采用外推内压法处理设施的干涉性，因此HTS算法得到的布局结构比HWL算法更加紧凑，其目标值更优，但梯度法比较耗时，HTS算法的运行时间s比HWL算法略有增加．图1和图2分别为HTS算法对于算例P6-6和P12-4在10次运行中找到的最优构型布局图．10.13245/j.hust.210206.F001图1算例P6-6最优构型布局图10.13245/j.hust.210206.F002图2算例P12-4最优构型布局图4.2　实际应用算例的测试对文献[5]中的实际应用算例(这里设施长宽固定)运行HTS算法10次，得到的最优目标费用为2.581 390 1×105元．虽然HTS算法所得的构型目标费用比HWL算法所得的构型目标费用2.570 855 0×105元略有增加，但仅增加了1 053.51元，对于数量级为十万的目标值而言可谓是无差别的．表4列出了HTS算法所得的构型中各个布局阶段的设施坐标(x，y)与方向d．10.13245/j.hust.210206.T004表4HTS算法所得的构型布局的设施坐标与方向设施编号第1阶段第2阶段第3阶段xydxydxyd111.007.39011.1910.00110.3013.590213.8010.8917.1813.40110.4910.591310.5010.6416.866.11014.057.431411.0114.3007.1810.00014.3013.590514.2513.4903.9510.00014.159.94167.1014.3013.406.50010.657.79177.0010.75110.115.85114.304.20087.207.00013.916.00010.554.2014.3　已知最优解算例的测试为了测试HTS算法在车间空间紧张的情况下的布局效果，进一步对文献[5]中两个已知最优解算例P3-3和P4-2进行测试，独立运行HTS算法10次．对于两个已知最优解算例，HTS算法与HWL算法的测试结果如表5所示．可以看出HTS算法得到了两个算例的已知最优解，而HWL算法只得到了算例P3-3的最优解，对于算例P4-2只得到了近似最优解．这主要是因为在HTS算法中，基于自适应步长的梯度算法对物料搬运费用进行最速下降搜索，在保证设施不干涉的同时，能尽可能降低物料搬运费用．在每一轮迭代直至到达局部极小值之后，算法采用启发式构型更新策略重新构造布局，有利于HTS算法跳出局部最优，从而找到最优布局．10.13245/j.hust.210206.T005表5已知最优解算例的测试结果算例HWLHTSE(X)/元s/minE(X)/元s/minP3-31 510.52.41 510.55.6P4-2848.53.5827.55.95 结语本研究针对在复杂市场环境下不等面积动态设施布局问题，首先基于混合整数规划建立了多阶段动态设施布置的连续性模型，然后基于拟物策略提出自适应梯度算法解决设施之间的干涉性约束问题，最后通过对传统禁忌搜索算法中禁忌对象和解的接受准则进行改进，并结合基于邻域构型的启发式布局更新策略，提出了一种基于改进的禁忌搜索策略的启发式算法．该算法通过禁忌搜索并结合启发式布局更新策略对解空间进行全局搜索，而自适应梯度算法则进行局部搜索．通过三组算例对算法性能进行测试，实验结果表明：将全局搜索与局部搜索相结合的搜索策略能有效地提高算法对全局最优解的搜索性能．