Solve time series model in several lengths

GAMS
1

Hello

I have a large MIP planning model.
The only set in the model ‘t’, describes time (and alias ‘k’ to help
formulate things like dynamic limits within the model).

I need to run the model through long input timeseries for analysis.
Because of solution performance I would like to cut it to
several shorter lenghts ‘automatically’ with code by defining a single run length,
e.g. 100 ‘t’ steps. The last one may be shorter, of course.

I’m thinking I’d need a set or sets to enumerate these calculation runs
and loop through them with multiple solves.

I haven’t been able to formulate this successfully, though.

Could you please nudge me to the right direction with how to create
these sub problems and how to loop solve over them.

Reading through documentation I’m guessing that by looping through the
whole of the set ‘t’ piecewise like this would produce the collective results
of single lengths in the resulting variables (t).

Thank you

2

This approach is often called a rolling horizon, often with some overlap between the periods (via lag operations) and is used a lot in energy system modeling. The general idea is to split your overall time horizon into chunks, iterate over the chunks and solve the model for the time steps in the chunks. Here is some sample code:

set tt 'time horizon' /t1*t8760/, t(tt) 'simulation periods';

parameter inflow(tt), outflow(tt), price(tt);
inflow(tt) = uniform(0,1); outflow(tt) = uniform(0,1); price(tt) = uniform(0,1);
scalar s_up 'storage capacity' / 3 /, s_init 'initial storage level' / [s_up/2] /; 

positive variable s(tt) 'storage', o(tt) 'overflow', u(tt) 'underflow', sell(tt) 'sales';
equation sbal(tt);
sbal(t(tt))..
  s(tt) =e= s(tt-1) + s_init$sameas(tt,'t1') + inflow(tt) - sell(tt) - outflow(tt) - o(tt) + u(tt);

variable obj; equation defobj;
defobj.. obj =e= 999*sum(t(tt), o(tt) + u(tt)) - sum(t(tt), price(tt)*sell(tt));

model storage /all/;
s.up(tt) = s_up;

option solveLink=%solveLink.loadLibrary%, limrow=0, limcol=0, solprint=off;

* Solve in one swoop
t(tt) = yes;
solve storage min obj us lp;
scalar best_obj; best_obj = obj.l;

parameter rep;
rep(tt,'best','sell') = sell.l(tt);
rep(tt,'best','price')$sell.l(tt) = price(tt);
rep(tt,'best','o') = o.l(tt);
rep(tt,'best','u') = u.l(tt);
rep(tt,'best','s') = s.l(tt);

* Solve with rolling horizon
set d 'days' / d1*d365 /, h 'hours' /h1*h24/;
set map(tt,d,h) / #tt:(#d.#h) /;

option solprint=silent;
scalar rolling_obj /0/;
loop(d,
   t(tt) = sum(map(tt,d,h), 1);
*  Simple forward looking strategy by fixing storage level to half of capacity
   s.fx(tt)$map(tt,d,'h24') = s_up/2;
   solve storage min obj us lp;
   rolling_obj = rolling_obj + obj.l;
   rep(t,'rolling','sell') = sell.l(t);
   rep(t,'rolling','price')$sell.l(t) = price(t);
   rep(t,'rolling','o') = o.l(t);
   rep(t,'rolling','u') = u.l(t);
   rep(t,'rolling','s') = s.l(t);
*  now fix this day's decision for variables with a lag (tt-1)
   s.fx(t) = s.l(t);   
);
option rep:4:1:2;
display best_obj, rolling_obj, rep;

This implements a simple storage facility with a know in and outflow. “Extra” commodity can be sold to the market at a varying price. Over- and under utilization is penalized. First, we solve the model (since this is simple) over the entire horizon. Next we iterate over the chunks and solve for the period in the chunk. Since the storage balance equation uses the storage level from the previous period, the last period from the last chunk n-1 is included in the model of chunk n. Hence this must be fixed in the model. I guess the rest of model is pretty much self-explanatory.

This model also shows the value of being able “to look into the future”. The model with a rolling horizon needs some tweaking not to deplete the storage at the end of the optimization period, even though that is the optimal solution for this chunk. I added the constraint that in the last period of the chunk the storage level needs to be at half of the capacity. With this the objective drops from 511352.801 to 308805.158. Still the solution over the entire period is significantly smaller, namely 199828.459.

Hope this helps,
-Michael

3

Thanks very much Michael,
I will try this out.

4

Hi,
I’m not sure I understand the details of this model.

  • Solve with rolling horizon
    set d ‘days’ / d1d365 /, h ‘hours’ /h1h24/;
    set map(tt,d,h) / #tt:(#d.#h) /;
    …
    loop(d,
    t(tt) = sum(map(tt,d,h), 1);
    …
    solve storage min obj us lp;

Could you please explain a bit more about how the loop works? I think I’m a bit lost with
the sets/indices and how the periods overlap in this example.

Also what is being iterated?

As a first step with my model I’m trying to calculate periods without overlap and with
no transfer of values between them.
Later I may need to transfer some (variable) values from the end of previous
period to the beginning of the next.

Thanks again,

5

This model solves batches of 24 hours, so it loops over days (d) and there no overlap between the periods of the batches. We just reach into the last hour of the previous day with s(tt-1) (as data). If you want overlap, you need to extend the set tt in each iteration (e.g. include the last 6 hours of the previous day as decision variables) and only fix up to the 7th hour before midnight rather all of the previous period (s.fx(t) = s.l(t);).

Good luck.

-Michael