Special Ordered Sets of type 2 (SOS2) for piecewise linearisation

[brynpickering]

With the new math formulation, it is possible to define simple piecewise linear constraints (e.g. for investment costs), possibly solving #107.

In these, a convex hull is created that approximates the non-linear curve and the structure of the model ensures that the variable will sit somewhere along the upper bound of that hull and therefore take on a value very similar to the nonlinear curve's values (see image, Fig 4.2 from my phd).

For example, efficiency might change with technology output. As the output increases on the Y axis, the curve lead to increasingly more input requirements, therefore simulating a reduced efficiency with greater input. This works because the model wants to maximise output for each unit input (otherwise it has to "spend" more on the input carrier, increasing the objective function value) and so will have to sit on the upper bounding curve.

The advantage of this approach is that the model can avoid lots of integer/binary variables. The disadvantage is that it can only be used to describe problems where things get worse with increased output/capacity of a technology. In reality, economies of scale lead to lower investment costs for increased capacity, and part-load efficiencies are often lower than full load efficiency.

To describe these problems, I believe we need to use Special Ordered Sets. in SOS of type 2 - SOS2 - each "joint" on the piecewise curve is described by a binary variable and only two adjacent variables may be non-zero at any time. This means only one section of a piecewise curve is "active" at a given moment and so input and output must be related according to the slope of the line between those two points. SOS2 come with the advantage of being able to describe most non-linear curve shapes in the form of a discrete set of straight lines. But they need many more binary variables so might make most problems intractable. Nevertheless, this should be an option for Calliope modellers.

Solvers have python APIs to define an SOS that we could tap into (Pyomo, Gurobi). They also have more generalised (and easier to use) piecewise linear constraint interfaces (Pyomo, Gurobi) but we still need to provide a way for users to construct them.

If we use the generalised piecewise linear constraint constructor, I would force it to use SOS2 with the "==" operator between x and y. This is partly because it is all that Gurobi can handle (Pyomo has many more options), but also because it is the most likely formulation of a piecewise constraint that someone will think of.

We can't use the standard math syntax we have. Instead we need something that allows reference to the x and y variables being compared and the values of x and y at each "breakpoint". Some ideas:

1

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: 
       variable: flow_cap
       values: [0, 1, 2, 3, 4]
    y: 
       variable: piecewise_cost_investment
       values: [0, 10, 15, 18, 20]

2

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: flow_cap
    x_values: [0, 1, 2, 3, 4]
    y: piecewise_cost_investment
    y_values: [0, 10, 15, 18, 20]

3

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: flow_cap
    y: piecewise_cost_investment
    breakpoints: 
      0: 0
      1: 10
      2: 15
      3: 18
      4: 20

4

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: flow_cap
    y: piecewise_cost_investment
    breakpoints: 
      - [0, 0]
      - [1, 10]
      - [2, 15]
      - [3, 18]
      - [4, 20]

5

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: flow_cap
    y: piecewise_cost_investment
    breakpoints:       
      - {x: 0, y: 0}
      - {x: 1, y: 10}
      - {x: 2, y: 15}
      - {x: 3, y: 18}
      - {x: 4, y: 20}

Option 1 and 2 are probably better as breakpoints will likely be defined per technology / node and that could be defined by an input parameter in those options. E.g.,:

1:

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: 
       variable: flow_cap
       values: my_constraint_breakpoint_data_x
    y: 
       variable: piecewise_cost_investment
       values: my_constraint_breakpoint_data_y

and...

techs:
  some_tech:
    my_constraint_breakpoint_data_x: 
      data: [0, 1, 2, 3, 4]
      index: [0, 1, 2, 3, 4]
      dims: breakpoints
    my_constraint_breakpoint_data_y:
      data: [0, 1, 2, 3, 4]
      index: [0, 10, 15, 18, 20]
      dims: breakpoints

However, there would need to be an additional model dimension for breakpoints that this parameter is indexed over. All the options would work if you have to define one piecewise constraint for each index item in foreach. Maybe that's OK because it's unlikely that lots of different piecewise constraints are going to be needed (one per tech for 5-15 techs?) but it is a break from the current math format.

Or, there could be one parameter per breakpoint, so the options would become, e.g.:

1:

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: 
       variable: flow_cap
       values: [my_constraint_breakpoint_data_1_x, my_constraint_breakpoint_data_2_x, ...]
    y: 
       variable: piecewise_cost_investment
       values: [my_constraint_breakpoint_data_1_y, my_constraint_breakpoint_data_2_y, ...]

4:

piecewise_constraints:
  my_constraint:
    foreach: [nodes, techs]
    x: flow_cap
    y: piecewise_cost_investment
    breakpoints:
      - [my_constraint_breakpoint_data_1_x, my_constraint_breakpoint_data_1_y]
      - [my_constraint_breakpoint_data_2_x, my_constraint_breakpoint_data_2_y]

[jmorrisnrel]

Has there been any consideration of perhaps farther down the line looking for areas where certain piecewise linear parameters could be linearized by splitting out the nodes/techs into separate linear technologies? A specific example would be increasing LCOE-based cost-capacity curves where economies of scale are outweighed by increasing LCOEs as less-desirable locations are included (assuming there is no base cost/y offset incurred). This would be highly dependent on other constraints fitting around splitting nodes within the model (i.e. a non-linear efficiency curve wouldn't work with a split technology without other processing done/a more complex constraint added), but many of our technical resource assessments output linearized costs for a number of points that would need to be aggregated into one region/node for a Calliope model, with build decisions reliably going from lowest to highest LCOE.

A modeler could preprocess that data to create linear sub-technologies to represent the different pieces of the curve, but if there are opportunities to do some of that within Calliope and increase model size/performance and ease of modeling that would be valuable to some users.

[brynpickering]

LCOE is an output of the model as it is the (total investment in capacity + the total operating costs ) / (total electricity output) that defines it. My interpretation of your case is that the closer you get to the use of all e.g., available area, the less desirable the remaining area is. This might be because the cost of land increases, access for labour and materials is restricted, or because the yield from these locations will be lower. All would adversely affect the LCOE as you have depicted in that curve.

Your figure looks like the combinatorial curve of economies of scale (the flattening off of the curve from zero to 20GW) and then pushing capacity into less desirable locations (the rising curve from 130GW). This could be depicted in Calliope using piecewise linearisation, but would require a different Y-axis. E.g., total investment cost. Perhaps two separate curves would be possible, one for cumulative capacity against investment cost and the other for cumulative capacity against maximum power output (i.e., capacity factor).

I'd be more inclined to say that defining different technologies (expensive_productive_pv, cheap_productive_pv, expensive_not_productive_pv, ...) would be better as it gives you the possibility of more fine-grain control and perhaps makes it easier to interpret the results.

Still, your example is a good one for using SOS2 since the curve changes the sign of its gradient. That's something that the convex hull approach can't handle (see case (d) in my image).

Let me know if I've got the wrong end of the stick as to the problem you're describing!