在 pyomo Python 中优化批处理 - 如何为每个批设置恒定的操作功率？

Question

我想在 Python 中使用 pyomo 优化工业批处理过程。目标函数是最小化系统成本。在我的简化示例中，我使用随机值指定生产成本。批处理的独特之处在于：

1. 系统一旦打开，必须运行四个或四个时间步的倍数（连续多个批次）。 --> 这已经实施了。
2. 一旦启动一个批次，该批次的运行功率应保持恒定。 --> 未解决

将间歇过程添加到具有恒定排放的存储中

我已经尝试过这个约束：

# def production_on_rule_3(model, t):
    if t>3:
        return model.production[t] >= model.production[t-1] * sum(model.start[tt] for tt in range(t-3, t) if tt > 0)
    else:
        return pyo.Constraint.Skip
model.production_on_con_3 = pyo.Constraint(model.timesteps, rule=production_on_rule_3)

但这样，我只保证，在这个过程中，运算能力不会变小。

有人知道如何解决这个问题吗？

代码：

import pyomo.environ as pyo
import numpy as np
import random 
import csv
import matplotlib.pyplot as plt

# create the model
model = pyo.ConcreteModel()

# Parameter EAF
# 1. Define the model, sets, and parameters.
T = 64  # number of timesteps
production_goal = 80
maxPower = 4
minPower = 0

# Parameter storage
max_capacity = 7.5  # Maximale Speicherkapazität des Batteriespeichers (z.B., in kWh)

# timesteps
model.timesteps = pyo.RangeSet(1, T) 

# Production cost function
production_cost = (np.sin(np.linspace(0, 8*np.pi, T)) + 1) * 50
production_cost = [_ for _ in range(T)]
production_cost = [random.randint(1,100) for _ in range(T)]



# 2. Define the decision variables.
# variables storage
model.charge = pyo.Var(model.timesteps, domain=pyo.NonNegativeReals, initialize=0.0)  # charge
model.discharge = pyo.Var(model.timesteps, domain=pyo.NonNegativeReals, initialize=0.0)  # discharge
model.energy = pyo.Var(model.timesteps, domain=pyo.NonNegativeReals, initialize=0.0)  # SOC
model.storage_status = pyo.Var(model.timesteps, domain=pyo.Binary, initialize=0.0)  # SOC > 0?

# variables EAF
model.production = pyo.Var(model.timesteps, within=pyo.Integers)  # production decision
model.running = pyo.Var(model.timesteps, within=pyo.Binary)  # plant status
model.start = pyo.Var(model.timesteps, within=pyo.Binary)  # start of the plant



# 3. Define the objective function.
def obj_rule(model):
    return sum(production_cost[t-1] * model.production[t] for t in model.timesteps)
model.obj = pyo.Objective(rule=obj_rule, sense=pyo.minimize)



# 4. Define the constraints.

# Constraints EAF
# Total production constraint
def total_production_rule(model):
    return sum(model.production[t] for t in model.timesteps) == production_goal
model.total_production_con = pyo.Constraint(rule=total_production_rule)

# Production is only possible if plant is on
def production_on_rule_1(model, t):
    return model.production[t] <= model.running[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.production[t] >= model.running[t] 
model.production_on_con_2 = pyo.Constraint(model.timesteps, rule=production_on_rule_2)

def production_on_rule_3(model, t):
    if t>3:
        return model.production[t] >= model.production[t-1] * sum(model.start[tt] for tt in range(t-3, t) if tt > 0)
    else:
        return pyo.Constraint.Skip
model.production_on_con_3 = pyo.Constraint(model.timesteps, rule=production_on_rule_3)

# Running if started in current or preceding 3 periods
def running_rule(model, t):
    return model.running[t] == sum(model.start[tt] for tt in range(t-3, t+1) if tt > 0)
model.running_con = pyo.Constraint(model.timesteps, rule=running_rule)

# Start if no start occurred in preceding 3 periods
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)

