Tutorial 3: Mixed Integer Linear Programming

This example is based on the urban scale example model, but with an override. In the model’s scenarios.yaml file overrides are defined which trigger binary and integer decision variables, creating a MILP model, rather than a conventional LP model.

Units

The capacity of a technology is usually a continuous decision variable, which can be within the range of 0 and energy_cap_max (the maximum capacity of a technology). In this model, we introduce a unit limit on the CHP instead:

    chp:
        constraints:
            units_max: 4
            energy_cap_per_unit: 300
            energy_cap_min_use: 0.2
        costs:
            monetary:
                energy_cap: 700
                purchase: 40000

A unit maximum allows a discrete, integer number of CHP to be purchased, each having a capacity of energy_cap_per_unit. Any of energy_cap_max, energy_cap_min, or energy_cap_equals are now ignored, in favour of units_max, units_min, or units_equals. A useful feature unlocked by introducing this is the ability to set a minimum operating capacity which is only enforced when the technology is operating. In the LP model, energy_cap_min_use would force the technology to operate at least at that proportion of its maximum capacity at each time step. In this model, the newly introduced energy_cap_min_use of 0.2 will ensure that the output of the CHP is 20% of its maximum capacity in any time step in which it has a non-zero output.

Purchase cost

The boiler does not have a unit limit, it still utilises the continuous variable for its capacity. However, we have introduced a purchase cost:

    boiler:
        costs:
            monetary:
                energy_cap: 35
                purchase: 2000

By introducing this, the boiler now has a binary decision variable associated with it, which is 1 if the boiler has a non-zero energy_cap (i.e. the optimisation results in investment in a boiler) and 0 if the capacity is 0. The purchase cost is applied to the binary result, providing a fixed cost on purchase of the technology, irrespective of the technology size. In physical terms, this may be associated with the cost of pipework, land purchase, etc. The purchase cost is also imposed on the CHP, which is applied to the number of integer CHP units in which the solver chooses to invest.

MILP functionality can be easily applied, but convergence is slower as a result of integer/binary variables. It is recommended to use a commercial solver (e.g. Gurobi, CPLEX) if you wish to utilise these variables outside this example model.

Asynchronous energy production/consumption

The heat pipes which distribute thermal energy in the network may be prone to dissipating heat in an unphysical way. I.e. given that they have distribution losses associated with them, in any given timestep, a link could produce and consume energy in the same timestep, losing energy to the atmosphere in both instances, but having a net energy transmission of zero. This allows e.g. a CHP facility to overproduce heat to produce more cheap electricity, and have some way of dumping that heat. The asynchronous_prod_con binary constraint ensures this phenomenon is avoided:

    heat_pipes:
        constraints:
            force_asynchronous_prod_con: true

Now, only one of carrier_prod and carrier_con can be non-zero in a given timestep. This constraint can also be applied to storage technologies, to similarly control charge/discharge.

Running the model

We now take you through running the model in a Jupyter notebook, which is included fully below. To download and run the notebook yourself, you can find it here. You will need to have Calliope installed.

Notebook

Calliope Urban Scale MILP Example Model

For more details on analysing input/output data, see the full urban scale example model

In [1]:
import calliope

# We increase logging verbosity
calliope.set_log_level('INFO')
In [2]:
model = calliope.examples.milp()

# Note, we see the overrides that we have applied printed here, thanks to inreasing logging verbosity
# We also see a warning that we are applying binary/integer decision variables in a model, to remind us that this
# model may take a while to run
[2019-05-27 15:05:00] INFO: Model: initialising
[2019-05-27 15:05:00] INFO: Applying the following overrides without a specific scenario name: ['milp']
[2019-05-27 15:05:00] INFO: Override applied to model.name: Urban-scale example model -> Urban-scale example model with MILP
`run.solver_options.mipgap`:0.05 applied from override as new configuration
`techs.boiler.costs.monetary.energy_cap`:35 applied from override as new configuration
`techs.boiler.costs.monetary.purchase`:2000 applied from override as new configuration
`techs.chp.constraints.energy_cap_min_use`:0.2 applied from override as new configuration
`techs.chp.constraints.energy_cap_per_unit`:300 applied from override as new configuration
`techs.chp.constraints.units_max`:4 applied from override as new configuration
Override applied to techs.chp.costs.monetary.energy_cap: 750 -> 700
`techs.chp.costs.monetary.purchase`:40000 applied from override as new configuration
`techs.heat_pipes.constraints.force_asynchronous_prod_con`:True applied from override as new configuration
[2019-05-27 15:05:00] INFO: Model: preprocessing stage 1 (model_run)
/Users/brynmorp/OneDrive/Projects/project-calliope/repos/calliope/calliope/core/preprocess/model_run.py:233: FutureWarning:

