2d_constraints.gms 115 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
$ontext
This file is part of Backbone.

Backbone is free software: you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.

Backbone is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License
along with Backbone.  If not, see <http://www.gnu.org/licenses/>.
$offtext

18

19
* =============================================================================
20
* --- Constraint Equation Definitions -----------------------------------------
21
22
23
24
* =============================================================================

* --- Energy Balance ----------------------------------------------------------

25
q_balance(gn(grid, node), msft(m, s, f, t))${   not p_gn(grid, node, 'boundAll')
26
                                                } .. // Energy/power balance dynamics solved using implicit Euler discretization
27
28
29
30

    // The left side of the equation is the change in the state (will be zero if the node doesn't have a state)
    + p_gn(grid, node, 'energyStoredPerUnitOfState')${gn_state(grid, node)} // Unit conversion between v_state of a particular node and energy variables (defaults to 1, but can have node based values if e.g. v_state is in Kelvins and each node has a different heat storage capacity)
        * [
31
32
            + v_state(grid, node, s, f+df_central(f,t), t)                   // The difference between current
            - v_state(grid, node, s+ds_state(grid,node,s,t), f+df(f,t+dt(t)), t+dt(t))       // ... and previous state of the node
33
34
35
36
37
38
39
40
            ]

    =E=

    // The right side of the equation contains all the changes converted to energy terms
    + p_stepLength(m, f, t) // Multiply with the length of the timestep to convert power into energy
        * (
            // Self discharge out of the model boundaries
41
            - p_gn(grid, node, 'selfDischargeLoss')${ gn_state(grid, node) }
42
                * v_state(grid, node, s, f+df_central(f,t), t) // The current state of the node
43
44

            // Energy diffusion from this node to neighbouring nodes
45
            - sum(gnn_state(grid, node, to_node),
46
                + p_gnn(grid, node, to_node, 'diffCoeff')
47
                    * v_state(grid, node, s, f+df_central(f,t), t)
48
49
50
                ) // END sum(to_node)

            // Energy diffusion from neighbouring nodes to this node
51
            + sum(gnn_state(grid, from_node, node),
52
                + p_gnn(grid, from_node, node, 'diffCoeff')
53
                    * v_state(grid, from_node, s, f+df_central(f,t), t) // Incoming diffusion based on the state of the neighbouring node
54
55
56
                ) // END sum(from_node)

            // Controlled energy transfer, applies when the current node is on the left side of the connection
57
            - sum(gn2n_directional(grid, node, node_),
58
                + (1 - p_gnn(grid, node, node_, 'transferLoss')) // Reduce transfer losses
59
                    * v_transfer(grid, node, node_, s, f, t)
60
                + p_gnn(grid, node, node_, 'transferLoss') // Add transfer losses back if transfer is from this node to another node
61
                    * v_transferRightward(grid, node, node_, s, f, t)
62
63
64
                ) // END sum(node_)

            // Controlled energy transfer, applies when the current node is on the right side of the connection
65
            + sum(gn2n_directional(grid, node_, node),
66
                + v_transfer(grid, node_, node, s, f, t)
67
                - p_gnn(grid, node_, node, 'transferLoss') // Reduce transfer losses if transfer is from another node to this node
68
                    * v_transferRightward(grid, node_, node, s, f, t)
69
70
71
72
                ) // END sum(node_)

            // Interactions between the node and its units
            + sum(gnuft(grid, node, unit, f, t),
73
                + v_gen(grid, node, unit, s, f, t) // Unit energy generation and consumption
74
                )
75
76

            // Spilling energy out of the endogenous grids in the model
77
            - v_spill(grid, node, s, f, t)${node_spill(node)}
78
79

            // Power inflow and outflow timeseries to/from the node
80
            + ts_influx_(grid, node, f, t, s)   // Incoming (positive) and outgoing (negative) absolute value time series
81
82

            // Dummy generation variables, for feasibility purposes
83
84
            + vq_gen('increase', grid, node, s, f, t) // Note! When stateSlack is permitted, have to take caution with the penalties so that it will be used first
            - vq_gen('decrease', grid, node, s, f, t) // Note! When stateSlack is permitted, have to take caution with the penalties so that it will be used first
85
    ) // END * p_stepLength
86
;
87
88

* --- Reserve Demand ----------------------------------------------------------
89
90
// NOTE! Currently, there are multiple identical instances of the reserve balance equation being generated for each forecast branch even when the reserves are committed and identical between the forecasts.
// NOTE! This could be solved by formulating a new "ft_reserves" set to cover only the relevant forecast-time steps, but it would possibly make the reserves even more confusing.
91

92
q_resDemand(restypeDirectionNode(restype, up_down, node), sft(s, f, t))
93
94
    ${  ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')
        and not [ restypeReleasedForRealization(restype)
95
                  and sft_realized(s, f, t)]
96
        } ..
97
98
    // Reserve provision by capable units on this node
    + sum(nuft(node, unit, f, t)${nuRescapable(restype, up_down, node, unit)},
99
        + v_reserve(restype, up_down, node, unit, s, f+df_reserves(node, restype, f, t), t)
100
101
102
103
            * [ // Account for reliability of reserves
                + 1${sft_realized(s, f+df_reserves(node, restype, f, t), t)} // reserveReliability limits the reliability of reserves locked ahead of time.
                + p_nuReserves(node, unit, restype, 'reserveReliability')${not sft_realized(s, f+df_reserves(node, restype, f, t), t)}
                ] // END * v_reserve
104
105
        ) // END sum(nuft)

106
    // Reserve provision from other reserve categories when they can be shared
107
    + sum((nuft(node, unit, f, t), restype_)${p_nuRes2Res(node, unit, restype_, up_down, restype)},
108
        + v_reserve(restype_, up_down, node, unit, s, f+df_reserves(node, restype_, f, t), t)
109
            * p_nuRes2Res(node, unit, restype_, up_down, restype)
110
111
112
113
114
            * [ // Account for reliability of reserves
                + 1${sft_realized(s, f+df_reserves(node, restype, f, t), t)} // reserveReliability limits the reliability of reserves locked ahead of time.
                + p_nuReserves(node, unit, restype, 'reserveReliability')${not sft_realized(s, f+df_reserves(node, restype, f, t), t)}
                    * p_nuReserves(node, unit, restype_, 'reserveReliability')
                ] // END * v_reserve
115
116
        ) // END sum(nuft)

117
    // Reserve provision to this node via transfer links
118
    + sum(gn2n_directional(grid, node_, node)${restypeDirectionNodeNode(restype, up_down, node_, node)},
119
        + (1 - p_gnn(grid, node_, node, 'transferLoss') )
120
            * v_resTransferRightward(restype, up_down, node_, node, s, f+df_reserves(node_, restype, f, t), t) // Reserves from another node - reduces the need for reserves in the node
121
        ) // END sum(gn2n_directional)
