Search code examples
pythonpython-3.xpyomo

Sparse sets in Pyomo Concrete Model


I want to create a sparse set in a ConcreteModel in pyomo. Here is a minimum example:

import pyomo
import pyomo.environ as pe

m = pe.ConcreteModel

m.P = pe.Set(initialize=['A', 'B', 'C'])
m.Q = pe.Set(initialize=[1, 2, 3, 4, 5, 6, 9, 10, 11])

D1 = {'A': [1, 2, 3], 'B': [9, 10, 11], 'C': [4, 5, 6]}
m.E = pe.Set(initialize=D1)

m.x = pe.Var(m.Q)

m.obj = pe.Objective(expr=1, sense=pe.minimize)


def constraint_rule(m: pe.ConcreteModel, p: str):
    return sum(m.x[i] for i in m.E[p]) <= 1


m.add_constraint = pe.Constraint(m.P, rule=constraint_rule)

opt = pe.SolverFactory('gurobi')
opt.solve(m)

When I run this model, I get the following message:

TypeError: valid_problem_types() missing 1 required positional argument: 'self'

Is this a problem in the construction of the set m.E?


Solution

  • You're getting that error because you're not creating an instance of the pyomo.environ.ConcreteModel but using an alias to it. You need to use parenthesis in m = pe.ConcreteModel()

    Now, I assume that you want to express your constraints something like this:

    x[1] + x[2] + x[3] <= 1
    x[9] + x[10] + x[11] <= 1
    x[4] + x[5] + x[6] <= 1
    

    using the actual sets created. Then you need to create model.E a a Subset of model.P, since the current values relays upong each value of model.P. Then you need to change it as follows:

    m.E = pe.Set(m.P, initialize=D1)
    

    The actual model will be something like this:

    import pyomo
    import pyomo.environ as pe
    
    m = pe.ConcreteModel()
    
    m.P = pe.Set(initialize=['A', 'B', 'C'])
    m.Q = pe.Set(initialize=[1, 2, 3, 4, 5, 6, 9, 10, 11])
    
    D1 = {'A': [1, 2, 3], 'B': [9, 10, 11], 'C': [4, 5, 6]}
    m.E = pe.Set(m.P, initialize=D1)
    
    m.x = pe.Var(m.Q)
    
    m.obj = pe.Objective(expr=1, sense=pe.minimize)
    
    
    def constraint_rule(m: pe.ConcreteModel, p: str):
        return sum(m.x[i] for i in m.E[p]) <= 1
    
    
    m.add_constraint = pe.Constraint(m.P, rule=constraint_rule)
    
    opt = pe.SolverFactory('gurobi')
    opt.solve(m)
    

    This will give you the following solution:

    >>>m.x.display()
    x : Size=9, Index=Q
        Key : Lower : Value : Upper : Fixed : Stale : Domain
          1 :  None :   1.0 :  None : False : False :  Reals
          2 :  None :   0.0 :  None : False : False :  Reals
          3 :  None :   0.0 :  None : False : False :  Reals
          4 :  None :   1.0 :  None : False : False :  Reals
          5 :  None :   0.0 :  None : False : False :  Reals
          6 :  None :   0.0 :  None : False : False :  Reals
          9 :  None :   1.0 :  None : False : False :  Reals
         10 :  None :   0.0 :  None : False : False :  Reals
         11 :  None :   0.0 :  None : False : False :  Reals
    

    I assume that you understand that this only generate a feasible solution, since the OF is a constant