2d_constraints.gms 98.1 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
27
28
29
30
                                            } .. // Energy/power balance dynamics solved using implicit Euler discretization

    // 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(to_node${ 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(from_node${ 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(node_${ 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(node_${ 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 -----------------------------------------------

213
214
215
q_maxDownward(m, s, gnuft(grid, node, unit, f, t))${msft(m, s, f, t)
                                                    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
262
    // 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
                + v_startup(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
263
                    * p_uCounter_runUpMin(unit, counter)
264
265
266
267
268
269
270
                ) // 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))
271
                * p_uCounter_shutdownMin(unit, counter)
272
            ) // END sum(shutdownCounter)
273

274
275
276
277
278
    // Consuming units, greater than maxCons
    // Available capacity restrictions
    - p_unit(unit, 'availability')
        * [
            // Capacity factors for flow units
279
            + sum(flowUnit(flow, unit),
280
                + ts_cf_(flow, node, f, t, s)
281
282
283
284
285
286
287
288
                ) // 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')
                * [
289
                    // Capacity online
290
291
                    + 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)}
292
293
294
295
296
297
298
299

                    // 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)
300
301
                    ] // END * p_gnu(unitSizeCons)
            ] // END * p_unit(availability)
302
;
303
304
305

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

306
307
308
q_maxUpward(m, s, gnuft(grid, node, unit, f, t))${msft(m, s, f, t)
                                                    and {
                                                 [   ord(t) < tSolveFirst + smax(restype, p_nReserves(node, restype, 'reserve_length')) // Unit is either providing
309
310
311
312
313
314
315
316
317
318
319
320
321
                                                    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))
                                                    ]
322
                                                }}..
323
    // Energy generation/consumption
324
    + v_gen(grid, node, unit, s, f, t)
325
326
327
328

    // Considering output constraints (e.g. cV line)
    + sum(gngnu_constrainedOutputRatio(grid, node, grid_output, node_, unit),
        + p_gnu(grid_output, node_, unit, 'cV')
329
            * v_gen(grid_output, node_, unit, s, f, t)
330
331
332
        ) // END sum(gngnu_constrainedOutputRatio)

    // Upwards reserve participation
333
    + sum(nuRescapable(restype, 'up', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
334
        + v_reserve(restype, 'up', node, unit, s, f+df_reserves(node, restype, f, t), t)
335
336
337
338
339
340
        ) // END sum(nuRescapable)

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

    // Consuming units
    + p_gnu(grid, node, unit, 'unitSizeCons')
341
        * sum(suft(effGroup, unit, f, t), // Uses the minimum 'lb' for the current efficiency approximation
342
343
344
345
            + p_effGroupUnit(effGroup, unit, 'lb')${not ts_effGroupUnit(effGroup, unit, 'lb', f, t)}
            + ts_effGroupUnit(effGroup, unit, 'lb', f, t)
            ) // END sum(effGroup)
        * [
346
347
            + 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)
348
349
350
351
352
353
354
            ] // 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
355
            + sum(flowUnit(flow, unit),
356
                + ts_cf_(flow, node, f, t, s)
357
358
359
360
361
362
363
364
                ) // 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')
                * [
365
                    // Capacity online
366
367
                    + 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)}
368
369
370
371
372
373
374
375

                    // 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)
376
377
                    ] // END * p_gnu(unitSizeGen)
            ] // END * p_unit(availability)
378

379
380
381
382
383
    // 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
                + v_startup(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
384
                    * p_uCounter_runUpMax(unit, counter)
385
386
387
388
389
390
391
                ) // 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))
392
                * p_uCounter_shutdownMax(unit, counter)
393
            ) // END sum(shutdownCounter)
394
;
395

396
397
* --- Reserve Provision of Units with Investments -----------------------------

398
399
400
401
402
403
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)
404
405
406
407
408
409
410
411
412
413
414
415
416

    =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_)
            ]