122
    + sum(gn2n_directional(grid, node, node_)${restypeDirectionNodeNode(restype, up_down, node_, node)},
123
        + (1 - p_gnn(grid, node, node_, 'transferLoss') )
124
            * v_resTransferLeftward(restype, up_down, node, node_, s, f+df_reserves(node_, restype, f, t), t) // Reserves from another node - reduces the need for reserves in the node
125
126
127
128
129
130
131
132
        ) // END sum(gn2n_directional)

    =G=

    // Demand for reserves
    + ts_reserveDemand_(restype, up_down, node, f, t)${p_nReserves(node, restype, 'use_time_series')}
    + p_nReserves(node, restype, up_down)${not p_nReserves(node, restype, 'use_time_series')}

133
134
    // Reserve demand increase because of units
    + sum(nuft(node, unit, f, t)${p_nuReserves(node, unit, restype, 'reserve_increase_ratio')}, // Could be better to have 'reserve_increase_ratio' separately for up and down directions
135
        + sum(gnu(grid, node, unit), v_gen(grid, node, unit, s, f, t)) // Reserve sets and variables are currently lacking the grid dimension...
136
137
138
            * p_nuReserves(node, unit, restype, 'reserve_increase_ratio')
        ) // END sum(nuft)

139
    // Reserve provisions to another nodes via transfer links
140
    + sum(gn2n_directional(grid, node, node_)${restypeDirectionNodeNode(restype, up_down, node, node_)},   // If trasferring reserves to another node, increase your own reserves by same amount
141
        + v_resTransferRightward(restype, up_down, node, node_, s, f+df_reserves(node, restype, f, t), t)
142
        ) // END sum(gn2n_directional)
143
    + sum(gn2n_directional(grid, node_, node)${restypeDirectionNodeNode(restype, up_down, node, node_)},   // If trasferring reserves to another node, increase your own reserves by same amount
144
        + v_resTransferLeftward(restype, up_down, node_, node, s, f+df_reserves(node, restype, f, t), t)
145
146
147
        ) // END sum(gn2n_directional)

    // Reserve demand feasibility dummy variables
148
149
    - vq_resDemand(restype, up_down, node, s, f+df_reserves(node, restype, f, t), t)
    - vq_resMissing(restype, up_down, node, s, f+df_reserves(node, restype, f, t), t)${ft_reservesFixed(node, restype, f+df_reserves(node, restype, f, t), t)}
150
;
151

152
153
154
155
* --- N-1 Reserve Demand ----------------------------------------------------------
// NOTE! Currently, there are multiple identical instances of the reserve balance equation being generated for each forecast branch even when the reserves are committed and identical between the forecasts.
// NOTE! This could be solved by formulating a new "ft_reserves" set to cover only the relevant forecast-time steps, but it would possibly make the reserves even more confusing.

156
q_resDemandLargestInfeedUnit(grid, restypeDirectionNode(restype, 'up', node), unit_fail(unit_), sft(s, f, t))
157
    ${  ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')
158
        and gn(grid, node)
159
160
161
        and not [ restypeReleasedForRealization(restype)
            and ft_realized(f, t)
            ]
162
        and p_nuReserves(node, unit_, restype, 'portion_of_infeed_to_reserve')
163
164
165
        } ..
    // Reserve provision by capable units on this node excluding the failing one
    + sum(nuft(node, unit, f, t)${nuRescapable(restype, 'up', node, unit) and (ord(unit_) ne ord(unit))},
166
        + v_reserve(restype, 'up', node, unit, s, f+df_reserves(node, restype, f, t), t)
167
168
169
170
            * [ // Account for reliability of reserves
                + 1${sft_realized(s, f+df_reserves(node, restype, f, t), t)} // reserveReliability limits the reliability of reserves locked ahead of time.
                + p_nuReserves(node, unit, restype, 'reserveReliability')${not sft_realized(s, f+df_reserves(node, restype, f, t), t)}
                ] // END * v_reserve
171
172
173
        ) // END sum(nuft)

    // Reserve provision from other reserve categories when they can be shared
174
    + sum((nuft(node, unit, f, t), restype_)${p_nuRes2Res(node, unit, restype_, 'up', restype)},
175
        + v_reserve(restype_, 'up', node, unit, s, f+df_reserves(node, restype_, f, t), t)
176
            * p_nuRes2Res(node, unit, restype_, 'up', restype)
177
178
179
180
181
            * [ // Account for reliability of reserves
                + 1${sft_realized(s, f+df_reserves(node, restype, f, t), t)} // reserveReliability limits the reliability of reserves locked ahead of time.
                + p_nuReserves(node, unit, restype, 'reserveReliability')${not sft_realized(s, f+df_reserves(node, restype, f, t), t)}
                    * p_nuReserves(node, unit, restype_, 'reserveReliability')
                ] // END * v_reserve
182
183
184
185
186
        ) // END sum(nuft)

    // Reserve provision to this node via transfer links
    + sum(gn2n_directional(grid, node_, node)${restypeDirectionNodeNode(restype, 'up', node_, node)},
        + (1 - p_gnn(grid, node_, node, 'transferLoss') )
187
            * v_resTransferRightward(restype, 'up', node_, node, s, f+df_reserves(node_, restype, f, t), t) // Reserves from another node - reduces the need for reserves in the node
188
189
190
        ) // END sum(gn2n_directional)
    + sum(gn2n_directional(grid, node, node_)${restypeDirectionNodeNode(restype, 'up', node_, node)},
        + (1 - p_gnn(grid, node, node_, 'transferLoss') )
191
            * v_resTransferLeftward(restype, 'up', node, node_, s, f+df_reserves(node_, restype, f, t), t) // Reserves from another node - reduces the need for reserves in the node
192
193
194
195
196
        ) // END sum(gn2n_directional)

    =G=

    // Demand for reserves of the failing one
197
    v_gen(grid,node,unit_,s,f,t) * p_nuReserves(node, unit_, restype, 'portion_of_infeed_to_reserve')
198
199

    // Reserve provisions to another nodes via transfer links
200
    + sum(gn2n_directional(grid, node, node_)${restypeDirectionNodeNode(restype, 'up', node, node_)},   // If trasferring reserves to another node, increase your own reserves by same amount
201
        + v_resTransferRightward(restype, 'up', node, node_, s, f+df_reserves(node, restype, f, t), t)
202
        ) // END sum(gn2n_directional)
203
    + sum(gn2n_directional(grid, node_, node)${restypeDirectionNodeNode(restype, 'up', node, node_)},   // If trasferring reserves to another node, increase your own reserves by same amount
204
        + v_resTransferLeftward(restype, 'up', node_, node, s, f+df_reserves(node, restype, f, t), t)
205
206
207
        ) // END sum(gn2n_directional)

    // Reserve demand feasibility dummy variables
