3

Below is a simple Pyomo script using the decorator syntax - I would like to understand how to use this syntax within a class - in this case inside Model.

None-class version

from pyomo.environ import *

import random

random.seed(1000)

model = AbstractModel()

model.N = Param(within=PositiveIntegers)
model.P = Param(within=RangeSet(1, model.N))
model.M = Param(within=PositiveIntegers)

model.Locations = RangeSet(1, model.N)
model.Customers = RangeSet(1, model.M)

model.d = Param(
    model.Locations,
    model.Customers,
    initialize=lambda n, m, model: random.uniform(1.0, 2.0),
    within=Reals,
)

model.x = Var(model.Locations, model.Customers, bounds=(0.0, 1.0))
model.y = Var(model.Locations, within=Binary)


@model.Objective()
def obj(model):
    return sum(
        model.d[n, m] * model.x[n, m] for n in model.Locations for m in model.Customers
    )


@model.Constraint(model.Customers)
def single_x(model, m):
    return (sum(model.x[n, m] for n in model.Locations), 1.0)


@model.Constraint(model.Locations, model.Customers)
def bound_y(model, n, m):
    return model.x[n, m] - model.y[n] <= 0.0


@model.Constraint()
def num_facilities(model):
    return sum(model.y[n] for n in model.Locations) == model.P

Decorator version within a class that doesn't work:

from pyomo.environ import *

import random

random.seed(1000)


class Model:
    def __init__(self):
        self.model = AbstractModel()

        self.model.N = Param(within=PositiveIntegers)
        self.model.P = Param(within=RangeSet(1, self.model.N))
        self.model.M = Param(within=PositiveIntegers)

        self.model.Locations = RangeSet(1, self.model.N)
        self.model.Customers = RangeSet(1, self.model.M)

        self.model.d = Param(
            self.model.Locations,
            self.model.Customers,
            initialize=lambda n, m, model: random.uniform(1.0, 2.0),
            within=Reals,
        )

        self.model.x = Var(
            self.model.Locations, self.model.Customers, bounds=(0.0, 1.0)
        )
        self.model.y = Var(self.model.Locations, within=Binary)

    @model.Objective()
    def obj(model):
        return sum(
            model.d[n, m] * model.x[n, m]
            for n in model.Locations
            for m in model.Customers
        )

    @model.Constraint(model.Customers)
    def single_x(model, m):
        return (sum(model.x[n, m] for n in model.Locations), 1.0)

    @model.Constraint(model.Locations, model.Customers)
    def bound_y(model, n, m):
        return model.x[n, m] - model.y[n] <= 0.0

    @model.Constraint()
    def num_facilities(model):
        return sum(model.y[n] for n in model.Locations) == model.P
Michael
  • 2,436
  • 1
  • 36
  • 57

2 Answers2

2

I'm not able to help you on this, I just have a few qustions:

  • do you know if the use of @model.Objective() (same for Constraint etc) is documented somewhere? I didn't know it existed, and it's awesome
  • why do you want your "function rules" to be methods of the class? couldn't you defined them as functions within the __init__ method?

I guess what I'm missing is the benefit of using a class in the first place. If you are just trying to wrap the model construction somehow, then a better approach is using a function:

def create_model():
    model = AbstractModel()

    ...

    @model.Constraint()
    def some_rule_function(model):
        ...

    ...

    return model

EDIT: if you really want to wrap everything into a class:

class Model:
    def __init__(self, model):
         self.model = model

    # alternative constructor:
    # def __init__(self):
    #     self.model = create_model()

    def construct(self, data):
        # get concrete model
        self.model = self.model.create_instance(data)

    def run(self, solver, **kwargs):
        with pe.SolverFactory(solver) as solver:
            solver.solve(self.model, **kwargs)

    def construct_and_run(self, data, solver, **kwargs):
        self.construct(data)
        self.data(solver, **kwargs)

    # other behavior you want to add to the class

example usage:

model = Model(create_model())
  • There is not much documentation around the use of decorator with pyomo however this example if coming from pyomo/examples/pyomo/p-median/decorated_pmedian.py on the pyomo github repository. The idea of using a class is to have a sort of wrapper for the model. I think my question goes beyond pyomo to more a general guidance on how to use a decorator within a class. Can the `create_model` method be part of the `Model` class itself in your example? – Michael Nov 29 '19 at 10:05
  • you have two choices: one is subclassing `AbstractModel` (which I don't recommend), the other is to have the `__init__` method of your class do what `create_model` does above, and still set `self.model = model` at the end. – Giorgio Balestrieri Nov 29 '19 at 10:20
  • Sorry, I don't really get the logic - My thinking was that I will wrap the full model, constraints, etc. into a class so I can instantiate it simply bypassing all the data and call a `run` method. This would create a level of abstraction. – Michael Nov 29 '19 at 15:49
  • I chanced slightly the example above, but I'm not 100% sure I understand what you want to achieve. The "class" that wraps everything is `AbstractModel` or `ConcreteModel`, that include all the components (you can access them through `model.component_name`). So you could just have a `run` or `main` function that creates the model first and then solves it. – Giorgio Balestrieri Nov 29 '19 at 17:42
2

Trying to answer your direct question, here's something that seems to work for me. My interpretation is that since your model is called self.model, the decorators should also match that.

Note that I used s as the first argument in the constraint method definitions just to see if it worked, but it could also be model or whatever you want to call it.

class Model:
    def __init__(self):
        self.model = pyo.AbstractModel()

        self.model.N = pyo.Param(initialize=5, within=pyo.PositiveIntegers)
        self.model.P = pyo.Param(initialize=3, within=pyo.RangeSet(1, self.model.N))
        self.model.M = pyo.Param(initialize=3, within=pyo.PositiveIntegers)

        self.model.Locations = pyo.RangeSet(1, self.model.N)
        self.model.Customers = pyo.RangeSet(1, self.model.M)

        self.model.d = pyo.Param(
            self.model.Locations,
            self.model.Customers,
            initialize=lambda n, m, model: random.uniform(1.0, 2.0),
            within=pyo.Reals,
        )

        self.model.x = pyo.Var(
            self.model.Locations, self.model.Customers, bounds=(0.0, 1.0)
        )
        self.model.y = pyo.Var(self.model.Locations, within=pyo.Binary)

        @self.model.Objective()
        def obj(s):
            return sum(
                s.d[n, m] * s.x[n, m]
                for n in s.Locations
                for m in s.Customers
            )

        @self.model.Constraint(self.model.Customers)
        def single_x(s, m):
            return (sum(s.x[n, m] for n in s.Locations), 1.0)

        @self.model.Constraint(self.model.Locations, self.model.Customers)
        def bound_y(s, n, m):
            return s.x[n, m] - s.y[n] <= 0.0

        @self.model.Constraint()
        def num_facilities(s):
            return sum(s.y[n] for n in s.Locations) == s.P

You would then be able to instantiate the model with model = Model(), though annoyingly (at least to me), all your Pyomo model components will be within the attribute model.model (e.g., model.model.P).

What I've done before to make the naming cleaner is to inherit from AbstractModel (though the other answer suggests that may not be good practice):

from pyomo.core.base.PyomoModel import AbstractModel

class Model(AbstractModel):
    def __init__(self):
        AbstractModel.__init__(self)

        self.N = pyo.Param(initialize=5, within=pyo.PositiveIntegers)
        self.P = pyo.Param(initialize=3, within=pyo.RangeSet(1, self.N))
        self.M = pyo.Param(initialize=3, within=pyo.PositiveIntegers)

        self.Locations = pyo.RangeSet(1, self.N)
        self.Customers = pyo.RangeSet(1, self.M)

        self.d = pyo.Param(
            self.Locations,
            self.Customers,
            initialize=lambda n, m, model: random.uniform(1.0, 2.0),
            within=pyo.Reals,
        )

        self.x = pyo.Var(
            self.Locations, self.Customers, bounds=(0.0, 1.0)
        )
        self.y = pyo.Var(self.Locations, within=pyo.Binary)

        @self.Objective()
        def obj(s):
            return sum(
                s.d[n, m] * s.x[n, m]
                for n in s.Locations
                for m in s.Customers
            )

        @self.Constraint(self.Customers)
        def single_x(s, m):
            return (sum(s.x[n, m] for n in s.Locations), 1.0)

        @self.Constraint(self.Locations, self.Customers)
        def bound_y(s, n, m):
            return s.x[n, m] - s.y[n] <= 0.0

        @self.Constraint()
        def num_facilities(s):
            return sum(s.y[n] for n in s.Locations) == s.P

In this case, you still instantiate as model = Model() but your Pyomo model components can be accessed as model.P.