基于改进模拟退火算法的非线性土方调配方案优化方法

时间：2024-07-28

范斌李迟典周诚

(1.武汉地铁集团有限公司质量安全监察部，武汉 430070；2.华中科技大学国家数字建造技术创新中心，武汉 430074；3.华中科技大学土木与水利工程学院，武汉 430074)

引言

土方工程是建设项目的重要一环，其顺利、稳定完工是整个工程按时完成的前提[1]。由于土方工程的复杂性，其成本通常占整个工程成本的较大部分[2]。土方调配方案是影响土方工程成本的重要因素，针对土方调配方案的优化对成本的控制有很大作用。

将土方调配过程简化并归纳出目标函数与约束函数，再利用数学模型对其求解是土方调配优化的一般思路。当前研究多集中于线性规划模型[3]。其中，Hare考虑路径可访问性并利用混合整数线性规划模型进行求解[4]；邢志国利用综合运距替换平均运距进行模型求解，使得土方调配优化结果与实际工程成本相适应[5]；闫若钰利用Revit对土方量进行计算并利用Dynamo实现土方调配线性规划自动求解[6]。除此之外，多目标规划与动态规划的研究也较多，Parente依据成本最低和持续时间最短构建多目标动态规划模型并利用多目标遗传算法进行土方调配方案求解[7]；黎天胜以成本最低和环保为目标构建土方调配模型来实现成本节约与环境保护[8]；常峻从实时监控系统获取信息来建立动态规划模型并据此对土方调配模型进行实时规划[9]。

但当前针对非线性土方调配的研究较少。Easa考虑实际土方施工中借方和弃方的单位成本随借方和弃方量的变化而变化，构建了一个非线性土方调配模型[10]。其他研究大多则是先利用启发式算法求解线性调配问题，再提出其用于求解非线性调配问题的可能性。

现实土方调配中存在借方阶梯费率、调配量过大导致交通拥堵、调配量不同而采用不同运输工具导致费率不同等多种非线性因素，单纯的线性规划模型已不能满足实际需要，亟需对非线性规划模型进行研究。

1 土方调配非线性规划模型

土方调配非线性规划模型最突出的特点是目标函数是一个非线性函数，即土方成本与土方量的函数关系是非线性的。土方成本的非线性主要是由于在考虑土方调配非线性因素的情况下，土方工程费率不再是一个固定值，而是会随着土方量的增加而发生变化的。分析各因素对土方成本的影响原理，在线性成本的基础上提出非线性土方调配的成本函数，如式1所示：

式1中，Vij、dTij为由挖方区i调配到填方区j的调配体积与转运运距，VFi、VCi、dCi、dDi为挖方区i的填方部分的体积、挖方部分体积、收集运距及摊铺运距；VFj、VCj、dCj、dDj为填方区j的填方部分的体积、挖方部分体积、收集运距及摊铺运距；fC、fT、fD、fW、fB为土方调配中的收集费率函数、转运费率函数、摊铺费率函数、弃方费率函数及借方费率函数；Nc、Nf则分别表示挖方区数量、填方区数量。费率函数是非线性土方成本的重要计算因素，一般为调配土方量的函数。如果费率函数是一个常数，则式1可以表示为线性调配成本；如果费率函数不是常数，则式1为非线性规划的目标函数。

对于土方调配问题，其约束一般为土方量约束。现假设土方施工分区整体土方量为填方大于挖方，需要借土且所有土方材料都能满足填方需求，不需要进行土方换填等操作。针对上述工程假设，土方量约束函数用式2到式6表示：

(2)

(3)

(4)

VF=0

(5)

Vij≥0

(6)

其中，VB表示借方体积，VF表示弃方体积，其余字母含义与式1相同。式2到式6的约束条件是假设填方量多余挖方量的情况下得到的，每个公式代表不同的含义。式2表示填方区收到的土方材料体积要小于或等于自己需求的材料体积，小于部分的体积由之后的借方体积补充；式3表示挖方区外运的体积要等于自己多余的体积；式4表示借方体积的计算，因为假设总挖方量要小于总填方量且所有挖方土壤都可以满足回填要求，因此借方体积为总填方体积与总挖方体积的差；式5表示弃方体积的计算，在假设条件下的弃方体积即为零；式6表示调配土方量要满足非负性约束。每个弃方区的弃方量约束和每个借方区的借方量约束在本文中不做详细考虑。

通过将约束条件整合到目标函数中来实现对约束的处理。等式约束通过消除一个参数来转化为不等式约束，例如对于式3所示的约束，需要将Vi1作为换元主体换到等式一边，如式7中的第二个式子所示；由于Vi1是一个非负参数，因此可以将其转换为一个不等式，如下列第三个式子所示：