417
418
419
420
        * p_unit(unit, 'availability') // Taking into account availability...
        * [
            // ... and capacity factor for flow units
            + sum(flowUnit(flow, unit),
421
                + ts_cf_(flow, node, f, t, s)
422
423
                ) // END sum(flow)
            + 1${not unit_flow(unit)}
424
425
426
            ] // How to consider reserveReliability in the case of investments when we typically only have "realized" time steps?
;

427
428
* --- Unit Startup and Shutdown -----------------------------------------------

429
q_startshut(m, s, uft_online(unit, f, t))$msft(m, s, f, t) ..
430
    // Units currently online
431
432
    + 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)}
433
434

    // Units previously online
435
436

    // The same units
437
    - 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))
438
                                                             and not uft_aggregator_first(unit, f, t) } // This reaches to tFirstSolve when dt = -1
439
    - 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))
440
441
442
443
                                                             and not uft_aggregator_first(unit, f, t) }

    // Aggregated units just before they are turned into aggregator units
    - sum(unit_${unitAggregator_unit(unit, unit_)},
444
445
        + 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))}
446
        )${uft_aggregator_first(unit, f, t)} // END sum(unit_)
447

448
449
    =E=

450
    // Unit startup and shutdown
451

452
    // Add startup of units dt_toStartup before the current t (no start-ups for aggregator units before they become active)
453
    + sum(unitStarttype(unit, starttype),
454
        + v_startup(unit, starttype, s, f+df(f,t+dt_toStartup(unit, t)), t+dt_toStartup(unit, t))
455
        )${not [unit_aggregator(unit) and ord(t) + dt_toStartup(unit, t) <= tSolveFirst + p_unit(unit, 'lastStepNotAggregated')]} // END sum(starttype)
456

457
458
459
460
    // 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
461

462
    // Shutdown of units at time t
463
    - v_shutdown(unit, s, f, t)
464
;
465

466
*--- Startup Type -------------------------------------------------------------
467
// !!! NOTE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
468
469
470
// 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.
471

472
473
q_startuptype(m, s, starttypeConstrained(starttype), uft_online(unit, f, t))
    ${msft(m, s, f, t) and unitStarttype(unit, starttype)} ..
474
475

    // Startup type
476
    + v_startup(unit, starttype, s, f, t)
477
478
479
480

    =L=

    // Subunit shutdowns within special startup timeframe
Topi Rasku's avatar
Topi Rasku committed
481
    + sum(unitCounter(unit, counter)${dt_starttypeUnitCounter(starttype, unit, counter)},
482
        + 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
483
            ${t_active(t+(dt_starttypeUnitCounter(starttype, unit, counter)+1))}
484
485
486
        ) // END sum(counter)

    // NOTE: for aggregator units, shutdowns for aggregated units are not considered
487
;
488

489

490
491
*--- Online Limits with Startup Type Constraints and Investments --------------

492
493
q_onlineLimit(m, s, uft_online(unit, f, t))${msft(m, s, f, t) and {
                                            p_unit(unit, 'minShutdownHours')
494
                                            or p_u_runUpTimeIntervals(unit)
495
496
                                            or unit_investLP(unit)
                                            or unit_investMIP(unit)
497
                                            }} ..
498
    // Online variables
499
500
    + 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)}
501
502
503
504
505
506

    =L=

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

507
    // Number of units unable to become online due to restrictions
Topi Rasku's avatar
Topi Rasku committed
508
    - sum(unitCounter(unit, counter)${dt_downtimeUnitCounter(unit, counter)},
509
        + 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
510
            ${t_active(t+(dt_downtimeUnitCounter(unit, counter) + 1))}
511
512
513
        ) // 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)