There will be no default cost class for the objective function in v0.7.0 (currently "monetary" with a weight of 1). Explicitly specify the cost class(es) you would like to use under `run.objective_options.cost_class`. E.g. `{"monetary": 1}` to replicate the current default.

[2019-05-27 15:05:01] INFO: Model: preprocessing stage 2 (model_data)
[2019-05-27 15:05:01] INFO: Model: preprocessing complete. Time since start: 0:00:01.353896
Warning: Possible issues found during model processing:
 * Integer and / or binary decision variables are included in this model. This may adversely affect solution time, particularly if you are using a non-commercial solver. To improve solution time, consider changing MILP related solver options (e.g. `mipgap`) or removing MILP constraints.
In [3]:
# Model inputs can be viewed at `model.inputs`. 
# Variables are indexed over any combination of `techs`, `locs`, `carriers`, `costs` and `timesteps`, 
# although `techs`, `locs`, and `carriers` are often concatenated. 
# e.g. `chp`, `X1`, `heat` -> `X1::chp::heat` 
model.inputs
Out[3]:
<xarray.Dataset>
Dimensions:                               (carrier_tiers: 3, carriers: 3, coordinates: 2, costs: 1, loc_carriers: 10, loc_tech_carriers_conversion_plus: 3, loc_techs: 26, loc_techs_area: 3, loc_techs_conversion: 2, loc_techs_conversion_plus: 1, loc_techs_export: 4, loc_techs_finite_resource: 9, loc_techs_investment_cost: 20, loc_techs_milp: 1, loc_techs_non_conversion: 23, loc_techs_om_cost: 9, loc_techs_supply_plus: 3, loc_techs_transmission: 10, locs: 4, techs: 9, timesteps: 48)
Coordinates:
  * loc_techs_transmission                (loc_techs_transmission) object 'N1::heat_pipes:X1' ... 'N1::heat_pipes:X2'
  * loc_techs_milp                        (loc_techs_milp) object 'X1::chp'
  * loc_techs_conversion_plus             (loc_techs_conversion_plus) object 'X1::chp'
  * loc_tech_carriers_conversion_plus     (loc_tech_carriers_conversion_plus) object 'X1::chp::gas' ... 'X1::chp::electricity'
  * loc_carriers                          (loc_carriers) object 'X1::heat' ... 'X1::electricity'
  * loc_techs_finite_resource             (loc_techs_finite_resource) object 'X3::demand_heat' ... 'X2::pv'
  * costs                                 (costs) object 'monetary'
  * loc_techs_investment_cost             (loc_techs_investment_cost) object 'X1::supply_gas' ... 'N1::heat_pipes:X2'
  * loc_techs                             (loc_techs) object 'X1::supply_gas' ... 'X2::pv'
  * locs                                  (locs) object 'N1' 'X2' 'X3' 'X1'
  * carrier_tiers                         (carrier_tiers) object 'out_2' ... 'out'
  * loc_techs_conversion                  (loc_techs_conversion) object 'X3::boiler' 'X2::boiler'
  * loc_techs_area                        (loc_techs_area) object 'X3::pv' ... 'X1::pv'
  * loc_techs_supply_plus                 (loc_techs_supply_plus) object 'X3::pv' ... 'X1::pv'
  * loc_techs_export                      (loc_techs_export) object 'X3::pv' ... 'X1::pv'
  * loc_techs_non_conversion              (loc_techs_non_conversion) object 'X1::supply_gas' ... 'X2::pv'
  * timesteps                             (timesteps) datetime64[ns] 2005-07-01 ... 2005-07-02T23:00:00
  * coordinates                           (coordinates) object 'x' 'y'
  * loc_techs_om_cost                     (loc_techs_om_cost) object 'X1::supply_gas' ... 'X2::boiler'
  * carriers                              (carriers) object 'gas' ... 'heat'
  * techs                                 (techs) object 'chp' ... 'demand_electricity'
