Piecewise model performance

Hello,
I wrote a piecewise object for oemof and I recognised some performance problems. The computer can solve this model really fast, but the result processing takes forever (multiple days). I saw that my piecewise object creates lots of indexed variables, which seems ok, but I think this is the reason why the processing takes so long.
Does anyone know how to speed up the processing?
This codesegment is way to slow for my model:

var_dict = {(str(bv).split('.')[0], str(bv).split('.')[-1], i): bv[i].value
            for bv in block_vars for i in getattr(bv, '_index')}

This is the piecewise object:
class GasBus(Bus):

def __init__(self, *args, **kwargs):
    super().__init__(*args, **kwargs)
    self.slack=kwargs.get('slack', False)
    self.p_max = kwargs.get('p_max', 1000)
    self.p_min = kwargs.get('p_min', -1000)

class GasLineBlock(SimpleBlock):

CONSTRAINT_GROUP = True

def __init__(self, *args, **kwargs):
    super().__init__(*args, **kwargs)

def _create(self, group=None):

    if group is None:
        return None

    m = self.parent_block()

    self.GASLINE_PL = Set(initialize=[n for n in group])

    self.GAS_BUSES = Set(initialize=[s for s in m.es.nodes if isinstance(s, GasBus)])



    self.pressure=Var(self.GAS_BUSES, m.TIMESTEPS, bounds=(0,1))

    self.delta_pressure = Var(self.GASLINE_PL, m.TIMESTEPS, bounds=(-1,1))

    self.energy_flow = Var(self.GASLINE_PL, m.TIMESTEPS)

    for n in self.GASLINE_PL:
        for t in m.TIMESTEPS:
            for ob in list(n.outputs.keys()):
                if ob.slack == True:
                    self.pressure[ob,t].value=1
                    self.pressure[ob,t].fix()

            for ob in list(n.inputs.keys()):
                if ob.slack == True:
                    self.pressure[ob,t].value=1
                    self.pressure[ob,t].fix()


    self.piecewise=Piecewise(self.GASLINE_PL,
                    m.TIMESTEPS,
                    self.energy_flow,
                    self.delta_pressure,
                    pw_pts=[-1, 0, 0.25, 0.50, 0.75, 1],
                    pw_constr_type='EQ',
                    pw_repn='CC',
                    f_rule=[0, 0, 0.50, 0.707107, 0.866025, 1, 3.1628])


    def flow_eq_pressure(block,n,t):
        expr = 0
        expr += (self.pressure[list(n.outputs.keys())[0],t]-self.pressure[list(n.inputs.keys())[0],t])
        expr += self.delta_pressure[n,t]

        return expr == 0

    self.flow_eq_pressure = Constraint(self.GASLINE_PL, m.TIMESTEPS, rule=flow_eq_pressure)

    def energy_flow_out(block, n, t):
        expr = 0
        expr += - m.flow[n, list(n.outputs.keys())[0], t]
        expr += self.energy_flow[n,t]*n.K_1

        return expr == 0

    self.energy_flow_out = Constraint(self.GASLINE_PL, m.TIMESTEPS, rule=energy_flow_out)

    def energy_flow_in(block, n, t):
        expr = 0
        expr += - m.flow[list(n.inputs.keys())[0],n, t]*n.conv_factor
        expr += self.energy_flow[n,t]*n.K_1

        return expr == 0

    self.energy_flow_in = Constraint(self.GASLINE_PL, m.TIMESTEPS, rule=energy_flow_in)

class GasLine(Transformer):

def __init__(self, *args, **kwargs):

    super().__init__(*args, **kwargs)
    self.conv_factor = kwargs.get('conv_factor', 1)
    self.K_1=kwargs.get('K_1')
    self.input_list = list(kwargs.get('input_list', []))
    self.output_list = list(kwargs.get('output_list', []))

    if len(self.inputs) > 1 or len(self.outputs) > 1:
        raise ValueError("Component GasLine must not have more than \
                         one input and one output!")

def constraint_group(self):
    return GasLineBlock

If anyone has an idea how to speed things up I would be very thankfull.
Best, Philipp

Is there nobody who has an idea?

Best, Philipp

You could use a code-profiler to find out which part of the code takes the time.

The profiler of the Spyder IDE is pretty easy to use. If it takes days you should reduce the number of time steps to speed up the problem for the analysis.

I still wonder why the processing could use days. You could also add messages between each line to make sure it is really the processing that takes so long.

from datetime import datetime
start = datetime.now()
[...]
print("START OPTIMISING:", datetime.now()-start)
model.solve(solver=my_solver)

print("START PROCESSING:", datetime.now()-start)
results = outputlib.processing.results(model)

After that you might be able to describe the problem more detailed.

I am not sure it might also be related to issue #630.

I have only the dimmest idea of how oemof works. But the OSeMOSYS folk had performance problems translating their model from MathProg into its solver‑oriented representation. @jonas.hoersch fixed this issue by recognizing that some of the deeply nested traversals were unnecessary and that the required functionality could be undertaken within earlier and/or parallel traversals. Nested traversals scale really badly as you are probably aware. This fix depends on the nature of your problem but if you can remove unneeded nested traversals then you should do so. This strategy is different from profiling (as indicated in the previous post) because it leverages off problem knowledge and not code enhancement. HTH, R.

Thank you @uwe.krien and @robbie.morrison for your answers. They gave us hints to identify the issue and find a solution.

We used datetime to identify the nested for-loop in oemof outputli - processing - create_dataframe()

var_dict = {(str(bv).split('.')[0], str(bv).split('.')[-1], i): bv[i].value
            for bv in block_vars for i in getattr(bv, '_index')}

as the cause of the poor performance, togehter with the piecewise model results as stated here https://github.com/oemof/oemof/pull/592#issuecomment-514685806:

The following problem arises in the outputlib when using pyomo’s Piecewise:
The blockvars for all other constraints have this form:

flow[gas,pwltf,t]
flow[electricity,demand,t]
PiecewiseLinearTransformerBlock.inflow[pwltf,t]

However, the blockvars produced by pyomo’s Piecewise look like this:

PiecewiseLinearTransformerBlock.piecewise[pwltf,t].CC_lambda[i]
PiecewiseLinearTransformerBlock.piecewise[pwltf,t].CC_bin_y[i]

, where i is the index for the different segments in the piecewise linearised function. Apparently, Piecewise creates child variables
that carry this new piecewise index.

This gives wrong results in outputlib.processing.create_dataframe() when creating the pyomo tuple/oemof tuple.
The index i is written instead of [pwlft,t] , which it should be according to what usually happens.

For large models it is the huge number of those pyomo piecewise blockvars that causes the extensive time consumption (and additionally, also causes wrong results). Since the results in the picewise blockvars are not important to our problem, we excluded them from processing. We only need those where a proper tuple can be created. Our var_dict is now created:

var_dict = {}
    for bv in block_vars:
        for i in getattr(bv, '_index'):
            if isinstance(i, tuple):
                var_dict[(str(bv).split('.')[0], str(bv).split('.')[-1], i)] = bv[i].value

Maybe this is not the best solution, but at the moment it works for our models.
Lukas

2 Likes

Great that you have a workable solution. By the way, here is the press release for the OSeMOSYS fix by @jonas.hoersch. His improvement decreased the model translation time “from around 15 hours to 9 minutes”. The was no change, of course, to the solution time.

Thanks from me too :smiley: