Defining piecewise linear constraints¶

In this tutorial, we use the national scale example model to implement a piecewise linear constraint. This constraint will represent a non-linear relationship between capacity and cost per unit capacity of Concentrating Solar Power (CSP).

In [1]:

import numpy as np
import plotly.express as px

import calliope
from calliope.io import read_rich_yaml

calliope.set_log_verbosity("INFO", include_solver_output=False)

import numpy as np import plotly.express as px import calliope from calliope.io import read_rich_yaml calliope.set_log_verbosity("INFO", include_solver_output=False)

Model setup¶

Defining our piecewise curve¶

In the base national scale model, the CSP has a maximum rated capacity of 10,000 kW and a cost to invest in that capacity of 1000 USD / kW.

In our updated model, the cost to invest in capacity will vary from 5000 USD / kW to 500 USD / kW as the CSP capacity increases:

In [2]:

capacity_steps = [0, 2500, 5000, 7500, 10000]
cost_steps = [0, 3.75e6, 6e6, 7.5e6, 8e6]

cost_per_cap = np.nan_to_num(np.divide(cost_steps, capacity_steps)).astype(int)

fig = px.line(
    x=capacity_steps,
    y=cost_steps,
    labels={"x": "Capacity (kW)", "y": "Investment cost (USD)"},
    markers="o",
    range_y=[0, 10e6],
    text=[f"{i} USD/kW" for i in cost_per_cap],
)
fig.update_traces(textposition="top center")
fig.show()

capacity_steps = [0, 2500, 5000, 7500, 10000] cost_steps = [0, 3.75e6, 6e6, 7.5e6, 8e6] cost_per_cap = np.nan_to_num(np.divide(cost_steps, capacity_steps)).astype(int) fig = px.line( x=capacity_steps, y=cost_steps, labels={"x": "Capacity (kW)", "y": "Investment cost (USD)"}, markers="o", range_y=[0, 10e6], text=[f"{i} USD/kW" for i in cost_per_cap], ) fig.update_traces(textposition="top center") fig.show()

[2026-04-22 17:27:56] WARNING  /tmp/ipykernel_3475/173276940.py:4: RuntimeWarning: invalid value encountered in divide
  cost_per_cap = np.nan_to_num(np.divide(cost_steps, capacity_steps)).astype(int)

We can then provide this data when we load our model:

Note

We must index our piecewise data over "breakpoints".

In [3]:

new_input_data = f"""
  data_definitions:
    capacity_steps:
      data: {capacity_steps}
      index: [0, 1, 2, 3, 4]
      dims: "breakpoints"
    cost_steps:
      data: {cost_steps}
      index: [0, 1, 2, 3, 4]
      dims: "breakpoints"
"""
print(new_input_data)

new_input_data = f""" data_definitions: capacity_steps: data: {capacity_steps} index: [0, 1, 2, 3, 4] dims: "breakpoints" cost_steps: data: {cost_steps} index: [0, 1, 2, 3, 4] dims: "breakpoints" """ print(new_input_data)

  data_definitions:
    capacity_steps:
      data: [0, 2500, 5000, 7500, 10000]
      index: [0, 1, 2, 3, 4]
      dims: "breakpoints"
    cost_steps:
      data: [0, 3750000.0, 6000000.0, 7500000.0, 8000000.0]
      index: [0, 1, 2, 3, 4]
      dims: "breakpoints"

Creating our piecewise constraint¶

We create the piecewise constraint by linking decision variables to the piecewise curve we have created. In this example, we need:

1. a new decision variable for investment costs that can take on the value defined by the curve at a given value of flow_cap;
2. to link that decision variable to our total cost calculation by injecting it into the cost_investment_flow_cap expression; and
3. to define the piecewise constraint.

Note that we replace the CSP investment cost related to its flow_cap here. So, although the technology still has a cost_flow_cap associated with it in the input data, this parameter will now have no effect.

In [4]:

new_math = """
  dimensions:
    breakpoints:
      dtype: integer
      ordered: true
      iterator: breakpoint
  parameters:
    capacity_steps:
      description: Capacity at each piecewise curve breakpoint
      unit: power
    cost_steps:
      description: Investment cost at each piecewise curve breakpoint
      unit: cost
  variables:
    piecewise_cost_investment:
      description: "Investment cost that increases monotonically"
      foreach: ["nodes", "techs", "carriers", "costs"]
      where: "[csp] in techs"
      bounds:
        min: 0
        max: .inf
      default: 0
      unit: cost
  global_expressions:
    cost_investment_flow_cap:
      where: (flow_cap AND (cost_flow_cap OR cost_flow_cap_per_distance)) OR piecewise_cost_investment
      equations:
        - expression: "$cost_sum * flow_cap"
          where: "NOT [csp] in techs"
        - expression: "piecewise_cost_investment"
          where: "[csp] in techs"
  piecewise_constraints:
    csp_piecewise_costs:
      description: "Set investment costs values along a piecewise curve using special ordered sets of type 2 (SOS2)."
      foreach: ["nodes", "techs", "carriers", "costs"]
      where: "piecewise_cost_investment"
      x_expression: "flow_cap"
      x_values: "capacity_steps"
      y_expression: "piecewise_cost_investment"
      y_values: "cost_steps"
"""

new_math = """ dimensions: breakpoints: dtype: integer ordered: true iterator: breakpoint parameters: capacity_steps: description: Capacity at each piecewise curve breakpoint unit: power cost_steps: description: Investment cost at each piecewise curve breakpoint unit: cost variables: piecewise_cost_investment: description: "Investment cost that increases monotonically" foreach: ["nodes", "techs", "carriers", "costs"] where: "[csp] in techs" bounds: min: 0 max: .inf default: 0 unit: cost global_expressions: cost_investment_flow_cap: where: (flow_cap AND (cost_flow_cap OR cost_flow_cap_per_distance)) OR piecewise_cost_investment equations: - expression: "$cost_sum * flow_cap" where: "NOT [csp] in techs" - expression: "piecewise_cost_investment" where: "[csp] in techs" piecewise_constraints: csp_piecewise_costs: description: "Set investment costs values along a piecewise curve using special ordered sets of type 2 (SOS2)." foreach: ["nodes", "techs", "carriers", "costs"] where: "piecewise_cost_investment" x_expression: "flow_cap" x_values: "capacity_steps" y_expression: "piecewise_cost_investment" y_values: "cost_steps" """

Building and checking the piecewise constraint¶

With our inputs and piecewise constraint defined, we can build our optimisation problem and inject this new math.

In [5]:

new_input_data_as_dict = read_rich_yaml(new_input_data)
new_math_as_dict = read_rich_yaml(new_math)

m = calliope.examples.national_scale(
    override_dict=new_input_data_as_dict,
    math_dict={"piecewise_math": new_math_as_dict},
    extra_math=["piecewise_math"],
)

new_input_data_as_dict = read_rich_yaml(new_input_data) new_math_as_dict = read_rich_yaml(new_math) m = calliope.examples.national_scale( override_dict=new_input_data_as_dict, math_dict={"piecewise_math": new_math_as_dict}, extra_math=["piecewise_math"], )

[2026-04-22 17:27:56] INFO     Math init | loading pre-defined math.

[2026-04-22 17:27:56] INFO     Math init | loading math files {'spores', 'operate', 'base', 'storage_inter_cluster', 'milp'}.

[2026-04-22 17:27:56] INFO     Model: preprocessing data

[2026-04-22 17:27:56] INFO     Math build | building applied math with ['base', 'piecewise_math'].

[2026-04-22 17:27:56] INFO     input data `color` not defined in model math; it will not be available in the optimisation problem.

[2026-04-22 17:27:56] INFO     input data `name` not defined in model math; it will not be available in the optimisation problem.

[2026-04-22 17:27:56] INFO     input data `color` not defined in model math; it will not be available in the optimisation problem.

[2026-04-22 17:27:56] INFO     input data `name` not defined in model math; it will not be available in the optimisation problem.

[2026-04-22 17:27:56] INFO     Model: initialisation complete

In [6]:

m.inputs.capacity_steps

m.inputs.capacity_steps

Out[6]:

<xarray.DataArray 'capacity_steps' (breakpoints: 5)> Size: 40B
array([    0.,  2500.,  5000.,  7500., 10000.])
Coordinates:
  * breakpoints  (breakpoints) int64 40B 0 1 2 3 4

xarray.DataArray

'capacity_steps'

• breakpoints: 5

• 0.0 2.5e+03 5e+03 7.5e+03 1e+04

  array([    0.,  2500.,  5000.,  7500., 10000.])

• Coordinates: (1)
  • breakpoints
    (breakpoints)
    int64
    0 1 2 3 4

    array([0, 1, 2, 3, 4])

• Indexes: (1)
  • breakpoints
    PandasIndex

    PandasIndex(Index([0, 1, 2, 3, 4], dtype='int64', name='breakpoints'))

• Attributes: (0)

In [7]:

m.inputs.cost_steps

m.build()

m.inputs.cost_steps m.build()

[2026-04-22 17:27:56] INFO     Model: backend build starting

[2026-04-22 17:27:57] INFO     Optimisation Model | parameters/lookups | Generated.

[2026-04-22 17:27:57] INFO     Optimisation Model | variables | Generated.

[2026-04-22 17:27:58] INFO     Optimisation Model | global_expressions | Generated.

[2026-04-22 17:27:59] INFO     Optimisation Model | constraints | Generated.

[2026-04-22 17:27:59] INFO     Optimisation Model | piecewise_constraints | Generated.

[2026-04-22 17:27:59] INFO     Optimisation Model | objectives | Generated.

[2026-04-22 17:27:59] INFO     Model: backend build complete

And we can see that our piecewise constraint exists in the built optimisation problem "backend"

In [8]:

m.backend.verbose_strings()
m.backend.get_piecewise_constraint("csp_piecewise_costs").to_series().dropna()

m.backend.verbose_strings() m.backend.get_piecewise_constraint("csp_piecewise_costs").to_series().dropna()

Out[8]:

nodes      techs  carriers  costs   
region1_1  csp    power     monetary    piecewise_constraints[csp_piecewise_costs][reg...
region1_2  csp    power     monetary    piecewise_constraints[csp_piecewise_costs][reg...
region1_3  csp    power     monetary    piecewise_constraints[csp_piecewise_costs][reg...
Name: csp_piecewise_costs, dtype: object

Solve the optimisation problem¶

Once we have all of our optimisation problem components set up as we desire, we can solve the problem.

In [9]:

m.solve()

m.solve()

[2026-04-22 17:27:59] INFO     Optimisation model | starting model in base mode.

[2026-04-22 17:28:01] INFO     Optimisation Model | postprocess | Generated.

[2026-04-22 17:28:01] INFO     Backend: solver finished running. Time since start of solving optimisation problem: 0:00:01.275987

[2026-04-22 17:28:01] INFO     Postprocessing: applied zero threshold 1e-10 to model results.

[2026-04-22 17:28:01] INFO     Postprocessing: ended. Time since start of solving optimisation problem: 0:00:01.310130

[2026-04-22 17:28:01] INFO     Backend: model solve completed. Time since start of solving optimisation problem: 0:00:01.310448

The results are stored in an array that can be accessed from m.results.

Analysing the outputs¶

In [10]:

# Absolute
csp_cost = m.results.cost_investment_flow_cap.sel(techs="csp")
csp_cost.to_series().dropna()

# Absolute csp_cost = m.results.cost_investment_flow_cap.sel(techs="csp") csp_cost.to_series().dropna()

Out[10]:

nodes      carriers  costs   
region1_1  power     monetary    8000000.0
region1_2  power     monetary          0.0
region1_3  power     monetary    3780385.9
Name: cost_investment_flow_cap, dtype: float64

In [11]:

# Relative to capacity
csp_cap = m.results.flow_cap.sel(techs="csp")
csp_cost_rel = csp_cost / csp_cap
csp_cost_rel.to_series().dropna()

# Relative to capacity csp_cap = m.results.flow_cap.sel(techs="csp") csp_cost_rel = csp_cost / csp_cap csp_cost_rel.to_series().dropna()

Out[11]:

nodes      carriers  costs   
region1_1  power     monetary     800.00000
region1_3  power     monetary    1492.00507
dtype: float64

In [12]:

# Plotted on our piecewise curve
fig.add_scatter(
    x=csp_cap.to_series().dropna().values,
    y=csp_cost.to_series().dropna().values,
    mode="markers",
    marker_symbol="cross",
    marker_size=10,
    marker_color="cyan",
    name="Installed capacity",
)
fig.show()

# Plotted on our piecewise curve fig.add_scatter( x=csp_cap.to_series().dropna().values, y=csp_cost.to_series().dropna().values, mode="markers", marker_symbol="cross", marker_size=10, marker_color="cyan", name="Installed capacity", ) fig.show()

Troubleshooting¶

If you are failing to load a piecewise constraint or it isn't working as expected, here are some common things to note:

1. The extent of your x_values and y_values will dictate the maximum values of your piecewise decision variables. In this example, we define capacity_steps over the full capacity range that we allow our CSP to cover in the model. However, if we set capacity_steps to [0, 2500, 5000, 7500, 9000] then flow_cap would never go above a value of 9000.

2. The x_values and y_values parameters must have the same number of breakpoints and be indexed over breakpoints. It is possible to extend these parameters to be indexed over other dimensions (e.g., different technologies with different piecewise curves) but it must always include the breakpoints dimension.

3. x_values must increase monotonically. That is, [0, 5000, 2500, 7500, 10000] is not valid for capacity_steps in this example. y_values, on the other hand, can vary any way you like; [0, 6e6, 3.75e6, 8e6, 7.5e6] is valid for cost_steps.

4. x_expression and y_expression must include reference to at least one decision variable. It can be a math expression, not only a single decision variable. flow_cap + storage_cap / 2 would be valid for x_expression in this example.

5. Piecewise constraints will make your problem more difficult to solve since each breakpoint adds a binary decision variable. Larger models with detailed piecewise constraints may not solve in a reasonable amount of time.