模拟退火求解作业车间调度问题：C语言实现与调参指南

简介作业车间调度问题是一类经典NP-Hard组合优化问题这份C语言资源包以最小化总完成时间为目标给出基于模拟退火算法的完整求解方案能够帮助计算机、工业工程专业学生及初级算法工程师理解元启发式方法在生产调度中的实际落地也适合作为课程设计或毕业设计的参考案例。压缩包采用zip格式共2个文件约210KB结构简洁cpp源码实现完整的退火求解流程docx文档则对JSP背景、模型构建、算法原理及实验结果进行说明便于逐行对照阅读。目前已有3498人下载学习具备一定的实践热度。源码中设计了作业工序、加工时间、设备需求等对象并覆盖初始调度生成、Metropolis接受准则、初始温度设定、冷却策略和停止条件等关键模块文档进一步补充了总完成时间等方案评价指标以及初始温度、冷却速率对搜索性能的影响分析同时介绍运行实例、参数设置与结果分析。二者结合可帮助读者从问题建模到算法编码形成完整闭环降低仿真实验和二次开发的上手门槛。1. 模拟退火解作业车间调度为什么我劝你别急着上遗传算法作业车间调度问题Job Shop Scheduling, JSP是排产场景里最常见的硬骨头n个工件的每道工序都要在指定机器上加工工序之间有先后约束机器在同一时刻只能干一件活目标是最小化最大完工时间makespan。这类问题规模稍大精确解就跑不动了于是很多人第一反应是上遗传算法。但我的血泪经验是如果只是想在C语言项目里快速落地一套能跑的调度求解器模拟退火反而更合适——它没有交叉变异那些复杂度只需要定义好邻域和温度表就能在代码量小一个量级的情况下达到接近遗传算法质量的解。这篇文章就围绕C语言实现模拟退火解JSP展开从数学模型到完整代码再到调参避坑适合有C基础、想做排产Demo或毕业设计的从业者。2. 从JSP数学模型到模拟退火四要素邻域、冷却、Metropolis准则与C语言落地2.1 JSP的工序约束与时间表建模JSP的标准定义是给定n个工件m台机器每个工件有m道工序第k道工序必须在指定机器上加工加工时间固定并且第k道工序必须等第k-1道工序完成才能开始。在大多数基准算例里每个工件恰好经过所有机器一次所以工序数等于机器数。这里有个初学者最容易忽略的点一台机器不能同时加工多个工件一个工件也不能同时上多台机器。所以解码时每个工序的开始时间不是简单地把上一个工序的结束时间相加而必须同时校验「该工件上一道工序完成时间」和「该机器可用时间」两条约束取两者中的最大值作为开工时间。数学上JSP可以建模为带优先约束的调度问题。设工件i第j道工序的完成时间为x则工序先后要求x_i,j x_i,j-1 p机器排他要求x_i,j 该机器上前一个工序的完成时间 p。这些约束在代码里最终体现为一个解码函数。我在做这个阶段时见过不少同学先把作业顺序排好然后简单用累计时间算完工时间结果所有解都是非法解跑出来的makespan比任何合法解都小其实是约束丢了。所以建模的关键不是背公式而是把解码函数写对。解码函数本质上是一个「离散事件推演」按给定的工序序列逐项推进每个工件和每台机器的状态最后取所有工件完成时间的最大值即为makespan。验证解码是否正确的办法也很直接手工把一个小算例的序列画成甘特图逐行对比程序输出。2.2 模拟退火的温度退火与邻域结构模拟退火源自金属退火固体加热到高温后缓慢冷却内能趋于最低。在JSP里状态就是一个合法的工序序列目标函数是makespan内能就是目标值。模拟退火四大要素初始解、初始温度、温度更新函数、邻域生成方式。其中邻域结构和温度表直接决定求解质量。对基于工序编码的解来说最常用的是交换两个位置swap。因为每个工件在序列里出现m次只要交换不改变每个工件的出现次数得到的新序列一定是合法解。随机交换两个位置天然保持这个性质所以实现简单。另一个常用邻域是逆序invert两个位置之间的子序列它能在一次操作内产生较大扰动适合跳出长周期局部陷阱。实际中我会混用内循环前半段用swap后半段用invert或者按概率以0.7/0.3选择。温度更新方面经典做法是几何退火T alpha * Talpha通常取0.85到0.95。alpha越接近1退火越慢解质量越高但耗时呈倍数增长。我一般先跑几次小规模算例观察目标值下降曲线再决定alpha。Metropolis准则是退火的核心如果新解比当前解好直接接受如果更差以概率exp(-delta / T)接受。这个「有概率接受坏解」的行为是模拟退火能跳出局部最优的关键。用一张表对比两种邻域的特点邻域方式扰动范围对目标值影响适用场景swap交换两点局部通常改变1~2道工序小容易拿到细粒度改进低温阶段精细搜索invert逆序子段全局可能改变多道工序顺序大有助于跳出局部陷阱高温阶段粗犷探索我在调参时遵循一个粗准则高温阶段多用invert低温阶段多用swap。代码里用随机概率控制比例比人为切换温度区间更简单也少一个需要调试的阈值。2.3 用C语言表示调度解编码、解码与合法性校验在C语言里实现JSP第一步要设计解的数据结构。我采用最主流的“基于工序的编码”用一个长度为n*m的整数数组表示解数组中的元素取值0到n-1表示工件编号同一个工件在数组中出现m次第k次出现表示该工件的第k道工序。例如2个工件各2道工序序列[0,1,0,1]表示工件0的第一道工序工件1的第一道工序工件0的第二道工序工件1的第二道工序。这种编码的好处是任何长度为n*m且满足每个工件恰好出现m次的整数序列都对应一个合法的工序先后关系。解码时按顺序扫描数组用两个数组记录状态job_done_time[i]表示工件i当前已完成工序的结束时间machine_done_time[j]表示机器j上一次加工完成时间。对当前工序开工时间为两者中的较大者。具体到C语言我用结构体保存实例数据#define MAX_JOBS 20 #define MAX_MACHINES 20 typedef struct { int num_jobs; // 工件数 int num_machines; // 机器数/每工件的工序数 int machine_seq[MAX_JOBS][MAX_MACHINES]; // 工件i第k道工序所用的机器 int time_seq[MAX_JOBS][MAX_MACHINES]; // 工件i第k道工序的加工时间 } JobShop; typedef struct { int *code; // 调度序列长度为 num_jobs * num_machines int makespan; } Solution;这里的machine_seq和time_seq是二维数组行号是工件编号列号是工序序号。例如machine_seq[2][1]4表示工件2的第二道工序使用机器4time_seq[2][1]6表示这道工序加工6个单位时间。项目输入数据通常以「每个工件的机器顺序和加工时间」给出读入后填这两个表即可。注意这里默认每个工件的工序数等于机器数这是经典JSP基准算例的设定。如果你的项目里有缺工序或返工的情况需要把二维数组改成结构体数组或链表代码复杂度会上去但核心解码逻辑不变。合法性校验的核心就是op_count计数器它保证同一工件的第k道工序不会在第k-1道工序之前被调度。3. 完整C语言实现从随机初解到邻域搜索的核心代码拆解3.1 数据结构与参数初始化下面给出一份可直接编译的C代码骨架。先看实例数据和全局声明#include stdio.h #include stdlib.h #include string.h #include math.h #include time.h #define MAX_JOBS 20 #define MAX_MACHINES 20 typedef struct { int num_jobs; int num_machines; int machine_seq[MAX_JOBS][MAX_MACHINES]; int time_seq[MAX_JOBS][MAX_MACHINES]; } JobShop; JobShop g_shop; double initial_temperature 50.0; double alpha 0.90; double min_temperature 0.1; int inner_iterations 200; int total_ops() { return g_shop.num_jobs * g_shop.num_machines; }逻辑说明g_shop 作为全局变量保存算例数据这样解码函数可以直接访问避免传指针带来的额外参数开销。alpha、初始温度、内循环次数都是全局变量方便在 main 里按算例规模临时调整。参数说明initial_temperature50 适合小算例如果目标值本身很大比如几百要按下一章的方法重新标定。inner_iterations200 对6x6算例够用对10x10算例需要提高到500以上。生成合法随机初始解void init_solution(int *code, int n) { int pos 0; for (int j 0; j g_shop.num_jobs; j) { for (int k 0; k g_shop.num_machines; k) { code[pos] j; } } // Fisher-Yates洗牌 for (int i n - 1; i 0; i--) { int j rand() % (i 1); int tmp code[i]; code[i] code[j]; code[j] tmp; } }逻辑说明先构造一个有序序列每个工件重复num_machines次再用Fisher-Yates洗牌随机打乱。洗牌比「每次随机选一个没满的工件」更稳定运行时间确定不会在最后一个位置陷入死循环。3.2 关键函数计算最大完工时间makespanmakespan是整个算法的目标函数必须保证正确且高效int calc_makespan(const int *code) { int job_done[MAX_JOBS] {0}; int machine_done[MAX_MACHINES] {0}; int op_count[MAX_JOBS] {0}; int total total_ops(); for (int idx 0; idx total; idx) { int job code[idx]; int op op_count[job]; int machine g_shop.machine_seq[job][op]; int proc g_shop.time_seq[job][op]; int start job_done[job]; if (machine_done[machine] start) { start machine_done[machine]; } int finish start proc; job_done[job] finish; machine_done[machine] finish; } int makespan 0; for (int j 0; j g_shop.num_jobs; j) { if (job_done[j] makespan) makespan job_done[j]; } return makespan; }逻辑说明线性扫描序列每次取一个工件编号通过 op_count[job] 知道这是该工件的第几道工序进而取出机器与加工时间。开工时间取两个约束的较晚者——该工件前序工序完成时间以及该机器释放时间。更新后把完成时间写回两个状态数组。时间复杂度O(n*m)在退火过程中可被反复调用。参数说明job_done 和 machine_done 初始为0表示时刻0所有工件和机器都可用。op_count[job] 从0递增保证工序严格按先后顺序执行。这个函数是只读的不修改 code所以可以安全传入新产生的邻居。3.3 邻域生成与Metropolis接受准则邻域生成采用 swap 和 invert 混合void generate_neighbor(int *code, int n) { int i rand() % n; int j rand() % n; while (i j) { j rand() % n; } // 保证交换的不是同一个工件位置否则序列不变 while (code[i] code[j]) { j rand() % n; } if ((rand() % 100) 30) { // 30%概率做逆序 int low (i j) ? i : j; int high (i j) ? i : j; while (low high) { int tmp code[low]; code[low] code[high]; code[high] tmp; low; high--; } } else { // 70%概率做交换 int tmp code[i]; code[i] code[j]; code[j] tmp; } }逻辑说明swap 适合小范围调整invert 适合大范围扰动。关键细节是第二个 while 循环——如果 code[i] 和 code[j] 是同一个工件编号交换后序列不变白白浪费一次目标函数评估。加这个判断能在小规模算例上显著避免无效迭代。Metropolis 接受函数int accept(int delta, double T) { if (delta 0) return 1; double r (double)rand() / RAND_MAX; return (r exp(-delta / T)); }delta 是新解makespan减去当前解makespan。delta 为负说明新解更好直接接受delta 为正时接受概率随温度升高而增大随差值的增大而指数下降。3.4 主函数把退火循环串起来并统计耗时有了以上几个函数主退火循环可以这样组织void run_simulated_annealing(int *best_code, int *best_energy) { int n total_ops(); int *cur_code malloc(n * sizeof(int)); int *new_code malloc(n * sizeof(int)); init_solution(cur_code, n); int cur_energy calc_makespan(cur_code); *best_energy cur_energy; memcpy(best_code, cur_code, n * sizeof(int)); double T initial_temperature; while (T min_temperature) { for (int iter 0; iter inner_iterations; iter) { memcpy(new_code, cur_code, n * sizeof(int)); generate_neighbor(new_code, n); int new_energy calc_makespan(new_code); int delta new_energy - cur_energy; if (accept(delta, T)) { memcpy(cur_code, new_code, n * sizeof(int)); cur_energy new_energy; if (new_energy *best_energy) { *best_energy new_energy; memcpy(best_code, new_code, n * sizeof(int)); } } } T * alpha; } free(cur_code); free(new_code); }逻辑说明每个温度下执行 inner_iterations 次内循环。每次生成新解根据 accept 决定是否替换当前解。即使接受坏解也不影响 best_code它始终保存历史最优解。这相当于“有记忆的退火”是我跑调度类问题的习惯。参数说明alpha 每轮乘一次总轮数约为 log(min_T / init_T) / log(alpha)。比如初温50终温0.1alpha0.9轮数约57轮。每轮200次内循环总共约1.1万次目标函数评估对6x6算例已经能稳定拿到接近最优的结果。main 函数里的调用int main() { srand(time(NULL)); // 这里读入算例并填入 g_shop // 例如 g_shop.num_jobs 6; g_shop.num_machines 6; // 逐行填写 machine_seq 和 time_seq int n total_ops(); int *best_code malloc(n * sizeof(int)); int best_energy 0; clock_t start clock(); run_simulated_annealing(best_code, best_energy); clock_t end clock(); printf(best makespan: %d\n, best_energy); printf(time: %.3f s\n, (double)(end - start) / CLOCKS_PER_SEC); free(best_code); return 0; }这里用 clock() 做耗时统计单线程调试足够。在多线程或长时间任务里建议换成高精度计时函数避免CLOCKS_PER_SEC的分辨率问题。4. 避坑指南模拟退火调参中的四个典型翻车现场模拟退火跑JSP代码本身不难难的是参数和邻域设计。以下四条是我在多个模拟项目里踩过并且事后觉得最典型的坑按现象、原因、解决三步记录。4.1 翻车一初温设太低搜索一开局就被冻结现象程序运行很快但结果很差甚至不如最简单的先到先得规则。观察接受率一开始就低于0.1后面几乎没有任何坏解被接受。原因初温T0太小Metropolis接受差解的概率exp(-delta/T)几乎为0。比如T01而邻域里makespan差值常见为2~3接受概率只有e^{-2}约0.13搜索很快退化成爬山失去了跳坑能力。解决初温应该与目标值的变动量级匹配。我一般跑10~20个随机初始解计算makespan的标准差sigma然后取初温为sigma的5~10倍。不想做统计的话可以直接用初始解makespan的1.5倍作为初温。判断标准是在第一轮内循环里坏解接受率应该在0.5~0.8之间。低于0.3说明初温太低高于0.95说明初温太高白白浪费时间。4.2 翻车二退火速度太快最后成了爬山算法现象初温合理但最终结果总比不上小规模算例已知最优解而且多跑几次结果波动很大。原因alpha设得太小比如0.8。22轮之后温度就从50降到0.1前面还没充分探索就“熄灭”了。几何退火的本质是在高温阶段多采样低温阶段精细挖优alpha直接在指数上影响总轮数。解决把alpha提到0.90~0.95并配合内循环次数一起调。我常用的经验配对alpha0.92内循环次数总工序数20。时间允许时alpha0.95配总工序数50能得到更稳的解。注意调大alpha后运行时间会显著增加先在6x6算例上对比再决定值不值得。4.3 翻车三邻域移动无效重复解导致无法逃逸现象makespan很早就停在某个值之后无论温度怎么降都不再更新。回看接受记录发现大部分新解的目标值都一样。原因swap邻域如果始终选到同一个工件编号的两个位置交换后序列不变。比如序列[0,1,0,1]交换第0和第2个位置两者都是0结果仍为[0,1,0,1]。这种无效移动在工件数少的算例里概率不低。另外如果两个位置相邻且属于同一工件的相同参数工序也会产生等价解。解决在generate_neighbor里强制code[i] ! code[j]否则重选j。更彻底的做法是当连续多轮比如20轮都没有改进时对当前解做一次“叛乱”随机取一段子序列整体逆序或者把当前最好解复制回来重新加热。这种重启操作能有效跳出平台期实现也简单。4.4 翻车四目标函数只算完工时间忽略了工序约束现象自己设计的小算例手工能算出最优解但程序运行得到的目标值比手工算还低明显不可能。原因解码实现有bug通常是没有按op_count[job]严格取工序而是按job的某道固定工序算导致某些工序的机器或前置约束被跳过。还有一种情况是忘掉把job_done[job]和machine_done[machine]都更新为finish只更新了机器没更新工件造成后一个工序提前开工。解决写一个打印调试函数对每个工序输出job、op_idx、machine、start、finish逐条检查哪些工序的start小于前序工序的finish哪些机器的start小于机器上上一工序的finish。特别要检查第一个工件的第一道工序它的start应为0。另一有效方式是把同一个解交给两个不同的解码函数对比结果我在某个图像处理Demo的调度子模块里就吃过这个亏后来把所有解码收敛到同一个实现彻底避免了口径不一致。5. 一个验证技巧用基准算例验收你的模拟退火求解器有了可以跑的代码下一步是确认它有没有写对。我的做法是先用小规模已知最优算例做硬验收再谈调参提升。5.1 小算例验收清单选择6工件6机器的小算例已知最优makespan为55。这个规模下正确实现的模拟退火即使参数不是最优跑几次也应该在55~58之间如果始终大于60优先检查解码函数而不是急着调参。验收流程如下生成随机初始解记录初始makespan通常应在80~100之间运行模拟退火输出best_energy重复运行20次统计最小、平均、波动范围把best_code打印成工序序列手工推演一遍约束确认合法。我常用的一个小数据样例机器顺序和加工时间见下表方便你对照。这张表是自构造的不代表任何标准算例最优值需要你自己穷举或程序跑出来。工件工序1机器/时间工序2机器/时间工序3机器/时间工序4机器/时间工序5机器/时间工序6机器/时间J0M2/3M1/7M5/4M3/6M4/5M0/2J1M1/6M4/2M0/5M5/8M3/4M2/3J2M3/5M0/4M4/7M2/6M1/3M5/1J3M5/2M2/5M1/3M4/6M0/7M3/4J4M0/4M3/3M2/6M1/5M5/2M4/8J5M4/5M5/1M3/7M0/2M2/6M1/4注意M0到M5只是机器编号与数组下标一一对应。填表时把M2/3写成machine_seq[0][0]2, time_seq[0][0]3即可。5.2 调参后的小动作验收通过后我还会做两个小动作第一把随机数种子固定下来用来做参数A/B对比第二记录每个温度下的接受率画出温度-接受率曲线。接受率从初温时的0.7~0.8逐渐下降到0.01左右才算是一次健康的退火。如果接受率掉得太快说明alpha太大或内循环不够如果一直不降说明初温过高或终温太高。在那以后我每次拿到一个新的JSP实例都先写解码函数然后拿已知最优的小算例跑通再进行模拟退火的调参不再像之前那样直接套模板。这套流程帮我在多个跨平台系统的排产模块里快速定位了问题也省了不少冤枉时间。希望帮到你。本文还有配套的精品资源点击获取