3c_inputsLoop.gms 27.2 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
* =============================================================================
* --- 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);
302
303
304
    // Make a temporary clone of tt_interval(t)
    option clear = tt_;
    tt_(tt_interval) = yes;
305
306
307
308

    // If stepsPerInterval equals one, simply use all the steps within the block
    if(mInterval(mSolve, 'stepsPerInterval', counter) = 1,

309
310
* --- Select time series data matching the intervals, for stepsPerInterval = 1, this is trivial.

311
312
313
        loop(ft(f_solve(f), tt_interval(t)),
            ts_unit_(unit_timeseries(unit), param_unit, f, t) // Only if time series enabled for the unit
                = ts_unit(unit, param_unit, f, t+dt_circular(t));
314
315
316
317
318
319
320
$ontext
* Should these be handled here at all? See above comment
            ts_effUnit_(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f_solve, t)
                = ts_effUnit(effSelector, unit, effSelector, param_eff, f_solve, t+dt_circular(t));
            ts_effGroupUnit_(effSelector, unit_timeseries(unit), param_eff, f_solve, t)
                = ts_effGroupUnit(effSelector, unit, param_eff, f_solve, t+dt_circular(t));
$offtext
321
            ts_cf_(flowNode(flow, node), fts(f, t, s))
322
323
324
                = ts_cf(flow, node, f + (df_scenario(f, t)$gn_scenarios(flow, node, 'ts_cf')),
                        t + (dt_scenarioOffset(flow, node, 'ts_cf', s)
                             + dt_circular(t)$(not gn_scenarios(flow, node, 'ts_cf'))));
325
            ts_influx_(gn(grid, node), fts(f, t, s))
326
327
328
                = ts_influx(grid, node, f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
                            t + (dt_scenarioOffset(grid, node, 'ts_influx', s)
                                 + dt_circular(t)$(not gn_scenarios(grid, node, 'ts_influx'))));
329
            // Reserve demand relevant only up until reserve_length
330
            ts_reserveDemand_(restypeDirectionNode(restype, up_down, node), f, t)
331
              ${ord(t) <= tSolveFirst + p_nReserves(node, restype, 'reserve_length')}
332
333
                = ts_reserveDemand(restype, up_down, node,
                                   f + (df_scenario(f,t)$gn_scenarios(restype, node, 'ts_reserveDemand')),
334
                                   t+dt_circular(t));
335
            ts_node_(gn_state(grid, node), param_gnBoundaryTypes, fts(f, t, s))
336
              $p_gnBoundaryPropertiesForStates(grid, node, param_gnBoundaryTypes, 'useTimeseries')
337
338
                = ts_node(grid, node, param_gnBoundaryTypes,
                          f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
339
340
                          t + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                               + dt_circular(t)$(not gn_scenarios(grid, node, 'ts_node'))));
341
342
            // Fuel price time series
            ts_fuelPrice_(fuel, t)
343
              ${ p_fuelPrice(fuel, 'useTimeSeries') }
344
345
346
                = ts_fuelPrice(fuel, t+dt_circular(t));
        ); // END loop(ft)

347
348
* --- If stepsPerInterval exceeds 1 (stepsPerInterval < 1 not defined) --------

349
350
351
352
    elseif mInterval(mSolve, 'stepsPerInterval', counter) > 1,

        // Select and average time series data matching the intervals, for stepsPerInterval > 1
        // Loop over the t:s of the interval
353
        loop(ft(f_solve(f), tt_interval(t)),
354
355
            // Select t:s within the interval
            Option clear = tt;
356
357
358
            tt(tt_(t_))
                ${  ord(t_) >= ord(t)
                    and ord(t_) < ord(t) + mInterval(mSolve, 'stepsPerInterval', counter)
359
360
                 }
                = yes;
361
362
363

            ts_unit_(unit_timeseries(unit), param_unit, f, t)
                = sum(tt(t_), ts_unit(unit, param_unit, f, t_+dt_circular(t_)))
364
365
366
367
368
369
370
371
372
373
                    / mInterval(mSolve, 'stepsPerInterval', counter);
$ontext
* Should these be handled here at all? See above comment
            ts_effUnit_(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f_solve, t)
                = sum(tt(t_), ts_effUnit(effSelector, unit, effSelector, param_eff, f_solve, t_+dt_circular(t_))))
                    / mInterval(mSolve, 'stepsPerInterval', counter);
            ts_effGroupUnit_(effSelector, unit_timeseries(unit), param_eff, f_solve, t)
                = sum(tt(t_), ts_effGroupUnit(effSelector, unit, param_eff, f_solve, t_+dt_circular(t_))))
                    / mInterval(mSolve, 'stepsPerInterval', counter);
