Device scheduler: mixed-integer linear model

Introduction

This generic device scheduler is able to handle a site with multiple devices, with various types of constraints on the site level and on the device level, and with multiple market commitments on the site level.

A typical example is a house with many devices. The commitments are assumed to be with regard to the flow of energy to the device (positive for consumption, negative for production). In practice, this generic scheduler is used in the StorageScheduler to schedule a storage device.

The solver minimizes the costs of deviating from the commitments. For a more detailed explanation of commitments in FlexMeasures, see Commitments.

The model is a mixed-integer linear program: binary variables model the sign of a device’s power, the direction of a commitment deviation (only when the cost curve is non-convex), and the choice of operation mode for devices with power bands. Without any of these, the model reduces to a plain linear program.

Note

The model below is built by device_scheduler(), using Pyomo. A second backend, highspy_optimization, builds the same model directly with the HiGHS Python API, which is much faster to construct. Both are fed by the same input preparation step, prepare_scheduling_problem(), and the symbols below refer to its output.

Notation

Indexes

Symbol

Variable in the Code

Description

\(c\)

c

Sub-commitments: one per commitment group, and one per deviation direction where a group spans several time steps (see Commitments and sub-commitments).

\(d\)

d

Devices, for example, a battery or a load.

\(j\)

j

0-indexed time dimension.

\(g\)

cg, cjg

Device groups within a device-scoped commitment \(c\).

\(s\)

sg

Stock groups: sets of devices that share one stock (e.g. one state-of-charge sensor).

\(e\)

eg

EMS constraint groups: sets of devices sharing one site-level capacity constraint (one per commodity).

\(b\)

db

Power bands (S2 operation modes) of a banded device \(d\).

\(k\)

coupling_group_range

Coupling groups: sets of devices whose flows are hard-coupled in fixed proportions.

\(n\)

balance_group_range

Balance groups: internal commodity nodes (e.g. a heat or steam network) whose flows must net to zero.

Note

The time index \(j\) has two interpretations: a time period or an instantaneous moment at the end of time period \(j\). For example, \(j\) in flow constraints correspond to time periods, whereas \(j\) used in a stock constraint refers to the end of time period \(j\).

Sets and mappings

The constraints below sum over the devices that a group covers. These sets name those memberships.

Symbol

Variable in the Code

Description

\(s(d)\)

device_to_group

The stock group that device \(d\) belongs to (its primary one, if it participates in several).

\(S(s)\)

group_to_devices

The devices that share stock group \(s\). Written as \(d \in s\) below.

\(E(e)\)

ems_constraint_device_groups

The devices covered by EMS constraint group \(e\).

\(\mathcal{D}(c)\)

commodity_devices

The devices covered by an EMS-level commitment \(c\): all devices, or those of the commitment’s commodity.

\(G(c,g)\)

device_group_lookup

The devices in device group \(g\) of a device-scoped commitment \(c\).

\(N(n)\)

balance_group_specs

The devices attached to internal commodity node \(n\).

Parameters

Symbol

Variable in the Code

Description

\(Price_{up}(c)\)

up_price

Price of incurring an upwards deviation in sub-commitment \(c\).

\(Price_{down}(c)\)

down_price

Price of incurring a downwards deviation in sub-commitment \(c\).

\(Commitment(c,j)\)

commitment_quantity

Committed quantity of sub-commitment \(c\) for time period \(j\) (a flow, or a stock for stock commitments).

\(\eta_{up}(d,j)\)

device_derivative_up_efficiency

Upwards conversion efficiency (stock increase : flow in).

\(\eta_{down}(d,j)\)

device_derivative_down_efficiency

Downwards conversion efficiency (flow out : stock decrease).

\(Stock_{min}(d,j)\)

device_min

Minimum stock of device \(d\) at the end of time period \(j\), relative to its initial stock.

\(Stock_{max}(d,j)\)

device_max

