Pyomo整数规划模型结果异常：员工排班求解值与预期不符

任务背景 #

Larry Edison是Buckly学院计算机中心主任，需要为中心排班。中心运营时间为早8点至午夜，各时段所需最低顾问数量如下：

可雇佣全职与兼职顾问：

• 全职顾问需连续工作8小时，可选早班（8AM-下午4点）、中班（正午-晚8点）、晚班（下午4点-午夜），时薪40欧元；
• 兼职顾问按上述4个4小时班次雇佣，时薪30欧元。

额外要求：每个时段内全职顾问数量至少为兼职的2倍。

要求 #

需构建带索引变量的线性规划模型，用Pyomo求解，且变量需为整数（无法雇佣部分员工）。需编写无参函数model_larry_edison，返回符合以下要求的Pyomo模型：

• 目标函数wage：最小化总薪资
• 变量：full_time（3个索引对应3个全职班次）、part_time（4个索引对应4个兼职班次）
• 8个约束：4个总人数约束（A_total、B_total、C_total、D_total）满足各时段最低需求；4个比例约束（A_ratio、B_ratio、C_ratio、D_ratio）满足全职与兼职人数比例要求。

问题 #

运行提供的代码后，得到的总薪资为4106.6666666，但预期结果为4160，需要排查模型问题并修正。

原代码 #

import pandas as pd
import pyomo.environ as pyo

def model_larry_edison():

  salary_df_ft = pd.DataFrame()
  salary_df_ft.at["FT1", "Salary"] = 40
  salary_df_ft.at["FT2", "Salary"] = 40
  salary_df_ft.at["FT3", "Salary"] = 40

  salary_df_pt = pd.DataFrame()
  salary_df_pt.at["PT1", "Salary"] = 30
  salary_df_pt.at["PT2", "Salary"] = 30
  salary_df_pt.at["PT3", "Salary"] = 30
  salary_df_pt.at["PT4", "Salary"] = 30

  salary_ft = salary_df_ft.astype(int)
  salary_pt = salary_df_pt.astype(int)

  model = pyo.ConcreteModel('ProductionPlan')

  model.salary_ft = pyo.Set(initialize = salary_ft.index)
  model.salary_pt = pyo.Set(initialize = salary_pt.index)

  model.full_time = pyo.Var(model.salary_ft, bounds=(0,None))
  model.part_time = pyo.Var(model.salary_pt, bounds=(0,None))

  model.wage = pyo.Objective(expr = 320*(model.full_time["FT1"]+model.full_time["FT2"]+model.full_time["FT3"]) + 120*(model.part_time["PT1"]+model.part_time["PT2"]+model.part_time["PT3"]+model.part_time["PT4"]), sense=pyo.minimize)

  model.A_total = pyo.Constraint(expr = (model.full_time["FT1"] + model.part_time["PT1"]) >= 4)
  model.B_total = pyo.Constraint(expr = (model.full_time["FT1"] + model.full_time["FT2"] + model.part_time["PT2"]) >= 8)
  model.C_total = pyo.Constraint(expr = (model.full_time["FT2"] + model.full_time["FT3"] + model.part_time["PT3"]) >= 10)
  model.D_total = pyo.Constraint(expr = (model.full_time["FT3"] + model.part_time["PT4"]) >= 6)
  model.A_ratio = pyo.Constraint(expr = model.full_time["FT1"] >= 2*model.part_time["PT1"])
  model.B_ratio = pyo.Constraint(expr = (model.full_time["FT1"] + model.full_time["FT2"]) >= 2*model.part_time["PT2"])
  model.C_ratio = pyo.Constraint(expr = (model.full_time["FT2"] + model.full_time["FT3"]) >= 2*model.part_time["PT3"])
  model.D_ratio = pyo.Constraint(expr = model.full_time["FT3"] >= 2*model.part_time["PT4"])

  return model
model = model_larry_edison()
results = pyo.SolverFactory('cbc').solve(model)
print(model.wage())

问题排查与修正 #

问题核心：变量未设置为整数类型。要求明确说明员工数量必须是整数，但原代码中full_time和part_time变量默认是连续型变量，导致求解器返回小数解，从而得到非整数的薪资结果。

修正方法：在定义变量时添加domain=pyo.Integers参数，强制变量为整数。

修正后的代码 #

import pandas as pd
import pyomo.environ as pyo

def model_larry_edison():

  salary_df_ft = pd.DataFrame()
  salary_df_ft.at["FT1", "Salary"] = 40
  salary_df_ft.at["FT2", "Salary"] = 40
  salary_df_ft.at["FT3", "Salary"] = 40

  salary_df_pt = pd.DataFrame()
  salary_df_pt.at["PT1", "Salary"] = 30
  salary_df_pt.at["PT2", "Salary"] = 30
  salary_df_pt.at["PT3", "Salary"] = 30
  salary_df_pt.at["PT4", "Salary"] = 30

  salary_ft = salary_df_ft.astype(int)
  salary_pt = salary_df_pt.astype(int)

  model = pyo.ConcreteModel('ProductionPlan')

  model.salary_ft = pyo.Set(initialize = salary_ft.index)
  model.salary_pt = pyo.Set(initialize = salary_pt.index)

  # 添加domain=pyo.Integers，强制变量为整数
  model.full_time = pyo.Var(model.salary_ft, bounds=(0,None), domain=pyo.Integers)
  model.part_time = pyo.Var(model.salary_pt, bounds=(0,None), domain=pyo.Integers)

  model.wage = pyo.Objective(expr = 320*(model.full_time["FT1"]+model.full_time["FT2"]+model.full_time["FT3"]) + 120*(model.part_time["PT1"]+model.part_time["PT2"]+model.part_time["PT3"]+model.part_time["PT4"]), sense=pyo.minimize)

  model.A_total = pyo.Constraint(expr = (model.full_time["FT1"] + model.part_time["PT1"]) >= 4)
  model.B_total = pyo.Constraint(expr = (model.full_time["FT1"] + model.full_time["FT2"] + model.part_time["PT2"]) >= 8)
  model.C_total = pyo.Constraint(expr = (model.full_time["FT2"] + model.full_time["FT3"] + model.part_time["PT3"]) >= 10)
  model.D_total = pyo.Constraint(expr = (model.full_time["FT3"] + model.part_time["PT4"]) >= 6)
  model.A_ratio = pyo.Constraint(expr = model.full_time["FT1"] >= 2*model.part_time["PT1"])
  model.B_ratio = pyo.Constraint(expr = (model.full_time["FT1"] + model.full_time["FT2"]) >= 2*model.part_time["PT2"])
  model.C_ratio = pyo.Constraint(expr = (model.full_time["FT2"] + model.full_time["FT3"]) >= 2*model.part_time["PT3"])
  model.D_ratio = pyo.Constraint(expr = model.full_time["FT3"] >= 2*model.part_time["PT4"])

  return model
model = model_larry_edison()
results = pyo.SolverFactory('cbc').solve(model)
print(model.wage())

运行修正后的代码，将得到预期的总薪资4160，符合整数规划的要求。

内容的提问来源于stack exchange，提问作者swnr31