208
209
    - vq_resDemand(restype, 'up', node, s, f+df_reserves(node, restype, f, t), t)
    - vq_resMissing(restype, 'up', node, s, f+df_reserves(node, restype, f, t), t)${ft_reservesFixed(node, restype, f+df_reserves(node, restype, f, t), t)}
210
;
211
212
* --- Maximum Downward Capacity -----------------------------------------------

Erkka Rinne's avatar
Erkka Rinne committed
213
q_maxDownward(gnu(grid, node, unit), msft(m, s, f, t))${gnuft(grid, node, unit, f, t)
214
215
                                                    and {
                                                    [   ord(t) < tSolveFirst + smax(restype, p_nReserves(node, restype, 'reserve_length')) // Unit is either providing
216
217
218
219
220
221
222
223
224
225
226
227
228
229
                                                        and sum(restype, nuRescapable(restype, 'down', node, unit)) // downward reserves
                                                        ]
                                                    // NOTE!!! Could be better to form a gnuft_reserves subset?
                                                    or [ // the unit has an online variable
                                                        uft_online(unit, f, t)
                                                        and [
                                                            (unit_minLoad(unit) and p_gnu(grid, node, unit, 'unitSizeGen')) // generators with a min. load
                                                            or p_gnu(grid, node, unit, 'maxCons') // or consuming units with an online variable
                                                            ]
                                                        ] // END or
                                                    or [ // consuming units with investment possibility
                                                        gnu_input(grid, node, unit)
                                                        and [unit_investLP(unit) or unit_investMIP(unit)]
                                                        ]
230
                                                    }} ..
231
    // Energy generation/consumption
232
    + v_gen(grid, node, unit, s, f, t)
233
234

    // Considering output constraints (e.g. cV line)
235
236
    + sum(gngnu_constrainedOutputRatio(grid, node, grid_output, node_, unit),
        + p_gnu(grid_output, node_, unit, 'cV')
237
            * v_gen(grid_output, node_, unit, s, f, t)
238
239
240
        ) // END sum(gngnu_constrainedOutputRatio)

    // Downward reserve participation
241
    - sum(nuRescapable(restype, 'down', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
242
        + v_reserve(restype, 'down', node, unit, s, f+df_reserves(node, restype, f, t), t) // (v_reserve can be used only if the unit is capable of providing a particular reserve)
243
244
245
246
247
248
        ) // END sum(nuRescapable)

    =G= // Must be greater than minimum load or maximum consumption  (units with min-load and both generation and consumption are not allowed)

    // Generation units, greater than minload
    + p_gnu(grid, node, unit, 'unitSizeGen')
249
        * sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
250
251
252
253
            + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
            + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
            ) // END sum(effGroup)
        * [ // Online variables should only be generated for units with restrictions
254
255
            + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f+df_central(f,t), t)} // LP online variant
            + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f+df_central(f,t), t)} // MIP online variant
256
257
            ] // END v_online

258
259
260
261
    // Units in run-up phase neet to keep up with the run-up rate
    + p_gnu(grid, node, unit, 'unitSizeGen')
        * sum(unitStarttype(unit, starttype)${uft_startupTrajectory(unit, f, t)},
            sum(runUpCounter(unit, counter)${t_active(t+dt_trajectory(counter))}, // Sum over the run-up intervals
262
263
264
265
266
267
                + [
                    + v_startup_LP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineLP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    + v_startup_MIP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineMIP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    ]
268
                    * p_uCounter_runUpMin(unit, counter)
269
270
271
272
273
274
275
                ) // END sum(runUpCounter)
            ) // END sum(unitStarttype)

    // Units in shutdown phase need to keep up with the shutdown rate
    + p_gnu(grid, node, unit, 'unitSizeGen')
        * sum(shutdownCounter(unit, counter)${t_active(t+dt_trajectory(counter)) and uft_shutdownTrajectory(unit, f, t)}, // Sum over the shutdown intervals
            + v_shutdown(unit, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
276
                * p_uCounter_shutdownMin(unit, counter)
277
            ) // END sum(shutdownCounter)
278

279
280
281
282
283
    // Consuming units, greater than maxCons
    // Available capacity restrictions
    - p_unit(unit, 'availability')
        * [
            // Capacity factors for flow units
284
            + sum(flowUnit(flow, unit),
285
                + ts_cf_(flow, node, f, t, s)
286
287
288
289
290
291
292
293
                ) // END sum(flow)
            + 1${not unit_flow(unit)}
            ] // END * p_unit(availability)
        * [
            // Online capacity restriction
            + p_gnu(grid, node, unit, 'maxCons')${not uft_online(unit, f, t)} // Use initial maximum if no online variables
            + p_gnu(grid, node, unit, 'unitSizeCons')
                * [
294
                    // Capacity online
295
296
                    + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
                    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
297
298
299
300
301
302
303
304

                    // Investments to additional non-online capacity
                    + sum(t_invest(t_)${    ord(t_)<=ord(t)
                                            and not uft_online(unit, f, t)
                                            },
                        + v_invest_LP(unit, t_)${unit_investLP(unit)} // NOTE! v_invest_LP also for consuming units is positive
                        + v_invest_MIP(unit, t_)${unit_investMIP(unit)} // NOTE! v_invest_MIP also for consuming units is positive
                        ) // END sum(t_invest)
305
306
                    ] // END * p_gnu(unitSizeCons)
            ] // END * p_unit(availability)
307
;
308
309
310

* --- Maximum Upwards Capacity ------------------------------------------------

Erkka Rinne's avatar
Erkka Rinne committed
311
q_maxUpward(gnu(grid, node, unit), msft(m, s, f, t))${gnuft(grid, node, unit, f, t)
312
313
                                                    and {
                                                 [   ord(t) < tSolveFirst + smax(restype, p_nReserves(node, restype, 'reserve_length')) // Unit is either providing
314
315
316
317
318
319
320
321
322
323
324
325
326
                                                    and sum(restype, nuRescapable(restype, 'up', node, unit)) // upward reserves
                                                    ]
                                                or [
                                                    uft_online(unit, f, t) // or the unit has an online variable
                                                        and [
                                                            [unit_minLoad(unit) and p_gnu(grid, node, unit, 'unitSizeCons')] // consuming units with min_load
                                                            or [p_gnu(grid, node, unit, 'maxGen')]                          // generators with an online variable
                                                            ]
                                                    ]
                                                or [
                                                    gnu_output(grid, node, unit) // generators with investment possibility
                                                    and (unit_investLP(unit) or unit_investMIP(unit))
                                                    ]
327
                                                }}..
328
    // Energy generation/consumption
329
    + v_gen(grid, node, unit, s, f, t)
330
331
332
333

    // Considering output constraints (e.g. cV line)
    + sum(gngnu_constrainedOutputRatio(grid, node, grid_output, node_, unit),
        + p_gnu(grid_output, node_, unit, 'cV')
334
            * v_gen(grid_output, node_, unit, s, f, t)
335
336
337
        ) // END sum(gngnu_constrainedOutputRatio)

    // Upwards reserve participation
338
    + sum(nuRescapable(restype, 'up', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
339
        + v_reserve(restype, 'up', node, unit, s, f+df_reserves(node, restype, f, t), t)
340
341
342
343
344
        ) // END sum(nuRescapable)

    =L= // must be less than available/online capacity

    // Consuming units
345
    - p_gnu(grid, node, unit, 'unitSizeCons')
346
        * sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
347
348
349
350
            + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
            + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
            ) // END sum(effGroup)
        * [
351
352
            + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)} // Consuming units are restricted by their min. load (consuming is negative)
            + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)} // Consuming units are restricted by their min. load (consuming is negative)