Maximum stock of device \(d\) at the end of time period \(j\), relative to its initial stock.

\(Stock_0(d)\)

initial_stock

Initial stock of device \(d\), shared by all devices in its stock group.

\(\epsilon(d,j)\)

device_efficiency

Storage efficiency (stock losses), shared by all devices in a stock group.

\(P_{max}(d,j)\)

device_derivative_max

Maximum flow of device \(d\) during time period \(j\).

\(P_{min}(d,j)\)

device_derivative_min

Minimum flow of device \(d\) during time period \(j\).

\(P^{ems}_{min}(e,j)\)

ems_derivative_min

Minimum aggregated flow of EMS constraint group \(e\) during time period \(j\).

\(P^{ems}_{max}(e,j)\)

ems_derivative_max

Maximum aggregated flow of EMS constraint group \(e\) during time period \(j\).

\(D(d,j)\)

stock_delta

Explicit stock gain or loss of device \(d\) during time period \(j\).

\(\gamma(k,d)\)

coupling_device_specs

Fixed proportion of device \(d\) within coupling group \(k\). Positive for inputs (consuming), negative for outputs (producing).

\(B_{min}(d,b)\), \(B_{max}(d,b)\)

band_lookup

Lower and upper flow bound of power band \(b\) of device \(d\).

\(M_d\)

Md

Big-M bounding device power: the largest absolute device flow limit (at least 1 MW).

\(M_c\)

Mc

Big-M bounding commitment deviations: the summed absolute device flow limits (at least 1 MW).

Variables

Symbol

Variable in the Code

Description

\(\Delta_{up}(c)\)

commitment_upwards_deviation

Upwards deviation from sub-commitment \(c\) (\(\geq 0\)). One variable per sub-commitment, not per time step.

\(\Delta_{down}(c)\)

commitment_downwards_deviation

Downwards deviation from sub-commitment \(c\) (\(\leq 0\)).

\(\sigma_c(c)\)

commitment_sign

Binary. Upwards deviation allowed if \(\sigma_c(c)=1\), downwards deviation otherwise.

\(P_{up}(d,j)\)

device_power_up

Upwards (consuming) power of device \(d\) during time period \(j\) (\(\geq 0\)).

\(P_{down}(d,j)\)

device_power_down

Downwards (producing) power of device \(d\) during time period \(j\) (\(\leq 0\)).

\(P^{ems}(d,j)\)

ems_power

Net flow of device \(d\) during time period \(j\).

\(\sigma(d,j)\)

device_power_sign

Binary. Upwards power activation if \(\sigma(d,j)=1\), downwards power activation otherwise.

\(Stock(s,j)\)

group_stock

Stock of stock group \(s\) at the end of time period \(j\).

\(\Delta Stock(s,j)\)

n/a

Auxiliary symbol used below: the stock change of stock group \(s\) during time period \(j\), before losses.

\(\alpha(k,j)\)

coupling_alpha

Common normalised flow level of coupling group \(k\) during time period \(j\).

\(y(d,b,j)\)

device_band

Binary. Device \(d\) operates in power band \(b\) during time period \(j\).

Commitments and sub-commitments

A commitment may declare a group per time step, which defines the set of time steps within which deviations are accounted for together (for example: only the highest breach per calendar month is penalised). Before the model is built, every commitment is split into sub-commitments: one per group, and — where a group spans multiple time steps — one per deviation direction. The index \(c\) therefore runs over sub-commitments, and each carries a single pair of deviation prices.

This is why the deviation variables \(\Delta_{up}(c)\) and \(\Delta_{down}(c)\) are not indexed by time: one deviation is paid for per sub-commitment, however many time steps it spans. For the common case of a per-time-step commitment (such as an energy tariff), each time step forms its own group, and the two coincide.

After solving, the costs of the sub-commitments are aggregated back to the original commitments (model.commitment_costs) and per commodity (model.commodity_costs).

