MixedInteger Quadratic Programming (MIQP) model does not converge
AnsweredHi,
I created a mixed integer quadratic optimization model in pyomo and want to solve it using gurobi.
My model uses a hourly timeseries of a single year, which describes the costs associated with the operation of a electricity grid.
My model now shall find a hourly timeseries of this year, which describes a grid tariff, which holds following three constraints:
1. only three tariff levels are allowed over the year
2. tariff time structure may differ between weekdays and weekend days
3. tariff time structure may differ between months but not within a single month.
I want to find the optimal tariff timeseries by minimizing the mean squared error of the hourly timeseries for costs and tariff.
Follwoing my model, the variable data holds a column called 'Netzentgelt' with the grid costs and has a datetimeindex in hourly resolution for a full year.
from pyomo.environ import ConcreteModel, Var, Objective, Constraint, SolverFactory, minimize, RangeSet, NonNegativeReals, Binary, quicksum
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
# Read data and preprocess
data = pd.read_csv('data/grid_cost_timeseries_h.csv', sep=';')
# Example costs time series
costs_time_series = data['Netzentgelt'].values
months = data.index.month
weekend_indicator = data.index.weekday >= 5 # 0=Monday, 6=Sunday
# Initialize model
model = ConcreteModel()
N = len(costs_time_series)
# Decision variables for tariffs
model.tariffs = Var(RangeSet(N), within=NonNegativeReals)
# Three potential group values
model.group_values = Var(RangeSet(3), within=NonNegativeReals)
model.group_membership = Var(RangeSet(N), RangeSet(3), within=Binary)
# Monthly and weekday/weekend groupings
model.monthly_group = Var(months.unique(), RangeSet(3), within=NonNegativeReals)
model.weekend_group = Var([0, 1], RangeSet(3), within=NonNegativeReals)
# Objective: Minimize squared error
def objective_rule(model):
return quicksum((model.tariffs[i]  costs_time_series[i1])**2 for i in model.tariffs)
model.objective = Objective(rule=objective_rule, sense=minimize)
# Constraint: Each tariff must be assigned to exactly one group
def one_group_constraint(model, i):
return sum(model.group_membership[i, j] for j in RangeSet(3)) == 1
model.one_group = Constraint(RangeSet(N), rule=one_group_constraint)
# Constraint: Link tariff values to group values
def link_tariff_group_constraint(model, i):
return model.tariffs[i] == sum(model.group_values[j] * model.group_membership[i, j] for j in RangeSet(3))
model.link_tariff_group = Constraint(RangeSet(N), rule=link_tariff_group_constraint)
# Monthly consistency constraint
def monthly_consistency_constraint(model, month, j):
return model.monthly_group[month, j] == model.group_values[j]
model.monthly_consistency = Constraint(months.unique(), RangeSet(3), rule=monthly_consistency_constraint)
# Weekday/Weekend constraint
def weekday_weekend_constraint(model, is_weekend, j):
return model.weekend_group[is_weekend, j] == model.group_values[j]
model.weekday_weekend = Constraint([0, 1], RangeSet(3), rule=weekday_weekend_constraint)
# Solve the model
solver = SolverFactory('gurobi') # Using Gurobi solver
result = solver.solve(model, tee=True)
# Extract solution
optimal_tariff_time_series = np.array([model.tariffs[i].value for i in RangeSet(N)])
group_values = np.array([model.group_values[j].value for j in RangeSet(3)])
group_membership = np.array([[model.group_membership[i,j].value for j in RangeSet(3)] for i in RangeSet(N)])
Does anyone of you has proposals how to simplify the model to improve runtime as the current model does not converge.
Are there also recommendations to finetune gurobi to improve solving time?
Following a screenshot from the gurobi log.

Hi,
Looking at the code, I'm wondering if (apart from performance) it does what you want it to do.
 The part about weekday/weekend seems disconnected from the rest of the model. The weekend_indicator is not being used elsewhere.
 Variables monthly_group and weekend_group seem overkill  irrespective of the month resp. [0, 1] value for the first index, you set them equal to the group value only looking at j. After that, you're not using those variables anymore either.