Vi1=Vci-VFi-(Vi2+Vi3+…+ViNf)

(7)

Vci-VFi-(Vi2+Vi3+…+ViNf)≥0

对于不等式约束则通过构建罚函数来进行转化。如果解满足不等式要求则罚函数应该为0，否则通过罚函数来对解的越界进行惩罚。本文以二次罚函数来构建罚函数：以式2为例，此约束是不等式约束，其罚函数构造如式8所示，其中σ为常数，控制罚函数F的放缩，对于不同的问题有不同的大小，需要根据实际情况考虑，本文通过多次尝试后确定σ取值为1 000。

(8)

为了利用启发式算法对有约束的非线性规划问题求解，需要将成本与约束整合为一个统一的目标函数，如式9所示。式9一共由两部分组成：第一部分即为土方成本，第二部分是罚函数。罚函数由式2、式3通过处理得到，式4和式5是借方和弃方体积，其成本不随调配方案改变而改变，不进行转换；式6为非负约束，不用转化，可以直接代入启发式算法进行计算。

通过上述论述，可以将非线性土方调配问题转化为式9所示的非负约束下的函数最小值问题，然后利用启发式算法进行求解。

2 改进模拟退火优化算法

模拟退火算法是一种通过模拟固体的退火过程进行优化的常用启发式算法，其算法流程可以用图1表示。首先，在解空间种随机产生一个新解，将这个新解作为一个初始解带入后续算法；通过一定的规律在旧解的基础上扰动产生新解并计算新旧两解目标函数值的变化量；如果变化量≤0，则接受新解；如果变化量>0，则按照式10所示的Metropolis准则判断是否接受，其中ΔE即为目标函数值的变化量；根据内循环次数要求执行内循环并判断当前温度下是否满足了终止条件，如果满足则输出最终解。

(10)

模拟退火算法可以跳出局部最优获得全局最优解，过程简单，鲁棒性强，在非线性优化问题的求解上有较大优势。但该算法也有一些缺点，例如算法运行时间较长，算法运行中可能会导致最优解丢失，不同参数对于算法结果影响较大等问题。

图1 模拟退火算法流程

针对算法运行时间较长的问题，本文参考L. Ingber提出的方法对降温函数与扰动方式进行改进[11]。其降温函数改进为式11所示的指数函数，其中，T0为初始温度，k为降温次数，α为降温系数，一般为[0.7， 1]之间取值，D一般取1，此降温函数可显著提高降温速度。改进的扰动方式如式12所示，旧解x的i维分量xi的上下限为[mini，maxi]，则其i维分量xi经过扰动后产生的新解为xi’，sgn为阶跃函数，T为当前情况下的温度，ui为[0， 1]上均匀分布的随机数。

T=T0αk1/D

(11)

针对算法运行中的最优解丢失问题，本文在算法中加入最优解保留策略，最优解仅记录整个过程中的最优目标函数值而不参与退火过程。

对于不同参数对算法结果影响较大，本文尝试通过对比调参来保证结果的可靠性。

3 算法验证与对比

3.1 案例简介

为了验证本文方法的效果，选用武汉某土方工程进行计算对比与验证。该工程填挖方量大，土方工期长，土方工程成本大，具有代表性。本文选择三个不同的算例进行算法分析：根据整体工程分包分区情况和施工段划分情况进行土方工程分区，具体分区数量为6，如图2所示；根据现场实际的土方工程分区进行分区调配，具体分区数量为13，如图3所示；根据区域联通性和区域面积进行分区，具体分区数量为29，如图4所示。其中13分区情况为实际的分区情况，故将其作为主要研究场景进行研究并进行算法比选。为了解不同问题规模下改进算法的求解性能，本文利用改进算法对6分区问题和29分区问题进行求解，分析不同的问题规模下算法的变化。参考实际情况，假设土方调配过程中土方转运量不大于20 000m3时，收集费率、转运费率和摊铺费率为6.65元/(km·m3)；如果土方转运量大于20 000m3，土方费率降低为原费率的0.9，即为5.985元/(km·m3)，不考虑其他影响土方成本的线性与非线性因素。

图2 6分区方案划分示意图

图3 13分区方案划分示意图

图4 29分区方案划分示意图

表1 各算法可选参数

3.2 算法调参及结果展示

分析前人的研究经验，对于本文的土方调配优化问题，拟采用改进模拟退火、遗传算法、粒子群算法和差分进化算法对比求解，各算法的可选参数组合如表1所示。其中改进模拟退火的冷却进度表还有一个参数为终止温度，这里将终止温度固定为10-7来减少调参。为了减小算法的随机性，每个参数组合均计算三次来减小误差，参数组合利用“算法种类(SA/GA/PSO/DE)-参数一—参数二—参数三—参数四”来表示。