353
354
355
356
357
358
359
            ] // END * p_gnu(unitSizeCons)

    // Generation units
    // Available capacity restrictions
    + p_unit(unit, 'availability') // Generation units are restricted by their (available) capacity
        * [
            // Capacity factor for flow units
360
            + sum(flowUnit(flow, unit),
361
                + ts_cf_(flow, node, f, t, s)
362
363
364
365
366
367
368
369
                ) // END sum(flow)
            + 1${not unit_flow(unit)}
            ] // END * p_unit(availability)
        * [
            // Online capacity restriction
            + p_gnu(grid, node, unit, 'maxGen')${not uft_online(unit, f, t)} // Use initial maxGen if no online variables
            + p_gnu(grid, node, unit, 'unitSizeGen')
                * [
370
                    // Capacity online
371
372
                    + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f ,t)}
                    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
373
374
375
376
377
378
379
380

                    // Investments to non-online capacity
                    + sum(t_invest(t_)${    ord(t_)<=ord(t)
                                            and not uft_online(unit, f ,t)
                                            },
                        + v_invest_LP(unit, t_)${unit_investLP(unit)}
                        + v_invest_MIP(unit, t_)${unit_investMIP(unit)}
                        ) // END sum(t_invest)
381
382
                    ] // END * p_gnu(unitSizeGen)
            ] // END * p_unit(availability)
383

384
385
386
387
    // Units in run-up phase neet to keep up with the run-up rate
    + p_gnu(grid, node, unit, 'unitSizeGen')
        * sum(unitStarttype(unit, starttype)${uft_startupTrajectory(unit, f, t)},
            sum(runUpCounter(unit, counter)${t_active(t+dt_trajectory(counter))}, // Sum over the run-up intervals
388
389
390
391
392
393
                + [
                    + v_startup_LP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineLP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    + v_startup_MIP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineMIP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    ]
394
                    * p_uCounter_runUpMax(unit, counter)
395
396
397
398
399
400
401
                ) // END sum(runUpCounter)
            ) // END sum(unitStarttype)

    // Units in shutdown phase need to keep up with the shutdown rate
    + p_gnu(grid, node, unit, 'unitSizeGen')
        * sum(shutdownCounter(unit, counter)${t_active(t+dt_trajectory(counter)) and uft_shutdownTrajectory(unit, f, t)}, // Sum over the shutdown intervals
            + v_shutdown(unit, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
402
                * p_uCounter_shutdownMax(unit, counter)
403
            ) // END sum(shutdownCounter)
404
;
405

406
407
* --- Reserve Provision of Units with Investments -----------------------------

408
409
410
411
412
413
q_reserveProvision(nuRescapable(restypeDirectionNode(restype, up_down, node), unit), sft(s, f, t))${ord(t) <= tSolveFirst + p_nReserves(node, restype, 'reserve_length')
                                                                                                    and nuft(node, unit, f, t)
                                                                                                    and (unit_investLP(unit) or unit_investMIP(unit))
                                                                                                    and not ft_reservesFixed(node, restype, f+df_reserves(node, restype, f, t), t)
                                                                                                   } ..
    + v_reserve(restype, up_down, node, unit, s, f+df_reserves(node, restype, f, t), t)
414
415
416
417
418
419
420
421
422
423
424
425
426

    =L=

    + p_nuReserves(node, unit, restype, up_down)
        * [
            + sum(grid, p_gnu(grid, node, unit, 'maxGen') + p_gnu(grid, node, unit, 'maxCons') )  // Reserve sets and variables are currently lacking the grid dimension...
            + sum(t_invest(t_)${ ord(t_)<=ord(t) },
                + v_invest_LP(unit, t_)${unit_investLP(unit)}
                    * sum(grid, p_gnu(grid, node, unit, 'unitSizeTot')) // Reserve sets and variables are currently lacking the grid dimension...
                + v_invest_MIP(unit, t_)${unit_investMIP(unit)}
                    * sum(grid, p_gnu(grid, node, unit, 'unitSizeTot')) // Reserve sets and variables are currently lacking the grid dimension...
                ) // END sum(t_)
            ]
427
428
429
430
        * p_unit(unit, 'availability') // Taking into account availability...
        * [
            // ... and capacity factor for flow units
            + sum(flowUnit(flow, unit),
431
                + ts_cf_(flow, node, f, t, s)
432
433
                ) // END sum(flow)
            + 1${not unit_flow(unit)}
434
435
436
            ] // How to consider reserveReliability in the case of investments when we typically only have "realized" time steps?
;

437
438
* --- Unit Startup and Shutdown -----------------------------------------------

439
q_startshut(m, s, uft_online(unit, f, t))$msft(m, s, f, t) ..
440
    // Units currently online
441
442
    + v_online_LP (unit, s, f+df_central(f,t), t)${uft_onlineLP (unit, f, t)}
    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
443
444

    // Units previously online
445
446

    // The same units
447
    - v_online_LP (unit, s+ds(s,t), f+df(f,t+dt(t)), t+dt(t))${ uft_onlineLP_withPrevious(unit, f+df(f,t+dt(t)), t+dt(t))
448
                                                             and not uft_aggregator_first(unit, f, t) } // This reaches to tFirstSolve when dt = -1
449
    - v_online_MIP(unit, s+ds(s,t), f+df(f,t+dt(t)), t+dt(t))${ uft_onlineMIP_withPrevious(unit, f+df(f,t+dt(t)), t+dt(t))
450
451
452
453
                                                             and not uft_aggregator_first(unit, f, t) }

    // Aggregated units just before they are turned into aggregator units
    - sum(unit_${unitAggregator_unit(unit, unit_)},
454
455
        + v_online_LP (unit_, s, f+df(f,t+dt(t)), t+dt(t))${uft_onlineLP_withPrevious(unit_, f+df(f,t+dt(t)), t+dt(t))}
        + v_online_MIP(unit_, s, f+df(f,t+dt(t)), t+dt(t))${uft_onlineMIP_withPrevious(unit_, f+df(f,t+dt(t)), t+dt(t))}
456
        )${uft_aggregator_first(unit, f, t)} // END sum(unit_)
457

458
459
    =E=

460
    // Unit startup and shutdown
461

462
    // Add startup of units dt_toStartup before the current t (no start-ups for aggregator units before they become active)
463
    + sum(unitStarttype(unit, starttype),
464
465
466
467
        + v_startup_LP(unit, starttype, s, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t))
            ${ uft_onlineLP(unit, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t)) }
        + v_startup_MIP(unit, starttype, s, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t))
            ${ uft_onlineMIP(unit, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t)) }