Do you see my concern? Have you tested the model on a smaller time horizon, to check its correctness?
Happy to think about improvements, but I fear that I don't fully understand the model yet. Most importantly: how many different tariff values can be chosen for the whole year? What does a "tariff level" in your bullet list refer to?
Kind regards,
Ronald0 
Hi Ronald,
thanks for your quick reply.
I did not tested the above model so far, as it did not yet converge for a timeframe which holds al least two months. But I created as a first step a simpler model, which holds only the first constraint:1. only three tariff levels are allowed over the year
# Read data and preprocess
data = pd.read_csv('data/grid_cost_timeseries_h.csv', sep=';')# Get only small subset
data = data.loc['20210101 00:00':'20210201 23:00']
# Example costs time series
costs_time_series = data['Netzentgelt'].values# Initialize model
model = ConcreteModel()# Define the length of the time series
N = len(costs_time_series)# Decision variables for tariffs
model.tariffs = Var(RangeSet(N), within=NonNegativeReals)# Three potential group values
model.group_values = Var(RangeSet(3), within=NonNegativeReals)
# Binary variables indicating group membership
model.group_membership = Var(RangeSet(N), RangeSet(3), within=Binary)# Objective: Minimize squared error
def objective_rule(model):
return quicksum((model.tariffs[i]  costs_time_series[i1])**2 for i in model.tariffs)
model.objective = Objective(rule=objective_rule, sense=minimize)# Constraint: Each tariff must be assigned to exactly one group
def one_group_constraint(model, i):
return sum(model.group_membership[i, j] for j in RangeSet(3)) == 1
model.one_group = Constraint(RangeSet(N), rule=one_group_constraint)# Constraint: Link tariff values to group values
def link_tariff_group_constraint(model, i):
return model.tariffs[i] == sum(model.group_values[j] * model.group_membership[i, j] for j in RangeSet(3))
model.link_tariff_group = Constraint(RangeSet(N), rule=link_tariff_group_constraint)# Solve the model
solver = SolverFactory('gurobi') # Using Gurobi solver
result = solver.solve(model, tee=True)# Extract solution
optimal_tariff_time_series = np.array([model.tariffs[i].value for i in RangeSet(N)])
group_values = np.array([model.group_values[j].value for j in RangeSet(3)])
group_membership = np.array([[model.group_membership[i,j].value for j in RangeSet(3)] for i in RangeSet(N)])print("Optimal Tariffs:", optimal_tariff_time_series)
print("Group Values:", group_values)This works fine but but also not converges for a full year.
I see your concerns for the integration of constraint 2 and constraint 3, this is not yet correctly implemented!
Regarding qour question: Over the whole year only three different tariff values are allowed to be used. Which does the model here already.
Constraint 2 shall guarantee that the time sequences of those three tariff values are only allowed to differ between weekdays and weekends. Constraint 3 shall guarantee that the time sequence within a month need to be the same but may differn between months.
I am happy about any recommendations how to improve the constraint 2 and 3 implementation and how to improve solver convergence.
I hope I could make my request a little bit more clear.
Thanks a lot!0 
Hi,
So if I understand correctly:
 We decide on three numbers  the three tariffs A, B, C
 Then, when looking at a specific month like January, then for each hour of a weekday (24 choices as the pattern applies to every weekday) we pick one of the three tariffs. So the weekday pattern might be 4x A, 4x B, 8x C, 4x B, 4x A and this would apply to every weekday in January.
 For the weekend days of January, you design a separate pattern. It applies to every weekend day in January.
 Then for February, you can construct a different pattern for the weekdays and another pattern for the weekend days.
Is that correct?  Cool puzzle by the way!
Kind regards,
Ronald0 
Hi Ronald,
you are totally correct that is the puzzle ; )
Thanks for your support!
0 
Hi Fabian,
Here's what I would probably do.
 First focus on adding the missing constraints. Performance comes later; the reason is that those missing constraints might greatly change the structure of the model.
 Regarding the formulation: You can simplify the problem by "merging" all weekdays of January. So, you only need one tariff variable for the price on a weekday in January 0:0000:59am. That gives you 12x2x24 variables. Based on those, you can create the objective function.
 Then, you need the three "main" variables for the three tariffs, and a matrix of booleans for the 12x2x24 by 3. You can link them in a way like you're doing it already in your code.
 Then, I would test for a smaller dataset. Don't just reduce the #months to 2, but also consider just looking at three or four hours a day. Make sure the results are really what you expect. Using a small dataset allows you to iterate much faster.
 Finally look at performance. One issue seems to be that the bound is 0 and doesn't really improve. If that's still the case when the model is complete, this might be something to focus on. There may be some tricks to provide a good lowerbound on the minimum deviation. For example, if you take the minimum and maximum cost for the 0:000:59am slot on a workday in January, you can already reason about the minimum deviation you will pay, knowing that you can only select a single tariff for all those slots.
Kind regards,
Ronald0 
Hi Ronald,
ok thanks I will try to further develop the model with these hints.
Fabian0
Please sign in to leave a comment.
Comments
6 comments