Data variables:
    energy_cap_min_use                    (loc_techs) float64 nan nan ... nan
    export_carrier                        (loc_techs_export) <U11 'electricity' ... 'electricity'
    resource_area_per_energy_cap          (loc_techs_area) int64 7 7 7
    units_max                             (loc_techs_milp) int64 4
    energy_prod                           (loc_techs) float64 1.0 1.0 ... 1.0
    resource_area_max                     (loc_techs_area) int64 1500 1500 1500
    parasitic_eff                         (loc_techs_supply_plus) float64 0.85 ... 0.85
    energy_con                            (loc_techs) float64 nan 1.0 ... nan
    lifetime                              (loc_techs) float64 25.0 25.0 ... 25.0
    energy_cap_per_unit                   (loc_techs) float64 nan nan ... nan
    energy_eff                            (loc_techs) float64 nan 0.9269 ... nan
    resource                              (loc_techs_finite_resource, timesteps) float64 -0.0156 ... 0.0
    resource_unit                         (loc_techs_finite_resource) <U15 'energy' ... 'energy_per_area'
    force_asynchronous_prod_con           (loc_techs) float64 nan 1.0 ... nan
    resource_eff                          (loc_techs_finite_resource) float64 nan ... 1.0
    energy_cap_max                        (loc_techs) float64 2e+03 ... 250.0
    force_resource                        (loc_techs_finite_resource) bool True ... True
    cost_depreciation_rate                (costs, loc_techs_investment_cost) float64 0.1102 ... 0.1102
    cost_om_annual                        (costs, loc_techs_om_cost) float64 nan ... nan
    cost_energy_cap                       (costs, loc_techs_investment_cost) float64 1.0 ... 0.9
    cost_om_prod                          (costs, loc_techs_om_cost) float64 nan ... nan
    cost_export                           (costs, loc_techs_om_cost, timesteps) float64 nan ... nan
    cost_purchase                         (costs, loc_techs_investment_cost) float64 nan ... nan
    cost_om_con                           (costs, loc_techs_om_cost) float64 0.025 ... 0.004
    distance                              (loc_techs_transmission) float64 3.0 ... 3.0
    lookup_remotes                        (loc_techs_transmission) <U18 'X1::heat_pipes:N1' ... 'X2::heat_pipes:N1'
    available_area                        (locs) float64 nan 1.3e+03 900.0 500.0
    loc_coordinates                       (coordinates, locs) int64 5 8 ... 3 7
    colors                                (techs) <U7 '#E4AB97' ... '#072486'
    inheritance                           (techs) <U29 'conversion_plus' ... 'demand'
    names                                 (techs) <U29 'Combined heat and power' ... 'Electrical demand'
    carrier_ratios                        (carrier_tiers, loc_tech_carriers_conversion_plus) float64 1.0 ... 1.0
    carrier_ratios_min                    (carrier_tiers, loc_techs_conversion_plus) float64 0.8 ... 0.8
    lookup_loc_carriers                   (loc_carriers) <U175 'X1::demand_heat::heat,X1::heat_pipes:N1::heat,X1::chp::heat' ... 'X1::power_lines:X2::electricity,X1::pv::electricity,X1::power_lines:X3::electricity,X1::chp::electricity,X1::demand_electricity::electricity,X1::supply_grid_power::electricity'
    lookup_loc_techs                      (loc_techs_non_conversion) <U35 'X1::supply_gas::gas' ... 'X2::pv::electricity'
    lookup_loc_techs_conversion           (carrier_tiers, loc_techs_conversion) object None ... 'X2::boiler::heat'
    lookup_primary_loc_tech_carriers_in   (loc_techs_conversion_plus) <U12 'X1::chp::gas'
    lookup_primary_loc_tech_carriers_out  (loc_techs_conversion_plus) <U20 'X1::chp::electricity'
    lookup_loc_techs_conversion_plus      (carrier_tiers, loc_techs_conversion_plus) object 'X1::chp::heat' ... 'X1::chp::electricity'
    lookup_loc_techs_export               (loc_techs_export) <U20 'X3::pv::electricity' ... 'X1::pv::electricity'
    lookup_loc_techs_area                 (locs) <U6 '' 'X2::pv' ... 'X1::pv'
    timestep_resolution                   (timesteps) float64 1.0 1.0 ... 1.0
    timestep_weights                      (timesteps) float64 1.0 1.0 ... 1.0
    max_demand_timesteps                  (carriers) datetime64[ns] 2005-07-01 ... 2005-07-01T07:00:00
