| |||||||||||
| |||||||||||
Python notebooks Description Python notebooks available in the GitHub repository python-notebooks.
Source Files By clicking on a file name, a preview is opened at the bottom of this page. unitcommitment_indicators.ipynb
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# **Solving an electricity generation problem using indicator constraints**"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"***unitcommitment_indicators.ipynb***\n",
"\n",
"This example shows how to model and solve an electricity generation problem typically found in power markets (see [Garver (1963)](https://ieeexplore.ieee.org/document/4501405)), showcasing the use of indicator constraints to model change state constraints when generators are turned on/off.\n",
"\n",
"© Copyright 2025 Fair Isaac Corporation\n",
"\n",
"Licensed under the Apache License, Version 2.0 (the \"License\"); you may not use this file except in compliance with the License. You may obtain a copy of the License at http://www.apache.org/licenses/LICENSE-2.0.\n",
" \n",
"Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an \"AS IS\" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.\n",
"\n",
"This example uses FICO® Xpress software. By running it, you agree to the Community License terms of the [Xpress Shrinkwrap License Agreement](https://community.fico.com/s/contentdocument/06980000002h0i5AAA) with respect to the FICO® Xpress software. See the [licensing options](https://www.fico.com/en/fico-xpress-trial-and-licensing-options) overview for additional details and information about obtaining a paid license."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Install the xpress package\n",
"%pip install -q xpress"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Problem description"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Four types of power generators are available to meet daily electricity demand and a security reserve of at least 20% above the estimated demand. Each type of generator has a set minimum and maximum power output.\n",
"A generator can only be started or stopped at the beginning of a time period. The objective is to determine which generators should be used, and at which power level, in each period so that total daily cost is minimized.\n",
"\n",
"The length and anticipated daily electricity demand for each of the seven planning periods of uneven length are as follows:\n",
"\n",
"| Time Period | Length (hours) | Demand (MW) |\n",
"| --- | --- | --- |\n",
"| 00h to 06h | 6 | 12000 |\n",
"| 06h to 09h | 3 | 32000 |\n",
"| 09h to 12h | 3 | 25000 |\n",
"| 12h to 14h | 2 | 36000 |\n",
"| 14h to 18h | 4 | 25000 |\n",
"| 18h to 22h | 4 | 30000 |\n",
"| 22h to 00h | 2 | 18000 |\n",
"\n",
"The characteristics of power generators are described in the table below, per generator type:\n",
"\n",
"| Type | No. units | Min. power | Max. power | Start-up cost | Hourly cost at min. output | Hourly cost per MW above min. output|\n",
"| --- | --- | --- | --- | --- | --- | --- |\n",
"| A | 10 | 750 | 1750 | 5000 | 2250 | 2.7 |\n",
"| B | 4 | 1000 | 1500 | 1600 | 1800 | 2.2 |\n",
"| C | 8 | 1200 | 2000 | 2400 | 3750 | 1.8 |\n",
"| D | 3 | 1800 | 3500 | 1200 | 4800 | 3.8 |"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Model formulation"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Three arrays of binary variables are required to determine when to 'start' and 'stop' generators, and decide which ones are set to 'work' in each time period. Another set of real variables are introduced to represent the additonal energy production ('padd') of each working unit type above its minimum output level. Variable 'start' is the difference between 'work' in a time period and 'work' in any preceding period, and 'stop' is the difference between 'work' in any previous time period and the current one.\n",
"\n",
"This unit commitment problem can be mathematically formulated as follows:\n",
"$$\\min \\sum_{u \\in UNITS} \\sum_{t \\in PERIODS} COST_{S_u}^{\\rm start} \\cdot start_{u,t} + LEN_t \\cdot (COST_{S_u}^{\\rm min} \\cdot work_{u,t} + COST_{S_u}^{\\rm add} \\cdot padd_{u,t}) \n",
"+ \\sum_{u \\in UNITS} \\sum_{t \\in PERIODS} stop_{u,t} \\cdot PEN $$\n",
"\n",
"Subject to:\n",
"\n",
"$$\n",
"\\begin{array}{llll}\n",
"& \\hbox{If generator starts in period:} \\\\\n",
"& \\qquad start_{u,t} \\geq work_{u,t} - work_{u,n} , \\qquad \\forall u \\in UNITS, \\forall t \\in PERIODS, n = (NT+t-1) \\hbox{ \\% } NT \\\\\n",
"& \\qquad start_{u,t} \\leq work_{u,t} , \\qquad \\forall u \\in UNITS, \\forall t \\in PERIODS \\\\\n",
"& \\hbox{If generator stops in period:} \\\\\n",
"& \\qquad stop_{u,t} \\geq work_{u,n} - work_{u,t} , \\qquad \\forall u \\in UNITS, \\forall t \\in PERIODS, n = (NT+t-1) \\hbox{ \\% } NT \\\\\n",
"& \\qquad stop_{u,t} \\leq 1 - work_{u,t}, \\qquad \\forall u \\in UNITS, \\forall t \\in PERIODS \\\\\n",
"& \\hbox{Limit on power production above minimum level:} \\\\\n",
"& \\qquad padd_{u,t}\\leq (P_{S_u}^{\\rm max}-P_{S_u}^{\\rm min}) \\cdot work_{u,t}, \\qquad \\forall u \\in UNITS, \\forall t \\in PERIODS \\\\\n",
"& \\hbox{Satisfy daily electricity demand:} \\\\\n",
"& \\qquad \\sum_{u \\in UNITS} P_{S_u}^{\\rm min} \\cdot work_{u,t} + padd_{u,t} \\geq D_t, \\qquad \\forall t \\in PERIODS \\\\\n",
"& \\hbox{Ensure security reserve of 20\\%:} \\\\\n",
"& \\qquad \\sum_{u \\in UNITS} P_{S_u}^{\\rm max} \\cdot work_{u,t} \\geq 1.20 \\cdot D_t, \\qquad \\forall t \\in PERIODS \\\\\n",
"\\end{array}\n",
"$$\n",
"\n",
"Where:\n",
"* $NT$ = number of time periods ($NT = |PERIODS|$)\n",
"* $PERIODS$ = set of time periods ($t \\in PERIODS$)\n",
"* $TYPES$ = set of generator types ($s \\in TYPES$)\n",
"* $UNITS$ = set of unit generators to be committed ($u \\in UNITS$)\n",
"* $AVAIL_s$ = available number of type $s$ generators ($|UNITS| = \\sum_{s \\in TYPES} AVAIL_s$)\n",
"* $LEN_t$ = length of time period $t$\n",
"* $DEM_t$ = demand in time period $t$\n",
"* $TYPE_u$ = type of generator $u$\n",
"* $COST_{s}^{\\rm start}$ = start-up cost of generator type $s$\n",
"* $COST_{s}^{\\rm min}$ = hourly cost of operating generator $s$ at minimum output\n",
"* $COST_{s}^{\\rm add}$ = cost/hour/MW of prod. above min. level\n",
"* $POW_{s}^{\\rm min}$,$POW_{s}^{\\rm max}$ = minimum and maximum power output of a generator type $s$\n",
"* $PEN$ = penalty to ensure that stop variables are 0 when not needed\n",
"\n",
"The decision variables are:\n",
"* $start_{u,t}$ = 1 if generator $u$ is started in period $t$, 0 otherwise\n",
"* $stop_{u,t}$ = 1 if generator $u$ is stopped in period $t$, 0 otherwise\n",
"* $work_{u,t}$ = 1 if generator $u$ is working during period $t$, 0 otherwise\n",
"* $padd_{u,t}$ = production level (MW) of generator $u$ above the minimum level $P_{S_u}^{\\rm min}$ in period $t$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Implementation of the basic model"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The code cell below demonstrates the implementation of the above model formulation and prints the results using the Xpress Python interface."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import xpress as xp\n",
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"\n",
"# Time periods\n",
"LEN = [6, 3, 3, 2, 4, 4, 2]\n",
"DEM = [12000, 32000, 25000, 36000, 25000, 30000, 18000]\n",
"\n",
"# Power plants\n",
"PMIN = [750, 1000, 1200, 1800] # minimum output (MW) per generator type\n",
"PMAX = [1750, 1500, 2000, 3500] # maximum output (MW) per generator type\n",
"CSTART = [5000, 1600, 2400, 1200] # start-up cost per generator type\n",
"CMIN = [2250, 1800, 3750, 4800] # hourly cost of operating generator type at minimum output (see PMIN)\n",
"CADD = [2.7, 2.2, 1.8, 3.8] # cost/hour/MW of prod. above min. level per generator type\n",
"AVAIL = [10, 4, 8, 3] # number of units per type\n",
"\n",
"NT = 7 # number of time periods\n",
"PERIODS = range(NT) # set of time periods\n",
"TYPES = range(4) # power generator types\n",
"UNITS = range(sum(AVAIL)) # power generation units\n",
"TYPE = [i for i in TYPES for p in range(AVAIL[i])] # associating units with types\n",
"PEN = 0.1 # penalty associated with stopping units\n",
"\n",
"# Create problem\n",
"p = xp.problem(\"Unit commitment\")\n",
"\n",
"# Create decision variables\n",
"start = p.addVariables(UNITS, PERIODS, vartype=xp.binary, name='start')\n",
"stop = p.addVariables(UNITS, PERIODS, vartype=xp.binary, name='stop')\n",
"work = p.addVariables(UNITS, PERIODS, vartype=xp.binary, name='work')\n",
"padd = p.addVariables(UNITS, PERIODS, name='padd')\n",
"\n",
"# If generator starts in period\n",
"p.addConstraint(start[u,t] >= work[u,t] - work[u,(NT+t-1) % NT] for u in UNITS for t in PERIODS)\n",
"p.addConstraint(start[u,t] <= work[u,t] for u in UNITS for t in PERIODS)\n",
"\n",
"# If generator stops before period\n",
"p.addConstraint(stop[u,t] >= work[u,(NT+t-1) % NT] - work[u,t] for u in UNITS for t in PERIODS)\n",
"p.addConstraint(stop[u,t] <= 1 - work[u,t] for u in UNITS for t in PERIODS)\n",
"\n",
"# Limit on power production above minimum level\n",
"p.addConstraint(padd[u,t] <= (PMAX[TYPE[u]]-PMIN[TYPE[u]]) * work[u,t] for u in UNITS for t in PERIODS)\n",
"\n",
"# Satisfy demands\n",
"p.addConstraint(xp.Sum(PMIN[TYPE[u]]*work[u,t] + padd[u,t] for u in UNITS) >= DEM[t] for t in PERIODS)\n",
"\n",
"# Security reserve of 20%\n",
"p.addConstraint(xp.Sum(PMAX[TYPE[u]]*work[u,t] for u in UNITS) >= 1.20 * DEM[t] for t in PERIODS)\n",
"\n",
"# Create and add the oObjective function of the problem (hint: compute 'daily cost' and 'penalty' separately)\n",
"Cost = xp.Sum(CSTART[TYPE[u]] * start[u,t] +\n",
" LEN[t] * (CMIN[TYPE[u]] * work[u,t] + CADD[TYPE[u]] * padd[u,t]) for u in UNITS for t in PERIODS)\n",
"\n",
"Penalty = PEN * xp.Sum(stop[u,t] for u in UNITS for t in PERIODS)\n",
"p.setObjective(Cost + Penalty)\n",
"\n",
"# Optimize the problem and print the daily cost, penalty and total objective value\n",
"p.optimize()\n",
"print(\"Daily cost:\", round(p.getSolution(Cost),2))\n",
"print(\"Penalty:\", round(p.getSolution(Penalty),2))\n",
"print(\"Objective value:\", round(p.attributes.objval,2))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next we define a function to print the solution in a user-friendly way."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Display solution\n",
"def print_sol(prob):\n",
" ct = 0\n",
" print(f\"Time period \", end=\"\")\n",
" for t in PERIODS:\n",
" print('{:5d}-{:2}'.format(ct,ct+LEN[t]), end=\"\")\n",
" ct += LEN[t]\n",
"\n",
" for u in UNITS:\n",
" print(\"\\n\",f\"\\nUnit {u+1} Working \", end=\"\")\n",
" for t in PERIODS:\n",
" print(\" off\" if prob.getSolution(work[u,t]) == 0 else \" on\", end=\"\")\n",
" print(\"\\n Status change\", end=\"\")\n",
" for t in PERIODS:\n",
" if prob.getSolution(stop[u,t]) > 0.5: print(\" stop\", end=\"\")\n",
" else:\n",
" if prob.getSolution(start[u,t]) > 0.5: print(\" start\", end=\"\")\n",
" else: print(\" -\", end=\"\")\n",
" print(\"\\n Total output \", end=\"\")\n",
" for t in PERIODS:\n",
" print('{:8d}'.format(round(prob.getSolution(PMIN[TYPE[u]]*work[u,t] + padd[u,t]))), end=\"\")\n",
" print(\"\\n of which add.\", end=\"\")\n",
" for t in PERIODS:\n",
" print('{:8d}'.format(round(prob.getSolution(padd[u,t]))), end=\"\")\n",
"\n",
"print_sol(p)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now we use the *matplotlib* pacakge to define a function that plots a bar chart with stacked bars of hourly power output per generator type per time period, the total power output and the total demand in each planning period. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def plot_sol(prob):\n",
"\n",
" # Get the total power output per generator type per planning period\n",
" outputs = [[0 for _ in PERIODS] for _ in TYPES]\n",
" reserve = [0 for _ in PERIODS]\n",
" for i in TYPES:\n",
" for t in PERIODS:\n",
" for u in UNITS:\n",
" if TYPE[u] == i and prob.getSolution(work[u,t]) > 0.5:\n",
" power_output = prob.getSolution(PMIN[TYPE[u]]*work[u,t] + padd[u,t])\n",
" outputs[i][t] += power_output\n",
" reserve[t] += PMAX[TYPE[u]] - power_output\n",
" \n",
" # Labels for the unit types\n",
" labels = [f\"Unit type {i + 1}\" for i in TYPES]\n",
"\n",
" # Create a stacked bar chart\n",
" fig, ax = plt.subplots()\n",
"\n",
" # Plot output for each unit type per planning period\n",
" for i in range(len(outputs)):\n",
" ax.bar(range(1,NT+1), outputs[i], label=labels[i], bottom=np.sum(outputs[:i], axis=0))\n",
"\n",
" # Plot the total reserve per planning period\n",
" ax.bar(range(1,NT+1), reserve, label='Total reserve', bottom=np.sum(outputs, axis=0), fill=False, edgecolor='black', linestyle='--')\n",
"\n",
" # Plot the demand data as a line\n",
" ax.plot(range(1,NT+1), DEM, color='black', marker='o', label='Demand (MW)')\n",
"\n",
" # Label and size the axes\n",
" ax.set_xlabel('Planning period', fontsize=13)\n",
" ax.set_ylabel('Total output (MW)', fontsize=13)\n",
"\n",
" # Add a legend\n",
" ax.legend(loc='upper left', bbox_to_anchor=(1,1), fontsize=13)\n",
"\n",
" # Show the plot\n",
" plt.show()\n",
"\n",
"plot_sol(p)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Adding indicator constraints"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In real-world applications, generators must often remain for some time in a certain state after they have been switched to that state, i.e. if a generator is turned ON/OFF, it must remain in that state for at least $ON_s^{\\rm min}$/$OFF_s^{\\rm min}$ periods, respectively. Such state change constraints can be formulated with the help of so-called **indicator constraints** as their enforcement depends on the value of a binary variable, called the 'indicator'.\n",
"\n",
"The indicator constraints for this problem can be formulated as follows:\n",
"\n",
"$$\n",
"\\begin{align*}\n",
"& \\hbox{Can only switch OFF at least $ON_s^{\\rm min}$ periods after it has been turned ON:} \\\\\n",
"& \\qquad \\sum_{j=t+1}^{t+ON_s^{\\rm min}-1} stop_{u,(j \\hbox{ \\% } NT)} \\leq 0, \\qquad \\forall u \\in \\mathcal{U}, \\forall t \\in \\mathcal{T}: start_{u,t} = 1 \\\\\n",
"& \\hbox{Can only switch ON at least $OFF_s^{\\rm min}$ periods after it has been turned OFF:} \\\\\n",
"& \\qquad \\sum_{j=t+1}^{t+OFF_s^{\\rm min}-1} start_{u,(j \\hbox{ \\% } NT)} \\leq 0, \\qquad \\forall u \\in \\mathcal{U}, \\forall t \\in \\mathcal{T}: stop_{u,t} = 1 \\\\\n",
"\\end{align*}\n",
"$$\n",
"\n",
"where $ON_s^{\\rm min}$,$OFF_s^{\\rm min}$ = minimum time intervals a generator type $s$ must be ON/OFF once it has switched to that state, respectively.\n",
"\n",
"Indicator constraints can conveniently be added by using the [problem.addIndicator()](https://www.fico.com/fico-xpress-optimization/docs/latest/solver/optimizer/python/HTML/chModeling.html?scroll=secModelingIndicator) method (note the \"==\" as the symbol for the equality on the indicator).\n",
"\n",
"The code cell below adds the indicator constraints to the model, triggers the optimization again and prints the solution. Comparing the outcomes with those from the previous code cells, we can verify that there is an extra cost by introducing the new constraints, as now units are not allowed to be turned ON and then OFF within three time periods. Also note that the indicators make the problem more difficult (hence slower) to solve."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Minimum time\n",
"ONMIN = [3, 3, 3, 3] # minimum time intervals a generator type $s$ must be ON once switched to that state\n",
"DWMIN = [3, 3, 3, 3] # minimum time intervals a generator type $s$ must be OFF once switched to that state\n",
"\n",
"# Indicator constraints\n",
"for u in UNITS:\n",
" for t in PERIODS:\n",
" # Can only switch off at least ONMIN periods later\n",
" p.addIndicator(start[u,t] == 1, xp.Sum(stop[u,j % NT] for j in range(t+1,t+ONMIN[TYPE[u]])) <= 0)\n",
" # Can only switch on at least DWMIN periods later\n",
" p.addIndicator(stop[u,t] == 1, xp.Sum(start[u,j % NT] for j in range(t+1,t+DWMIN[TYPE[u]])) <= 0)\n",
"\n",
"# Re-optimize the problem and print the daily cost, penalty and total objective value\n",
"p.controls.outputlog = 0\n",
"p.optimize()\n",
"print(\"Daily cost:\", round(p.getSolution(Cost),2))\n",
"print(\"Penalty:\", round(p.getSolution(Penalty),2))\n",
"print(\"Objective value:\", round(p.attributes.objval,2))\n",
"\n",
"print_sol(p)\n",
"plot_sol(p)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.12.9"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
| |||||||||||
| © Copyright 2025 Fair Isaac Corporation. |