514
    - sum(unitAggregator_unit(unit, unit_),
Topi Rasku's avatar
Topi Rasku committed
515
        + sum(unitCounter(unit, counter)${dt_downtimeUnitCounter(unit, counter)},
516
            + 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
517
                ${t_active(t+(dt_downtimeUnitCounter(unit, counter) + 1))}
518
519
            ) // END sum(counter)
        )${unit_aggregator(unit)} // END sum(unit_)
520
521
522

    // Investments into units
    + sum(t_invest(t_)${ord(t_)<=ord(t)},
523
524
        + v_invest_LP(unit, t_)${unit_investLP(unit)}
        + v_invest_MIP(unit, t_)${unit_investMIP(unit)}
525
526
527
        ) // END sum(t_invest)
;

528
529
530
531
*--- 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).
532
533
q_onlineOnStartUp(s, uft_online(unit, f, t))
    ${sft(s, f, t) and sum(starttype, unitStarttype(unit, starttype))}..
534
535

    // Units currently online
536
537
    + 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)}
538
539
540
541

    =G=

    + sum(unitStarttype(unit, starttype),
542
        + v_startup(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
543
544
545
      ) // END sum(starttype)
;

546
547
q_offlineAfterShutdown(s, uft_online(unit, f, t))
    ${sft(s, f, t) and sum(starttype, unitStarttype(unit, starttype))}..
548

549
550
551
552
553
    // Number of existing units
    + p_unit(unit, 'unitCount')

    // Investments into units
    + sum(t_invest(t_)${ord(t_)<=ord(t)},
554
555
        + v_invest_LP(unit, t_)${unit_investLP(unit)}
        + v_invest_MIP(unit, t_)${unit_investMIP(unit)}
556
557
        ) // END sum(t_invest)

558
    // Units currently online
559
560
    - 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)}
561
562
563

    =G=

564
    + v_shutdown(unit, s, f, t)
565
566
;

567
568
*--- Minimum Unit Uptime ------------------------------------------------------

569
570
q_onlineMinUptime(m, s, uft_online(unit, f, t))
    ${msft(m, s, f, t) and  p_unit(unit, 'minOperationHours')} ..
571
572

    // Units currently online
573
574
    + 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)}
575
576
577
578

    =G=

    // Units that have minimum operation time requirements active
Topi Rasku's avatar
Topi Rasku committed
579
    + sum(unitCounter(unit, counter)${dt_uptimeUnitCounter(unit, counter)},
580
        + sum(unitStarttype(unit, starttype),
581
            + v_startup(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))
Topi Rasku's avatar
Topi Rasku committed
582
                ${t_active(t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1))}
583
            ) // END sum(starttype)
584
585
586
        ) // 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
587
588
    + sum(unitAggregator_unit(unit, unit_),
        + sum(unitCounter(unit, counter)${dt_uptimeUnitCounter(unit, counter)},
589
            + sum(unitStarttype(unit, starttype),
590
                + v_startup(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))
Topi Rasku's avatar
Topi Rasku committed
591
                    ${t_active(t+(dt_uptimeUnitCounter(unit, counter)+dt_toStartup(unit, t) + 1))}
592
593
594
                ) // END sum(starttype)
            ) // END sum(counter)
        )${unit_aggregator(unit)} // END sum(unit_)
595
596
;

597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
* --- 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)
;

623
* --- Ramp Constraints --------------------------------------------------------
624
625
626
627

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

629
630
    + v_genRamp(grid, node, unit, s, f, t)
        * p_stepLength(m, f, t)
631

632
    =E=
633

634
    // Change in generation over the interval: v_gen(t) - v_gen(t-1)
635
    + v_gen(grid, node, unit, s, f, t)
636

637
    // Unit generation at t-1 (except aggregator units right before the aggregation threshold, see next term)
638
    - 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)}
639
640
    // Unit generation at t-1, aggregator units right before the aggregation threshold
    + sum(unit_${unitAggregator_unit(unit, unit_)},
641
        - v_gen(grid, node, unit_, s+ds(s,t), f+df(f,t+dt(t)), t+dt(t))
642
      )${uft_aggregator_first(unit, f, t)}
