零部件循环取货越库物流系统仿真优化 施 文，刘志学，杨 威 （华中科技大学 管理学院，湖北 武汉 430084） 摘要：在众多汽车制造商产能扩大的背景下，以零部件平均流程时间为绩效指标，深入研究了基于第三方物流零部件循环取货越库配送物流模式的优化设计问题，建立符合实际物流运作的仿真模型。由于存在较多的仿真输入(因子)，提出首先应用序贯分支筛选关键因子，再针对关键因子建立响应面与Kriging元模型优化的研究步骤，解决了由于“高维因子”而导致的较难优化的问题。运用调研数据进行仿真实验计算，结果表明关键因子的最优配置可使物流绩效达到较为理想的水平；Kriging元模型的预测结果略优于响应面的最佳绩效值，该模型与方法为产能扩大后TPL物流改进方案的制定提供了可靠的解决思路，具有良好的应用价值。 关键词：零部件物流；第三方物流；越库；仿真；序贯分支；响应面元模型；Kriging元模型 中图分类号：F274 文献标志码：A Screening and optimizing the key factors of milk-run and cross-docking parts logistics system: Sequential bifurcation, response surface and kriging methods SHI Wen，LIU Zhi-xue，YANG Wei （School of Management, Huazhong University of Science and Technology, Wuhan 430074, China） Abstract: Automobile and logistics are the two critical restructuring and invigorating industries in China. In the context of assembly plant capacity expansion, with the purpose of improving system performance, this paper elaborates on the optimal design of key operations in the logistics mode of auto parts milk run and cross docking dealt with by a TPL company. We build a fresh simulation model which integrated the operations of milk-run and cross-docking based on practical research data, and propose a hybrid methods that integrate Sequential Bifurcation, Response Surface Methodology (RSM) and Kriging. The results of specified simulation design indicate that the logistics performances are greatly improved when the TPL simply optimizes the key operations, resulting in a great reduction of TPL’s workload. In addition, the optimal prediction of Kriging method shows a better result than RSM. This study makes the contribution on the decision support automobile parts TPL Company when automobile assembly plant capacity expands, and guides to solve the optimization issue of complex simulated system. Key words：auto parts logistics; third-party logistics; cross-docking; simulation; sequential bifurcation; response surface meta-model; kriging meta-model 1 引言 国内汽车制造企业零部件入厂物流(零部件由供应商至主机厂)大多采用基于第三方物流（Third-Party Logistics, TPL）的循环取货(Milk-run)和越库配送(Cross-docking) 模式(以下简称TPL-MRCD)，其中TPL配合主机厂JIT生产实施JIT配送，对上游零部件供应商循环取货后，经越库配送将零部件送抵主机厂生产线线边仓库。在TPL-MRCD系统中，存在众多影响JIT成功实施的因子(仿真实验中的术语)，如循环取货车辆数、越库配送车辆数、越库卸载区/转载区库台数、越库 收稿日期：2011-11-07；修订日期：2011-12-28。Received 07 Nov. 2011;accepted 28 . 基金项目：国家自然科学基金资助项目(71072034)。Foundation item: Project supported by the National Natural Science Foundation, China(No. 71072034). 1
内装卸载设备数等。因子对系统物流绩效的影响在主机厂产能扩大的背景下尤为凸显。据《中国汽车行业发展报告2010》预测，未来10年国内汽车需求增长率将保持在10%左右，包括武汉东风乘用车﹑神龙汽车﹑东风本田在内的众多整车厂纷纷设立新厂。产能的扩大势必带来的零部件物流量的增加，这对TPL的JIT实施提出更高的要求，如何优化配置系统影响因子是TPL的重要任务之一。 由于所建模型涵盖了零部件Milk-run和越库配送两部分物流运作，故首先总结和分析了与此对应的两部分文献。此外，由于所建模型为一类随机仿真模型（同一仿真因子组合下系统重复运行所得仿真输出不同），笔者综述了现有文献针对该类模型的优化方法。 循环取货与越库配送研究 有关汽车零部件物流的研究多数集中于上游零部件Milk-run问题，文献[1-2]对日本汽车制造厂的零部件Milk-run展开研究，他们认为基于TPL的越库物流模式在小批量、高频率的零部件配送模式中极为重要，但未提出任何数学模型。文献[3]探讨了国内基于第三方物流的汽车零部件物流Milk-run模式的实施情况。有关越库的技术性研究始于20世纪90年代，文献[4]研究了越库中心(Cross-Docking Center, CDC)作为供应链某一层级的网络选址问题。文献[5-6]探讨了货物直接运输或经CDC转运的网络设计问题。文献[7]将CDC与车辆路径问题联合考虑，采用禁忌搜索算法找出最优的车辆调度路径。文献[8-10]研究了越库卸载区与装载区各一个库台(1对1)的单库台车辆调度模型，在此基础上，文献[11]拓展到1对2或2对1的多库台模型。 综上可知，现有文献多以解析法建立Milk-run或越库内部运作的模型，较少将两者结合起来分析整体运作对系统绩效的影响，尤其缺乏针对汽车行业零部件Milk-run取货和CDC配送运作的集成模型。这主要因为两方面单独的模型已较为复杂，多数为NP难或NP完全。然而许多模型的假设与实际情形差距较远，如文献[8-11]假设较少的装卸载库台，实际中库台数一般为实际中库台数一般为[12]6～200个，笔者调研的东本储运、长安民生不低于10个。仿真是解决复杂问题强有力的工具，本文基于零部件入厂物流构建了TPL-MRCD离散事件(随机)仿真模型，该模型更符合零部件物流的实际运作流程。 (随机)仿真优化研究 在供应链(随机)仿真系统中，多数学者基于仿真实验的输入输出数据直接建立响应面(RSM)元模型并通过非线性规划找到了输入因子的最佳组合，如文献[13-17]。元模型是仿真模型输入(因子)输出[18](系统绩效或响应)关系的数学接近式，最终目的用于验证仿真模型，进行灵敏度分析和优化。上述文献中的因子数量较少，一般为5个，其所建元模型需要的仿真次数小于100。当仿真系统存在[19]大量因子时，直接建立元模型所需的仿真实验次数将成指数倍增加，这将大大提高计算机与仿真分析者的工作量。本文提出采用序贯分支(SB)与元模型相结合的方法筛选并优化模型中的关键因子，即先筛选出众多因子中的关键因子，再基于关键因子建立元模型。相比其他筛选技术，在符合假设的前提[18]下SB法更为高效。在现有供应链仿真领域的文献中，未发现因子筛选与元模型优化相结合的仿真案例。另一方面，RSM也称低阶多项式回归元模型，它适用于拟合局部输入输出数据。相比RSM，2
[20-21][22]Kriging从更大的实验区域内拟合数据，具有较高的柔性。Kriging源于地质统计学，文献首次将Kriging引入到随机仿真领域，目前有关Kriging的优化研究仍以小型随机仿真模型为主，如(s, S)库存模[23-24][20]型和M/M/1排队模型，缺乏基于实际背景的复杂仿真案例的优化研究，并缺乏两类元模型间的比较。 鉴于以上模型与方法论两方面内容，本文首先建立TPL-MRCD模型，通过灵敏度分析探讨产能扩大后物流绩效的变化规律；将该模型与序贯分支筛选设计及两类元模型相结合，深入分析TPL-MRCD物流系统的最佳改进方案，为TPL运作决策的调整提供可靠的解决方案，为两类元模型的比较提供一定的理论参考。 2 问题描述 国内绝大多数汽车零部件供应商分布于珠三角、长三角、环渤海湾三大汽车产业带，而某些主机厂(Assemble plant, AP)位于距离较远的中西部地区，如东风本田、东风神龙、长安福特等，并且通过合资的方式成立了为其服务的第三方物流企业(TPL)，如：东本储运、捷富凯—鸿泰、长安民生等。从规模经济角度考虑，TPL在各汽车产业集聚区设立CDC。汽车制造商、零部供应商和TPL事先签订合作协议，委托TPL实施零部件物流Milk-run，即根据零部件供应商的分布、零部件的补货频率、补货数量、需求特征和对生产的重要性等因素，划分若干Milk-run区域。图1展示了整个TPL-MRCD物流流程，通常情况下，从某个Milk-run取货的零部件数量不够整车(Less-than-truckload, LTL)，因此该车将LTL零部件先运输至区域CDC(支线运输)，经过一定时间的集并，再整车(Truckload，TL)运输至厂边仓库最终到达主机厂(干线运输)。 图1 TPL-MRCD零部件物流 图1中间部分展示了CDC内部的基本结构。库台(Doors)是卸载区(Receiving dock)和装载区(Shipping dock)最主要的组成部分，一般地，单个库台仅容纳一辆货车。一旦入库车辆(Inbound trucks, IT, 图1中红色卡车)抵达CDC，调度室视库台的繁忙程度决定该车辆下一位置，若卸载区有空闲库台，入库车辆则随机分配至其中之一；相反，若卸载区库台已全部占用，入库车辆则先进入车辆等候区等待，直到某一库台空闲。入库车辆停靠在库台时，操作工人分派至此卸货扫描。扫描完成后，3
零部件根据装载区出库车辆(Outbound trucks, OT)的到达情况决定下一步的路线：路线1，若装载区库台停靠有OT，则零部件通过传送带直接运至装载区，其中，卸载区库台与传送带间及传送带与装载区间的两段运输由叉车完成；路线2，若装载区库台未停有OT，零部件则由叉车先运至暂存区暂存，待OT到达后，再由叉车运至该OT停靠的库台。由于路线2的零部件有一定存储时间，故装载区优先选择路线2的零部件。零部件在装载区装至满载后再干线运输至厂边仓库。此外，零部件不允许在CDC暂存区存储太长的时间，一般不超过24小时，否则由额外的货车将其运至厂边仓库。 远距离的零部件供应商需要配合下游主机厂JIT生产，实施JIT配送。产量的扩大对零部件的流程时间提出了更高的要求。因此，本文以零部件平均流程时间(Cycle time, CT)作为绩效指标，CT是[8-11]衡量系统作业效率高低的标准，部分越库文献也将此作为研究目标。由于集成了CDC的外部运作，本文中CT包括零部件的Milk-run时间、运输时间和CDC内的运作时间三部分。 2 仿真模型 逻辑模型 [25]本文利用Rockwell公司开发的面向对象的仿真平台ARENA 作为仿真工具。为了建模和测试的顺利进行，本文将TPL-MRCD模型分割为八个相互联系的子模型(Submodel)，图2不同颜色的模块标示了各个子模型(逻辑控制模型未在图中显示)，子模型间以粗线条连接。鉴于篇幅，本文介绍循环取货及支线运输、卸载区和装载区三个子模型： (1)子模型1 循环取货及支线运输模型。 支线运输车辆或称入库车辆(IT)为模型中支线运输资源，在主机厂一个生产周期内数量为TQ，下标M和C分别表示Milk-run和CDC。IT有忙和闲两种状态，忙态从IT接到调度指令开M,C始，到CDC卸载区零部件卸载完毕结束；闲态的IT载着空容器自动返回至上次服务的Milk-run区。CDC对第i个Milk-run发出循环取货指令(主机厂的订货指令)的平均间隔时间BT服从均值为M(i)μ的指数分布；若此时支线运输有空闲的IT，则立即将距该Milk-run区最近的IT调度至此；否BTM(i)则等待空闲的IT，其中，第和个Milk-run间的平均距离为ID，，为TPL划M(i),M(j)分的Milk-run数。第个Milk-run供应商数量为，代表供应商总数。第个Milk-run单次循环取货的供应商所占比例为，即每次循环取货的供应商数为。单个供应商零部件取货时间服从参数为的均匀分布，供应商单位零部件装载时间服从参数的三角分布。取完货后IT零部件装载量=, ，其中为IT的容量，为单次循环取货应达到的最小装载量。因此，第个Milk-run的取货时间。第个Milk-run与CDC的距离为，支线运输车辆行驶的平均速度为固定值。 （2）子模型2 卸载区模型。 IT到达CDC触发卸载区事件，该事件包括库台分配、零部件卸货和扫描三个子事件。卸载区4
的库台数为(下标表示CDC)，若IT到达时处于忙态的库台超过，则该IT先到车辆等候区排队等待直到空闲库台出现，默认情况下等待区服务规则为先到先服务()。每一库台配备足够数量的作业工人，单位数量零部件卸货时间服从参数为的三角分布。因此，某辆IT的卸货时间为。扫描与卸货在卸载区库台同步进行，包装的朝向错误或标签损坏时有发生，需要工人二次扫描或手动输入条码，单位数量零部件扫描时间服从参数为的均匀分布，一次扫描通过率为，第二次操作时间与第一次相同。零部件下一步路线取决于卸载区是否停靠有OT，若有，则零部件直接进入传送区由传送带运至该OT，否则进入暂存区存储。 （3）子模型3 装载区模型。 装载区库台数等于，若OT到达时，装载区无空闲库台，则车辆到卸货区车辆等待区排队等待，默认情况下等候区排队规则为。库台配备相应数量的操作工人，整车OT的零部件装载量为，单位数量零部件装载时间服从参数为的三角分布，因此，每辆OT的装载时间为LT*Q。 CTL 图2 TPL-MRCD物流逻辑流程 算例分析 本文模型为终态(Terminating)仿真模型，数据来源为笔者所调研的武汉某家汽车TPL企业。它5
所服务的汽车制造企业以30 d为一个计划生产周期，因此，共有30个仿真日(每日工作12 h，共360 h)，另有5个仿真日的系统Warm-up时间。由于是随机性仿真，为了提高仿真输出的精确度和实验设计的效率，每个仿真场景的重复次数为5。 本模型的基本假设包括：①每日工作12 h，8:00-20:00；②运输车辆IT和OT仅有忙和闲两种状态；③支线运输与干线运输交通状况良好；④CDC暂存区容量无限；⑤任一卸载区库台至装载区库台距离为一固定值；⑥CDC内部操作工人数足以满足库内零部件操作需求；⑦卸载区(装载区)库台同质，即卸载区(装载区)库台的卸载时间(装载时间)服从相同的参数分布；⑧汽车零部件在供应商处已完成包装，CDC内无包装流程；⑨需求在该计划周期内处于恒定状态；⑩初始时刻IT与OT停靠在CDC；⑾整个供应链完全信息共享，因此主机厂发送采购指令后，TPL立即发车循环取货。 3 产能扩大对系统绩效的影响 图3a～图3d展示了产能扩大后TPL企业部分物流设备的使用情况(横坐标为仿真时间，纵坐标为部分物流设备的占用数，如运输车辆、库台数、叉车、拖车等)。可以看出，产能越小，企业设备利用率越低，设备状态的变动幅度越大；随着产能的增大，企业资源利用率也随之提高，资源越来越接近于瓶颈状态。 表1列出了不同组合成比例取值时重复进行10次仿真后得到的系统绩效平均值，结合表1中的部分与图3a～图3b可看出，随着的增加，CT的期望与标准差有一定的下降，下降原因是在其他条件不变的情况下，零部件在CDC内能更短快地形成，即减少了零部件在CDC中的等待时间；此时和的期望成比例增加，总体而言系统绩效未向不利方向转变。由上述分析可知，在短时间内，产能的增加对TPL-MRCD系统绩效无不利影响；相反，系统绩效可获得一定的改善。从长期来看，由于受到资源瓶颈的约束，产能的持续增加对系统绩效产生不利影响，图3c～图3d中部分资源长期处于瓶颈状态(表现为一条直线)，与之对应的表1中零部件平均物流时间急剧增加()，供应商发货量(19 )虽有所增加，但收货量维持在一个较为固定的数值(18 500)，说明TPL运能已达到极限，系统内积压了较多的在制品库存，TPL需改进物流运作以满足新的物流需求。 表1 不同发货频率下的系统绩效 绩效指标 (供应商发货量)(厂边仓库收货量) 期望标准差期望 标准差期望 标准差 1 1 3 3 7 7 9 9 15 15 19 18 26 18 6
(c) 为{, , , , } (a) 为{, , , , } (d) 为{, ,1, , } (b) 为{, , 3, 5, } 红线代表，蓝线代表，黑线代表，粉红线代表，绿线代表 图3 不同发货频率下系统部分资源利用率 由以上分析可知，随着产能的持续增大，企业需要采取有效的改进措施应对这一变化，使物流绩效达到较为理想的水平；否则，系统将长期处于超负荷状态，这将严重影响主机厂的正常生产。 4 筛选设计与元模型优化 SB筛选重要因子 本模型中影响系统绩效因子共有58个(见附件)，其中可变因子21个(不可变因子包括：、、、、等供应链结构既定值)，具有“高维性”。21个可控因子高低水平见表2所示，单个因子从低水平变为高水平不会减少仿真输出的期望，如第9号因子(收货区库台)个数由20(低水平)变为10(高水平)的过程中，CT(非严格)单调上升，这符合现实情形。高低水平的判断源于TPL企业和分析员的先验知识，对无法判定的因子需个别处理。 表2 TPL-MRCD系统影响因子 序号 参数 低水平 高水平 1 UNIF(10, 15) UNIF(15,20) 2 UNIF(1,2,3) UNIF(2,3,4) 3 TRIA(1,2,3) TRIA(3,4,5) 4 EXPO(20) EXPO(30) 5 TRIA(2,3,4) TRIA(4,5,6) 6 TRIA(,,)TRIA(,,) 7 40% 60% 8 1%2%9 20 10 10 2010 11 20 10 12 402013 70 50 14 243615 IT排队规则SPT FIFO 16 30207
17 100 75 18 10075 19 5 10 20 2010 21 24 12 步骤1 将21个因子聚合为一个因子组(Group)，检测该因子组是否有重要效应。如果有，转入第2步；否则，表明21个因子对系统绩效无显著影响，停止筛选。第2步将21个因子分为两个子组(Subgroup)，即分支，检测两个子组的重要性。按照类似的步骤，将重要的因子组分为更小的子组，放弃不重要的子组，最终筛选出重要因子。 由于模型因子间有交互作用，如对CT的影响大小与其他因子取值有关，故一个有效元[26]模型是具有二阶交互作用和噪声的一阶多项式： (1) 。 。 (2) 其中，表示元模型的响应，21为因子的个数，为第个因子的主效应()，为因子与的交互作用()，噪声是随机仿真产生的近似误差，通过(2)式标准化，即[－1,1]；为的原始值，和为的上下限，如第9号因子标准化后。标准化意味着因子可以按主效应大小排序，即最重要的因子主效应值最大，最不重要的因子效应值接近于0。 为了方便起见，采用以下符号： ，第次仿真重复中第1到号因子高水平，第到21号因子低水平所获得的仿真输出值； ，第次仿真重复中第1到号因子低水平，到21号因子高水平所获得的仿真输出值； ，第到号的主效应之和，即：=； 根据式(1)，基于第次重复的无偏群组效应估计量和单个无偏因子效应估计量分别为式(3)和式(4)，进行m次重复仿真，从而得到每个估计效应的均值和标准差，见式(5)。同理亦可得出群组效应的均值和标准差，最后通过t-test检验效应的显著性。 ， (3) ， (4) 。 (5) 由于TPL-MRCD是随机性仿真系统，故进行m=5次重复仿真实验。表3给出了前四个场景的重复输出及群组效应值，即：场景①，所有因子都置于低水平的极端场景(此时CT最小)；场景②，8
所有因子都处于高水平的极端场景(此时CT最大)；场景③，第一次分支后，1～16号因子处于低水平，17～21号因子置于高水平的场景；场景④，第一次分支后，1～16号因子处于高水平，17～21号因子置于低水平的场景。表中数据显示：场景①时，平均输出为= 1；场景②时，= 3。因此，1～21号因子群组效应估计量由式(3)获得，=(－ 1)/2= 1，该群组效应的标准差=，可见群组效应显著。 4步骤2 将当前群组的21个因子分成2个子群，第一个子群放入16(2)个因子，第2个子群放入剩余的5个因子。由表4可见，场景③时，= 4；场景④时，= 6，得到(1－16)和(17－21)两个子群效应估计值及标准差为= 6，=，= 5，= 0，说明(1－16)为较重要的子群，重要因子的筛选需进一步分支。 表3 TPL-MRCD仿真中前四个场景重复5次的仿真实验结果 重复 1 6 4 9 0 2 5 8 9 6 3 8 5 9 4 1 7 4 9 5 7 1 0 均值 1 4 5 1 标准差 1 0 8 图4展示了SB的后续步骤及其计算结果，当估计得出的子群主效应不显著为正时，停止进一步的分割。该模型SB过程停止在第13步，上限= 3。获得的因子清单中有4个因子，图4中黑色边框标识，重要因子占总因子数的19%。最重要的因子为第13号()因子，次重要因子包括第12号()、第9号()以及第18号()因子。 图4 TPL-MRCD模型SB步骤 9
由SB分析的结果可知，干线运输车辆是影响CT的关键因子，随着主机厂产能扩大，干线运输将出现运力不足的问题，这将导致零部件在CDC内的过度积压，延长零部件在暂存区内的等待时间，从而增加在制品库存，这不利于主机厂JIT生产的实施。这一结论也解释了第3节中产能扩大后实际收货量维持在一个较为固定的数值的现象(18500)。因此，TPL应首先增加干线运输的车辆数，确保干线运输的运能。此外，适当增加支线运输车辆数、CDC卸载区库台数和更换性能更好的车辆对系统绩效改进也有较大的帮助。 元模型优化重要因子 企业不可盲目地调整关键因子，本模型中大幅度增加，和将造成资源的浪费，其最优取值仍需由进一步确定。 由于重要因子中低水平值已为其边界值，故优化中仅分析，和三个因子。通过进一步的小规模仿真实验发现，参数组合{70,40,20}未使CT达到最优，表5中参数组合{80,70,30}和{100,100,50}的绩效指标值优于{70,40,20}的结果(= 6< 4)；而{80,70,30}，{100,100,50}仿真结果一致说明此时系统中重要因子已有空闲资源，该组合下的绩效值(= 6)为重要因子能使CT达到的最小值，因此重要因子最佳的资源配置取值界于[{70,40,20},{80,70,30}]之间。 表4 重要因子三组参数组合仿真数据 参数组合 仿真次数 1 7 2 2 5 8 3 4 3 4 3 0 5 3 6 平均值 4 6 标准差 1 1 响应面元模型及优化 对三因子[{70,40,20},{80,70,30}]建立RSM元模型，它是一种通过估计参数组合，使给定目标函数达到最小的启发式优化方法。与其他仿真优化方法相比，RSM需要相对较少的仿真重复和因子组合，并且与昂贵仿真(单次仿真运行需要大量时间)运行所需的时间相比，RSM启发式搜索所需的数[18]学和统计计算可以忽略。 响应(CT)和自变量(重要因子)之间的关系形式是未知的，RSM的第一步就是寻求响应和自变量集合{x,x,x}(代表)之间真实函数的一个合适的逼近式。当远离最优点时，123通常假定在的一个小区域范围内是真实模型曲面的合适近似，一般采用分辨率Ⅲ设计方法，若CT适合重要因子的线性函数建模，则近似函数是一阶元模型： 。 (6) 式中：，和是未知参数，通过仿真实验的数据由最小二乘法估计，为模型的随机误差项。为确定下一步因子的方向，RSM将式中拟合的梯度用到最速下降法之中到达最优点区域，若最10
优点位置方差分析显示系统有曲度出现，则改用更高阶的多项式，一般为中心复合设计(CCD)拟合的二阶元模型： 。 (7) 最后，对(7)式应用规范分析检查最优区域的外形并用数学规划求得极值，见式(8): 。 (8) 采用分辨率Ⅲ设计方法拟合三因子一阶响应面元模型，拟合结果见式，(9)。由于模型曲度显著，如图5所示，故将一阶模型(6)替换为二阶模型(7)。 26 25 24 23 A: TQc,f B: TQm,c 图5 二阶响应面元模型 对该区域采用3因子中心复合设计，包括组全析因设计，6组轴点设计及1组中心点设计，其中，非中心点仿真重复次数为5，中心点仿真重复次数为10，总的实验次数为8*5+6*5+10=80，最终拟合的二阶响应面方程为式。(10)。决定系数=和调整的= 5说明式。(10)能够较好地代表TPL-MRCD仿真系统，由式(8)和(10)求得最优的因子组合{80,62,23}，此时=。 ，(9) 。(10) Kriging元模型及优化 Kriging元模型一般形式为： 11 CT
·y(x)=μ+z(x) 。 (11) μ为仿真输出在整个实验区域的均值，单个μ并不代表响应面是平坦的，也可用回归模型替换μ； z(x)为仿真输出与均值的偏差构成的一个随机过程，是具有均值为0的协方差平稳过程。 一般地，Kriging使用如下线性预测量： n·y(x)λY(x)=λ'(x,X)Y(X) 。 (12) 0∑ii0i=1权重λ'(x,X)随着预测的输入x和旧的设计点X=(x,x,x)之间的距离的减小而增加，X决定00123了仿真输出向量Y(X)。Kriging选择最佳线性无偏估计量(Best linear unbiased predictor)作为权重λ'(x,X)的选择标准。 0[18, 24, 27]最佳权重可证明取值如下： -11−1'Γγ-1λ=Γ[γ+1] 。 (13) 0-11'Γ1Γ=Cov(y,y)是旧输出(已经仿真过的因子组合的输出)之间协方差的n×n维正半定对称矩阵，ii'i,i'=1,...,n，n为仿真过的因子组合数。 γ=Cov(y,y)是n个旧输出y与需要预测的输出y之间的协方差的n维向量。 i0i0最后，式(11)～式(13)表明： ·-1−1-1'-1-1−1-1y(x)=(1'Γ1)1'ΓY(X)+γΓ(Y(X)(1'Γ1)1'ΓY(X)1) 。 (14) 0这样，式(14)取决于协方差(相关系数)Γ与γ，式(11)表明协方差Γ与γ由输入组合间的距离决定，Kriging假设k维输入变量的相关函数为k个一维函数乘积： kkppjjρ=exp(−θ|x−x|)=exp(−θ|x−x|) 。 (15) ∑jj;ij;c∏jj;ij;cj=1j=1式中：|x−x|表示仿真输入j在因子组合i与c间的距离；θ为因子j的重要性，θ越大，j;ij;cjj因子j的影响越小；P表示相关函数的平滑度。 j[18]采用拉丁超立方体采样(Latin Hypercube Sampling, LHS)生成重要因子与CT仿真实验样本点。[24]采样个数为10×n，其中，n=3为因子组合数，共30个采样点，如图6所示。 30个仿真组合(x,x,x)重复3次的平均运行结果见附件，根据此值拟合的主要参数为：θ=( 1233, 0, 0)，μ=，其他参数结果见附件。利用式(14)预测该实验区域的CT值，得到1 ¶000个因子组合预测结果(见附件)，argminCT={80,63,23}，此时CT=。 12
图6 重要因子仿真实验点LHS采样 不同运作方案及元模型的比较 图7展示了产能扩大后因子在不同水平下进行25次独立仿真实验后得到的系统绩效指标值。“方案一”曲线表示主机厂产能增大后，TPL未做任何调整时零部件平均物流时间呈现出的变化趋势，由图中可看出，CT的期望值较大，标准差波动也较大，系统处于极不稳定的状态，说明此时零部件物流绩效难以达到汽车行业现行的JIT生产要求，不利于主机厂的正常生产，TPL物流运作亟待改善。“方案二”曲线表示仅调整重要因子时的绩效，图7中系统绩效明显改善，其期望值和标准差均大幅度减小，系统状态相对稳定。这说明在主机厂产能增大的过程中，找出影响系统绩效的关键物流运作是TPL工作的重点。“方案三”与“方案四”时，系统绩效进一步改进，CT期望值与标准差进一步缩小，此时为重要因子所能发挥的最大效用。“方案一”与“方案三、四”曲线间的差值即为重要因子的绩效改进贡献值，本模型中重要因子贡献值达到%。由此可见，TPL只优化少数的关键物流运作，便可使绩效恢复到较为合理的水平，这将大大减少TPL的工作量。“方案五”曲线表示改进所有物流运作得到的绩效，该组绩效值为五种组合下的最优值，但TPL所需改进的运作内容也最多。“方案三、四”与“方案五”曲线间的差值为非重要因子的绩效改进贡献值。本模型中，非重要因子贡献值仅为%。由上述分析可以看出，随着主机厂产能扩大，TPL企业开始无需调整所有运作，盲目地扩大所有资源未必有理想的效果，如本文中增加装载区库台数(第10号因子)对绩效改进贡献不大，而首先应重点改进关键物流运作，使其效用最大化；若关键运作最优水平达到的绩效仍无法满足汽车制造商的JIT生产需求，则再按因子重要性逐步改进其他的物流运作，使系统逐步达到最佳绩效，如本文中第2、17和1号因子，如图5所示。 图7显示，Kriging最优因子组合绩效值略优于RSM最优绩效值，这一结果源于Kriging初始样本点LHS的空间填充(Space-filling)特性，而RSM中的中心复合设计实验点则较为稀疏；此外，Kriging的空间关联性的考虑也使得预测更为准确。以上计算分析结果充分显示，提出的数学模型和方法是可靠且实用的，具有较高的操作性和可行性。13
方案一：TPL无调整任何运作(所有因子高水平)，CT=； 方案二：TPL仅调整重要运作(重要因子低水平，其余因子高水平)，CT=； 方案三：TPL重要运作调整至RSM最优(重要因子最优水平，其余因子高水平)，CT=；120 方案四：TPL重要运作调整至Kriging最优(重要因子最优水平，其余因子高水平)，CT=； 方案五：TPL调整所有运作(重要因子最优水平，其余因子低水平)，CT=仿真次数0510152025 图3 不同因子水平的仿真结果比较 5 结束语 为分析主机厂产能扩大对汽车零部件物流产生的影响以及TPL企业如何识别与优化影响其系统绩效的其关键因子以满足新的物流需求，构建了从供应商到越库再到主机厂场边仓库的零部件离散仿真模型。运用序贯分支筛选出影响系统绩效的关键因子，并通过建立关键因子的响应面与Kriging元模型，得到系统绩效最优时的关键因子配置。研究结果表明：TPL集中调整系统的关键因子可极大改进系统绩效；Kriging的预测值略优于响应面的最佳预测结果。本文建立的仿真模型及提出的筛选与优化方法能辅助汽车行业TPL企业管理者设置运作参数，为主机厂产能扩大后企业决策的制定提供一定的理论基础和决策依据。本研究也尚存不足，拟从以下诸方面进一步改进：进一步研究仿真系统的参数可靠性和系统稳健性；考虑更多的随机因素，如操作工人的调度、轮班制度和资源故障等；综合考虑多种绩效评价指标，如物流量、成本等；建立加入供给、需求等不确定性因素的进行系统的稳健设计；深入挖掘一下日本、美国汽车物流技术发展动态，并与国内的技术研究比较，以利于提高我国的汽车物流技术的应用水平。 参考文献： [1] NEMOTO T, HAYASHI K, HASHIMOTO M. Milk-run logistics by Japanese automobile manufacturers in Thailand[J]. Procedia - Social and Behavioral Sciences，2010, 2(3): 5980-5989. [2] KANEKO J, NOJIRI W. The logistics of Just-in-Time between parts suppliers and car assemblers in Japan[J]. Journal of Transport Geography，2008, 16(3): 155-173. [3] ZUO Xiaolu, LIU Zhixue, ZHENG Changzheng. Analysis and optimization of Milk-run mode for auto parts logistics[J]. Automotive Engineering，2011, 33(1): 79-84(in Chinese).[左晓露, 刘志学, 郑长征. 汽车零部件循环取货物流模式的分析与优化[J]. 汽车工程，2011, 33(1): 79-84.] [4] ROSS A, JAYARAMAN V. An evaluation of new heuristics for the location of cross-docks distribution centers in supply chain network design[J]. Computers & Industrial Engineering，2008, 55(1): 64-79. [5] MA H, MIAO Z, LIM A, et al. Crossdocking distribution networks with setup cost and time window 14 CT(小时)
constraint[J]. Omega，2011, 39(1): 64-72. [6] MUSA R, ARNAOUT J P, JUNG H. Ant colony optimization algorithm to solve for the transportation problem of cross-docking network[J]. Computers & Industrial Engineering，2010, 59(1): 85-92. [7] LEE Y H, JUNG J W, LEE K M. Vehicle routing scheduling for cross-docking in the supply chain[J]. Computers & Industrial Engineering，2006, 51(2): 247-256. [8] YU W, EGBELU P J. Scheduling of inbound and outbound trucks in cross docking systems with temporary storage[J]. European Journal of Operational Research， 2008, 184(1): 377-396. [9] CHEN F, LEE C Y. Minimizing the makespan in a two-machine cross-docking flow shop problem[J]. European Journal of Operational Research，2009, 193(1): 59-72. [10] VAHDANI B, ZANDIEH M. Scheduling trucks in cross-docking systems: robust meta-heuristics[J]. Computers & Industrial Engineering，2010, 58(1): 12-24. [11] CHEN F, SONG K. Minimizing makespan in two-stage hybrid cross docking scheduling problem[J]. Computers & Operations Research，2009, 36(6): 2066-2073. [12] GUE K R. The effects of trailer scheduling on the layout of freight terminals[J]. Transportation Science， 1999, 33(4): 419-428. [13] YAL INKAYA O, MIRAC BAYHAN G. Modelling and optimization of average travel time for a metro line by simulation and response surface methodology[J]. European Journal of Operational Research，2009, 196(1): 225-233. [14] SHANG J, LI S, TADIKAMALLA P. Operational design of a supply chain system using the Taguchi method, response surface methodology, simulation, and optimization[J]. International Journal of Production Research， 2004, 42(18): 3823-3849. [15] KUMAR S, NOTTESTAD D A. Capacity design: an application using discrete-event simulation and designed experiments[J]. IIE Transactions， 2006, 38(9): 729-736. [16] NOGUERA J H, WATSON E F. Response surface analysis of a multi-product batch processing facility using a simulation metamodel[J]. International Journal of Production Economics，2006, 102(2): 333-343. [17] REIS DOS SANTOS M I, PORTA NOVA A M O. Statistical fitting and validation of non-linear simulation metamodels: a case study[J]. European Journal of Operational Research，2006, 171(1): 53-63. [18] KLEIJNEN J P C. Design and analysis of simulation experiments[M]. New York，.,USA: Springer, 2008. [19] KLEIJNEN J P C, SANCHEZ S M, LUCAS T W, et al. State-of-the-art review: a user’s guide to the brave new world of designing simulation experiments[J]. Informs Journal on Computing，2005, 17(3): 263-289. [20] ANKENMAN B, NELSON B L, STAUM J. Stochastic Kriging for Simulation Metamodeling[J]. Operations Research，2010, 58(2): 371-382. [21] KLEIJNEN J P C, VAN BEERS W C M. Robustness of Kriging when interpolating in random simulation with heterogeneous variances: Some experiments[J]. European Journal of Operational Research，2005, 165(3): 826-834. [22] KLEIJNEN J P C, VAN BEERS W C M. Application-driven sequential designs for simulation experiments: Kriging metamodelling[J]. Journal of the Operational Research Society，2004, 55(8): 876-883. 15
[23] KLEIJNEN J P C, BEERS W V, NIEUWENHUYSE I V. Constrained optimization in expensive simulation: Novel approach[J]. European Journal of Operational Research，2010, 202(1): 164-174. [24] ZAKERIFAR M, BILES W E, EVANS G W. Kriging metamodeling in multiple-objective simulation optimization[J]. Simulation，2011, 87(10): 843-856. [25] KELTON W D, SADOWSKI R P, STURROCK D T. Simulation with ARENA [M]. 4th York .,USA: McGraw-Hill，2007. [26] WAN H, ANKENMAN B E, NELSON B L. Controlled sequential bifurcation: A new factor-screening method for discrete-event simulation[J]. Operations Research，2006, 54(4): 743-755. [27] KLEIJNEN J P C. Kriging metamodeling in simulation: a review[J]. European Journal of Operational Research，2009, 192(3): 707-716. 作者简介： 施 文（1981-），男，湖北武汉人，博士研究生，研究方向：物流系统仿真及仿真优化，E-mail：shi3wen@； 刘志学（1963-），男，湖南武冈人，教授，博士生导师，研究方向：物流与供应链管理； 杨 威（1987-），男，河南信阳人，博士研究生，研究方向：物流与供应链管理。 16