468
        )${not [unit_aggregator(unit) and ord(t) + dt_toStartup(unit, t) <= tSolveFirst + p_unit(unit, 'lastStepNotAggregated')]} // END sum(starttype)
469

470
471
472
473
    // NOTE! According to 3d_setVariableLimits,
    // cannot start a unit if the time when the unit would become online is outside
    // the horizon when the unit has an online variable
    // --> no need to add start-ups of aggregated units to aggregator units
474

475
    // Shutdown of units at time t
476
    - v_shutdown(unit, s, f, t)
477
;
478

479
*--- Startup Type -------------------------------------------------------------
480
// !!! NOTE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
481
482
483
// This formulation doesn't work as intended when unitCount > 1, as one recent
// shutdown allows for multiple hot/warm startups on subsequent time steps.
// Pending changes.
484

485
486
q_startuptype(m, s, starttypeConstrained(starttype), uft_online(unit, f, t))
    ${msft(m, s, f, t) and unitStarttype(unit, starttype)} ..
487
488

    // Startup type
489
490
    + v_startup_LP(unit, starttype, s, f, t)${ uft_onlineLP(unit, f, t) }
    + v_startup_MIP(unit, starttype, s, f, t)${ uft_onlineMIP(unit, f, t) }
491
492
493
494

    =L=

    // Subunit shutdowns within special startup timeframe
Topi Rasku's avatar
Topi Rasku committed
495
    + sum(unitCounter(unit, counter)${dt_starttypeUnitCounter(starttype, unit, counter)},
496
        + v_shutdown(unit, s, f+df(f,t+(dt_starttypeUnitCounter(starttype, unit, counter)+1)), t+(dt_starttypeUnitCounter(starttype, unit, counter)+1))
Topi Rasku's avatar
Topi Rasku committed
497
            ${t_active(t+(dt_starttypeUnitCounter(starttype, unit, counter)+1))}
498
499
500
        ) // END sum(counter)

    // NOTE: for aggregator units, shutdowns for aggregated units are not considered
501
;
502

503

504
505
*--- Online Limits with Startup Type Constraints and Investments --------------

506
507
q_onlineLimit(m, s, uft_online(unit, f, t))${msft(m, s, f, t) and {
                                            p_unit(unit, 'minShutdownHours')
508
                                            or p_u_runUpTimeIntervals(unit)
509
510
                                            or unit_investLP(unit)
                                            or unit_investMIP(unit)
511
                                            }} ..
512
    // Online variables
513
514
    + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f ,t)}
515
516
517
518
519
520

    =L=

    // Number of existing units
    + p_unit(unit, 'unitCount')

521
    // Number of units unable to become online due to restrictions
Topi Rasku's avatar
Topi Rasku committed
522
    - sum(unitCounter(unit, counter)${dt_downtimeUnitCounter(unit, counter)},
523
        + v_shutdown(unit, s, f+df(f,t+(dt_downtimeUnitCounter(unit, counter) + 1)), t+(dt_downtimeUnitCounter(unit, counter) + 1))
Topi Rasku's avatar
Topi Rasku committed
524
            ${t_active(t+(dt_downtimeUnitCounter(unit, counter) + 1))}
525
526
527
        ) // END sum(counter)

    // Number of units unable to become online due to restrictions (aggregated units in the past horizon or if they have an online variable)
528
    - sum(unitAggregator_unit(unit, unit_),
Topi Rasku's avatar
Topi Rasku committed
529
        + sum(unitCounter(unit, counter)${dt_downtimeUnitCounter(unit, counter)},
530
            + v_shutdown(unit_, s, f+df(f,t+(dt_downtimeUnitCounter(unit, counter) + 1)), t+(dt_downtimeUnitCounter(unit, counter) + 1))
Topi Rasku's avatar
Topi Rasku committed
531
                ${t_active(t+(dt_downtimeUnitCounter(unit, counter) + 1))}
532
533
            ) // END sum(counter)
        )${unit_aggregator(unit)} // END sum(unit_)
534
535
536

    // Investments into units
    + sum(t_invest(t_)${ord(t_)<=ord(t)},
537
538
        + v_invest_LP(unit, t_)${unit_investLP(unit)}
        + v_invest_MIP(unit, t_)${unit_investMIP(unit)}
539
540
541
        ) // END sum(t_invest)
;

542
543
544
545
*--- Both q_offlineAfterShutdown and q_onlineOnStartup work when there is only one unit.
*    These equations prohibit single units turning on and off at the same time step.
*    Unfortunately there seems to be no way to prohibit this when unit count is > 1.
*    (it shouldn't be worthwhile anyway if there is a startup cost, but it can fall within the solution gap).
546
547
q_onlineOnStartUp(s, uft_online(unit, f, t))
    ${sft(s, f, t) and sum(starttype, unitStarttype(unit, starttype))}..
548
549

    // Units currently online
550
551
    + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
552
553
554
555

    =G=

    + sum(unitStarttype(unit, starttype),
556
557
558
559
        + v_startup_LP(unit, starttype, s, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t)) //dt_toStartup displaces the time step to the one where the unit would be started up in order to reach online at t
            ${ uft_onlineLP(unit, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t)) }
        + v_startup_MIP(unit, starttype, s, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t)) //dt_toStartup displaces the time step to the one where the unit would be started up in order to reach online at t
            ${ uft_onlineMIP(unit, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t)) }
560
561
562
      ) // END sum(starttype)
;

563
564
q_offlineAfterShutdown(s, uft_online(unit, f, t))
    ${sft(s, f, t) and sum(starttype, unitStarttype(unit, starttype))}..
565

566
567
568
569
570
    // Number of existing units
    + p_unit(unit, 'unitCount')

    // Investments into units
    + sum(t_invest(t_)${ord(t_)<=ord(t)},
571
572
        + v_invest_LP(unit, t_)${unit_investLP(unit)}
        + v_invest_MIP(unit, t_)${unit_investMIP(unit)}
573
574
        ) // END sum(t_invest)

575
    // Units currently online
576
577
    - v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
    - v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
578
579
580

    =G=

581
    + v_shutdown(unit, s, f, t)
582
583
;

584
585
*--- Minimum Unit Uptime ------------------------------------------------------