Cost function

The cost function quantifies the total cost of upwards and downwards deviations from the different sub-commitments.

(1)\[\min \sum_{c} [\Delta_{up}(c) \cdot Price_{up}(c) + \Delta_{down}(c) \cdot Price_{down}(c)]\]

Note that \(\Delta_{down}(c) \leq 0\), so a positive downwards deviation price yields a negative contribution, and vice versa. This mirrors the sign convention described under Commitments.

State dynamics

Stock is tracked per stock group \(s\): a set of devices that share one stock, such as several converters filling one buffer. Devices that are not part of a declared stock group form a stock group of their own. Because the stock is a property of the group, all its devices must declare the same storage efficiency and the same initial stock.

To simplify the description of the model, the auxiliary variable \(\Delta Stock(s,j)\) is introduced in the documentation. It represents the change of \(Stock(s,j)\), taking into account conversion efficiencies but not considering the storage losses.

(2)\[ \Delta Stock(s,j) = \sum_{d \in s} \left[ \frac{P_{down}(d,j)}{\eta_{down}(d,j)} + P_{up}(d,j) \cdot \eta_{up}(d,j) + D(d,j) \right]\]

The stock is then defined recursively, rather than as a running sum over all preceding time steps, which keeps the number of nonzeros in the model linear (rather than quadratic) in the scheduling horizon:

(3)\[ Stock(s, j) = a(s,j) \cdot Stock(s, j-1) + b(s,j) \cdot \Delta Stock(s,j)\]

with \(Stock(s, -1) = Stock_0(d)\) for any device \(d \in s\) (the group’s devices share one initial stock), and with the loss coefficients

(4)\[\begin{split} (a(s,j), b(s,j)) = \begin{cases} (1, 1) & \text{if } \epsilon(d,j) = 1 \\ \left(\epsilon(d,j), \frac{\epsilon(d,j) - 1}{\log(\epsilon(d,j))}\right) & \text{otherwise} \end{cases}\end{split}\]

for any device \(d \in s\) (they all share the group’s storage efficiency).

Note

This is the linear treatment of storage losses: the stock is assumed to change at a constant rate, while losses decay exponentially, within each time step. The scheduler models losses this way exclusively. The alternative treatments (perfect, left and right) still exist in apply_stock_changes_and_losses(), which reconstructs a stock series from a known power series.

Constraints

Device bounds

Stock bounds are expressed relative to the device’s initial stock, and refer to the stock of the group the device belongs to:

(5)\[ Stock_{min}(d,j) \leq Stock(s(d), j) - Stock_0(d) \leq Stock_{max}(d,j)\]
(6)\[ P_{min}(d,j) \leq P_{up}(d,j) + P_{down}(d,j)\leq P_{max}(d,j)\]
(7)\[ min(P_{min}(d,j),0) \leq P_{down}(d,j)\leq 0\]
(8)\[ 0 \leq P_{up}(d,j)\leq max(P_{max}(d,j),0)\]

Upwards/Downwards activation selection

Avoid simultaneous upwards and downwards activation during the same time period.

(9)\[ P_{up}(d,j) \leq M_d \cdot \sigma(d,j)\]
(10)\[ -P_{down}(d,j) \leq M_d \cdot (1-\sigma(d,j))\]

The same trick prevents a sub-commitment from deviating in both directions at once:

(11)\[ \Delta_{up}(c) \leq M_c \cdot \sigma_c(c)\]
(12)\[ -\Delta_{down}(c) \leq M_c \cdot (1-\sigma_c(c))\]

These two constraints are only added when the summed deviation prices do not describe a convex cost curve. For a convex cost curve, deviating in both directions is never optimal anyway, so dropping them leaves the problem a pure LP. The direct HiGHS backend then omits \(\sigma_c(c)\) as well, while the Pyomo backend still declares the variable, leaving it unreferenced.