将6分区调配问题带入改进模拟退火计算得到1 050个运行结果。因为数据量过多，将所有数据按照马尔可夫链长度数N进行分类，以初始温度和降温系数的参数组合为横轴，土方成本为纵轴绘制图形，如图5所示。由图可知，最优值为1 960 892元，当N大于500 000时，基本所有参数组合均接近这个最优值。例如起始温度=100，马尔可夫链长度=500 000，降温系数=0.75得到的最优值就是1 960 892元，此参数下目标函数随时间变化如图6所示。

图5 6分区改进模拟退火算法结果汇总

按照同样的方法对分区数量为13、29的两种情况求解，得到如下结果：对于13分区问题，调参结果如图7所示，改进算法得到的最小成本为1 765 672元，参数组合SA-1000000-300000-0.85，此参数组合下的算法随时间收敛如图8所示；对于29分区，结果如图9所示，成本的最小值为1 922 424元，是在参数组合SA-1000000-600000-0.8处得到的，此参数组合下的算法随时间收敛如图10所示。

图6 6分区改进模拟退火最优解之一的收敛情况

图7 13分区改进模拟退火算法结果汇总

图8 13分区下改进模拟退火最优解收敛过程

图9 29分区改进模拟退火算法结果汇总

图10 29分区下改进模拟退火最优解变化过程

图11 13分区下遗传算法种群数量和终止代数对比图

对于遗传算法，种群数量和终止代数是影响算法运行时间的主要参数，因此根据种群数量和终止代数的组合可以将算法结果分为18类，选取每类最优的参数组合，得到表3所示的表格。将表3得到的18个数据按照种群数量N分类汇总可以得到如图11所示的不同参数对比图。从图11可以得出针对13分区问题，遗传算法最优的参数组合为GA-3000-3000-0.99-0.001，得到的最优解为1 772 355元，最优参数组合的目标函数随算法运行时间的图像，如图12所示。

图12 13分区下遗传算法最优解收敛过程

表3 13分区遗传算法交叉概率和变异概率对比表

对于粒子群算法，种群数量和终止代数是影响算法运行时间的主要参数，根据种群数量和终止代数的组合可以将算法结果分为9类，得到表4和图13所示的不同参数对比图。从图13可以得出，针对13分区问题，粒子群算法最优的参数组合为GA-3000-12000-1.1-1，得到的最优解为1 858 936元，最优参数组合的目标函数随算法运行时间的图像如图14所示。

表4 13分区粒子群算法惯性权重和学习因子对比表

图13 13分区下粒子群算法种群数量和终止代数对比图

图14 13分区粒子群算法最优解收敛过程

同理，对于差分进化算法，可以得到表5和图15。最优参数组合的目标函数随算法运行时间变化的图像如图16所示。

表5 13分区差分进化算法交叉概率和放缩因子对比表

图15 13分区下差分进化算法种群数量和终止代数对比图

图16 13分区差分进化算法最优解收敛过程

3.3 算法结果分析

以13分区的不同算法收敛情况进行对比分析，得到算法对非线性调配的模型的求解效果，汇总四种算法的最优参数组合，最优目标函数和算法运行时间，得到表6所示的各种算法效果对比表。

表6 13分区非线性调配的各算法效果对比

由表6可以得出，对于13分区的情况，改进模拟退火得到最优解；遗传算法与改进模拟退火的计算时间大致相同，得到的解比改进算法的最优解大0.4%；粒子群算法计算时间最长，得到成本解比改进算法大5.3%；差分进化算法的计算时间最短，但是其成本相较于改进算法求得的最优解大1.2%。综上，改进模拟退火的解最优且运行时间较短，选取改进模拟退火作为非线性土方调配的最优求解方法。

在此基础上，为了解不同问题规模对改进模拟退火的求解性能和参数组合的影响，利用改进模拟退火对该项目6分区、13分区和29分区情况下的非线性土方调配求解，得到不同分区数量下最优参数组合、最优目标函数和最优组合的函数图像，如表7所示。

表7 不同分区下改进模拟退火算法效果对比

由前文可知，分区数量为6时，有多个参数组合都能得到最优解，此情况下的求解难度较小；分区数量为13时，仅有一个最优参数组合且得到的成本解是收敛的，比6分区得到的成本小，求解难度适中；分区数量为29的情况下也仅有一个最优参数组合，此参数组合下的成本仍有下降可能，求解难度最大。通过上述分析可知，随着分区数量的增加，求解难度逐步增加，对于不同规模的问题，最优参数组合也不同。