586
587
q_onlineMinUptime(m, s, uft_online(unit, f, t))
    ${msft(m, s, f, t) and  p_unit(unit, 'minOperationHours')} ..
588
589

    // Units currently online
590
591
    + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
592
593
594
595

    =G=

    // Units that have minimum operation time requirements active
596
597
598
    + sum(unitCounter(unit, counter)${  dt_uptimeUnitCounter(unit, counter)
                                        and t_active(t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)) // Don't sum over counters that don't point to an active time step
                                        },
599
        + sum(unitStarttype(unit, starttype),
600
601
602
603
            + v_startup_LP(unit, starttype, s, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1))
                ${ uft_onlineLP(unit, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)) }
            + v_startup_MIP(unit, starttype, s, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1))
                ${ uft_onlineMIP(unit, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)) }
604
            ) // END sum(starttype)
605
606
607
        ) // END sum(counter)

    // Units that have minimum operation time requirements active (aggregated units in the past horizon or if they have an online variable)
Topi Rasku's avatar
Topi Rasku committed
608
    + sum(unitAggregator_unit(unit, unit_),
609
610
611
        + sum(unitCounter(unit, counter)${  dt_uptimeUnitCounter(unit, counter)
                                            and t_active(t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)) // Don't sum over counters that don't point to an active time step
                                            },
612
            + sum(unitStarttype(unit, starttype),
613
614
615
616
                + v_startup_LP(unit, starttype, s, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1))
                    ${ uft_onlineLP(unit, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)) }
                + v_startup_MIP(unit, starttype, s, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1))
                    ${ uft_onlineMIP(unit, f+df(f,t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)), t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1)) }
617
618
619
                ) // END sum(starttype)
            ) // END sum(counter)
        )${unit_aggregator(unit)} // END sum(unit_)
620
621
;

622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
* --- Cyclic Boundary Conditions for Online State -----------------------------

q_onlineCyclic(uss_bound(unit, s_, s), m)${ ms(m, s_)
                                            and ms(m, s)
                                            and tSolveFirst = mSettings(m, 't_start')
                                            }..

    // Initial value of the state of the unit at the start of the sample
    + sum(mst_start(m, s, t),
        + sum(sft(s, f, t),
            + v_online_LP(unit, s, f+df(f,t+dt(t)), t+dt(t))${uft_onlineLP_withPrevious(unit, f+df(f,t+dt(t)), t+dt(t))}
            + v_online_MIP(unit, s, f+df(f,t+dt(t)), t+dt(t))${uft_onlineMIP_withPrevious(unit, f+df(f,t+dt(t)), t+dt(t))}
            ) // END sum(ft)
        ) // END sum(mst_start)

    =E=

    // State of the unit at the end of the sample
    + sum(mst_end(m, s_, t_),
        + sum(sft(s_, f_, t_),
            + v_online_LP(unit, s_, f_, t_)${uft_onlineLP(unit, f_, t_)}
            + v_online_MIP(unit, s_, f_, t_)${uft_onlineMIP(unit, f_, t_)}
            ) // END sum(ft)
        ) // END sum(mst_end)
;

648
* --- Ramp Constraints --------------------------------------------------------
649
650
651
652

q_genRamp(m, s, gnuft_ramp(grid, node, unit, f, t))${  ord(t) > msStart(m, s) + 1
                                                       and msft(m, s, f, t)
                                                       } ..
653

654
655
    + v_genRamp(grid, node, unit, s, f, t)
        * p_stepLength(m, f, t)
656

657
    =E=
658

659
    // Change in generation over the interval: v_gen(t) - v_gen(t-1)
660
    + v_gen(grid, node, unit, s, f, t)
661

662
    // Unit generation at t-1 (except aggregator units right before the aggregation threshold, see next term)
663
    - v_gen(grid, node, unit, s+ds(s,t), f+df(f,t+dt(t)), t+dt(t))${not uft_aggregator_first(unit, f, t)}
664
665
    // Unit generation at t-1, aggregator units right before the aggregation threshold
    + sum(unit_${unitAggregator_unit(unit, unit_)},
666
        - v_gen(grid, node, unit_, s+ds(s,t), f+df(f,t+dt(t)), t+dt(t))
667
      )${uft_aggregator_first(unit, f, t)}
668
;
669

670
* --- Ramp Up Limits ----------------------------------------------------------
671
672
673
674

q_rampUpLimit(m, s, gnuft_ramp(grid, node, unit, f, t))${  ord(t) > msStart(m, s) + 1
                                                           and msft(m, s, f, t)
                                                           and p_gnu(grid, node, unit, 'maxRampUp')
675
676
677
678
679
                                                           and [ sum(restype, nuRescapable(restype, 'up', node, unit))
                                                                 or uft_online(unit, f, t)
                                                                 or unit_investLP(unit)
                                                                 or unit_investMIP(unit)
                                                                 ]
680
                                                           } ..
681
    + v_genRamp(grid, node, unit, s, f, t)
682
    + sum(nuRescapable(restype, 'up', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
683
        + v_reserve(restype, 'up', node, unit, s, f+df_reserves(node, restype, f, t), t) // (v_reserve can be used only if the unit is capable of providing a particular reserve)
684
685
686
687
688
        ) // END sum(nuRescapable)
        / p_stepLength(m, f, t)

    =L=

689
    // Ramping capability of units without an online variable
690
691
692
693
694
695
696
697
698
699
700
701
    + (
        + ( p_gnu(grid, node, unit, 'maxGen') + p_gnu(grid, node, unit, 'maxCons') )${not uft_online(unit, f, t)}
        + sum(t_invest(t_)${ ord(t_)<=ord(t) },
            + v_invest_LP(unit, t_)${not uft_onlineLP(unit, f, t) and unit_investLP(unit)}
                * p_gnu(grid, node, unit, 'unitSizeTot')
            + v_invest_MIP(unit, t_)${not uft_onlineMIP(unit, f, t) and unit_investMIP(unit)}
                * p_gnu(grid, node, unit, 'unitSizeTot')
          )
      )
        * p_gnu(grid, node, unit, 'maxRampUp')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]

702
    // Ramping capability of units with an online variable
703
    + (
704
705
        + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
        + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
706
707
708
709
710
      )
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnu(grid, node, unit, 'maxRampUp')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]

711
712
713
714
715
716
717
718
719
720
721
722
723
    // Generation units not be able to ramp from zero to min. load within one time interval according to their maxRampUp
    + sum(unitStarttype(unit, starttype)${   uft_online(unit, f, t)
                                             and gnu_output(grid, node, unit)
                                             and not uft_startupTrajectory(unit, f, t)
                                             and ( + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                                                       + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                                                       + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
                                                     ) // END sum(effGroup)
                                                       / p_stepLength(m, f, t)
                                                   - p_gnu(grid, node, unit, 'maxRampUp')
                                                       * 60 > 0
                                                   )
                                             },