Grid constraints

(13)\[P^{ems}(d,j) = P_{up}(d,j) + P_{down}(d,j)\]

Site-level capacity is enforced per EMS constraint group \(e\), over the devices \(E(e)\) it covers. The StorageScheduler uses one group per commodity, so each commodity gets its own site-level capacity constraint. A single group covering all devices is the default (and the historical behaviour).

(14)\[ P^{ems}_{min}(e,j) \leq \sum_{d \in E(e)} P^{ems}(d,j) \leq P^{ems}_{max}(e,j)\]

EMS-level commitments

Commitments that do not name a device apply to the site as a whole. Writing \(\mathcal{D}(c)\) for the devices such a commitment covers — all devices, or, if the commitment names a commodity, the devices of that commodity — and

(15)\[\Xi(c,j) = Commitment(c,j) + \Delta_{down}(c) + \Delta_{up}(c) - \sum_{d \in \mathcal{D}(c)} P^{ems}(d,j)\]

the constraint is bounded on the side(s) for which the sub-commitment carries a deviation price:

Sub-commitment prices

Constraint

Meaning

Both prices given

\(0 \leq \Xi \leq 0\)

The commitment is met exactly, or paid for in both directions.

Only an upwards deviation price

\(\Xi \geq 0\)

Flow above the committed quantity is a breach; staying below is free.

Only a downwards deviation price

\(\Xi \leq 0\)

Flow below the committed quantity is a breach; staying above is free.

Device-scoped commitments

A commitment may instead name devices, optionally organised in device groups \(g\) (the device_group column; by default each device forms its own group). It is then bound once per device group, against the aggregate of that group’s devices \(G(c,g)\). For a flow commitment, the aggregate is a flow:

(16)\[\Xi(c,j,g) = Commitment(c,j) + \Delta_{down}(c) + \Delta_{up}(c) - \sum_{d \in G(c,g)} P^{ems}(d,j)\]

and for a stock commitment, it is a stock change since the start of the schedule:

(17)\[\Xi(c,j,g) = Commitment(c,j) + \Delta_{down}(c) + \Delta_{up}(c) - \sum_{d \in G(c,g)} [Stock(s(d), j) - Stock_0(d)]\]

Both are bounded exactly as in the table above. A stock commitment that names a stock group couples to that group as a whole, through a single representative device, so a shared stock is not counted more than once.

Hard flow coupling

Devices in a coupling group \(k\) are forced to operate in fixed proportion to one another — for example a CHP unit whose gas input, heat output and power output move together. A free variable \(\alpha(k,j)\) represents the group’s common normalised flow level:

(18)\[P^{ems}(d,j) = \gamma(k,d) \cdot \alpha(k,j) \qquad \forall d \in k\]

Balance groups (internal commodity nodes)

An internal commodity node — a heat or steam network without its own grid connection — stores nothing itself: everything produced into the node must be consumed from it within the same time step.

(19)\[\sum_{d \in N(n)} P^{ems}(d,j) = 0\]

Derivative efficiencies and stock deltas describe each device’s own stock-side conversion, and do not enter this commodity-side balance. To give a node storage, include a storage device in the group: its flow absorbs the imbalance, while its stock is bounded by its own device constraints.

Power bands (S2 operation modes)

A device may declare a list of signed power bands, and must then operate within exactly one of them at every time step. Which band it runs in is a free binary decision:

(20)\[\sum_{b} y(d,b,j) = 1\]
(21)\[\sum_{b} y(d,b,j) \cdot B_{min}(d,b) \leq P_{up}(d,j) + P_{down}(d,j) \leq \sum_{b} y(d,b,j) \cdot B_{max}(d,b)\]

Because exactly one band binary is 1, each sum selects the bound of the chosen band. Bands are what make a device’s feasible operating region non-convex (for instance, “off, or between 40% and 100%”), at the cost of one binary variable per device per band per time step.