643
;
644

645
* --- Ramp Up Limits ----------------------------------------------------------
646
647
648
649

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')
650
651
652
653
654
                                                           and [ sum(restype, nuRescapable(restype, 'up', node, unit))
                                                                 or uft_online(unit, f, t)
                                                                 or unit_investLP(unit)
                                                                 or unit_investMIP(unit)
                                                                 ]
655
                                                           } ..
656
    + v_genRamp(grid, node, unit, s, f, t)
657
    + sum(nuRescapable(restype, 'up', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
658
        + 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)
659
660
661
662
663
        ) // END sum(nuRescapable)
        / p_stepLength(m, f, t)

    =L=

664
    // Ramping capability of units without an online variable
665
666
667
668
669
670
671
672
673
674
675
676
    + (
        + ( 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]

677
    // Ramping capability of units with an online variable
678
    + (
679
680
        + 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)}
681
682
683
684
685
      )
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnu(grid, node, unit, 'maxRampUp')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]

686
687
688
689
690
691
692
693
694
695
696
697
698
    // 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
                + v_startup(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                    * [
                        + 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)
699

700
    // Shutdown of consumption units from full load
701
    + v_shutdown(unit, s, f, t)${uft_online(unit, f, t) and gnu_input(grid, node, unit)}
702
        * p_gnu(grid, node, unit, 'unitSizeTot')
703
;
704

705
* --- Ramp Down Limits --------------------------------------------------------
706
707
708
709

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')
710
711
712
713
714
                                                             and [ sum(restype, nuRescapable(restype, 'down', node, unit))
                                                                   or uft_online(unit, f, t)
                                                                   or unit_investLP(unit)
                                                                   or unit_investMIP(unit)
                                                                   ]
715
                                                             } ..
716
    + v_genRamp(grid, node, unit, s, f, t)
717
    - sum(nuRescapable(restype, 'down', node, unit)${ord(t) < tSolveFirst + p_nReserves(node, restype, 'reserve_length')},
718
        + 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)
719
720
721
722
723
        ) // END sum(nuRescapable)
        / p_stepLength(m, f, t)

    =G=

724
    // Ramping capability of units without online variable
725
726
727
728
729
730
731
732
733
734
735
736
    - (
        + ( 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]

737
    // Ramping capability of units that are online
738
    - (
739
740
        + 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)}
741
742
743
744
745
      )
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnu(grid, node, unit, 'maxRampDown')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]

746
    // Shutdown of generation units from full load
747
    - v_shutdown(unit, s, f, t)${   uft_online(unit, f, t)
748
749
750
                                    and gnu_output(grid, node, unit)
                                    and not uft_shutdownTrajectory(unit, f, t)
                                    }
751
        * p_gnu(grid, node, unit, 'unitSizeTot')
752

753
754
    // Units in shutdown phase need to keep up with the shutdown ramp rate
    - p_gnu(grid, node, unit, 'unitSizeGen')
755
756
757
758
759
760
761
762
763
        * [
            + 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)
764
765
            // 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))
766
767
768
                * p_unit(unit, 'rampSpeedFromMinload')
            ]
        * 60 // Unit conversion from [p.u./min] to [p.u./h]
769
770
;

771
772
773
774
775
776
777
* --- 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))
                                                          } ..

778
    + v_genRamp(grid, node, unit, s, f, t)
779

780
    =E=
781

782
783
    // Upward and downward ramp categories
    + sum(slack${ gnuft_rampCost(grid, node, unit, slack, f, t) },
784
785
        + v_genRampUpDown(grid, node, unit, slack, s, f, t)$upwardSlack(slack)
        - v_genRampUpDown(grid, node, unit, slack, s, f, t)$downwardSlack(slack)
786
      ) // END sum(slack)
787
788
;

