Skip to content
GitLab
Commit b227d09c authored by Erkka Rinne's avatar Erkka Rinne
Browse files

Use actual stochastic programming for long-term horizon 

Also enables the use of scenario tree construction with SCENRED2
parent e65277f8
Hide whitespace changes
Inline Side-by-side
inc/1a_definitions.gms
View file @ b227d09c
......@@ -50,7 +50,8 @@ Sets
t_forecastStart, // Time step for first reading the forecasts (not necessarily t_start)
t_forecastJump, // Number of time steps between each update of the forecasts
t_improveForecast "Number of time steps ahead of time on which the forecast is improved on each solve"
sampleLength "Length of sample in time steps for creating stocahstic scenarios from time series data"
scenarios "Number of long-term scenarios used"
scenarioLength "Length of scenario in time steps for creating stocahstic scenarios from time series data"
// Features
t_trajectoryHorizon, // Length of the horizon when start-up and shutdown trajectories are considered (in time steps)
......
inc/1b_sets.gms
View file @ b227d09c
......@@ -119,11 +119,12 @@ Sets
t_latestForecast(t) "t for the latest forecast that is available"
gnss_bound(grid, node, s, s) "Bound the samples so that the node state at the last interval of the first sample equals the state at the first interval of the second sample"
uss_bound(unit, s, s) "Bound the samples so that the unit online state at the last interval of the first sample equals the state at the first interval of the second sample"
s_parallel(s) "Samples which are treated as parallel"
s_active(s) "Samples with non-zero probability in the current model solve"
s_stochastic(s) "Samples used for stochastic programming"
ss(s, s) "Previous sample of sample"
s_prev(s) "Temporary set for previous sample"
longtermSamples(*, node, *) "Which grid/flow, node and timeseries/param have data for long-term scenarios"
s_scenario(s, scenario) "Which samples belong to which scenarios"
gn_scenarios(*, node, *) "Which grid/flow, node and timeseries/param have data for long-term scenarios"
* --- Sets used for the changing unit aggregation and efficiency approximations
uft(unit, f, t) "Active units on intervals, enables aggregation of units for later intervals"
......@@ -165,7 +166,7 @@ Option clear = modelSolves;
Option clear = ms;
Option clear = mf;
mf_realization(mType, 'f00') = yes;
Option clear = longtermSamples;
Option clear = gn_scenarios;
Option clear = mTimeseries_loop_read;
alias(m, mSolve);
......
inc/1c_parameters.gms
View file @ b227d09c
......@@ -27,6 +27,7 @@ Scalars
currentForecastLength "Length of the forecast in the curren solve, minimum of unchanging and decreasing forecast lengths"
count "General counter"
count_lambda, count_lambda2 "Counter for lambdas"
count_sample "Counter for samples"
cum_slope "Cumulative for slope"
cum_lambda "Cumulative for lambda"
heat_rate "Heat rate temporary parameter"
......@@ -90,10 +91,10 @@ Parameters
Parameters
p_msWeight(mType, s) "Weight of sample"
p_msProbability(mType, s) "Probability to reach sample conditioned on anchestor samples"
p_msProbability_orig(mType, s) "Original probabilities of samples in model"
p_mfProbability(mType, f) "Probability of forecast"
p_msft_probability(mType, s, f, t) "Probability of forecast"
p_sProbability(s) "Probability of sample"
p_scenProbability(scenario) "Probability of scenarios"
;
Scalar p_sWeightSum "Sum of sample weights";
......@@ -111,7 +112,7 @@ Parameters
dt_downtimeUnitCounter(unit, counter) "Displacement needed to account for downtime constraints (in time steps)"
dt_uptimeUnitCounter(unit, counter) "Displacement needed to account for uptime constraints (in time steps)"
dt_trajectory(counter) "Run-up/shutdown trajectory time index displacement"
dt_sampleOffset(*, node, *, s) "Time offset to make periodic time series data (for grid/flow, unit, label) to go into different samples"
dt_scenarioOffset(*, node, *, s) "Time offset to make periodic time series data (for grid/flow, unit, label) to go into different scenarios"
// Forecast displacement arrays
df(f, t) "Displacement needed to reach the realized forecast on the current time step"
......
inc/3a_periodicInit.gms
View file @ b227d09c
......@@ -56,28 +56,8 @@ $offtext
// Select samples for the model
if (not sum(s, ms(m, s)), // unless they have been provided as input
ms(m, s)$(ord(s) <= mSettings(m, 'samples')) = yes;
if (mSettings(m, 'samples') = 0, // Use all samples if mSettings/samples is 0
ms(m, s) = p_msProbability(m, s);
);
);
// Calculate which samples are treated as parallel and the previous samples
loop(ms_initial(m, s_), // Select the root sample
loop(ms(m, s)$(not sameas(s, s_)), // Select other samples than root
// If two samples share same starting time, treat them as parallel
if(msStart(m, s) = msStart(m, s - 1),
s_parallel(s) = yes;
s_parallel(s - 1) = yes;
);
// Set previous samples for samples
if(msEnd(m, s_) = msStart(m, s), ss(s, s_) = yes);
if(msEnd(m, s - 1) = msStart(m, s), ss(s, s - 1) = yes);
);
);
// Store original probabilities
p_msProbability_orig(m, s) = p_msProbability(m, s);
// Select forecasts in use for the models
if (not sum(f, mf(m, f)), // unless they have been provided as input
mf(m, f)$(ord(f) <= 1 + mSettings(m, 'forecasts')) = yes; // realization needs one f, therefore 1 + number of forecasts
......@@ -89,8 +69,6 @@ $offtext
// Select combinations of models, samples and forecasts to be solved
msf(m, s, f_solve(f))$(ms(m, s) and mf(m, f)) = yes;
msf(m, s_parallel(s), f_solve(f)) = mf_central(m, f); // Parallel samples only have central forecast
// Check the modelSolves for preset patterns for model solve timings
// If not found, then use mSettings to set the model solve timings
......
inc/3b_periodicLoop.gms
View file @ b227d09c
......@@ -21,105 +21,107 @@ $offtext
* --- Variables ---------------------------------------------------------------
// Free Variables
Option clear = v_gen;
Option clear = v_state;
Option clear = v_genRamp;
Option clear = v_transfer;
// Integer Variables
Option clear = v_online_MIP;
Option clear = v_invest_MIP;
Option clear = v_investTransfer_MIP;
// SOS2 Variables
Option clear = v_sos2;
// Positive Variables
Option clear = v_fuelUse;
Option clear = v_startup;
Option clear = v_shutdown;
Option clear = v_genRampUpDown;
Option clear = v_spill;
Option clear = v_transferRightward;
Option clear = v_transferLeftward;
Option clear = v_resTransferRightward;
Option clear = v_resTransferLeftward;
Option clear = v_reserve;
Option clear = v_investTransfer_LP;
Option clear = v_online_LP;
Option clear = v_invest_LP;
// Feasibility control
Option clear = v_stateSlack;
Option clear = vq_gen;
Option clear = vq_resDemand;
Option clear = vq_resMissing;
// Free Variables
Option clear = v_gen;
Option clear = v_state;
Option clear = v_genRamp;
Option clear = v_transfer;
// Integer Variables
Option clear = v_online_MIP;
Option clear = v_invest_MIP;
Option clear = v_investTransfer_MIP;
// SOS2 Variables
Option clear = v_sos2;
// Positive Variables
Option clear = v_fuelUse;
Option clear = v_startup;
Option clear = v_shutdown;
Option clear = v_genRampUpDown;
Option clear = v_spill;
Option clear = v_transferRightward;
Option clear = v_transferLeftward;
Option clear = v_resTransferRightward;
Option clear = v_resTransferLeftward;
Option clear = v_reserve;
Option clear = v_investTransfer_LP;
Option clear = v_online_LP;
Option clear = v_invest_LP;
// Feasibility control
Option clear = v_stateSlack;
Option clear = vq_gen;
Option clear = vq_resDemand;
Option clear = vq_resMissing;
* --- Equations ---------------------------------------------------------------
// Objective Function, Energy Balance, and Reserve demand
Option clear = q_obj;
Option clear = q_balance;
Option clear = q_resDemand;
// Unit Operation
Option clear = q_maxDownward;
Option clear = q_maxUpward;
Option clear = q_reserveProvision;
Option clear = q_startshut;
Option clear = q_startuptype;
Option clear = q_onlineLimit;
Option clear = q_onlineMinUptime;
Option clear = q_onlineCyclic;
Option clear = q_onlineOnStartUp;
Option clear = q_offlineAfterShutdown;
Option clear = q_genRamp;
Option clear = q_rampUpLimit;
Option clear = q_rampDownLimit;
Option clear = q_rampUpDown;
Option clear = q_rampSlack;
Option clear = q_outputRatioFixed;
Option clear = q_outputRatioConstrained;
Option clear = q_conversionDirectInputOutput;
Option clear = q_conversionSOS2InputIntermediate;
Option clear = q_conversionSOS2Constraint;
Option clear = q_conversionSOS2IntermediateOutput;
Option clear = q_fuelUseLimit;
// Energy Transfer
Option clear = q_transfer;
Option clear = q_transferRightwardLimit;
Option clear = q_transferLeftwardLimit;
Option clear = q_resTransferLimitRightward;
Option clear = q_resTransferLimitLeftward;
Option clear = q_reserveProvisionRightward;
Option clear = q_reserveProvisionLeftward;
// State Variables
Option clear = q_stateSlack;
Option clear = q_stateUpwardLimit;
Option clear = q_stateDownwardLimit;
Option clear = q_boundStateMaxDiff;
Option clear = q_boundCyclic;
// Policy
Option clear = q_inertiaMin;
Option clear = q_instantaneousShareMax;
Option clear = q_constrainedOnlineMultiUnit;
Option clear = q_capacityMargin;
Option clear = q_constrainedCapMultiUnit;
Option clear = q_emissioncap;
Option clear = q_energyShareMax;
Option clear = q_energyShareMin;
// Objective Function, Energy Balance, and Reserve demand
Option clear = q_obj;
Option clear = q_balance;
Option clear = q_resDemand;
// Unit Operation
Option clear = q_maxDownward;
Option clear = q_maxUpward;
Option clear = q_reserveProvision;
Option clear = q_startshut;
Option clear = q_startuptype;
Option clear = q_onlineLimit;
Option clear = q_onlineMinUptime;
Option clear = q_onlineCyclic;
Option clear = q_onlineOnStartUp;
Option clear = q_offlineAfterShutdown;
Option clear = q_genRamp;
Option clear = q_rampUpLimit;
Option clear = q_rampDownLimit;
Option clear = q_rampUpDown;
Option clear = q_rampSlack;
Option clear = q_outputRatioFixed;
Option clear = q_outputRatioConstrained;
Option clear = q_conversionDirectInputOutput;
Option clear = q_conversionSOS2InputIntermediate;
Option clear = q_conversionSOS2Constraint;
Option clear = q_conversionSOS2IntermediateOutput;
Option clear = q_fuelUseLimit;
// Energy Transfer
Option clear = q_transfer;
Option clear = q_transferRightwardLimit;
Option clear = q_transferLeftwardLimit;
Option clear = q_resTransferLimitRightward;
Option clear = q_resTransferLimitLeftward;
Option clear = q_reserveProvisionRightward;
Option clear = q_reserveProvisionLeftward;
// State Variables
Option clear = q_stateSlack;
Option clear = q_stateUpwardLimit;
Option clear = q_stateDownwardLimit;
Option clear = q_boundStateMaxDiff;
Option clear = q_boundCyclic;
// Policy
Option clear = q_inertiaMin;
Option clear = q_instantaneousShareMax;
Option clear = q_constrainedOnlineMultiUnit;
Option clear = q_capacityMargin;
Option clear = q_constrainedCapMultiUnit;
Option clear = q_emissioncap;
Option clear = q_energyShareMax;
Option clear = q_energyShareMin;
* --- Temporary Time Series ---------------------------------------------------
// Initialize temporary time series
Option clear = ts_unit_;
* Option clear = ts_effUnit_;
* Option clear = ts_effGroupUnit_;
Option clear = ts_influx_;
Option clear = ts_cf_;
Option clear = ts_unit_;
Option clear = ts_reserveDemand_;
Option clear = ts_node_;
// Initialize temporary time series
Option clear = ts_unit_;
*Option clear = ts_effUnit_;
*Option clear = ts_effGroupUnit_;
Option clear = ts_influx_;
Option clear = ts_cf_;
Option clear = ts_unit_;
Option clear = ts_reserveDemand_;
Option clear = ts_node_;
* =============================================================================
* --- Determine the forecast-intervals included in the current solve ----------
......@@ -128,26 +130,6 @@ $offtext
// Determine the time steps of the current solve
tSolveFirst = ord(tSolve); // tSolveFirst: the start of the current solve, t0 used only for initial values
* --- Define sample offsets for creating stochastic scenarios -----------------
Option clear = dt_sampleOffset;
loop(gn(grid, node),
loop(longtermSamples(grid, node, timeseries),
dt_sampleOffset(grid, node, timeseries, s)$(ord(s) > 1)
= (ord(s) - 2) * mSettings(mSolve, 'sampleLength');
);
loop(longtermSamples(grid, node, param_gnBoundaryTypes),
dt_sampleOffset(grid, node, param_gnBoundaryTypes, s)$(ord(s) > 1)
= (ord(s) - 2) * mSettings(mSolve, 'sampleLength');
);
);
loop((flowNode(flow, node), longtermSamples(flow, node, timeseries)),
dt_sampleOffset(flow, node, timeseries, s)$(ord(s) > 1)
= (ord(s) - 2) * mSettings(mSolve, 'sampleLength');
);
* --- Build the forecast-time structure using the intervals -------------------
// Initializing forecast-time structure sets
......@@ -156,12 +138,15 @@ Option clear = msft;
Option clear = mft;
Option clear = ft;
Option clear = sft;
Option clear = s_active, clear = s_stochastic, clear = s_scenario, clear = ss;
Option clear = mst_start, clear = mst_end;
// Initialize the set of active t:s, counters and interval time steps
Option clear = t_active;
Option clear = tt_block;
Option clear = cc;
tCounter = 1;
count_sample = 1;
// Determine the set of active interval counters (or blocks of intervals)
cc(counter)${ mInterval(mSolve, 'stepsPerInterval', counter) }
......@@ -181,7 +166,7 @@ currentForecastLength
// Is there any case where t_forecastLength should be larger than t_horizon? Could happen if one doesn't want to join forecasts at the end of the solve horizon.
// If not, add a check for currentForecastLength <= mSettings(mSolve, 't_horizon')
// and change the line below to 'tSolveLast = ord(tSolve) + mSettings(mSolve, 't_horizon');'
tSolveLast = ord(tSolve) + max(currentForecastLength, min(mSettings(mSolve, 't_horizon'), smax(s, msEnd(mSolve, s)))); // tSolveLast: the end of the current solve
tSolveLast = ord(tSolve) + mSettings(mSolve, 't_horizon'); // tSolveLast: the end of the current solve
Option clear = t_current;
t_current(t_full(t))
${ ord(t) >= tSolveFirst
......@@ -201,7 +186,8 @@ loop(cc(counter),
option clear = tt;
tt(t_current(t))
${ord(t) >= tSolveFirst + tCounter
and ord(t) <= min(tSolveFirst + mInterval(mSolve, 'lastStepInIntervalBlock', counter),
and ord(t) <= min(tSolveFirst
+ mInterval(mSolve, 'lastStepInIntervalBlock', counter),
tSolveLast)
} = yes;
......@@ -258,9 +244,6 @@ loop(cc(counter),
// Set of locked combinations of forecasts and intervals for the reserves?
// Reduce the model dimension
ft(f_solve, tt_interval(t)) = mft(mSolve, f_solve, t);
// Update tActive
t_active(tt_interval) = yes;
......@@ -270,6 +253,36 @@ loop(cc(counter),
and ord(t) < msEnd(mSolve, s) + tSolveFirst // Move the samples along with the dispatch
} = mft(mSolve, f_solve, t);
// Create stochastic programming scenarios
if(mSettings(mSolve, 'scenarios'),
// Select root sample and central forecast
loop((ms_initial(mSolve, s_), mf_central(mSolve, f_)),
s_active(s_) = yes;
loop(scenario$(ord(scenario) <= mSettings(mSolve, 'scenarios')),
s_scenario(s_, scenario) = yes;
loop(tt_interval(t)$(ord(t) >= msEnd(mSolve, s_) + tSolveFirst),
mft(mSolve, f_, t) = yes;
loop(s$(ord(s) = mSettings(mSolve, 'samples') + count_sample),
s_active(s) = yes;
s_stochastic(s) = yes;
msft(mSolve, s, f_, t) = yes;
s_scenario(s, scenario) = yes;
p_msProbability(mSolve, s) = p_scenProbability(scenario);
msStart(mSolve, s) = ord(t) - tSolveFirst;
msEnd(mSolve, s)
= ord(t) - tSolveFirst
+ mInterval(mSolve, 'stepsPerInterval', counter);
);
count_sample = count_sample + 1;
);
);
msf(mSolve, s, f_) = s_active(s);
);
ms(mSolve, s) = s_active(s);
);
// Reduce the model dimension
ft(f_solve, tt_interval(t)) = mft(mSolve, f_solve, t);
// Reduce the sample dimension
sft(s, f, t)$msft(mSolve, s, f, t) = ft(f, t);
......@@ -278,6 +291,37 @@ loop(cc(counter),
); // END loop(counter)
* Build stochastic tree by definfing previous samples
Option clear = s_prev;
loop(scenario$(ord(scenario) <= mSettings(mSolve, 'scenarios')),
loop(s_scenario(s, scenario),
if(not ms_initial(mSolve, s), ss(s, s_prev) = yes);
Option clear = s_prev; s_prev(s) = yes;
);
);
* --- Define sample offsets for creating stochastic scenarios -----------------
Option clear = dt_scenarioOffset;
loop(s_scenario(s, scenario)$(ord(s) > 1 and ord(scenario) > 1),
loop(gn_scenarios(grid, node, timeseries),
dt_scenarioOffset(grid, node, timeseries, s)
= (ord(scenario) - 1) * mSettings(mSolve, 'scenarioLength');
);
loop(gn_scenarios(grid, node, param_gnBoundaryTypes),
dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
= (ord(scenario) - 1) * mSettings(mSolve, 'scenarioLength');
);
loop(gn_scenarios(flow, node, timeseries),
dt_scenarioOffset(flow, node, timeseries, s)
= (ord(scenario) - 1) * mSettings(mSolve, 'scenarioLength');
);
);
* --- Determine various other forecast-time sets required for the model -------
// Set of realized intervals in the current solve
......@@ -328,30 +372,28 @@ mft_lastSteps(mSolve, ft(f,t))
= yes
;
// If this is the very first solve
if(tSolveFirst = mSettings(mSolve, 't_start'),
// Sample start and end intervals
loop(ms(mSolve, s),
tmp = 1;
tmp_ = 1;
loop(t_active(t),
if(tmp and ord(t) - tSolveFirst + 1 > msStart(mSolve, s),
mst_start(mSolve, s, t) = yes;
tmp = 0;
);
if(tmp_ and ord(t) - tSolveFirst + 1 > msEnd(mSolve, s),
mst_end(mSolve, s, t+dt(t)) = yes;
tmp_ = 0;
);
); // END loop(t_active)
// If the last interval of a sample is in mft_lastSteps, the method above does not work
if(tmp_,
mst_end(mSolve, s, t)${sum(f_solve, mft_lastSteps(mSolve, f_solve, t))} = yes;
// Sample start and end intervals
loop(ms(mSolve, s),
tmp = 1;
tmp_ = 1;
loop(t_active(t),
if(tmp and ord(t) - tSolveFirst + 1 > msStart(mSolve, s),
mst_start(mSolve, s, t) = yes;
tmp = 0;
);
); // END loop(ms)
// Displacement from the first interval of a sample to the previous interval is always -1
dt(t)${sum(ms(mSolve, s)$(not s_parallel(s)), mst_start(mSolve, s, t))} = -1;
); // END if(tSolveFirst)
if(tmp_ and ord(t) - tSolveFirst + 1 > msEnd(mSolve, s),
mst_end(mSolve, s, t+dt(t)) = yes;
tmp_ = 0;
);
); // END loop(t_active)
// If the last interval of a sample is in mft_lastSteps, the method above does not work
if(tmp_,
mst_end(mSolve, s, t)${sum(f_solve, mft_lastSteps(mSolve, f_solve, t))} = yes;
);
); // END loop(ms)
// Displacement from the first interval of a sample to the previous interval is always -1,
// except for stochastic samples
dt(t)${sum(ms(mSolve, s)$(not s_stochastic(s)), mst_start(mSolve, s, t))} = -1;
// Forecast index displacement between realized and forecasted intervals
// NOTE! This set cannot be reset without references to previously solved time steps in the stochastic tree becoming ill-defined!
......@@ -389,17 +431,6 @@ ft_reservesFixed(node, restype, f_solve(f), t_active(t))
}
= yes;
// Calculate sample displacements
Options clear = ds, clear = ds_state;
loop(ms(mSolve, s)$msStart(msolve, s), // Get all samples with defined start
loop(s_$(msEnd(mSolve, s_) = msStart(mSolve, s)), // Get the previous sample
ds(s, t)$(ord(t) = tSolveFirst + msStart(mSolve, s)) = -(ord(s) - ord(s_));
);
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);
);
* =============================================================================
* --- Defining unit aggregations and ramps ------------------------------------
* =============================================================================
......
inc/3c_inputsLoop.gms
View file @ b227d09c
......@@ -303,9 +303,9 @@ loop(cc(counter),
* --- Select time series data matching the intervals, for stepsPerInterval = 1, this is trivial.
loop(ft(f_solve, tt_interval(t)),
ts_unit_(unit_timeseries(unit), param_unit, f_solve, t) // Only if time series enabled for the unit
= ts_unit(unit, param_unit, f_solve, t+dt_circular(t));
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));
$ontext
* Should these be handled here at all? See above comment
ts_effUnit_(effGroupSelectorUnit(effSelector, unit_timeseries(unit), effSelector), param_eff, f_solve, t)
......@@ -313,21 +313,21 @@ $ontext
ts_effGroupUnit_(effSelector, unit_timeseries(unit), param_eff, f_solve, t)
= ts_effGroupUnit(effSelector, unit, param_eff, f_solve, t+dt_circular(t));
$offtext
ts_cf_(flowNode(flow, node), f_solve, t, s)$msf(mSolve, s, f_solve)
= ts_cf(flow, node, f_solve, t + (dt_sampleOffset(flow, node, 'ts_cf', s)
+ dt_circular(t)$(not longtermSamples(flow, node, 'ts_cf'))));
ts_influx_(gn(grid, node), f_solve, t, s)$msf(mSolve, s, f_solve)
= ts_influx(grid, node, f_solve, t + (dt_sampleOffset(grid, node, 'ts_influx', s)
+ dt_circular(t)$(not longtermSamples(grid, node, 'ts_influx'))));
ts_cf_(flowNode(flow, node), f, t, s)$sft(s, f, t)
= ts_cf(flow, node, f_solve, t + (dt_scenarioOffset(flow, node, 'ts_cf', s)
+ dt_circular(t)$(not gn_scenarios(flow, node, 'ts_cf'))));
ts_influx_(gn(grid, node), f, t, s)$sft(s, f, t)
= ts_influx(grid, node, f, t + (dt_scenarioOffset(grid, node, 'ts_influx', s)
+ dt_circular(t)$(not gn_scenarios(grid, node, 'ts_influx'))));
// Reserve demand relevant only up until reserve_length
ts_reserveDemand_(restypeDirectionNode(restype, up_down, node), f_solve, t)
${ord(t) <= tSolveFirst + p_nReserves(node, restype, 'reserve_length')}
= ts_reserveDemand(restype, up_down, node, f_solve, t+dt_circular(t));
ts_node_(gn_state(grid, node), param_gnBoundaryTypes, f_solve, t, s)
${p_gnBoundaryPropertiesForStates(grid, node, param_gnBoundaryTypes, 'useTimeseries')
and msf(mSolve, s, f_solve)}
= ts_node(grid, node, param_gnBoundaryTypes, f_solve, t + (dt_sampleOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t)$(not longtermSamples(grid, node, 'ts_node'))));
$p_gnBoundaryPropertiesForStates(grid, node, param_gnBoundaryTypes, 'useTimeseries')
= ts_node(grid, node, param_gnBoundaryTypes, f_solve, t + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t)$(not gn_scenarios(grid, node, 'ts_node'))));
// Fuel price time series
ts_fuelPrice_(fuel, t)
= ts_fuelPrice(fuel, t+dt_circular(t));
......@@ -339,7 +339,7 @@ $offtext
// Select and average time series data matching the intervals, for stepsPerInterval > 1
// Loop over the t:s of the interval
loop(ft(f_solve, tt_interval(t)),
loop(ft(f_solve(f), tt_interval(t)),
// Select t:s within the interval
Option clear = tt;
tt(tt_(t_))
......@@ -347,8 +347,9 @@ $offtext
and ord(t_) < ord(t) + mInterval(mSolve, 'stepsPerInterval', counter)
}
= yes;
ts_unit_(unit_timeseries(unit), param_unit, f_solve, t)
= sum(tt(t_), ts_unit(unit, param_unit, f_solve, t_+dt_circular(t_)))
ts_unit_(unit_timeseries(unit), param_unit, f, t)
= sum(tt(t_), ts_unit(unit, param_unit, f, t_+dt_circular(t_)))
/ mInterval(mSolve, 'stepsPerInterval', counter);
$ontext
* Should these be handled here at all? See above comment
......@@ -359,33 +360,33 @@ $ontext
= sum(tt(t_), ts_effGroupUnit(effSelector, unit, param_eff, f_solve, t_+dt_circular(t_))))
/ mInterval(mSolve, 'stepsPerInterval', counter);
$offtext
ts_influx_(gn(grid, node), f_solve, t, s)$msf(mSolve, s, f_solve)
= sum(tt(t_), ts_influx(grid, node, f_solve, t_ + (dt_sampleOffset(grid, node, 'ts_influx', s)
+ dt_circular(t_)$(not longtermSamples(grid, node, 'ts_influx')))))
ts_influx_(gn(grid, node), f_solve, t, s)$sft(s, f, t)
= sum(tt(t_), ts_influx(grid, node, f_solve, t_ + (dt_scenarioOffset(grid, node, 'ts_influx', s)
+ dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_influx')))))
/ mInterval(mSolve, 'stepsPerInterval', counter);
ts_cf_(flowNode(flow, node), f_solve, t, s)$msf(mSolve, s, f_solve)
= sum(tt(t_), ts_cf(flow, node, f_solve, t_ + (dt_sampleOffset(flow, node, 'ts_cf', s)
+ dt_circular(t_)$(not longtermSamples(flow, node, 'ts_cf')))))
ts_cf_(flowNode(flow, node), f_solve, t, s)$sft(s, f, t)
= sum(tt(t_), ts_cf(flow, node, f_solve, t_ + (dt_scenarioOffset(flow, node, 'ts_cf', s)
+ dt_circular(t_)$(not gn_scenarios(flow, node, 'ts_cf')))))
/ mInterval(mSolve, 'stepsPerInterval', counter);
// Reserves relevant only until reserve_length
ts_reserveDemand_(restypeDirectionNode(restype, up_down, node), f_solve, t)
ts_reserveDemand_(restypeDirectionNode(restype, up_down, node), f, t)
${ord(t) <= tSolveFirst + p_nReserves(node, restype, 'reserve_length') }
= sum(tt(t_), ts_reserveDemand(restype, up_down, node, f_solve, t_+dt_circular(t_)))
/ mInterval(mSolve, 'stepsPerInterval', counter);
ts_node_(gn_state(grid, node), param_gnBoundaryTypes, f_solve, t, s)
ts_node_(gn_state(grid, node), param_gnBoundaryTypes, f, t, s)
${p_gnBoundaryPropertiesForStates(grid, node, param_gnBoundaryTypes, 'useTimeseries')
and msf(mSolve, s, f_solve)}
and sft(s, f, t)}
// Take average if not a limit type
= (sum(tt(t_), ts_node(grid, node, param_gnBoundaryTypes, f_solve, t_ + (dt_sampleOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t_)$(not longtermSamples(grid, node, 'ts_node')))))
= (sum(tt(t_), ts_node(grid, node, param_gnBoundaryTypes, f_solve, 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')))
// Maximum lower limit
+ smax(tt(t_), ts_node(grid, node, param_gnBoundaryTypes, f_solve, t_ + (dt_sampleOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t_)$(not longtermSamples(grid, node, 'ts_node')))))
+ smax(tt(t_), ts_node(grid, node, param_gnBoundaryTypes, f_solve, t_ + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node')))))
$sameas(param_gnBoundaryTypes, 'downwardLimit')
// Minimum upper limit
+ smin(tt(t_), ts_node(grid, node, param_gnBoundaryTypes, f_solve, t_ + (dt_sampleOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t_)$(not longtermSamples(grid, node, 'ts_node')))))
+ smin(tt(t_), ts_node(grid, node, param_gnBoundaryTypes, f_solve, t_ + (dt_scenarioOffset(grid, node, param_gnBoundaryTypes, s)
+ dt_circular(t_)$(not gn_scenarios(grid, node, 'ts_node')))))
$sameas(param_gnBoundaryTypes, 'upwardLimit');
// Fuel price time series
ts_fuelPrice_(fuel, t)
......@@ -425,8 +426,6 @@ loop(effLevelGroupUnit(effLevel, effGroup, unit)${ mSettingsEff(mSolve, effLeve
* =============================================================================
* --- Scenario reduction ------------------------------------------------------
s_active(s) = ms(mSolve, s);
if(active(mSolve, 'scenred'),
$$include 'inc/scenred.gms'
);
......@@ -439,6 +438,16 @@ p_msft_probability(msft(mSolve, s, f, t))
p_mfProbability(mSolve, f_)) * p_msProbability(mSolve, s);
* --- Calculate sample displacements ------------------------------------------
Options clear = ds, clear = ds_state;
loop(ss(s, s_),
ds(s, t)$(ord(t) = tSolveFirst + msStart(mSolve, s)) = -(ord(s) - ord(s_));
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);
);
* --- Smooting of stochastic scenarios ----------------------------------------
$ontext