# Shut-down window constraint
M_starts = 5 # max possible starts in 20 periods
def shutdown_window_rule(model, t):
    if t > 5:  # Limiting timesteps to avoid out-of-bounds errors
        return sum(model.start[tt] for tt in range(t+3, t+21) if tt <= T) <= (1 - (model.start[t-4] - model.start[t])) * M_starts
    # else:
    #     return model.start[t] == 0
    return pyo.Constraint.Skip
model.shutdown_window_con = pyo.Constraint(model.timesteps, rule=shutdown_window_rule)


# Constraints storage
def capacity_rule1(model, t):
    return model.energy[t] <= max_capacity
model.capacity_con1 = pyo.Constraint(model.timesteps, rule=capacity_rule1)

def capacity_rule2(model, t):
    return model.energy[t] >= 0
model.capacity_con2 = pyo.Constraint(model.timesteps, rule=capacity_rule2)

def storage_status_rule1(model, t):
    if t > 1:
        return model.storage_status[t] <= model.energy[t-1]
    else:
        return pyo.Constraint.Skip
model.storage_status_con1 = pyo.Constraint(model.timesteps, rule=storage_status_rule1)

def storage_status_rule2(model, t):
    if t > 1:
        return model.storage_status[t] * max_capacity >= model.energy[t-1]
    else:
        return pyo.Constraint.Skip
model.storage_status_con2 = pyo.Constraint(model.timesteps, rule=storage_status_rule2)

def energy_balance_rule(model, t):
    if t == 1:
        return model.energy[t] == model.charge[t] - model.discharge[t]
    else:
        return model.energy[t] == model.energy[t - 1] + model.charge[t] - model.discharge[t]
model.energy_balance_con = pyo.Constraint(model.timesteps, rule=energy_balance_rule)

def charging_rule(model,t):
    return model.charge[t] == model.production[t]
model.charging_con = pyo.Constraint(model.timesteps, rule=charging_rule)

def discharging_rule(model,t):
    return model.discharge[t] == maxPower * (1 - 0.5/8) * model.storage_status[t]
model.discharging_con = pyo.Constraint(model.timesteps, rule=discharging_rule)

Answer 1

我假设“功率”与生产水平同义，因此如果在时间 5 开始运行，功率为 3，它将在第 5、6、7 和 8 期生产 3 个小部件，总共 12 个小部件您的生产目标...

所以这里有一个想法：重新表述你的问题。它可能变得更容易，因为生产是持续运行的，对吧？因此，如上所述，如果我告诉您它从第 5 期开始，幂为 3，您就知道结果 ==> 您可能不需要运行/生产变量。

因此，尝试将您的

start

变量重新表示为二进制，以指示开始发生（以帮助覆盖非重叠开始的约束），并将另一个

run_power

重新表示为整数，即任何时间段开始的功率级别。看看您是否可以仅根据这两个变量来考虑约束和成本函数，这两个变量仍然按时间索引。

继续上面的例子：

• 如果我告诉你一次运行从 t=5 的 3 次方开始，你能使用

  run_power

  变量计算出该特定运行的成本吗？ （是的。） 类似：

---

 cost[5] = run_power[5] * sum(price[5-8])

• 您能否将所有可能具有

  run_power

  的期间的数据进行汇总，以通过一点创造力获得总成本？ （是的。）

• 您可以重新使用约束来防止多次重叠运行和二进制

  start

  变量吗？ （是）

• 您可以使用大 M 约束在同一时间步将

  start

  变量链接到

  run_power

  变量吗？

• 您能否（类似于成本）汇总生产的小部件总数以进行约束，使其大于或等于需求？ ...

Answer 2

也许我自己解决了。我再次删除了运行变量并实现了 run_power 变量。 Run_power 描述了开始时的功率（整个批次也是如此）。我调整了目标函数：

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)

此外，我将生产变量更改为：

def production_rule(model, t):
return model.production[t] == sum(model.start[tt] for tt in range(t-3, t+1) if tt > 0) * sum(model.run_power[tt] for tt in range(t-3, t+1) if tt > 0)
model.running_con = pyo.Constraint(model.timesteps, rule=production_rule)

通过这种方法，我可以拥有连续的充电变量（生产），同时使用 run_power 来确保块生产。

看起来很简单，所以我不确定我是否错过了什么。