3c_inputsLoop.gms 24 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
$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
17

18
* =============================================================================
19
* --- Update the Forecast Data ------------------------------------------------
20
* =============================================================================
21

22
put log 'ord tSolve: ';
23
put log ord(tSolve):0:0 /;
24
putclose log;
25

26
27
// Determine the necessary horizon for updating data
option clear = tmp;
28
29
tmp = max(  mSettings(mSolve, 't_forecastLengthUnchanging'),
            mSettings(mSolve, 't_forecastLengthDecreasesFrom')
30
31
32
33
34
35
36
            );

// Find time steps until the forecast horizon
option clear = tt_forecast;
tt_forecast(t_current(t))
    ${ ord(t) <= tSolveFirst + tmp }
    = yes;
37

38
if (ord(tSolve) = tForecastNext(mSolve) - mSettings(mSolve, 't_forecastJump'), // tForecastNext updated already in periodicLoop!
39

40
41
42
    // Update ts_unit
    if (mTimeseries_loop_read(mSolve, 'ts_unit'),
        put_utility 'gdxin' / '%input_dir%/ts_unit/' tSolve.tl:0 '.gdx';
43
        execute_load ts_unit_update=ts_unit;
44
        ts_unit(unit_timeseries(unit), param_unit, f_solve(f), tt_forecast(t)) // Only update if time series enabled for the unit
45
            ${not mf_realization(mSolve, f) // Realization not updated
46
              and (mSettings(mSolve, 'onlyExistingForecasts')
47
                   -> ts_unit_update(unit, param_unit, f, t)) // Update only existing values (zeroes need to be EPS)
48
                }
49
            = ts_unit_update(unit, param_unit, f, t);
50
    ); // END if('ts_unit')
51
52
53
54
$ontext
* !!! NOTE !!!
* These probably shouldn't be read at all, as p_effUnit and p_effGroupUnit are
* not input data, but calculated based on p_unit
55
    // Update ts_effUnit
56
    if (mTimeseries_loop_read(mSolve, 'ts_effUnit'),
57
        put_utility 'gdxin' / '%input_dir%/ts_effUnit/' tSolve.tl:0 '.gdx';
58
        execute_load ts_effUnit_update=ts_effUnit;
59
        ts_effUnit(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f_solve(f), tt_forecast(t)) // Only update if time series enabled for the unit
60
            ${  not mf_realization(mSolve, f) // Realization not updated
61
*               and ts_effUnit_update(effSelector, unit, effSelector, param_eff, f, t) // Update only existing values (zeroes need to be EPS)
62
                }
63
            = ts_effUnit_update(effSelector, unit, effSelector, param_eff, f, t);
64
    ); // END if('ts_effUnit')
65
66

    // Update ts_effGroupUnit
67
    if (mTimeseries_loop_read(mSolve, 'ts_effGroupUnit'),
68
        put_utility 'gdxin' / '%input_dir%/ts_effGroupUnit/' tSolve.tl:0 '.gdx';
69
        execute_load ts_effGroupUnit_update=ts_effGroupUnit;
70
        ts_effGroupUnit(effSelector, unit_timeseries(unit), param_eff, f_solve(f), tt_forecast(t)) // Only update if time series enabled for the unit
71
            ${  not mf_realization(mSolve, f) // Realization not updated
72
*               and ts_effGroupUnit_update(effSelector, unit, param_eff, f, t) // Update only existing values (zeroes need to be EPS)
73
                }
74
            = ts_effGroupUnit_update(effSelector, unit, param_eff, f, t);
75
76
    ); // END if('ts_effGroupUnit')
$offtext
77
    // Update ts_influx
78
    if (mTimeseries_loop_read(mSolve, 'ts_influx'),
79
        put_utility 'gdxin' / '%input_dir%/ts_influx/' tSolve.tl:0 '.gdx';
80
81
        execute_load ts_influx_update=ts_influx;
        ts_influx(gn(grid, node), f_solve(f), tt_forecast(t))
82
            ${  not mf_realization(mSolve, f) // Realization not updated
83
                and (mSettings(mSolve, 'onlyExistingForecasts')
84
                     -> ts_influx_update(grid, node, f, t)) // Update only existing values (zeroes need to be EPS)
85
                }
86
87
            = ts_influx_update(grid, node, f, t);
    ); // END if('ts_influx')
88
89

    // Update ts_cf
90
    if (mTimeseries_loop_read(mSolve, 'ts_cf'),
91
        put_utility 'gdxin' / '%input_dir%/ts_cf/' tSolve.tl:0 '.gdx';
92
93
        execute_load ts_cf_update=ts_cf;
        ts_cf(flowNode(flow, node), f_solve(f), tt_forecast(t))
94
            ${  not mf_realization(mSolve, f) // Realization not updated
95
                and (mSettings(mSolve, 'onlyExistingForecasts')
96
                     -> ts_cf_update(flow, node, f, t)) // Update only existing values (zeroes need to be EPS)
97
                }
98
99
            = ts_cf_update(flow, node, f, t);
    ); // END if('ts_cf')
100
101

    // Update ts_reserveDemand
102
    if (mTimeseries_loop_read(mSolve, 'ts_reserveDemand'),
103
        put_utility 'gdxin' / '%input_dir%/ts_reserveDemand/' tSolve.tl:0 '.gdx';
104
105
        execute_load ts_reserveDemand_update=ts_reserveDemand;
        ts_reserveDemand(restypeDirectionNode(restype, up_down, node), f_solve(f), tt_forecast(t))
106
            ${  not mf_realization(mSolve, f) // Realization not updated
107
                and (mSettings(mSolve, 'onlyExistingForecasts')
108
                     -> ts_reserveDemand_update(restype, up_down, node, f, t)) // Update only existing values (zeroes need to be EPS)
109
                }
110
111
            = ts_reserveDemand_update(restype, up_down, node, f, t);
    ); // END if('ts_reserveDemand')
112
113

    // Update ts_node
114
    if (mTimeseries_loop_read(mSolve, 'ts_node'),
115
        put_utility 'gdxin' / '%input_dir%/ts_node/' tSolve.tl:0 '.gdx';
116
117
        execute_load ts_node_update=ts_node;
        ts_node(gn(grid, node), param_gnBoundaryTypes, f_solve(f), tt_forecast(t))
118
            ${  not mf_realization(mSolve, f) // Realization not updated
119
                and (mSettings(mSolve, 'onlyExistingForecasts')
120
                     -> ts_node_update(grid, node, param_gnBoundaryTypes, f ,t)) // Update only existing values (zeroes need to be EPS)
121
                }
122
123
124
125
126
            = ts_node_update(grid, node, param_gnBoundaryTypes, f, t);
    ); // END if('ts_node')

* --- NO FORECAST DIMENSION, SHOULD THESE BE HANDLED SEPARATELY? --------------
// Currently, only updated until the forecast horizon, but is this correct?
127
128

    // Update ts_fuelPriceChange
129
    if (mTimeseries_loop_read(mSolve, 'ts_fuelPriceChange'),
130
        put_utility 'gdxin' / '%input_dir%/ts_fuelPriceChange/' tSolve.tl:0 '.gdx';
131
132
        execute_load ts_fuelPriceChange_update=ts_fuelPriceChange;
        ts_fuelPriceChange(fuel, tt_forecast(t))
133
*            ${ ts_fuelPriceChange_update(fuel, t) } // Update only existing values (zeroes need to be EPS)
134
135
            = ts_fuelPriceChange_update(fuel, t);
    ); // END if('ts_fuelPriceChange')
136
137

    // Update ts_unavailability
138
    if (mTimeseries_loop_read(mSolve, 'ts_unavailability'),
139
        put_utility 'gdxin' / '%input_dir%/ts_unavailability/' tSolve.tl:0 '.gdx';
140
141
        execute_load ts_unavailability_update=ts_unavailability;
        ts_unavailability(unit, tt_forecast(t))
142
*            ${ ts_unavailability_update(unit, t) } // Update only existing values (zeroes need to be EPS)
143
144
            = ts_unavailability_update(unit, t);
    ); // END if('ts_unavailability')
145

146
); // END if(tForecastNext)
147