Attributes:
    calliope_version:    0.6.4-dev
    applied_overrides:   milp
    scenario:            milp
    defaults:            available_area: null\ncarrier_ratios: {}\ncharge_rat...
    allow_operate_mode:  1
In [4]:
# Individual data variables can be accessed easily, `to_pandas()` reformats the data to look nicer
# Here we look at one of the MILP overrides that we have added, the fixed `purchase` cost
model.inputs.cost_purchase.to_pandas().dropna(axis=1)
Out[4]:
loc_techs_investment_cost X3::boiler X1::chp X2::boiler
costs
monetary 2000.0 40000.0 2000.0
In [5]:
# Solve the model. Results are loaded into `model.results`. 
# By including logging (see package importing), we can see the timing of parts of the run, as well as the solver's log
model.run()
[2019-05-27 15:05:01] INFO: Backend: starting model run
[2019-05-27 15:05:02] INFO: Backend: model generated. Time since start: 0:00:02.023389
[2019-05-27 15:05:02] INFO: Backend: sending model to solver
[2019-05-27 15:05:04] INFO: Backend: solver finished running. Time since start: 0:00:04.462383
[2019-05-27 15:05:04] INFO: Backend: loaded results
[2019-05-27 15:05:04] INFO: Backend: generated solution array. Time since start: 0:00:04.551883
[2019-05-27 15:05:04] INFO: Postprocessing: started
[2019-05-27 15:05:05] INFO: Postprocessing: zero threshold of 1e-10 not required
[2019-05-27 15:05:05] INFO: Postprocessing: Model was feasible, deleting unmet_demand variable
[2019-05-27 15:05:05] INFO: Postprocessing: ended. Time since start: 0:00:04.944535
In [6]:
# Model results are held in the same structure as model inputs. 
# The results consist of the optimal values for all decision variables, including capacities and carrier flow
# There are also results, like system capacity factor and levelised costs, which are calculated in postprocessing
# before being added to the results Dataset

model.results
Out[6]:
<xarray.Dataset>
Dimensions:                          (carriers: 3, costs: 1, loc_tech_carriers_con: 19, loc_tech_carriers_export: 4, loc_tech_carriers_prod: 21, loc_techs: 26, loc_techs_area: 3, loc_techs_asynchronous_prod_con: 6, loc_techs_cost: 20, loc_techs_investment_cost: 20, loc_techs_milp: 1, loc_techs_om_cost: 9, loc_techs_purchase: 2, loc_techs_supply_plus: 3, techs: 9, timesteps: 48)
Coordinates:
  * loc_techs_milp                   (loc_techs_milp) object 'X1::chp'
  * loc_techs_asynchronous_prod_con  (loc_techs_asynchronous_prod_con) object 'N1::heat_pipes:X1' ... 'N1::heat_pipes:X2'
  * costs                            (costs) object 'monetary'
  * loc_techs_investment_cost        (loc_techs_investment_cost) object 'X1::supply_gas' ... 'N1::heat_pipes:X2'
  * loc_tech_carriers_con            (loc_tech_carriers_con) object 'X1::power_lines:X2::electricity' ... 'N1::heat_pipes:X1::heat'
  * loc_tech_carriers_prod           (loc_tech_carriers_prod) object 'X3::pv::electricity' ... 'X3::heat_pipes:N1::heat'
  * loc_techs                        (loc_techs) object 'X1::supply_gas' ... 'X2::pv'
  * loc_techs_area                   (loc_techs_area) object 'X3::pv' ... 'X1::pv'
  * loc_techs_supply_plus            (loc_techs_supply_plus) object 'X3::pv' ... 'X1::pv'
  * timesteps                        (timesteps) datetime64[ns] 2005-07-01 ... 2005-07-02T23:00:00
  * loc_tech_carriers_export         (loc_tech_carriers_export) object 'X3::pv::electricity' ... 'X1::chp::electricity'
  * loc_techs_om_cost                (loc_techs_om_cost) object 'X1::supply_gas' ... 'X2::boiler'
  * loc_techs_purchase               (loc_techs_purchase) object 'X3::boiler' 'X2::boiler'
  * carriers                         (carriers) object 'gas' ... 'heat'
  * loc_techs_cost                   (loc_techs_cost) object 'X1::supply_gas' ... 'X2::heat_pipes:N1'
  * techs                            (techs) object 'chp' ... 'demand_electricity'