Niina Helistö's avatar
Niina Helistö committed
789
* --- Upward and downward ramps constrained by slack boundaries ---------------
790
791
792
793
794

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

795
    + v_genRampUpDown(grid, node, unit, slack, s, f, t)
796

797
    =L=
798
799

    // Ramping capability of units without an online variable
800
801
802
803
804
805
806
807
808
809
810
    + (
        + ( 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]
811
812

    // Ramping capability of units with an online variable
813
    + (
814
815
        + 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)}
816
817
818
819
      )
        * p_gnu(grid, node, unit, 'unitSizeTot')
        * p_gnuBoundaryProperties(grid, node, unit, slack, 'rampLimit')
        * 60   // Unit conversion from [p.u./min] to [p.u./h]
820

821
822
823
824
825
826
827
828
829
830
831
832
833
    // 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
                + v_startup(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
                    * [
                        + 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)
834
835

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

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

847
    // Units in shutdown phase need to keep up with the shutdown ramp rate
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
    + 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')
                        + ( 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)
            // 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))
                * p_unit(unit, 'rampSpeedFromMinload')
            ]
        * 60 // Unit conversion from [p.u./min] to [p.u./h]
863
;
864

865
866
* --- Fixed Output Ratio ------------------------------------------------------

867
q_outputRatioFixed(gngnu_fixedOutputRatio(grid, node, grid_, node_, unit), sft(s, f, t))${  uft(unit, f, t)
868
869
870
                                                                                        } ..

    // Generation in grid
871
    + v_gen(grid, node, unit, s, f, t)
872
        / p_gnu(grid, node, unit, 'conversionFactor')
873
874
875
876

    =E=

    // Generation in grid_
877
    + v_gen(grid_, node_, unit, s, f, t)
878
        / p_gnu(grid_, node_, unit, 'conversionFactor')
879
;
880
881
882

* --- Constrained Output Ratio ------------------------------------------------

883
q_outputRatioConstrained(gngnu_constrainedOutputRatio(grid, node, grid_, node_, unit), sft(s, f, t))${  uft(unit, f, t)
884
885
886
                                                                                                    } ..

    // Generation in grid
887
    + v_gen(grid, node, unit, s, f, t)
888
        / p_gnu(grid, node, unit, 'conversionFactor')
889
890
891
892

    =G=

    // Generation in grid_
893
    + v_gen(grid_, node_, unit, s, f, t)
894
        / p_gnu(grid_, node_, unit, 'conversionFactor')
Juha Kiviluoma's avatar
Juha Kiviluoma committed
895
;
896
897
898

* --- Direct Input-Output Conversion ------------------------------------------

899
q_conversionDirectInputOutput(s, suft(effDirect(effGroup), unit, f, t))$sft(s, f, t) ..
900
901

    // Sum over endogenous energy inputs
902
    - sum(gnu_input(grid, node, unit)${not p_gnu(grid, node, unit, 'doNotOutput')},
903
        + v_gen(grid, node, unit, s, f, t)
904
905
906
907
        ) // END sum(gnu_input)

    // Sum over fuel energy inputs
    + sum(uFuel(unit, 'main', fuel),
908
        + v_fuelUse(fuel, unit, s, f, t)
909
910
911
912
913
914
        ) // END sum(uFuel)

    =E=

    // Sum over energy outputs
    + sum(gnu_output(grid, node, unit),
915
        + v_gen(grid, node, unit, s, f, t)
916
            * [ // efficiency rate
917
                + p_effUnit(effGroup, unit, effGroup, 'slope')${ not ts_effUnit(effGroup, unit, effGroup, 'slope', f, t) }
918
                + ts_effUnit(effGroup, unit, effGroup, 'slope', f, t)
919
920
921
                ] // END * v_gen
        ) // END sum(gnu_output)

922
    // Consumption of keeping units online (no-load fuel use)