$offtext
374
            ts_influx_(gn(grid, node), fts(f, t, s))
375
376
                = sum(tt(t_), ts_influx(grid, node,
                                        f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_influx')),
377
378
                                        t_ + (dt_scenarioOffset(grid, node, 'ts_influx', s)
                                              + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_influx')))))
379
                    / mInterval(mSolve, 'stepsPerInterval', counter);
380
            ts_cf_(flowNode(flow, node), fts(f, t, s))
381
382
                = sum(tt(t_), ts_cf(flow, node,
                                    f + (df_scenario(f, t)$gn_scenarios(flow, node, 'ts_cf')),
383
384
                                    t_ + (dt_scenarioOffset(flow, node, 'ts_cf', s)
                                          + dt_circular(t_)$(not gn_scenarios(flow, node, 'ts_cf')))))
385
386
                    / mInterval(mSolve, 'stepsPerInterval', counter);
            // Reserves relevant only until reserve_length
387
            ts_reserveDemand_(restypeDirectionNode(restype, up_down, node), f, t)
388
              ${ord(t) <= tSolveFirst + p_nReserves(node, restype, 'reserve_length')  }
389
                = sum(tt(t_), ts_reserveDemand(restype, up_down, node,
390
391
                                               f + (df_scenario(f, t)$gn_scenarios(restype, node, 'ts_reserveDemand')),
                                               t_+ dt_circular(t_)))
392
                    / mInterval(mSolve, 'stepsPerInterval', counter);
393
394
            ts_node_(gn_state(grid, node), param_gnBoundaryTypes, fts(f, t, s))
              $p_gnBoundaryPropertiesForStates(grid, node, param_gnBoundaryTypes, 'useTimeseries')
395
                   // Take average if not a limit type
396
                = (sum(tt(t_), ts_node(grid, node, param_gnBoundaryTypes,
397
398
399
400
401
402
                                       f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
                                       t_ + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                                             + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node')))))
                    / mInterval(mSolve, 'stepsPerInterval', counter))$(not (sameas(param_gnBoundaryTypes, 'upwardLimit')
                                                                            or sameas(param_gnBoundaryTypes, 'downwardLimit')
                                                                            or slack(param_gnBoundaryTypes)))
403
                  // Maximum lower limit
404
                  + smax(tt(t_), ts_node(grid, node, param_gnBoundaryTypes,
405
406
407
408
                                         f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
                                         t_ + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                                               + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node')))))
                      $(sameas(param_gnBoundaryTypes, 'downwardLimit') or downwardSlack(param_gnBoundaryTypes))
409
                  // Minimum upper limit
410
                  + smin(tt(t_), ts_node(grid, node, param_gnBoundaryTypes,
411
                                         f + (df_scenario(f, t)$gn_scenarios(grid, node, 'ts_node')),
412
413
                                         t_ + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
                                               + dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node')))))
414
                       $(sameas(param_gnBoundaryTypes, 'upwardLimit') or upwardSlack(param_gnBoundaryTypes));
415
416
            // Fuel price time series
            ts_fuelPrice_(fuel, t)
417
                ${ p_fuelPrice(fuel, 'useTimeSeries') }
418
419
420
421
422
423
424
                = sum(tt(t_), ts_fuelPrice(fuel, t_+dt_circular(t_)))
                    / mInterval(mSolve, 'stepsPerInterval', counter);
            ); // END loop(ft)

    ); // END if(stepsPerInterval)
); // END loop(counter)

425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
* --- 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
449
450
451
452
453
454
455
456
457
* =============================================================================
* --- Input data processing ---------------------------------------------------
* =============================================================================

* --- Scenario reduction ------------------------------------------------------
if(active(mSolve, 'scenred'),
    $$include 'inc/scenred.gms'
);

458
459
460
461
462
463
464
465
* --- 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);


466
467
* --- Calculate sample displacements ------------------------------------------
Options clear = ds, clear = ds_state;
Erkka Rinne's avatar
Erkka Rinne committed
468
469
loop((mst_start(mSolve, s, t), ss(s, s_)),
    ds(s, t) = -(ord(s) - ord(s_));
470
471
472
473
474
475
    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
476
477
* --- Smooting of stochastic scenarios ----------------------------------------
$ontext
478
479
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
480
481
482
483
484
485

[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

486
487
488
489
490
// Do smoothing
loop((ms_initial(mSolve, s_), mf_central(mSolve, f_), t_)
    $(ord(t_) = msEnd(mSolve, s_)),
    $$batinclude 'inc/smoothing.gms' ts_influx
    $$batinclude 'inc/smoothing.gms' ts_cf
Erkka Rinne's avatar
Erkka Rinne committed
491
);