148
* =============================================================================
149
* --- Optional forecast improvement -------------------------------------------
150
151
* =============================================================================

152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
if(mSettings(mSolve, 't_improveForecast'),
* Linear improvement of the central forecast towards the realized forecast,
* while preserving the difference between the central forecast and the
* remaining forecasts

    // Determine the set of improved time steps
    option clear = tt;
    tt(tt_forecast(t))
        ${ ord(t) <= tSolveFirst + mSettings(mSolve, 't_improveForecast') }
        = yes;

    // Temporary forecast displacement to reach the central forecast
    option clear = ddf;
    ddf(f_solve(f))
        ${ not mf_central(mSolve, f) }
        = sum(mf_central(mSolve, f_), ord(f_) - ord(f));

    // Temporary forecast displacement to reach the realized forecast
    option clear = ddf_;
    ddf_(f_solve(f))
        ${ not mf_realization(mSolve, f) }
        = sum(mf_realization(mSolve, f_), ord(f_) - ord(f));

* --- Calculate the other forecasts relative to the central one ---------------

    loop(f_solve(f)${ not mf_realization(mSolve, f) and not mf_central(mSolve, f) },
        // ts_unit
179
180
181
182
        ts_unit(unit_timeseries(unit), param_unit, f, tt(t))// Only update for units with time series enabled
            = ts_unit(unit, param_unit, f, t) - ts_unit(unit, param_unit, f+ddf(f), t);
$ontext
* Should these be handled here at all? See above note
183
        // ts_effUnit
184
185
        ts_effUnit(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f, tt(t)) // Only update for units with time series enabled
            = ts_effUnit(effSelector, unit, effSelector, param_eff, f, t) - ts_effUnit(effSelector, unit, effSelector, param_eff, f+ddf(f), t);
186
        // ts_effGroupUnit
187
188
        ts_effGroupUnit(effSelector, unit_timeseries(unit), param_eff, f, tt(t)) // Only update for units with time series enabled
            = ts_effGroupUnit(effSelector, unit, param_eff, f, t) - ts_effGroupUnit(effSelector, unit, param_eff, f+ddf(f), t);
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
$offtext
        // ts_influx
        ts_influx(gn(grid, node), f, tt(t))
            = ts_influx(grid, node, f, t) - ts_influx(grid, node, f+ddf(f), t);
        // ts_cf
        ts_cf(flowNode(flow, node), f, tt(t))
            = ts_cf(flow, node, f, t) - ts_cf(flow, node, f+ddf(f), t);
        // ts_reserveDemand
        ts_reserveDemand(restypeDirectionNode(restype, up_down, node), f, tt(t))
            = ts_reserveDemand(restype, up_down, node, f, t) - ts_reserveDemand(restype, up_down, node, f+ddf(f), t);
        // ts_node
        ts_node(gn(grid, node), param_gnBoundaryTypes, f, tt(t))
            = ts_node(grid, node, param_gnBoundaryTypes, f, t) - ts_node(grid, node, param_gnBoundaryTypes, f+ddf(f), t);
    ); // END loop(f_solve)

* --- Linear improvement of the central forecast ------------------------------

    loop(mf_central(mSolve, f),
        // ts_unit
208
        ts_unit(unit_timeseries(unit), param_unit, f, tt(t)) // Only update for units with time series enabled
209
            = [ + (ord(t) - tSolveFirst)
210
                    * ts_unit(unit, param_unit, f, t)
211
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
212
                    * ts_unit(unit, param_unit, f+ddf_(f), t)
213
                ] / mSettings(mSolve, 't_improveForecast');
214
215
$ontext
* Should these be handled here at all? See above note
216
        // ts_effUnit
217
        ts_effUnit(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f, tt(t)) // Only update for units with time series enabled
218
            = [ + (ord(t) - tSolveFirst)
219
                    * ts_effUnit(effSelector, unit, effSelector, param_eff, f, t)
220
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
221
                    * ts_effUnit(effSelector, unit, effSelector, param_eff, f+ddf_(f), t)
222
223
                ] / mSettings(mSolve, 't_improveForecast');
        // ts_effGroupUnit
224
        ts_effGroupUnit(effSelector, unit_timeseries(unit), param_eff, f, tt(t)) // Only update for units with time series enabled
225
            = [ + (ord(t) - tSolveFirst)
226
                    * ts_effGroupUnit(effSelector, unit, param_eff, f, t)
227
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
228
                    * ts_effGroupUnit(effSelector, unit, param_eff, f+ddf_(f), t)
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
                ] / mSettings(mSolve, 't_improveForecast');
$offtext
        // ts_influx
        ts_influx(gn(grid, node), f, tt(t))
            = [ + (ord(t) - tSolveFirst)
                    * ts_influx(grid, node, f, t)
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
                    * ts_influx(grid, node, f+ddf_(f), t)
                ] / mSettings(mSolve, 't_improveForecast');
        // ts_cf
        ts_cf(flowNode(flow, node), f, tt(t))
            = [ + (ord(t) - tSolveFirst)
                    * ts_cf(flow, node, f, t)
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
                    * ts_cf(flow, node, f+ddf_(f), t)
                ] / mSettings(mSolve, 't_improveForecast');
        // ts_reserveDemand
        ts_reserveDemand(restypeDirectionNode(restype, up_down, node), f, tt(t))
            = [ + (ord(t) - tSolveFirst)
                    * ts_reserveDemand(restype, up_down, node, f, t)
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
                    * ts_reserveDemand(restype, up_down, node, f+ddf_(f), t)
                ] / mSettings(mSolve, 't_improveForecast');
        // ts_node
        ts_node(gn(grid, node), param_gnBoundaryTypes, f, tt(t))
            = [ + (ord(t) - tSolveFirst)
                    * ts_node(grid, node, param_gnBoundaryTypes, f, t)
                + (tSolveFirst - ord(t) + mSettings(mSolve, 't_improveForecast'))
                    * ts_node(grid, node, param_gnBoundaryTypes, f+ddf_(f), t)
                ] / mSettings(mSolve, 't_improveForecast');
    ); // END loop(mf_central)

* --- Recalculate the other forecasts based on the improved central forecast --

    loop(f_solve(f)${ not mf_realization(mSolve, f) and not mf_central(mSolve, f) },
        // ts_unit
265
266
267
268
        ts_unit(unit_timeseries(unit), param_unit, f, tt(t)) // Only update for units with time series enabled
            = ts_unit(unit, param_unit, f, t) + ts_unit(unit, param_unit, f+ddf(f), t);
$ontext
* Should these be handled here at all? See above note
269
        // ts_effUnit
270
271
        ts_effUnit(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f, tt(t)) // Only update for units with time series enabled
            = ts_effUnit(effSelector, unit, effSelector, param_eff, f, t) + ts_effUnit(effSelector, unit, effSelector, param_eff, f+ddf(f), t);
272
        // ts_effGroupUnit
273
274
        ts_effGroupUnit(effSelector, unit_timeseries(unit), param_eff, f, tt(t)) // Only update for units with time series enabled
            = ts_effGroupUnit(effSelector, unit, param_eff, f, t) + ts_effGroupUnit(effSelector, unit, param_eff, f+ddf(f), t);
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
$offtext
        // ts_influx
        ts_influx(gn(grid, node), f, tt(t))
            = ts_influx(grid, node, f, t) + ts_influx(grid, node, f+ddf(f), t);
        // ts_cf
        ts_cf(flowNode(flow, node), f, tt(t))
            = max(min(ts_cf(flow, node, f, t) + ts_cf(flow, node, f+ddf(f), t), 1), 0); // Ensure that capacity factor forecasts remain between 0-1
        // ts_reserveDemand
        ts_reserveDemand(restypeDirectionNode(restype, up_down, node), f, tt(t))
            = max(ts_reserveDemand(restype, up_down, node, f, t) + ts_reserveDemand(restype, up_down, node, f+ddf(f), t), 0); // Ensure that reserve demand forecasts remains positive
        // ts_node
        ts_node(gn(grid, node), param_gnBoundaryTypes, f, tt(t))
            = ts_node(grid, node, param_gnBoundaryTypes, f, t) + ts_node(grid, node, param_gnBoundaryTypes, f+ddf(f), t);
    ); // END loop(f_solve)

); // END if(t_improveForecast)
291