724
725
726
727
        + v_startup_LP(unit, starttype, s, f, t)
            ${ uft_onlineLP(unit, f, t) }
        + v_startup_MIP(unit, starttype, s, f, t)
            ${ uft_onlineMIP(unit, f, t) }
728
729
730
731
732
733
734
735
736
737
738
739
      ) // END sum(starttype)
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * (
            + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
              ) // END sum(effGroup)
                / p_stepLength(m, f, t)
            - p_gnu(grid, node, unit, 'maxRampUp')
                * 60   // Unit conversion from [p.u./min] to [p.u./h]
          ) // END * v_startup

740
741
742
743
    // Units in the run-up phase need to keep up with the run-up rate
    + p_gnu(grid, node, unit, 'unitSizeTot')
        * sum(unitStarttype(unit, starttype)${uft_startupTrajectory(unit, f, t)},
            sum(runUpCounter(unit, counter)${t_active(t+dt_trajectory(counter))}, // Sum over the run-up intervals
744
745
746
747
748
749
                + [
                    + v_startup_LP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineLP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    + v_startup_MIP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineMIP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    ]
750
751
752
753
754
755
756
757
                    * [
                        + p_unit(unit, 'rampSpeedToMinLoad')
                        + ( p_gnu(grid, node, unit, 'maxRampUp') - p_unit(unit, 'rampSpeedToMinLoad') )${ not runUpCounter(unit, counter+1) } // Ramp speed adjusted for the last run-up interval
                            * ( p_u_runUpTimeIntervalsCeil(unit) - p_u_runUpTimeIntervals(unit) )
                        ]
                    * 60 // Unit conversion from [p.u./min] into [p.u./h]
                ) // END sum(runUpCounter)
            ) // END sum(unitStarttype)
758

759
    // Shutdown of consumption units according to maxRampUp
760
    + v_shutdown(unit, s, f, t)${uft_online(unit, f, t) and gnu_input(grid, node, unit)}
761
        * p_gnu(grid, node, unit, 'unitSizeTot')
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
        * p_gnu(grid, node, unit, 'maxRampUp')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
    // Consumption units not be able to ramp from min. load to zero within one time interval according to their maxRampUp
    + v_shutdown(unit, s, f, t)${   uft_online(unit, f, t)
                                    and gnu_input(grid, node, unit)
                                    and ( + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                                              + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                                              + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
                                              ) // END sum(effGroup)
                                              / p_stepLength(m, f, t)
                                          - p_gnu(grid, node, unit, 'maxRampUp')
                                              * 60 > 0
                                          )
                                    }
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * (
            + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
                ) // END sum(effGroup)
                / p_stepLength(m, f, t)
            - p_gnu(grid, node, unit, 'maxRampUp')
                * 60   // Unit conversion from [p.u./min] to [p.u./h]
          ) // END * v_startup
786
;
787

788
* --- Ramp Down Limits --------------------------------------------------------
789
790
791
792

q_rampDownLimit(m, s, gnuft_ramp(grid, node, unit, f, t))${  ord(t) > msStart(m, s) + 1
                                                             and msft(m, s, f, t)
                                                             and p_gnu(grid, node, unit, 'maxRampDown')
793
794
795
796
797
                                                             and [ sum(restype, nuRescapable(restype, 'down', node, unit))
                                                                   or uft_online(unit, f, t)
                                                                   or unit_investLP(unit)
                                                                   or unit_investMIP(unit)
                                                                   ]
798
                                                             } ..
799
    + v_genRamp(grid, node, unit, s, f, t)
800
    - sum(nuRescapable(restype, 'down', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
801
        + v_reserve(restype, 'down', node, unit, s, f+df_reserves(node, restype, f, t), t) // (v_reserve can be used only if the unit is capable of providing a particular reserve)
802
803
804
805
806
        ) // END sum(nuRescapable)
        / p_stepLength(m, f, t)

    =G=

807
    // Ramping capability of units without online variable
808
809
810
811
812
813
814
815
816
817
818
819
    - (
        + ( p_gnu(grid, node, unit, 'maxGen') + p_gnu(grid, node, unit, 'maxCons') )${not uft_online(unit, f, t)}
        + sum(t_invest(t_)${ ord(t_)<=ord(t) },
            + v_invest_LP(unit, t_)${not uft_onlineLP(unit, f, t) and unit_investLP(unit)}
                * p_gnu(grid, node, unit, 'unitSizeTot')
            + v_invest_MIP(unit, t_)${not uft_onlineMIP(unit, f, t) and unit_investMIP(unit)}
                * p_gnu(grid, node, unit, 'unitSizeTot')
          )
      )
        * p_gnu(grid, node, unit, 'maxRampDown')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]

820
    // Ramping capability of units that are online
821
    - (
822
823
        + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
        + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
824
825
826
827
828
      )
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnu(grid, node, unit, 'maxRampDown')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]

829
    // Shutdown of generation units according to maxRampDown
830
    - v_shutdown(unit, s, f, t)${   uft_online(unit, f, t)
831
832
833
                                    and gnu_output(grid, node, unit)
                                    and not uft_shutdownTrajectory(unit, f, t)
                                    }
834
        * p_gnu(grid, node, unit, 'unitSizeTot')
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
        * p_gnu(grid, node, unit, 'maxRampDown')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
    // Generation units not be able to ramp from min. load to zero within one time interval according to their maxRampDown
    - v_shutdown(unit, s, f, t)${   uft_online(unit, f, t)
                                    and gnu_output(grid, node, unit)
                                    and not uft_shutdownTrajectory(unit, f, t)
                                    and ( + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                                              + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                                              + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
                                            ) // END sum(effGroup)
                                            / p_stepLength(m, f, t)
                                          - p_gnu(grid, node, unit, 'maxRampDown')
                                              * 60 > 0
                                        )
                                }
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * (
            + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
                ) // END sum(effGroup)
                / p_stepLength(m, f, t)
            - p_gnu(grid, node, unit, 'maxRampDown')
                * 60   // Unit conversion from [p.u./min] to [p.u./h]
          ) // END * v_shutdown
860

861
862
    // Units in shutdown phase need to keep up with the shutdown ramp rate
    - p_gnu(grid, node, unit, 'unitSizeGen')
