基于Pyomo的Python批量过程优化：如何设置批次恒定运行功率？

我正在用Python的Pyomo工具优化工业批量生产过程，目标是最小化系统成本，生产成本采用随机值定义。该批量过程有两个核心特性：

• 系统启动后必须连续运行4个或4的倍数个时间步长（此功能已实现）
• 批次启动后，该批次的运行功率需保持恒定（这是目前未解决的问题）

尝试过的约束代码 #

我曾尝试添加以下约束，但这些约束只能保证运行功率不会在过程中降低，无法实现功率恒定的要求：

# 尝试确保工厂启动后，当前及后续三个时段功率相等
def block_production_rule_1(model, t):
    if t <= T-3:
        return model.production[t+1]  >= model.production[t] * model.start[t]
    return pyo.Constraint.Skip
model.block_production_rule_1 = pyo.Constraint(model.timesteps, rule=block_production_rule_1)

def block_production_rule_2(model,t):
    if t <= T-3:
        return model.production[t+2] >= model.production[t] * model.start[t]
    return pyo.Constraint.Skip
model.block_production_rule_2 = pyo.Constraint(model.timesteps, rule=block_production_rule_2)

def block_production_rule_3(model,t):
    if t <= T-3:
        return model.production[t+3] >= model.production[t] * model.start[t]
    return pyo.Constraint.Skip
model.block_production_rule_3 = pyo.Constraint(model.timesteps, rule=block_production_rule_3)

更新：问题解决 #

我修改了目标函数，同时调整了代码其他部分，目前功能已正常运行。

修改后的目标函数 #

def obj_rule(model):
    return sum((production_cost[t-1] + production_cost[t] + production_cost[t+1] + production_cost[t+2]) * model.run_power[t] for t in model.timesteps)
model.obj = pyo.Objective(rule=obj_rule, sense=pyo.minimize)

完整代码 #

import pyomo.environ as pyo
import numpy as np
import random 


# 定义模型
model = pyo.ConcreteModel()

# 1. 定义模型集合与参数
T = 56  # 时间步长总数
production_goal = 32
maxPower = 2
minPower = 0
model.timesteps = pyo.RangeSet(1, T)

# 生产成本函数
# production_cost = [1, 2, 0, 3, 4, 1, 0, 2, 4, 1, 3, 2, 0, 4, 2, 1, 0, 4, 3, 2, 1, 0, 4, 3, 2, 1, 0, 4, 3, 2, 1,0,3,3,3]
production_cost = (np.sin(np.linspace(0, 8*np.pi, T)) + 1) * 50
production_cost = [_ for _ in range(T+3)]
production_cost = [random.randint(0,4) for _ in range(T+3)]
production_cost_sum = []
for i in range(len(production_cost)):
    total = 0
    for j in range(i, min(i + 4, len(production_cost))):  # 处理末尾边界问题
        total += production_cost[j]
    production_cost_sum.append(total)

# 2. 定义决策变量
model.start = pyo.Var(model.timesteps, within=pyo.Binary)  # 工厂启动标志
model.run_power = pyo.Var(model.timesteps, within=pyo.Integers)  # 工厂运行功率

# 3. 定义目标函数
def obj_rule(model):
    return sum((production_cost[t-1] + production_cost[t] + production_cost[t+1] + production_cost[t+2]) * model.run_power[t] for t in model.timesteps)
model.obj = pyo.Objective(rule=obj_rule, sense=pyo.minimize)

# 4. 定义约束条件
# 总产量约束
def total_production_rule(model):
    return sum(model.run_power[t] for t in model.timesteps) == production_goal/4
model.total_production_con = pyo.Constraint(rule=total_production_rule)

# 仅当工厂启动时才能生产
def production_on_rule_1(model, t):
    return model.run_power[t] <= model.start[t] * maxPower
model.production_on_con_1 = pyo.Constraint(model.timesteps, rule=production_on_rule_1)

def production_on_rule_2(model, t):
    return model.run_power[t] >= model.start[t] 
model.production_on_con_2 = pyo.Constraint(model.timesteps, rule=production_on_rule_2)


# 启动约束：前3个时段未启动才能启动
def start_rule(model, t):
    return model.start[t] <= 1 - sum(model.start[tt] for tt in range(t-3, t) if tt > 0)
model.start_con = pyo.Constraint(model.timesteps, rule=start_rule)




# 5. 使用Gurobi求解模型
solver = pyo.SolverFactory('gurobi')
results = solver.solve(model)






import matplotlib.pyplot as plt
# 6. 绘制结果
x_values = [model.run_power[t]() for t in model.timesteps]
i = 0
while i < len(x_values):
    if x_values[i] > 0:
        value_to_propagate = x_values[i]
        for j in range(i, min(i + 4, len(x_values))):
            x_values[j] = value_to_propagate
        i += 4
    else:
        i += 1
y_values = [model.run_power[t]() for t in model.timesteps]
batch_over_values = [model.start[t]() for t in model.timesteps]

# 创建3个子图
fig, (ax1, ax2, ax3) = plt.subplots(3, 1, sharex=True, figsize=(15, 10))

# 绘制生产决策
ax1.step(range(1, T+1), x_values, where='post', label="Production (x)", color='blue')
ax1.set_ylabel("数值")
ax1.legend()
ax1.grid(True)

# 绘制工厂状态
ax2.step(range(1, T+1), y_values, where='post', label="工厂运行中", color='red', linestyle='--')
ax2.step(range(1, T+1), batch_over_values, where='post', label="工厂启动", color='blue', linestyle=':')
ax2.set_ylabel("数值")
ax2.legend()
ax2.grid(True)

# 统一x轴范围
ax3.step(range(1, T+4), production_cost, where='post', label="生产成本", color='green', linestyle=':')
ax3.set_ylabel("数值")
ax3.set_xlabel("时间步长")
ax3.legend()
ax3.grid(True)
ax3.set_xlim(1, T+4)  # 设置x轴范围

# 添加每个时间步长的刻度
ax3.set_xticks(range(1, T+5))  # 用5而非4，确保包含最后一个时间步长

# 绘制累计生产成本
ax3.step(range(1, T+4), production_cost_sum, where='post', label="累计生产成本", color='red', linestyle=':')
ax3.set_ylabel("数值")
ax3.set_xlabel("时间步长")
ax3.legend()
ax3.grid(True)

for t in range(1, T+4):
    value = production_cost_sum[t-1]
    ax3.text(t, value, str(value), ha='center', va='bottom')

# 显示图形
plt.show()

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