Data variables:
    energy_cap                       (loc_techs) float64 740.7 238.7 ... 142.0
    carrier_prod                     (loc_tech_carriers_prod, timesteps) float64 0.0 ... 0.8603
    carrier_con                      (loc_tech_carriers_con, timesteps) float64 -96.48 ... 0.0
    cost                             (costs, loc_techs_cost) float64 529.0 ... 0.05836
    resource_area                    (loc_techs_area) float64 350.0 994.1 0.0
    resource_con                     (loc_techs_supply_plus, timesteps) float64 0.0 ... 0.0
    resource_cap                     (loc_techs_supply_plus) float64 38.95 ... 0.0
    carrier_export                   (loc_tech_carriers_export, timesteps) float64 0.0 ... 0.0
    cost_var                         (costs, loc_techs_om_cost, timesteps) float64 5.832 ... 0.004252
    cost_investment                  (costs, loc_techs_investment_cost) float64 0.4472 ... 0.05836
    purchased                        (loc_techs_purchase) float64 0.0 1.0
    units                            (loc_techs_milp) float64 1.0
    operating_units                  (loc_techs_milp, timesteps) float64 1.0 ... 1.0
    prod_con_switch                  (loc_techs_asynchronous_prod_con, timesteps) float64 1.0 ... 0.0
    capacity_factor                  (loc_tech_carriers_prod, timesteps) float64 0.0 ... 0.08272
    systemwide_capacity_factor       (carriers, techs) float64 0.0 0.0 ... nan
    systemwide_levelised_cost        (carriers, costs, techs) float64 inf ... nan
    total_levelised_cost             (carriers, costs) float64 0.03962 ... 0.1099
Attributes:
    calliope_version:       0.6.4-dev
    applied_overrides:      milp
    scenario:               milp
    defaults:               available_area: null\ncarrier_ratios: {}\ncharge_...
    allow_operate_mode:     1
    model_config:           calliope_version: 0.6.4\nname: Urban-scale exampl...
    run_config:             backend: pyomo\nbigM: 1000000.0\ncyclic_storage: ...
    termination_condition:  optimal
    solution_time:          3.038068
    time_finished:          2019-05-27 15:05:04
In [7]:
# We can sum operating units of CHP over all locations and turn the result into a pandas DataFrame
df_units = model.get_formatted_array('operating_units').sum('locs').to_pandas().T

#The information about the dataframe tells us about the amount of data it holds in the index and in each column
df_units.info()
<class 'pandas.core.frame.DataFrame'>
DatetimeIndex: 48 entries, 2005-07-01 00:00:00 to 2005-07-02 23:00:00
Data columns (total 1 columns):
chp    48 non-null float64
dtypes: float64(1)
memory usage: 768.0 bytes
In [8]:
# Using .head() to see the first few rows of operating units

df_units.head()
Out[8]:
techs chp
timesteps
2005-07-01 00:00:00 1.0
2005-07-01 01:00:00 1.0
2005-07-01 02:00:00 1.0
2005-07-01 03:00:00 1.0
2005-07-01 04:00:00 1.0
In [9]:
# We can plot this by using the timeseries plotting functionality.
# The top-left dropdown gives us the chance to scroll through other timeseries data too.

model.plot.timeseries()
In [10]:
# plot.capacities gives a graphical view of the non-timeseries variables, both input and output

# Note, because we fix unit size, CHP now has a maximum capacity of 300kW, 
# compared to 260kW in the non-MILP case

model.plot.capacity()

See the Calliope documentation for more details on setting up and running a Calliope model.

ript src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.1.10/require.min.js">