英文:
how to write a pyomo optimization to select optimal volume of renewables?
问题
Here's the translation of the code part you provided:
import os
import pandas as pd
from pyomo.environ import *
import datetime
def model_to_df(model, first_period, last_period):
# Need to increase the first & last hour by 1 because of pyomo indexing
periods = range(model.T[first_period + 1], model.T[last_period + 1] + 1)
spot = [value(model.spot_v[i]) for i in periods]
load = [value(model.load_v[i]) for i in periods]
slr1 = [value(model.slr1_size[i]) for i in periods]
slr2 = [value(model.slr2_size[i]) for i in periods]
slr3 = [value(model.slr3_size[i]) for i in periods]
wnd1 = [value(model.wnd1_size[i]) for i in periods]
wnd2 = [value(model.wnd2_size[i]) for i in periods]
wnd3 = [value(model.wnd3_size[i]) for i in periods]
slr1_gen = [value(model.slr1_gen[i]) for i in periods]
slr2_gen = [value(model.slr2_gen[i]) for i in periods]
slr3_gen = [value(model.slr3_gen[i]) for i in periods]
wnd1_gen = [value(model.wnd1_gen[i]) for i in periods]
wnd2_gen = [value(model.wnd2_gen[i]) for i in periods]
wnd3_gen = [value(model.wnd3_gen[i]) for i in periods]
total_gen = [value(model.total_gen[i]) for i in periods]
spill_gen = [value(model.spill_gen[i]) for i in periods]
spill_gen_sum = [value(model.spill_gen_sum[i]) for i in periods]
spill_binary = [value(model.spill_binary[i]) for i in periods]
self_cons = [value(model.self_cons[i]) for i in periods]
df_dict = {
'Period': periods,
'spot': spot,
'load': load,
'slr1': slr1,
'slr2': slr2,
'slr3': slr3,
'wnd1': wnd1,
'wnd2': wnd2,
'wnd3': wnd3,
'slr1_gen': slr1_gen,
'slr2_gen': slr2_gen,
'slr3_gen': slr3_gen,
'wnd1_gen': wnd1_gen,
'wnd2_gen': wnd2_gen,
'wnd3_gen': wnd3_gen,
'total_gen': total_gen,
'spill_gen': spill_gen,
'self_cons': self_cons,
'spill_gen_sum': spill_gen_sum,
'spill_binary': spill_binary
}
df = pd.DataFrame(df_dict)
return df
LOCATION = r"C:\cbc-win64"
os.environ["PATH"] = LOCATION + ";" + os.environ["PATH"]
df = pd.DataFrame({
'hour': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
'load': [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100],
'spot': [65.4, 62.7, 60.9, 60.3, 61.8, 64.5, 65.9, 57.9, 39.7, 28.3, 20.9, 16.3, 18.1, 23.9, 32.3, 43.2, 59.3, 76.3, 80.5, 72.5, 73.1, 69.0, 67.9, 67.7],
'slr1': [0.00, 0.00, 0.00, 0.00, 0.00, 0.04, 0.20, 0.44, 0.60, 0.69, 0.71, 0.99, 1.00, 0.66, 0.75, 0.63, 0.52, 0.34, 0.14, 0.02, 0.00, 0.00, 0.00, 0.00],
<details>
<summary>英文:</summary>
**Background**
I am trying to write a pyomo optimization which takes in a customer's electricity load and the generation data of several renewable projects, then optimally solves for the lowest cost selection of renewable projects to minimize electricity consumption, subject to a few constraints.
**What I've tried**
Using the pyomo readthedocs and stackoverflow. I've written my first attempt (below), but I am having two issues.
**Problem**
1. ERROR: Rule failed for Expression 'd_spill_var' with index 0: PyomoException:
Cannot convert non-constant Pyomo expression
I think this is because I am trying to return a max(expr, 0) value for one of my dependent Expresions. However even if I change this, I still get issue 2 below;
2. RuntimeError: Cannot write legal LP file. Objective 'objective' has nonlinear terms that are not quadratic.
**Help Requested**
Could someone please point me in the right direction to solving the above two issues? Any help would be greatly appreciated!
**Code OLD v1.0**
import os
import pandas as pd
from pyomo.environ import *
import datetime
def model_to_df(model, first_period, last_period):
# Need to increase the first & last hour by 1 because of pyomo indexing
periods = range(model.T[first_period + 1], model.T[last_period + 1] + 1)
spot = [value(model.spot[i]) for i in periods]
load = [value(model.load[i]) for i in periods]
slr1 = [value(model.slr1_size[i]) for i in periods]
slr2 = [value(model.slr2_size[i]) for i in periods]
slr3 = [value(model.slr3_size[i]) for i in periods]
wnd1 = [value(model.wnd1_size[i]) for i in periods]
wnd2 = [value(model.wnd2_size[i]) for i in periods]
wnd3 = [value(model.wnd3_size[i]) for i in periods]
d_slrgen_var = [value(model.d_slrgen_var[i]) for i in periods]
d_wndgen_var = [value(model.d_wndgen_var[i]) for i in periods]
d_spill_var = [value(model.d_spill_var[i]) for i in periods]
d_selfcons_var = [value(model.d_selfcons_var[i]) for i in periods]
df_dict = {
'Period': periods,
'spot': spot,
'load': load,
'slr1': slr1,
'slr2': slr2,
'slr3': slr3,
'wnd1': wnd1,
'wnd2': wnd2,
'wnd3': wnd3,
'd_slrgen_var': d_slrgen_var,
'd_wndgen_var': d_wndgen_var,
'd_spill_var': d_spill_var,
'd_selfcons_var': d_selfcons_var
}
df = pd.DataFrame(df_dict)
return df
LOCATION = r"C:\cbc-win64"
os.environ["PATH"] = LOCATION + ";" + os.environ["PATH"]
df = pd.DataFrame({
'hour': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
'load': [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100],
'spot': [65.4, 62.7, 60.9, 60.3, 61.8, 64.5, 65.9, 57.9, 39.7, 28.3, 20.9, 16.3, 18.1, 23.9, 32.3, 43.2, 59.3, 76.3, 80.5, 72.5, 73.1, 69.0, 67.9, 67.7],
'slr1': [0.00, 0.00, 0.00, 0.00, 0.00, 0.04, 0.20, 0.44, 0.60, 0.69, 0.71, 0.99, 1.00, 0.66, 0.75, 0.63, 0.52, 0.34, 0.14, 0.02, 0.00, 0.00, 0.00, 0.00],
'slr2': [0.00, 0.00, 0.00, 0.00, 0.03, 0.19, 0.44, 0.68, 1.00, 0.83, 0.90, 0.88, 0.98, 0.94, 0.83, 0.70, 0.36, 0.11, 0.02, 0.00, 0.00, 0.00, 0.00, 0.00],
'slr3': [0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.03, 0.17, 0.39, 0.87, 0.91, 1.00, 0.89, 0.71, 0.71, 0.85, 0.63, 0.52, 0.32, 0.12, 0.02, 0.00, 0.00, 0.00],
'wnd1': [1.00, 0.72, 0.74, 0.94, 0.69, 0.90, 0.92, 0.76, 0.51, 0.35, 0.31, 0.34, 0.37, 0.28, 0.35, 0.40, 0.39, 0.32, 0.42, 0.48, 0.74, 0.63, 0.80, 0.97],
'wnd2': [0.95, 0.67, 0.82, 0.48, 0.51, 0.41, 0.33, 0.42, 0.34, 0.30, 0.39, 0.29, 0.34, 0.55, 0.67, 0.78, 0.84, 0.73, 0.77, 0.89, 0.76, 0.97, 1.00, 0.91],
'wnd3': [0.32, 0.35, 0.38, 0.28, 0.33, 0.38, 0.41, 0.38, 0.51, 0.65, 0.54, 0.88, 0.93, 0.89, 0.90, 1.00, 0.90, 0.76, 0.76, 0.92, 0.71, 0.56, 0.52, 0.40]
})
first_model_period = df['hour'].iloc[0]
last_model_period = df['hour'].iloc[-1]
# **********************
# Build Model
# **********************
model = ConcreteModel()
# Fixed Paramaters
model.T = Set(initialize=df.index.tolist(), doc='hourly intervals', ordered=True)
model.load_v = Param(model.T, initialize=df.load, doc='customers load', within=Any)
model.spot_v = Param(model.T, initialize=df.spot, doc='spot price for each interval', within=Any)
model.slr1 = Param(model.T, initialize=df.slr1, doc='1MW output solar farm 1', within=Any)
model.slr2 = Param(model.T, initialize=df.slr2, doc='1MW output solar farm 2', within=Any)
model.slr3 = Param(model.T, initialize=df.slr3, doc='1MW output solar farm 3', within=Any)
model.wnd1 = Param(model.T, initialize=df.wnd1, doc='1MW output wind farm 1', within=Any)
model.wnd2 = Param(model.T, initialize=df.wnd2, doc='1MW output wind farm 2', within=Any)
model.wnd3 = Param(model.T, initialize=df.wnd3, doc='1MW output wind farm 3', within=Any)
# Variable Parameters
model.slr1_flag = Var(model.T, doc='slr 1 on / off', within=Binary, initialize=0)
model.slr2_flag = Var(model.T, doc='slr 2 on / off', within=Binary, initialize=0)
model.slr3_flag = Var(model.T, doc='slr 3 on / off', within=Binary, initialize=0)
model.wnd1_flag = Var(model.T, doc='wnd 1 on / off', within=Binary, initialize=0)
model.wnd2_flag = Var(model.T, doc='wnd 2 on / off', within=Binary, initialize=0)
model.wnd3_flag = Var(model.T, doc='wnd 3 on / off', within=Binary, initialize=0)
model.slr1_size = Var(model.T, bounds=(0, 1500), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.slr2_size = Var(model.T, bounds=(0, 1500), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.slr3_size = Var(model.T, bounds=(0, 1500), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.wnd1_size = Var(model.T, bounds=(0, 1500), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.wnd2_size = Var(model.T, bounds=(0, 1500), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.wnd3_size = Var(model.T, bounds=(0, 1500), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.total_gen = Var(model.T, initialize=0, within=NonNegativeReals)
# Dependent Expression Parameters
def dependent_solar_gen(model, t):
"Total selected solar Generation"
return (model.slr1[t] * model.slr1_flag[t] * model.slr1_size[t]) + \
(model.slr2[t] * model.slr2_flag[t] * model.slr2_size[t]) + \
(model.slr3[t] * model.slr3_flag[t] * model.slr3_size[t])
model.d_slrgen_var = Expression(model.T, rule=dependent_solar_gen)
def dependent_wind_gen(model, t):
"Total selected wind Generation"
return (model.wnd1[t] * model.wnd1_flag[t] * model.wnd1_size[t]) + \
(model.wnd2[t] * model.wnd2_flag[t] * model.wnd2_size[t]) + \
(model.wnd3[t] * model.wnd3_flag[t] * model.wnd3_size[t])
model.d_wndgen_var = Expression(model.T, rule=dependent_wind_gen)
def dependent_spill(model, t):
"Volume of energy not consumed by customer (spilled into grid)"
expr = (model.d_slrgen_var[t] + model.d_wndgen_var[t]) - model.load_v[t]
return max(0, expr)
model.d_spill_var = Expression(model.T, rule=dependent_spill)
def dependent_self_cons(model, t):
"Volume of energy consumed by customer"
expr = (model.d_slrgen_var[t] + model.d_wndgen_var[t]) - model.d_spill_var[t]
return expr
model.d_selfcons_var = Expression(model.T, rule=dependent_self_cons)
# -----------------------
# Constraints
# -----------------------
def min_spill(model, t):
"Limit spill renewables to 10% of total"
return model.d_spill_var[t] <= 0.1 * (model.d_slrgen_var[t] + model.d_wndgen_var[t])
model.min_spill_c = Constraint(model.T, rule=min_spill)
def load_match(model, t):
"contract enough renewables to offset 100% load, even if its not time matched"
return (model.d_slrgen_var[t] + model.d_wndgen_var[t]) >= model.load_v[t]
model.load_match_c = Constraint(model.T, rule=load_match)
# **********************
# Define the income, expenses, and profit
# **********************
green_income = sum(model.spot_v[t] * model.d_spill_var[t] for t in model.T)
black_cost = sum(model.spot_v[t] * (model.load_v[t] - model.d_selfcons_var[t]) for t in model.T)
slr_cost = sum(40 * model.d_slrgen_var[t] for t in model.T)
wnd_cost = sum(70 * model.d_wndgen_var[t] for t in model.T)
profit = green_income - black_cost - slr_cost - wnd_cost
model.objective = Objective(expr=profit, sense=maximize)
# Solve the model
# solver = SolverFactory('glpk')
solver = SolverFactory('cbc')
solver.solve(model, timelimit=10)
results_df = model_to_df(model, first_period=first_model_period, last_period=last_model_period)
print(results_df)
**Solved**
Thanks @airliquid, I managed to solve the issue thanks to your advice.
What I had to do was linearize the max constraint, redefine dependent expressions as constraints, remove some redundant variables, and change the last two constraints to summations
It's not the prettiest answer, but it works!
**UPDATED CODE (V2)**
import os
import pandas as pd
from pyomo.environ import *
import datetime
def model_to_df(model, first_period, last_period):
# Need to increase the first & last hour by 1 because of pyomo indexing
periods = range(model.T[first_period + 1], model.T[last_period + 1] + 1)
spot = [value(model.spot_v[i]) for i in periods]
load = [value(model.load_v[i]) for i in periods]
slr1 = [value(model.slr1_size[i]) for i in periods]
slr2 = [value(model.slr2_size[i]) for i in periods]
slr3 = [value(model.slr3_size[i]) for i in periods]
wnd1 = [value(model.wnd1_size[i]) for i in periods]
wnd2 = [value(model.wnd2_size[i]) for i in periods]
wnd3 = [value(model.wnd3_size[i]) for i in periods]
slr1_gen = [value(model.slr1_gen[i]) for i in periods]
slr2_gen = [value(model.slr2_gen[i]) for i in periods]
slr3_gen = [value(model.slr3_gen[i]) for i in periods]
wnd1_gen = [value(model.wnd1_gen[i]) for i in periods]
wnd2_gen = [value(model.wnd2_gen[i]) for i in periods]
wnd3_gen = [value(model.wnd3_gen[i]) for i in periods]
total_gen = [value(model.total_gen[i]) for i in periods]
spill_gen = [value(model.spill_gen[i]) for i in periods]
spill_gen_sum = [value(model.spill_gen_sum[i]) for i in periods]
spill_binary = [value(model.spill_binary[i]) for i in periods]
self_cons = [value(model.self_cons[i]) for i in periods]
df_dict = {
'Period': periods,
'spot': spot,
'load': load,
'slr1': slr1,
'slr2': slr2,
'slr3': slr3,
'wnd1': wnd1,
'wnd2': wnd2,
'wnd3': wnd3,
'slr1_gen': slr1_gen,
'slr2_gen': slr2_gen,
'slr3_gen': slr3_gen,
'wnd1_gen': wnd1_gen,
'wnd2_gen': wnd2_gen,
'wnd3_gen': wnd3_gen,
'total_gen': total_gen,
'spill_gen': spill_gen,
'self_cons': self_cons,
'spill_gen_sum': spill_gen_sum,
'spill_binary': spill_binary
}
df = pd.DataFrame(df_dict)
return df
LOCATION = r"C:\cbc-win64"
os.environ["PATH"] = LOCATION + ";" + os.environ["PATH"]
df = pd.DataFrame({
'hour': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
'load': [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100],
'spot': [65.4, 62.7, 60.9, 60.3, 61.8, 64.5, 65.9, 57.9, 39.7, 28.3, 20.9, 16.3, 18.1, 23.9, 32.3, 43.2, 59.3, 76.3, 80.5, 72.5, 73.1, 69.0, 67.9, 67.7],
'slr1': [0.00, 0.00, 0.00, 0.00, 0.00, 0.04, 0.20, 0.44, 0.60, 0.69, 0.71, 0.99, 1.00, 0.66, 0.75, 0.63, 0.52, 0.34, 0.14, 0.02, 0.00, 0.00, 0.00, 0.00],
'slr2': [0.00, 0.00, 0.00, 0.00, 0.03, 0.19, 0.44, 0.68, 1.00, 0.83, 0.90, 0.88, 0.98, 0.94, 0.83, 0.70, 0.36, 0.11, 0.02, 0.00, 0.00, 0.00, 0.00, 0.00],
'slr3': [0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.03, 0.17, 0.39, 0.87, 0.91, 1.00, 0.89, 0.71, 0.71, 0.85, 0.63, 0.52, 0.32, 0.12, 0.02, 0.00, 0.00, 0.00],
'wnd1': [1.00, 0.72, 0.74, 0.94, 0.69, 0.90, 0.92, 0.76, 0.51, 0.35, 0.31, 0.34, 0.37, 0.28, 0.35, 0.40, 0.39, 0.32, 0.42, 0.48, 0.74, 0.63, 0.80, 0.97],
'wnd2': [0.95, 0.67, 0.82, 0.48, 0.51, 0.41, 0.33, 0.42, 0.34, 0.30, 0.39, 0.29, 0.34, 0.55, 0.67, 0.78, 0.84, 0.73, 0.77, 0.89, 0.76, 0.97, 1.00, 0.91],
'wnd3': [0.32, 0.35, 0.38, 0.28, 0.33, 0.38, 0.41, 0.38, 0.51, 0.65, 0.54, 0.88, 0.93, 0.89, 0.90, 1.00, 0.90, 0.76, 0.76, 0.92, 0.71, 0.56, 0.52, 0.40]
})
# df.to_csv('example.csv', index=False)
first_model_period = df['hour'].iloc[0]
last_model_period = df['hour'].iloc[-1]
# **********************
# Build Model
# **********************
model = ConcreteModel()
# Fixed Paramaters
model.T = Set(initialize=df.index.tolist(), doc='hourly intervals', ordered=True)
model.load_v = Param(model.T, initialize=df.load, doc='customers load', within=Any)
model.spot_v = Param(model.T, initialize=df.spot, doc='spot price for each interval', within=Any)
model.slr1 = Param(model.T, initialize=df.slr1, doc='1MW output solar farm 1', within=Any)
model.slr2 = Param(model.T, initialize=df.slr2, doc='1MW output solar farm 2', within=Any)
model.slr3 = Param(model.T, initialize=df.slr3, doc='1MW output solar farm 3', within=Any)
model.wnd1 = Param(model.T, initialize=df.wnd1, doc='1MW output wind farm 1', within=Any)
model.wnd2 = Param(model.T, initialize=df.wnd2, doc='1MW output wind farm 2', within=Any)
model.wnd3 = Param(model.T, initialize=df.wnd3, doc='1MW output wind farm 3', within=Any)
# Variable Parameters
model.slr1_size = Var(model.T, bounds=(0, 1000), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.slr2_size = Var(model.T, bounds=(0, 1000), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.slr3_size = Var(model.T, bounds=(0, 1000), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.wnd1_size = Var(model.T, bounds=(0, 1000), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.wnd2_size = Var(model.T, bounds=(0, 1000), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.wnd3_size = Var(model.T, bounds=(0, 1000), doc='selected size in MWs', initialize=0, within=NonNegativeIntegers)
model.slr1_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.slr2_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.slr3_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.wnd1_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.wnd2_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.wnd3_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.total_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.spill_gen_sum = Var(model.T, initialize=0, within=Reals)
model.spill_binary = Var(model.T, doc='to get max', within=Binary, initialize=0)
model.spill_gen = Var(model.T, initialize=0, within=NonNegativeReals)
model.self_cons = Var(model.T, initialize=0, within=NonNegativeReals)
# -----------------------
# Constraints
# -----------------------
# SIZE CONSTRAINTS
def slr1_size(model, t):
"slr1 size"
return model.slr1_size[t] == model.slr1_size[1]
model.slr_size1_c = Constraint(model.T, rule=slr1_size)
def slr2_size(model, t):
"slr2 size"
return model.slr2_size[t] == model.slr2_size[1]
model.slr_size2_c = Constraint(model.T, rule=slr2_size)
def slr3_size(model, t):
"slr3 size"
return model.slr3_size[t] == model.slr3_size[1]
model.slr_size3_c = Constraint(model.T, rule=slr3_size)
def wnd1_size(model, t):
"wnd1 size"
return model.wnd1_size[t] == model.wnd1_size[1]
model.wnd_size1_c = Constraint(model.T, rule=wnd1_size)
def wnd2_size(model, t):
"wnd2 size"
return model.wnd2_size[t] == model.wnd2_size[1]
model.wnd_size2_c = Constraint(model.T, rule=wnd2_size)
def wnd3_size(model, t):
"wnd3 size"
return model.wnd3_size[t] == model.wnd3_size[1]
model.wnd_size3_c = Constraint(model.T, rule=wnd3_size)
# GENERATION EXPRESSIONS / CONSTRAINTS
def slr1_gen(model, t):
"solar 1 generation"
return model.slr1_gen[t] == model.slr1[t] * model.slr1_size[t]
model.slr_gen1_c = Constraint(model.T, rule=slr1_gen)
def slr2_gen(model, t):
"solar 2 generation"
return model.slr2_gen[t] == model.slr2[t] * model.slr2_size[t]
model.slr_gen2_c = Constraint(model.T, rule=slr2_gen)
def slr3_gen(model, t):
"solar 3 generation"
return model.slr3_gen[t] == model.slr3[t] * model.slr3_size[t]
model.slr_gen3_c = Constraint(model.T, rule=slr3_gen)
def wnd1_gen(model, t):
"wind 1 generation"
return model.wnd1_gen[t] == model.wnd1[t] * model.wnd1_size[t]
model.wnd_gen1_c = Constraint(model.T, rule=wnd1_gen)
def wnd2_gen(model, t):
"wind 2 generation"
return model.wnd2_gen[t] == model.wnd2[t] * model.wnd2_size[t]
model.wnd_gen2_c = Constraint(model.T, rule=wnd2_gen)
def wnd3_gen(model, t):
"wind 3 generation"
return model.wnd3_gen[t] == model.wnd3[t] * model.wnd3_size[t]
model.wnd_gen3_c = Constraint(model.T, rule=wnd3_gen)
# TOTAL GENERATION
def total_gen(model, t):
"sum of generation"
return model.total_gen[t] == model.slr1_gen[t] + model.slr2_gen[t] + model.slr3_gen[t] + \
model.wnd1_gen[t] + model.wnd2_gen[t] + model.wnd3_gen[t]
model.total_gen_c = Constraint(model.T, rule=total_gen)
# SPILL GENERATION
def spill_gen_sum(model, t):
"X >= x1"
return model.spill_gen_sum[t] == model.total_gen[t] - model.load_v[t]
model.spill_gen_sum_c = Constraint(model.T, rule=spill_gen_sum)
def spill_check_one(model, t):
"X >= x1"
return model.spill_gen[t] >= model.spill_gen_sum[t]
model.spill_check_one_c = Constraint(model.T, rule=spill_check_one)
def spill_check_two(model, t):
"X >= x2"
return model.spill_gen[t] >= 0
model.spill_check_two_c = Constraint(model.T, rule=spill_check_two)
def spill_binary_one(model, t):
"X <= x1 + M(1-y)"
return model.spill_gen[t] <= model.spill_gen_sum[t] + 9999999*(1-model.spill_binary[t])
model.spill_binary_one_c = Constraint(model.T, rule=spill_binary_one)
def spill_binary_two(model, t):
"X <= x2 + My"
return model.spill_gen[t] <= 9999999*model.spill_binary[t]
model.spill_binary_two_c = Constraint(model.T, rule=spill_binary_two)
# SELF CONS
def self_cons(model, t):
"X <= x2 + My"
return model.self_cons[t] == model.total_gen[t] - model.spill_gen[t]
model.self_cons_c = Constraint(model.T, rule=self_cons)
# ACTUAL CONSTRAINTS
def min_spill(model, t):
"Limit spill renewables to 10% of total"
return sum(model.spill_gen[t] for t in model.T) <= 0.2 * sum(model.total_gen[t] for t in model.T)
model.min_spill_c = Constraint(model.T, rule=min_spill)
def load_match(model, t):
"contract enough renewables to offset 100% load, even if its not time matched"
return sum(model.total_gen[t] for t in model.T) >= sum(model.load_v[t] for t in model.T)
model.load_match_c = Constraint(model.T, rule=load_match)
# **********************
# Define the battery income, expenses, and profit
# **********************
green_income = sum(model.spot_v[t] * model.spill_gen[t] for t in model.T)
black_cost = sum(model.spot_v[t] * (model.load_v[t] - model.self_cons[t]) for t in model.T)
slr_cost = sum(40 * (model.slr1_gen[t] + model.slr2_gen[t] + model.slr3_gen[t]) for t in model.T)
wnd_cost = sum(70 * (model.wnd1_gen[t] + model.wnd2_gen[t] + model.wnd3_gen[t]) for t in model.T)
cost = black_cost + slr_cost + wnd_cost - green_income
model.objective = Objective(expr=cost, sense=minimize)
# Solve the model
# solver = SolverFactory('glpk')
solver = SolverFactory('cbc')
solver.solve(model)
# , timelimit=10
results_df = model_to_df(model, first_period=first_model_period, last_period=last_model_period)
print(results_df)
results_df.to_csv('temp.csv', index=False)
</details>
# 答案1
**得分**: 2
以下是翻译的内容：
1. 您的模型写得很好（清晰而有条理），但存在一些重要的线性规划问题。以下是一些建议，以帮助您解决这些问题....
您通过将变量相乘来使此模型非线性化，这会给求解器带来巨大的复杂性，在这种情况下是不必要的。您可以使用一些约束，包括一个"big-M"约束，以消除您在将标志乘以容量（两个变量...因此非线性性）部分的非线性性。您实际上正在执行以下操作：
output[t] = energy[t] * size[t] * running[t]
其中energy是一个参数，其他的是变量，而`output`则被合并到一个表达式中。您可以通过将output作为一个新的变量，并使用3个约束来线性化它：
output[t] >= energy[t] * size[t] - (1 - running[t]) * M
output[t] <= energy[t] * size[t]
output[t] <= running[t] * M
这里的M是能量*大小的逻辑上限。有时候上限和下限都是不需要的，但是由于超额约束，您可能需要这个来强制执行您想要的"相等"约束。我建议将`output`设为一个新变量，以减轻您的复杂表达式。
2. 在尝试创建"合并表达式"时，您实际上正在创建一组"索引表达式"。只需在模型的T上使用求和语句：
model.d_slrgen_var = Expression(model.T, rule=dependent_solar_gen)
创建了一组表达式，而不是一个求和。
3. `max()`是一个非线性运算符，您不能在LP表达式中使用它或返回它。因此，引入另一个类似于以下的变量：
spill[t] ∈ {非负实数 }
spill[t] >= ... 对于模型T中的每个t
当您解决了所有这些问题并使其正常工作后，您可以考虑重新组织模型，并通过`[source, time]`双重索引主要变量，而不是按源命名它们，这样会更清晰一些，但可以留到以后再考虑。如果您遇到困难，请回复评论...
<details>
<summary>英文:</summary>
Your model is well written (clear and organized), but shows some significant Linear Programming issues. A couple pointers to help you along the way....
1. You are making this model non-linear by multiplying variables together, which adds HUGE complication for the solver and is not necessary in this case. You can use a couple of constraints, including a "big-M" constraint to remove the portion where you are multiplying your flag times capacity (both variables...hence the nonlinearity). You are essentially doing this:
###
output[t] = energy[t] * size[t] * running[t]
where energy is a parameter, the others are variables, and `output` is rolled up into an expression. You can linearize that by just making output a variable and using 3 constraints:
output[t] >= energy[t] * size[t] - (1 - running[t]) * M
output[t] <= energy[t] * size[t]
output[t] <= running[t] * M
Where M is some logical upper bound on energy*size. Sometimes the upper + the lower are not needed, but you probably need this to enforce the "equality" constraint you want because of the overage constraint.
I would make `output` a new variable to alleviate the complex expressions you have.
2. In your attempt to make the "rolled up expressions" you are instead making a set of "indexed expressions". Just use a summation statement over model.T. This:
###
model.d_slrgen_var = Expression(model.T, rule=dependent_solar_gen)
makes a **set** of expressions, not a summation
3. `max()` is a nonlinear operator, and you cannot return or use that in a LP expression. So.... Introduce another variable similar to:
###
spill[t] ∈ {non-negative REALS }
spill[t] >= ... for each t in model.T
After you manage to get all that ironed out ( ;) ) and working, you might want to re-organize the model and double-index your main variables by `[source, time]` instead of naming them all by source, which is a bit cleaner, but window dressing for later.
Comment back if you are stuck...
</details>
通过集体智慧和协作来改善编程学习和解决问题的方式。致力于成为全球开发者共同参与的知识库，让每个人都能够通过互相帮助和分享经验来进步。
评论