PYOMO: LP Heat storage optimalization problem, I want to define the domain of a variable with discrete floats

Question

I am working on optimizing a heat dispatch model for a District Heating (DH) network. In this model, I have different heat assets that can produce heat at varying prices. The challenge arises from the fact that the price of heat produced by the Combined Heat and Power (CHP) plant depends on its heat load.

I want to include this dependency in the optimization process. Specifically, I aim to constrain the CHP heat load (model.CHP[t]) and its corresponding price ($model.price_CHP[t]$) to discrete points. At each time step, I have predefined arrays containing different heat load points ($chp.Line_PQ$) and their corresponding prices ($chp.Price_PQ[0]$).

Here are the relevant variables and parameters from my Pyomo model:

Variablesyour text

    model.CHP = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.EB = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.HP = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.BO = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.price_CHP = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)```

    model.CHP = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.EB = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.HP = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.BO = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)
    model.price_CHP = pyomo.Var(model.T, domain=pyomo.NonNegativeReals)

Parameters

    model.price_EB = pyomo.Param(model.T, initialize=Day_price_EB)
    model.price_HP = pyomo.Param(model.T, initialize=Day_price_HP)
    model.price_BO = pyomo.Param(model.T, initialize=Day_price_BO)

And here's the objective function I'm using:

    def objective_func(model):
        return sum([model.CHP[t]*model.price_CHP[t] + model.HP[t]*model.price_HP[t] + model.EB[t]*model.price_EB[t] +      model.BO[t]*model.price_BO[t] for t in model.T])

For the constraints, I've attempted to implement them as follows:

    discrete_points_chp = chp.Line_PQ
    discrete_prices_chp = chp.Price_PQ[0]

    def constrain_chp_to_discrete_points(model, t):
        return model.CHP[t] in set(discrete_points_chp)

    model.constr_chp_discrete = pyomo.Constraint(model.T, rule=constrain_chp_to_discrete_points)

    def constrain_price_chp_to_discrete_points(model, t):
    return model.price_CHP[t] in set(discrete_prices_chp)

However, when I run the optimization, I encounter issues with these constraints. Could someone please provide guidance on how to correctly implement such constraints in Pyomo?

SO I'm trying to use a set to define the discrete options for my variables but this does not work....

Answer 1

The constraint

 x ∈ {p[0],p[1],...p[n-1]}   (p[i] are constants)

can be implemented as:

 x = sum(p[i]*y[i])
 sum(y[i]) = 1
 y[i] ∈ {0,1}  (binary variable)

Answer 2

The aproach by Erwin worked for me, I had some struggles to implement it so therefore I want to share my example. I make use of the binary variable y(t,n) this will shift over all my time points and work load points and corresponding prices.

Code:

model.y = pyomo.Var(model.T, model.n , domain=pyomo.Binary)

def exactly_one_rule(model, t):
return sum(model.y[t, i] for i in model.n) == 1

def compute_Load_rule(model,t):
return model.CHP_Load[t] == sum(model.set_Loads[i+1] * model.y[t,i] for i in range(n))

model.compute_load_constraint = pyomo.Constraint(model.T, rule=compute_Load_rule)
model.compute_load_constraint.pprint()

# Compute x as the sum of p[i] * y[i]
def compute_Price_rule(model,t):
    return model.Price_CHP[t] == sum(model.Price_CHP_Set[i+1] * model.y[t,i] for i in range(n))

model.compute_price_constraint = pyomo.Constraint(model.T, rule=compute_Price_rule)
model.compute_price_constraint.pprint()

if you need any further help, feel free to reach out!