863
864
865
866
867
868
869
870
871
        * [
            + sum(shutdownCounter(unit, counter)${t_active(t+dt_trajectory(counter)) and uft_shutdownTrajectory(unit, f, t)}, // Sum over the shutdown intervals
                + v_shutdown(unit, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                    * [
                        + p_unit(unit, 'rampSpeedFromMinLoad')
                        + ( p_gnu(grid, node, unit, 'maxRampDown') - p_unit(unit, 'rampSpeedFromMinLoad') )${ not shutdownCounter(unit, counter-1) } // Ramp speed adjusted for the first shutdown interval
                            * ( p_u_shutdownTimeIntervalsCeil(unit) - p_u_shutdownTimeIntervals(unit) )
                        ]
                ) // END sum(shutdownCounter)
872
            // Units need to be able to shut down after shut down trajectory
873
            + v_shutdown(unit, s, f+df(f, t+dt_toShutdown(unit, t)), t+dt_toShutdown(unit, t))${uft_shutdownTrajectory(unit, f, t)}
874
875
876
                * p_unit(unit, 'rampSpeedFromMinload')
            ]
        * 60 // Unit conversion from [p.u./min] to [p.u./h]
877
878
879
880
881
882
883
884
885
886
887
888
889

    // Consumption units not be able to ramp from zero to min. load within one time interval according to their maxRampDown
    - sum(unitStarttype(unit, starttype)${   uft_online(unit, f, t)
                                             and gnu_input(grid, node, unit)
                                             and ( + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                                                       + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                                                       + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
                                                     ) // END sum(effGroup)
                                                       / p_stepLength(m, f, t)
                                                   - p_gnu(grid, node, unit, 'maxRampDown')
                                                       * 60 > 0
                                                   )
                                             },
890
891
892
893
        + v_startup_LP(unit, starttype, s, f, t)
            ${ uft_onlineLP(unit, f, t) }
        + v_startup_MIP(unit, starttype, s, f, t)
            ${ uft_onlineMIP(unit, f, t) }
894
895
896
897
898
899
900
901
902
903
904
      ) // END sum(starttype)
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * (
            + sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
                + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
                + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
              ) // END sum(effGroup)
                / p_stepLength(m, f, t)
            - p_gnu(grid, node, unit, 'maxRampDown')
                * 60   // Unit conversion from [p.u./min] to [p.u./h]
          ) // END * v_startup
905
906
;

907
908
909
910
911
912
913
* --- Ramps separated into upward and downward ramps --------------------------

q_rampUpDown(m, s, gnuft_ramp(grid, node, unit, f, t))${  ord(t) > msStart(m, s) + 1
                                                          and msft(m, s, f, t)
                                                          and sum(slack, gnuft_rampCost(grid, node, unit, slack, f, t))
                                                          } ..

914
    + v_genRamp(grid, node, unit, s, f, t)
915

916
    =E=
917

918
919
    // Upward and downward ramp categories
    + sum(slack${ gnuft_rampCost(grid, node, unit, slack, f, t) },
920
921
        + v_genRampUpDown(grid, node, unit, slack, s, f, t)$upwardSlack(slack)
        - v_genRampUpDown(grid, node, unit, slack, s, f, t)$downwardSlack(slack)
922
      ) // END sum(slack)
923
924
;

Niina Helistö's avatar
Niina Helistö committed
925
* --- Upward and downward ramps constrained by slack boundaries ---------------
926
927
928
929
930

q_rampSlack(m, s, gnuft_rampCost(grid, node, unit, slack, f, t))${  ord(t) > msStart(m, s) + 1
                                                                    and msft(m, s, f, t)
                                                                    } ..

931
    + v_genRampUpDown(grid, node, unit, slack, s, f, t)
932

933
    =L=
934
935

    // Ramping capability of units without an online variable
936
937
938
939
940
941
942
943
944
945
946
    + (
        + ( p_gnu(grid, node, unit, 'maxGen') + p_gnu(grid, node, unit, 'maxCons') )${not uft_online(unit, f, t)}
        + sum(t_invest(t_)${ ord(t_)<=ord(t) },
            + v_invest_LP(unit, t_)${not uft_onlineLP(unit, f, t) and unit_investLP(unit)}
                * p_gnu(grid, node, unit, 'unitSizeTot')
            + v_invest_MIP(unit, t_)${not uft_onlineMIP(unit, f, t) and unit_investMIP(unit)}
                * p_gnu(grid, node, unit, 'unitSizeTot')
          )
      )
        * p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
947
948

    // Ramping capability of units with an online variable
949
    + (
950
951
        + v_online_LP(unit, s, f+df_central(f,t), t)${uft_onlineLP(unit, f, t)}
        + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
952
953
954
955
      )
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
956

957
958
959
960
    // Units in the run-up phase need to keep up with the run-up rate
    + p_gnu(grid, node, unit, 'unitSizeTot')
        * sum(unitStarttype(unit, starttype)${uft_startupTrajectory(unit, f, t)},
            sum(runUpCounter(unit, counter)${t_active(t+dt_trajectory(counter))}, // Sum over the run-up intervals
961
962
963
964
965
966
                + [
                    + v_startup_LP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineLP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    + v_startup_MIP(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                        ${ uft_onlineMIP(unit, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter)) }
                    ]
967
968
969
970
971
972
973
974
                    * [
                        + p_unit(unit, 'rampSpeedToMinLoad')
                        + ( p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit') - p_unit(unit, 'rampSpeedToMinLoad') )${ not runUpCounter(unit, counter+1) } // Ramp speed adjusted for the last run-up interval
                            * ( p_u_runUpTimeIntervalsCeil(unit) - p_u_runUpTimeIntervals(unit) )
                        ]
                    * 60 // Unit conversion from [p.u./min] into [p.u./h]
                ) // END sum(runUpCounter)
            ) // END sum(unitStarttype)
975
976

    // Shutdown of consumption units from full load
977
    + v_shutdown(unit, s, f, t)${uft_online(unit, f, t) and gnu_input(grid, node, unit)}
978
979
980
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
981

982
    // Shutdown of generation units from full load and ramping of units in the beginning of the shutdown phase
983
    + v_shutdown(unit, s, f, t)${uft_online(unit, f, t) and gnu_output(grid, node, unit)}
984
985
986
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
987

988
    // Units in shutdown phase need to keep up with the shutdown ramp rate
989
990
991
992
993
994
    + p_gnu(grid, node, unit, 'unitSizeGen')
        * [
            + sum(shutdownCounter(unit, counter)${t_active(t+dt_trajectory(counter)) and uft_shutdownTrajectory(unit, f, t)}, // Sum over the shutdown intervals
                + v_shutdown(unit, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                    * [
                        + p_unit(unit, 'rampSpeedFromMinLoad')
995
                        + ( p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit') - p_unit(unit, 'rampSpeedFromMinLoad') )${ not shutdownCounter(unit, counter-1) } // Ramp speed adjusted for the first shutdown interval
996
997
998
999
1000
                            * ( p_u_shutdownTimeIntervalsCeil(unit) - p_u_shutdownTimeIntervals(unit) )
                        ]
                ) // END sum(shutdownCounter)
            // Units need to be able to shut down after shut down trajectory
            + v_shutdown(unit, s, f+df(f, t+dt_toShutdown(unit, t)), t+dt_toShutdown(unit, t))