923
924
925
    + sum(gnu_output(grid, node, unit),
        + p_gnu(grid, node, unit, 'unitSizeGen')
        ) // END sum(gnu_output)
926
        * [ // Unit online state
927
928
            + 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)}
929
930
931
932
933
934

            // Run-up and shutdown phase efficiency correction
            // Run-up 'online state'
            + 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
                    + v_startup(unit, starttype, s, f+df(f, t+dt_trajectory(counter)), t+dt_trajectory(counter))
935
                        * p_uCounter_runUpMin(unit, counter)
936
937
938
939
940
941
                        / p_unit(unit, 'op00') // Scaling the p_uCounter_runUp using minload
                    ) // END sum(runUpCounter)
                ) // END sum(unitStarttype)
            // Shutdown 'online state'
            + 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))
942
                    * p_uCounter_shutdownMin(unit, counter)
943
944
                        / p_unit(unit, 'op00') // Scaling the p_uCounter_shutdown using minload
                ) // END sum(shutdownCounter)
945
946
            ] // END * sum(gnu_output)
        * [
947
948
            + p_effGroupUnit(effGroup, unit, 'section')${not ts_effUnit(effGroup, unit, effDirect, 'section', f, t)}
            + ts_effUnit(effGroup, unit, effGroup, 'section', f, t)
949
            ] // END * sum(gnu_output)
950
;
951
952
953

* --- SOS2 Efficiency Approximation -------------------------------------------

954
q_conversionSOS2InputIntermediate(s, suft(effLambda(effGroup), unit, f, t))$sft(s, f, t) ..
955
956

    // Sum over endogenous energy inputs
957
    - sum(gnu_input(grid, node, unit)${not p_gnu(grid, node, unit, 'doNotOutput')},
958
        + v_gen(grid, node, unit, s, f, t)
959
960
961
962
        ) // END sum(gnu_input)

    // Sum over fuel energy inputs
    + sum(uFuel(unit, 'main', fuel),
963
        + v_fuelUse(fuel, unit, s, f, t)
964
965
        ) // END sum(uFuel)

966
    =E=
967
968
969
970
971

    // Sum over the endogenous outputs of the unit
    + sum(gnu_output(grid, node, unit), p_gnu(grid, node, unit, 'unitSizeGen'))
        * [
            // Consumption of generation
972
            + sum(effGroupSelectorUnit(effGroup, unit, effSelector),
973
                + v_sos2(unit, s, f, t, effSelector)
974
975
976
977
978
979
980
981
982
                    * [ // Operation points convert the v_sos2 variables into share of capacity used for generation
                        + p_effUnit(effGroup, unit, effSelector, 'op')${not ts_effUnit(effGroup, unit, effSelector, 'op', f, t)}
                        + ts_effUnit(effGroup, unit, effSelector, 'op', f, t)
                        ] // END * v_sos2
                    * [ // Heat rate
                        + p_effUnit(effGroup, unit, effSelector, 'slope')${not ts_effUnit(effGroup, unit, effSelector, 'slope', f, t)}
                        + ts_effUnit(effGroup, unit, effSelector, 'slope', f, t)
                        ] // END * v_sos2
                ) // END sum(effSelector)
983
           ]
984
;
985
986
987

* --- SOS 2 Efficiency Approximation Online Variables -------------------------

988
q_conversionSOS2Constraint(s, suft(effLambda(effGroup), unit, f, t))$sft(s, f, t) ..
989
990

    // Total value of the v_sos2 equals the number of online units
991
    + sum(effGroupSelectorUnit(effGroup, unit, effSelector),
992
        + v_sos2(unit, s, f, t, effSelector)
993
994
995
996
997
        ) // END sum(effSelector)

    =E=

    // Number of units online
998
    + v_online_MIP(unit, s, f+df_central(f,t), t)${uft_onlineMIP(unit, f, t)}
999
1000

    // Run-up and shutdown phase efficiency approximation