GAMS二次开发教程17生产计划优化实战——多期多工厂产能-库存-运输 MILP版本与许可声明工具/软件GAMS版本号以官方发行说明为准本篇以gams.com/latest/docs指向的当前分发为口径求解器侧示例取自 CPLEX 手册CPLEX 22.1以官方发行说明为准。语言/环境GAMS 语言 GDX模型类型为混合整数线性规划Mixed Integer Linear Program, MIP。本文目标把三个月产量怎么排、货怎么调、库存留多少、缺货认不认这一整类产销协同问题翻译成一份可求解、可报表、可回归的 MIP。许可边界求解 MIP 需要支持 MIP 的求解器许可demo license 可用规模上限以官方许可说明为准本篇的相关性重跑会连续求解三次 MIP注意许可类型是否限制实际使用次数见第 15 篇提及的许可会话问题。一句话结论多期多工厂计划的骨架只有四块——三张清单定边界决策变量 / 目标 / 约束、库存平衡方程定时间轴inv(p,t) e inv0(p)$(ord(t)1) inv(p,t-1)$(ord(t)1) y(p,t) - sum(m$link(p,m), x(p,m,t))两个互斥的$条件把首期与后续期一次写对、binary 最紧 big-M 定开停capMaxLim(p,t) .. y(p,t) l capMax(p) * run(p,t) ;这里的M就直接取该装置的产能上界永远不要用一个够大就行的数、稀疏运输变量定网络运输变量与求和都带同一个$link(p,m)而缺货必须用重惩罚项而不是硬约束否则需求高峰会让整个 MIP 变成不可行。〇、本篇要解决的认知问题Q1一个业务问题产销协同怎么变成模型决策变量 / 目标 / 约束三张清单各自怎么列Q2多期模型的库存平衡方程怎么写才不会在首期越界、又能在末期正确结转Q3开工与否为什么必须引入binary变量big-M 该取多大开了就必须产够的最小运行量怎么加Q4运输/调拨的网络流结构怎么表达为什么运输变量的声明与求和要带同一个$link(p,m)Q5目标函数的成本项该怎么构成缺货为什么用惩罚项需求 ±10%的敏感性重跑怎么组织、结果怎么收集一、机制解析1.1 从业务到模型先写三张清单业务侧的说法通常是这三个月各厂该产多少、货往哪调、库存留多少、实在供不上就缺货。要落成模型必须先把它拆成三张清单这三张清单必须在动笔写 GAMS 之前先用中文写出来清单本问题的答案落到 GAMS 的语句决策变量要决定什么各厂各期产量、各连线各期运量、各厂各期末库存、各市场各期缺货量、各厂各期是否开工VariablePositive VariableBinary Variable目标要最小化/最大化什么总成本最小固定开工成本 可变生产成本 运输成本 库存持有成本 缺货惩罚cost方程 Solve ... using mip minimizing z约束什么不能被违反产能上界与最小运行量、库存平衡、需求满足可缺货capMaxLim/capMinLim/invBal/demSat这张表的价值在于没有被列进约束清单的业务规则就不该出现在方程里反过来凡是业务上绝对不允许的就必须列进来。所有我们没有约束它它怎么自己跑偏了的困惑答案都在这里。1.2 库存平衡方程多期模型的时间轴多期模型的时间耦合完全由库存平衡方程承担。它是本期与上期之间唯一的桥期末库存(t) 期初库存(t) 本期产量 - 本期发货量其中期初库存在首期是给定的初始库存inv0(p)在之后各期就是上一期的期末库存inv(p,t-1)。在 GAMS 里这两种情形用两个互斥的$条件一次写完invBal(p,t) .. inv(p,t) e inv0(p)$(ord(t) 1) inv(p,t-1)$(ord(t) 1) y(p,t) - sum(m$link(p,m), x(p,m,t)) ;两条必须记住的细节inv(p,t-1)在首期ord(t)1是越界的靠$(ord(t) 1)把它关掉GAMS 不会去访问不存在的t0e是等式注意不是那属于 Python/GAMSpy 侧写法。末期结转则是免费的——只要让inv(p,t12)留在模型里并进入目标函数holdCost * inv(p,t)模型自己就会权衡提前生产备货与当期生产的成本不需要额外的终端条件。1.3binary big-M 最小运行量开停逻辑的标准三件套产能上限这句话在业务上有两种截然不同的含义必须分清连续上限装置在运行状态下能产多少y(p,t) l capMax(p)开停上限装置不开工就一点也不能产必须引入Binary Variable run(p,t)并用 big-M 把它与产能挂钩。GAMS 里有一个直接支持带条件的上界的写法y.up(p,t)$(...) ...但它只能表达静态上界表达不了由变量run决定的耦合。因此开停耦合的标准写法就是capMaxLim(p,t) .. y(p,t) l capMax(p) * run(p,t) ; capMinLim(p,t) .. y(p,t) g capMin(p) * run(p,t) ;两条一起看才完整第一条说不开工产量为零第二条说开了工产量必须达到最小经济运行负荷。没有第二条模型的最优解会变成每期都开工、每期产一点点——因为固定成本被分摊到极小批量上就会出现业务上根本做不到的排产方案。big-M 取值的铁律经验法则big-M只能取该约束下最紧的合法上界。本问题里就是capMax(p)——它是y(p,t)自身允许的最大值再大就是纯粹的松弛会让分支定界更慢、线性松弛更弱这一条不构成性能承诺只是建模的正确性要求。永远不要写100000 * run(p,t)这种够大就行的数。1.4 运输/调拨稀疏网络流运输变量天然是一个网络流结构x(p,m,t)从工厂p流向市场m。它有两条必须同时满足的纪律声明稀疏只声明实际存在的连线Set link(p,m)而不是p × m的全笛卡尔积第 16 篇 1.3 节使用端也带同一条件方程里写sum((p,m,t)$link(p,m), shipCost(p,m) * x(p,m,t))与sum(p$link(p,m), x(p,m,t))。声明与使用两端带同一个条件是稀疏建模最容易漏的一条。如果只在方程标号处带条件但求和里引入了一个不属于link的(p,m)组合那个组合的x会因为出现在方程里而被自动隐含扩展进定义域——稀疏声明就白做了。1.5 目标函数五项成本与缺货惩罚本篇的目标函数由五项构成每一项都对应一条业务口径成本项表达式业务含义固定开工成本sum((p,t), fixCost(p) * run(p,t))只要开工就发生与产量无关可变生产成本sum((p,t), varCost(p) * y(p,t))与产量成正比运输成本sum((p,m,t)$link(p,m), shipCost(p,m) * x(p,m,t))按连线、按量计费库存持有成本sum((p,t), holdCost * inv(p,t))占用期末库存推动按需生产缺货惩罚sum((m,t), shortCost * short(m,t))惩罚未满足的需求缺货为什么用惩罚项而不是硬约束sum(p$link(p,m), x(p,m,t)) g dem(m,t)因为硬约束意味着需求必须被满足一旦产能不足需求高峰、检修期整个 MIP 直接不可行你得到的是没有方案而不是代价最小的方案。改成需求满足方程带一个缺口变量short(m,t)、缺口在目标里被重罚以后模型永远有可行解且在产能充裕时short会被压到零。业务上这叫软约束 惩罚是产能规划类模型的标配。1.6 求解、报表与敏感性重跑模型类型由于存在run(p,t)这样的binary变量模型类型是mip。求解语句为Solve plan using mip minimizing z ;求解器用option MIP CPLEX ;指定MIP 的容差plan.optCR是相对最优性容差第 11 篇它决定搜索到多大的相对 gap 就可以停。教学示例取0.01它的语义与各求解器的实现差异以官方文档为准状态卡口MIP 的状态码与 LP 不同modelStat可能出现整数解这类取值完整码表以官方文档为准。工程上的做法是卡住plan.solveStat 1再加一条plan.modelStat的检查敏感性重跑用第 9 篇的loop范式把需求整体缩放 0.9 / 1.0 / 1.1每个场景重求解一次把 KPI 收进一个二维参数kpiS(s,*)。注意binary变量在多次solve之间会带着上一次的解第 13 篇的多 solve 自动传递初始点这对 MIP 来说是冷却起点还是热起点取决于求解器是否使用铁律 6因此重跑结论只能写成本算例上如此不做性能承诺铁律 7。二、完整代码与逐行剖析代码 17-1planning_core.gms—— 数据层 模型层不含Solve$title planning_core - multi period multi plant capacity-inventory-transport MILP core * 适用分发口径Variable/Equation/Model 语法以官方 Model and Solve 章节为准 * 求解许可本文件不含 SOLVE无需求解器许可 * 平台GAMS 语句通用 * 数值口径以下全部为教学示例值真实项目必须替换为数据源并标注来源铁律 10 Set t planning periods (months) / t1*t12 / p production plants / p1*p3 / m markets / m1*m4 / ; * 稀疏网络只声明实际存在的 plant-market 连线而不是 p x m 的全笛卡尔积 Set link(p,m) existing plant-market lanes / p1.(m1,m2) p2.(m2,m3) p3.(m3,m4) / ; Parameter dem(m,t) demand at market m in period t (cases) capMax(p) max production rate at plant p per period (cases) capMin(p) min production rate at plant p per period when running (cases) fixCost(p) fixed cost at plant p per period when running (USD) varCost(p) variable production cost at plant p (USD per case) shipCost(p,m) shipping cost on lane (p,m) (USD per case) inv0(p) initial inventory at plant p (cases) holdCost inventory holding cost (USD per case per period) / 0.8 / shortCost shortage penalty (USD per case) / 25.0 / ; dem(m,t) 120 12 * ord(t) 10 * ord(m) ; capMax(p) 900 ; capMin(p) 100 ; fixCost(p) 4000 ; varCost(p) 2.0 0.4 * ord(p) ; shipCost(p,m)$link(p,m) 1.2 0.3 * ord(p) 0.2 * ord(m) ; inv0(p) 50 ; Variable y(p,t) production quantity at plant p in period t (cases) x(p,m,t) shipment from plant p to market m in period t (cases) inv(p,t) end-of-period inventory at plant p in period t (cases) short(m,t) unmet demand at market m in period t (cases) run(p,t) binary: plant p is running in period t z total cost (USD) ; Positive Variable y, x, inv, short ; Binary Variable run ; Equation cost total fixed, variable, shipping, holding and shortage cost capMaxLim(p,t) max production rate, coupled to the run binary (big-M) capMinLim(p,t) min production rate when running invBal(p,t) inventory balance at plant p in period t demSat(m,t) demand satisfaction at market m in period t ; cost .. z e sum((p,t), fixCost(p) * run(p,t)) sum((p,t), varCost(p) * y(p,t)) sum((p,m,t)$link(p,m), shipCost(p,m) * x(p,m,t)) sum((p,t), holdCost * inv(p,t)) sum((m,t), shortCost * short(m,t)) ; * big-M 取该装置产能上界本身这是最紧的合法取值不用够大就行的数 capMaxLim(p,t) .. y(p,t) l capMax(p) * run(p,t) ; * 最小运行量不开工则产量为零开工则必须达到最小经济运行负荷 capMinLim(p,t) .. y(p,t) g capMin(p) * run(p,t) ; * 首期接初始库存之后各期接上期期末库存两个互斥的 $ 条件 invBal(p,t) .. inv(p,t) e inv0(p)$(ord(t) 1) inv(p,t-1)$(ord(t) 1) y(p,t) - sum(m$link(p,m), x(p,m,t)) ; * 软约束需求满足方程带缺口变量缺口在目标里被重罚 demSat(m,t) .. sum(p$link(p,m), x(p,m,t)) short(m,t) e dem(m,t) ; Model plan / all / ;逐行剖析Set link(p,m)用p1.(m1,m2)圆括号简写声明一条工厂的全部下线把网络结构写成数据而不是维度积shipCost(p,m)$link(p,m) ...只在连线上赋成本其余格保持 0 且不会被引用Binary Variable run ;声明开停指示变量——注意它是单独的Binary语句与Positive Variable y, x, inv, short ;分开写因为二值变量不具备非负连续的语义cost方程用五行sum把五项成本逐项列出顺序与 1.5 节的表完全一致这样目标里有哪几块钱一眼可查capMaxLim(p,t) .. y(p,t) l capMax(p) * run(p,t) ;是 big-M 耦合M就是capMax(p)capMinLim与它成对缺了这条就会出现每期都开工、每期产一点点的假最优invBal里inv0(p)$(ord(t) 1)与inv(p,t-1)$(ord(t) 1)是两个互斥条件首期只接初始库存、其余期只接上期期末inv(p,t-1)在首期不会越界demSat用e把发货 缺口 需求写成恒等式于是short(m,t)的含义精确等于缺口量缺货 KPI 可以直接用它的levelModel plan / all / ;把五条方程一次装配。代码 17-2planning_run.gms—— 单次求解 期间 × 工厂报表 GDX 导出$title planning_run - single MIP solve plus period by plant report and gdx export * 适用分发口径option / 模型属性optCR / resLim / iterLim以官方 Option Statement 页为准 * 求解许可本文件执行 SOLVE需要支持 MIP 的求解器许可以官方许可说明为准 * 平台GDX 跨平台通用 * 前置条件与 planning_core.gms 放在同一目录 $include planning_core.gms option MIP CPLEX ; * 求解器选择MIP 需要支持整数规划的求解器 plan.optCR 0.01 ; * 相对最优性容差搜索到 1% 相对 gap 即可停 plan.resLim 300 ; * 秒级资源上限先排除被时间截断这一可能 plan.iterLim 2000000 ; * 迭代上限 Solve plan using mip minimizing z ; * ---- 状态卡口MIP 的 modelStat 码表以官方文档为准 ---- abort$(plan.solveStat 1) SOLVE did not terminate normally: check solvestat ; abort$(plan.modelStat 1) MODEL is not optimal: check modelstat ; * ---- 报表参数期间 x 工厂的产量与期末库存以及按连线的运量 ---- Parameter prod(p,t) production by plant and period (cases) invEnd(p,t) end-of-period inventory by plant and period (cases) shipLane(p,m,t) shipped quantity by lane and period (cases) kpi(*) key performance indicators ; prod(p,t) y.l(p,t) ; invEnd(p,t) inv.l(p,t) ; shipLane(p,m,t)$link(p,m) x.l(p,m,t) ; kpi(totalCost) z.l ; kpi(fixCostPart) sum((p,t), fixCost(p) * run.l(p,t)) ; kpi(varCostPart) sum((p,t), varCost(p) * y.l(p,t)) ; kpi(shipCostPart) sum((p,m,t)$link(p,m), shipCost(p,m) * x.l(p,m,t)) ; kpi(holdCostPart) sum((p,t), holdCost * inv.l(p,t)) ; kpi(shortCostPart) sum((m,t), shortCost * short.l(m,t)) ; kpi(totalShort) sum((m,t), short.l(m,t)) ; kpi(plantPeriods) sum((p,t), run.l(p,t)) ; option decimals 2 ; option dispWidth 14 ; Display kpi, prod, invEnd, plan.modelStat, plan.solveStat ; * ---- 落一份文本报表便于归档、diff 与人工抽查 ---- file rep / planning_report.txt / ; put rep ; put period plant production inventory / ; loop((p,t), put t.tl:8, p.tl:8, y.l(p,t):14:1, inv.l(p,t):14:1 /) ; putclose rep ; * ---- 层间契约把解与 KPI 写成 GDX供 Python / Engine 侧消费 ---- execute_unload planning_solution.gdx, prod, invEnd, shipLane, kpi, run ;逐行剖析$include planning_core.gms让本文件只负责求解 输出两件事这正是第 16 篇的四层分离option MIP CPLEX ;指定求解器第 12 篇的选型矩阵MIP 与 LP 的求解器选项是分开的不要用option LP去指定 MIP 求解器plan.optCR 0.01是模型属性层的设置优先级高于option语句含义是相对最优性容差abort$(plan.solveStat 1)卡是否正常终止abort$(plan.modelStat 1)卡是否拿到最优两道闸门缺一不可——MIP 最危险的失败模式是被resLim截断却拿了一个漂亮的整数解prod(p,t) y.l(p,t)一行就把期间 × 工厂的产量矩阵整理好.l取变量的 levelKPI 里把成本按项拆开固定 / 可变 / 运输 / 库存 / 缺货这是产能规划类项目最常被追问的一张表——“这块钱花在哪了”put rep ;与loop((p,t), put ...)生成一份可归档的文本报表t.tl/p.tl取集合元素的文本标签左对齐输出最后execute_unload planning_solution.gdx, ...把解与 KPI 一起写成 GDX成为下游的接口。代码 17-3planning_sens.gms—— 需求 ±10% 的敏感性重跑与回归闸门$title planning_sens - demand plus or minus ten percent sensitivity rerun * 适用分发口径loop 与 Set/Parameter 语法以官方文档为准 * 求解许可本文件对每个场景执行一次 SOLVE共三次均需要 MIP 求解器许可 * 平台GAMS 语句通用 * 前置条件与 planning_core.gms 放在同一目录 $include planning_core.gms Set s demand scenarios / low, base, high / ; Parameter scale(s) demand multiplier of each scenario kpiS(s,*) key performance indicators by scenario demBase(m,t) baseline demand kept for the reruns ; scale(low) 0.9 ; scale(base) 1.0 ; scale(high) 1.1 ; demBase(m,t) dem(m,t) ; option MIP CPLEX ; plan.optCR 0.01 ; plan.resLim 300 ; loop(s, dem(m,t) demBase(m,t) * scale(s) ; Solve plan using mip minimizing z ; abort$(plan.solveStat 1) a scenario did not terminate normally ; kpiS(s,totalDemand) sum((m,t), dem(m,t)) ; kpiS(s,totalCost) z.l ; kpiS(s,plantPeriods) sum((p,t), run.l(p,t)) ; kpiS(s,totalShort) sum((m,t), short.l(m,t)) ; kpiS(s,totalInv) sum((p,t), inv.l(p,t)) ; ) ; option decimals 2 ; option dispWidth 16 ; Display kpiS ; * ---- 回归闸门需求上升时总成本不应低于需求下降时的总成本 ---- abort$(kpiS(high,totalCost) kpiS(low,totalCost)) REGRESSION: cost must not fall when demand rises ; execute_unload planning_sensitivity.gdx, kpiS ;逐行剖析$include planning_core.gms让数据 模型只写一次敏感性文件只写驱动逻辑demBase(m,t) dem(m,t) ;先留一份基线需求因为循环里dem会被反复覆盖——不先备份第二次循环就会拿上一次缩放后的需求再缩放一次这是这类脚本最经典的静默 bugloop(s, ...)是第 9 篇的参数化重求解范式循环体内第一件事是重算需求第二件事就是Solveabort$(plan.solveStat 1)放在循环内任何一个场景失败都立刻中止避免两个场景成功、一个场景挂了却还在往下汇总kpiS(s,...)是二维参数第一维是场景、第二维是自由指标用*声明这张表就是敏感性分析的最终交付物末尾的abort$(kpiS(high,totalCost) kpiS(low,totalCost))是业务级回归闸门需求整体上升而总成本反而更低只能说明数据、缩放或模型被改坏了——把它写成硬门禁比人眼扫表可靠得多最后execute_unload把敏感性表也导出成 GDX。三、常见报错与排查排查 3-1模型不可行modelStat显示不可行。现象需求高峰月份无解。根因把需求必须满足写成了硬约束sum(p$link(p,m), x(p,m,t)) g dem(m,t)而该期产能含capMin/capMax与检修根本不够。解法改成软约束demSat(m,t) .. sum(p$link(p,m), x(p,m,t)) short(m,t) e dem(m,t) ;把缺口交给short(m,t)并给一个足够高的shortCost这样模型永远可行缺货量本身就是一条业务指标。别用把shortCost设成很小的数来粉饰——那等于告诉模型缺货很便宜最优解会真的去缺货。排查 3-2最优解出现每期都开工、每期只产一点点。现象run.l(p,t)几乎全为 1产量却极低。根因只写了 big-M 上界y l capMax * run漏了最小运行量下界于是固定成本被摊到极小批量上反而更便宜。解法补上capMinLim(p,t) .. y(p,t) g capMin(p) * run(p,t) ;让开工必须至少产到最小经济运行负荷。改完要复核kpi(plantPeriods)开工的厂·期数是否回到业务上可执行的量级。排查 3-3稀疏网络漏了条件模型规模暴涨。现象明明只声明了几条连线方程项数却按全笛卡尔积增长。根因某处求和写了sum((p,m,t), ...)而漏了$link(p,m)或某个方程标号处漏了条件导致不属于link的组合被自动隐含扩展进定义域。解法“声明与使用两端带同一个条件”——所有涉及运输的项统一写成sum((p,m,t)$link(p,m), ...)与sum(p$link(p,m), ...)改完用.lst的Equation Listing复核项数第 16 篇 1.3 节。排查 3-4MIP 目标值看起来合理但其实是被时限截断的解。现象z.l与预期量级相符solveStat却不是 1达到resLim。根因MIP 达到资源上限时会返回当前最好的整数解而modelStat可能指示整数解而非最优完整码表以官方文档为准。解法把abort$(plan.solveStat 1)与abort$(plan.modelStat 1)固定成求解层的尾巴判断算不完还是模型难先看是不是被resLim/iterLim卡住第 16 篇 1.6 节再考虑放宽optCR或调求解器参数。排查 3-5敏感性循环里第二个场景的数据被二次缩放。现象low场景的成本反而比重跑前更高/更低得离谱。根因循环里直接写dem(m,t) dem(m,t) * scale(s)用的是上一次循环已经缩放过的需求缩放被累乘。解法像代码 17-3 那样先备份基线demBase(m,t) dem(m,t) ;循环里从demBase重新计算。另一个同源隐患是inv/short等变量的level在多个solve之间自动传递第 13 篇——它不改变最优解但会让你以为新一轮是冷启动判断收敛行为时要显式清零或明确接受残留起点并按铁律 6 区分GAMS 传了和求解器用了。四、动手练习练习 1基线跑通依次运行gams planning_run.gms。判定(a) 生成planning_report.txt与planning_solution.gdx(b)modelStat与solveStat都是 1©kpi(totalCost)等于fixCostPart varCostPart shipCostPart holdCostPart shortCostPart用Display核对到小数点后两位。练习 2缺货取舍把shortCost从 25 改成 60重跑。判定(a)kpi(totalShort)不会上升(b) 能解释缺货惩罚提高后模型是否愿意多开工、多留库存。练习 3最小运行量的必要性把capMinLim方程从Model plan / all / ;里去掉或把capMin(p)设为 0重跑。判定能指出kpi(plantPeriods)与逐期产量分布的变化并说明为什么去掉下界会得到业务上不可执行的点点开工方案。练习 4big-M 松紧对照把capMaxLim里的capMax(p) * run(p,t)改成100000 * run(p,t)重跑。判定(a) 最优目标值不变(b).lst里求解器报告的节点数/迭代量发生变化© 能写出结论——big-M 取最紧上界是建模正确性要求不是性能调优。练习 5敏感性重跑运行gams planning_sens.gms。判定(a) 生成planning_sensitivity.gdx(b)kpiS(high,totalCost)不小于kpiS(low,totalCost)© 能把demBase那一行删掉复现二次缩放的静默错误并说清原因。五、小结与下一篇预告本篇把产销协同问题落成了可交付的 MIP核心是四件事三张清单把业务翻译成模型决策变量 / 目标 / 约束没被列进约束清单的规则不该出现在方程里库存平衡方程承担全部时间耦合inv0(p)$(ord(t)1)与inv(p,t-1)$(ord(t)1)两个互斥条件一次写对首期与后续期末期结转由目标函数免费完成binary 最紧 big-M 最小运行量是开停逻辑的标准三件套只写上界会出现每期都开工、每期产一点点的假最优big-M 永远取该约束最紧的合法上界稀疏运输变量搭网络流声明与使用两端都带$link(p,m)而缺货必须用重惩罚的软约束否则需求高峰会让模型直接不可行。最后用loop KPI 表做 ±10% 敏感性重跑并用一条业务级abort$回归闸门把它锁住。四条铁律在这里落地许可先讲清2MIP 需要支持整数规划的求解器许可、初值机制说明白6多solve之间的level传递与求解器是否采纳、调优不承诺性能7、数值不臆造10全部教学示例值。第 18 篇《优化调度实战离散/连续时间调度模型与 MIP 建模技巧》把时间从期间细化到时间格离散时间与连续时间建模的取舍精度 vs 规模、用binary big-M 或 SOS 表达顺序与重叠、启动/停机成本与最小运行/停机时间、时间窗与优先级约束以及用mipTrace与初值给 MIP 搜索加速的工程手法最后用put输出甘特式文本报表——那是在同一个时间轴上把计划的月换成调度的小时。本篇认知问题回显FAQQ1一个业务问题产销协同怎么变成模型决策变量 / 目标 / 约束三张清单各自怎么列A做法是在动笔写 GAMS 之前先用中文把三张清单写出来。决策变量清单回答要决定什么本问题里是各厂各期产量y(p,t)、各连线各期运量x(p,m,t)、各厂各期末库存inv(p,t)、各市场各期缺货量short(m,t)、各厂各期是否开工run(p,t)分别用VariablePositive VariableBinary Variable落地。目标清单回答要最小化什么总成本最小由固定开工成本、可变生产成本、运输成本、库存持有成本、缺货惩罚五项构成写在cost方程里、由Solve plan using mip minimizing z驱动。约束清单回答什么不能被违反产能上界与最小运行量、库存平衡、需求满足可缺货分别落在capMaxLim/capMinLim/invBal/demSat四条方程上。这张表的关键价值是双向的没被列进约束清单的业务规则不该出现在方程里业务上绝对不允许的就必须列进来——所有没约束它、它怎么跑偏了的困惑答案都在这一步。Q2多期模型的库存平衡方程怎么写才不会在首期越界、又能在末期正确结转A多期模型的时间耦合完全由库存平衡方程承担语义是期末库存(t) 期初库存(t) 本期产量 − 本期发货量。首期的期初库存是给定初始库存之后各期的期初库存是上期期末库存两种情形用两个互斥的$条件一次写完inv(p,t) e inv0(p)$(ord(t) 1) inv(p,t-1)$(ord(t) 1) y(p,t) - sum(m$link(p,m), x(p,m,t)) ;。要点有三inv(p,t-1)在首期是越界的靠$(ord(t) 1)关掉它GAMS 就不会去访问不存在的t0e是 GAMS 的等式运算符不是后者属于 Python/GAMSpy 侧写法末期结转不需要额外终端条件——只要让inv(p,t12)留在模型里并进入目标函数holdCost * inv(p,t)模型自己就会权衡提前备货与当期生产的成本。Q3开工与否为什么必须引入binary变量big-M 该取多大开了就必须产够的最小运行量怎么加A因为产能上限有两种含义必须分清连续上限装置运行状态下能产多少y(p,t) l capMax(p)与开停上限不开工就一点也不能产。后者是一个由变量决定的上界用静态的y.up(p,t)$(...)表达不了必须引入Binary Variable run(p,t)并用 big-M 耦合capMaxLim(p,t) .. y(p,t) l capMax(p) * run(p,t) ;。big-M 的取值铁律是只能取该约束下最紧的合法上界——本问题里就是capMax(p)本身y(p,t)允许的最大值再大就是纯粹松弛永远不要写100000 * run(p,t)这种够大就行的数。最小运行量用下界形式加capMinLim(p,t) .. y(p,t) g capMin(p) * run(p,t) ;。只写上界、不写下界最优解会变成每期都开工、每期产一点点固定成本被摊到极小批量上这在业务上根本做不到是产能模型最典型的假最优。Q4运输/调拨的网络流结构怎么表达为什么运输变量的声明与求和要带同一个$link(p,m)A运输变量天然是网络流结构x(p,m,t)从工厂流向市场。两条纪律必须同时满足。第一是声明稀疏只声明实际存在的连线Set link(p,m) existing plant-market lanes / p1.(m1,m2) p2.(m2,m3) p3.(m3,m4) / ;而不是p × m的全笛卡尔积。第二是使用端也带同一条件目标里写sum((p,m,t)$link(p,m), shipCost(p,m) * x(p,m,t))需求满足里写sum(p$link(p,m), x(p,m,t))。为什么两端都要带因为 GAMS 会把出现在方程里的索引组合隐含扩展进符号的定义域——如果只在方程标号处带条件求和里引入了一个不属于link的(p,m)组合稀疏声明就白做了模型规模会按全笛卡尔积增长。改完后用.lst的Equation Listing复核项数是验证手段。Q5目标函数的成本项该怎么构成缺货为什么用惩罚项需求 ±10%的敏感性重跑怎么组织、结果怎么收集A目标由五项构成每项对应一条业务口径固定开工成本sum((p,t), fixCost(p) * run(p,t))、可变生产成本sum((p,t), varCost(p) * y(p,t))、运输成本sum((p,m,t)$link(p,m), shipCost(p,m) * x(p,m,t))、库存持有成本sum((p,t), holdCost * inv(p,t))、缺货惩罚sum((m,t), shortCost * short(m,t))。缺货用惩罚项而不用硬约束是因为硬约束意味着需求必须满足一旦产能不足需求高峰、检修整个 MIP 直接不可行你得到的是没有方案而不是代价最小的方案改成需求满足方程带一个缺口变量short(m,t)、缺口在目标里被重罚后模型永远可行产能充裕时缺口会被压到零。敏感性重跑用loop范式定义Set s / low, base, high /与scale(s) 0.9 / 1.0 / 1.1先把基线需求备份到demBase(m,t)否则循环里会对已缩放的需求二次缩放这是最经典的静默 bug循环体内重算需求、Solve、把结果收进二维参数kpiS(s,*)最后可用一条业务级回归闸门abort$(kpiS(high,totalCost) kpiS(low,totalCost))把需求升而成本降这类数据/模型损坏挡在交付之前。注意 MIP 在多次solve之间会带着上一次解的level但求解器是否采纳为起点取决于它的选项铁律 6因此重跑结论只写本算例上如此不作性能承诺铁律 7。