292
293
294
295
296
297
298
299
300
301
302
* =============================================================================
* --- Aggregate time series data for the time intervals -----------------------
* =============================================================================

// Loop over the defined blocks of intervals
loop(cc(counter),

    // Retrieve interval block time steps
    option clear = tt_interval;
    tt_interval(t) = tt_block(counter, t);

303
304
    // Select and average time series data matching the intervals
    ts_unit_(unit_timeseries(unit), param_unit, ft(f, tt_interval(t)))
305
306
        = sum(tt_aggregate(t, t_),
            ts_unit(unit, param_unit, f, t_+dt_circular(t_))
307
308
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);
309
310
$ontext
* Should these be handled here at all? See above comment
311
    ts_effUnit_(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, ft(f, tt_interval(t)))
312
313
        = sum(tt_aggregate(t, t_),
            ts_effUnit(effSelector, unit, effSelector, param_eff, f, t_+dt_circular(t_))
314
315
316
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);
    ts_effGroupUnit_(effSelector, unit_timeseries(unit), param_eff, ft(f, tt_interval(t)))
317
318
        = sum(tt_aggregate(t, t_),
            ts_effGroupUnit(effSelector, unit, param_eff, f, t_+dt_circular(t_))
319
320
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);
321
$offtext
322
    ts_influx_(gn(grid, node), fts(f, tt_interval(t), s))
323
        = sum(tt_aggregate(t, t_),
324
325
            ts_influx(grid, node,
                f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_influx')),
326
327
                t_+ (   + dt_scenarioOffset(grid, node, 'ts_influx', s)
                        + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_influx'))))
328
329
330
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);
    ts_cf_(flowNode(flow, node), fts(f, tt_interval(t), s))
331
        = sum(tt_aggregate(t, t_),
332
333
            ts_cf(flow, node,
                f + (df_scenario(f, t)$gn_scenarios(flow, node, 'ts_cf')),
334
335
                t_+ (   + dt_scenarioOffset(flow, node, 'ts_cf', s)
                        + dt_circular(t_)$(not gn_scenarios(flow, node, 'ts_cf'))))
336
337
338
339
340
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);
    // Reserves relevant only until reserve_length
    ts_reserveDemand_(restypeDirectionNode(restype, up_down, node), ft(f, tt_interval(t)))
      ${ord(t) <= tSolveFirst + p_nReserves(node, restype, 'reserve_length')  }
341
        = sum(tt_aggregate(t, t_),
342
343
            ts_reserveDemand(restype, up_down, node,
                f + (df_scenario(f, t)$gn_scenarios(restype, node, 'ts_reserveDemand')),
344
                t_+ dt_circular(t_))
345
346
347
348
349
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);
    ts_node_(gn_state(grid, node), param_gnBoundaryTypes, fts(f, tt_interval(t), s))
      $p_gnBoundaryPropertiesForStates(grid, node, param_gnBoundaryTypes, 'useTimeseries')
           // Take average if not a limit type
350
        = (sum(tt_aggregate(t, t_),
351
352
                ts_node(grid, node, param_gnBoundaryTypes,
                    f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
353
354
                    t_+ (   + dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                            + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node'))))
355
356
357
358
359
            )
            / mInterval(mSolve, 'stepsPerInterval', counter))$( not (sameas(param_gnBoundaryTypes, 'upwardLimit')
                                                                or sameas(param_gnBoundaryTypes, 'downwardLimit')
                                                                or slack(param_gnBoundaryTypes)))
          // Maximum lower limit
360
          + smax(tt_aggregate(t, t_),
361
362
                ts_node(grid, node, param_gnBoundaryTypes,
                    f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
363
364
                    t_+ (   + dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                            + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node'))))
365
366
367
                )
                $(sameas(param_gnBoundaryTypes, 'downwardLimit') or downwardSlack(param_gnBoundaryTypes))
          // Minimum upper limit
368
          + smin(tt_aggregate(t, t_),
369
370
                ts_node(grid, node, param_gnBoundaryTypes,
                    f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
371
372
                    t_+ (   + dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                            + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node'))))
373
374
375
376
377
                )
                $(sameas(param_gnBoundaryTypes, 'upwardLimit') or upwardSlack(param_gnBoundaryTypes));
    // Fuel price time series
    ts_fuelPrice_(fuel, tt_interval(t))
        ${ p_fuelPrice(fuel, 'useTimeSeries') }
378
379
        = sum(tt_aggregate(t, t_),
            + ts_fuelPrice(fuel, t_+dt_circular(t_))
380
381
382
            )
            / mInterval(mSolve, 'stepsPerInterval', counter);

383
384
); // END loop(counter)

385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
* --- Process unit time series data -------------------------------------------

// Calculate time series form parameters for units using direct input output conversion without online variable
// Always constant 'lb', 'rb', and 'section', so need only to define 'slope'.
loop(effGroupSelectorUnit(effDirectOff, unit_timeseries(unit), effDirectOff_),
    ts_effUnit(effDirectOff, unit, effDirectOff_, 'slope', ft(f, t))
        ${ sum(eff, ts_unit(unit, eff, f, t)) } // NOTE!!! Averages the slope over all available data.
        = sum(eff${ts_unit(unit, eff, f, t)}, 1 / ts_unit(unit, eff, f, t))
            / sum(eff${ts_unit(unit, eff, f, t)}, 1);
); // END loop(effGroupSelectorUnit)

// NOTE! Using the same methodology for the directOn and lambda approximations in time series form might require looping over ft(f,t) to find the min and max 'eff' and 'rb'
// Alternatively, one might require that the 'rb' is defined in a similar structure, so that the max 'rb' is located in the same index for all ft(f,t)

// Calculate unit wide parameters for each efficiency group
loop(effLevelGroupUnit(effLevel, effGroup, unit)${  mSettingsEff(mSolve, effLevel)
                                                    and p_unit(unit, 'useTimeseries')
                                                    },
    ts_effGroupUnit(effGroup, unit, 'lb', ft(f, t))${   sum(effSelector, ts_effUnit(effGroup, unit, effSelector, 'lb', f, t))}
        = smin(effSelector${effGroupSelectorUnit(effGroup, unit, effSelector)}, ts_effUnit(effGroup, unit, effSelector, 'lb', f, t));
    ts_effGroupUnit(effGroup, unit, 'slope', ft(f, t))${sum(effSelector, ts_effUnit(effGroup, unit, effSelector, 'slope', f, t))}
        = smin(effSelector$effGroupSelectorUnit(effGroup, unit, effSelector), ts_effUnit(effGroup, unit, effSelector, 'slope', f, t)); // Uses maximum efficiency for the group
); // END loop(effLevelGroupUnit)

Erkka Rinne's avatar
Erkka Rinne committed
409
410
411
412
* =============================================================================
* --- Input data processing ---------------------------------------------------
* =============================================================================

413
$ifthen defined scenario
Erkka Rinne's avatar
Erkka Rinne committed
414
* --- Scenario reduction ------------------------------------------------------
415
if(active(mSolve, 'scenred') and mSettings('schedule', 'scenarios') > 1,
Erkka Rinne's avatar
Erkka Rinne committed
416
417
    $$include 'inc/scenred.gms'
);
418
$endif
Erkka Rinne's avatar
Erkka Rinne committed
419

420
421
422
423
424
425
426
427
* --- Update probabilities ----------------------------------------------------
Option clear = p_msft_probability;
p_msft_probability(msft(mSolve, s, f, t))
    = p_mfProbability(mSolve, f)
        / sum(f_${ft(f_, t)},
              p_mfProbability(mSolve, f_)) * p_msProbability(mSolve, s);


428
429
* --- Calculate sample displacements ------------------------------------------
Options clear = ds, clear = ds_state;
Erkka Rinne's avatar
Erkka Rinne committed
430
431
loop((mst_start(mSolve, s, t), ss(s, s_)),
    ds(s, t) = -(ord(s) - ord(s_));
432
433
434
435
436
437
    ds_state(gn_state(grid, node), s, t)
      ${not sum(s__, gnss_bound(grid, node, s__, s))
        and not sum(s__, gnss_bound(grid, node, s, s__))} = ds(s, t);
);


Erkka Rinne's avatar
Erkka Rinne committed
438
439
* --- Smooting of stochastic scenarios ----------------------------------------
$ontext
440
441
Smoothen the scenarios following the methodology presented in [1, p. 443].
This avoids a discontinuity `jump' after the initial sample.
Erkka Rinne's avatar
Erkka Rinne committed
442
443
444
445
446
447

[1] A. Helseth, B. Mo, A. Lote Henden, and G. Warland, "Detailed long-term hydro-
    thermal scheduling for expansion planning in the Nordic power system," IET Gener.
    Transm. Distrib., vol. 12, no. 2, pp. 441 - 447, 2018.
$offtext

448
// Do smoothing
449
450
451
452
453
454
455
456
if(mSettings(mSolve, 'scenarios'),  // Only do smooting if using long-term scenarios
    // Select the initial sample, first `t` not in it and `f` of the last `t` in it
    loop((ms_initial(mSolve, s_), t_, ft(f_, t__))
        $[ord(t_) = msEnd(mSolve, s_) and mst_end(mSolve, s_, t__)
          and (mf_realization(mSolve, f_) xor mf_central(mSolve, f_))],
        $$batinclude 'inc/smoothing.gms' ts_influx
        $$batinclude 'inc/smoothing.gms' ts_cf
    );
457
); // END if('scenarios')