[Feature Request] Enable n-Step Ahead Forecasts in FromNPY Chronics Handler

[RalfCortes]

Environment

• Grid2op version: 1.12.1
• System: linux

Summary

Extend the FromNPY chronics handler to support multi-step (n-th step ahead) forecasts, instead of being limited to 1-step ahead forecasts.

This would significantly improve compatibility with modern forecasting pipelines and reinforcement learning setups that rely on multi-horizon predictions.

Current Behavior

Using grid2op, it is currently possible to inject custom time series dynamically through the FromNPY chronics handler:

import grid2op
import numpy as np
from lightsim2grid import LightSimBackend
from datetime import timedelta
from grid2op.Chronics import FromNPY

env = grid2op.make(
    "l2rpn_case14_sandbox",
    backend=LightSimBackend(),
)

init_obs = env.reset()

# Artificially create a timeseries :
load_p, load_q, prod_p, prod_v = (
    init_obs.load_p,
    init_obs.load_q,
    init_obs.prod_p,
    init_obs.prod_v,
)

load_p = np.concatenate([[load_p] * 48], axis=0)
load_q = np.concatenate([[load_q] * 48], axis=0)
prod_p = np.concatenate([[prod_p] * 48], axis=0)
prod_v = np.concatenate([[prod_v] * 48], axis=0)

load_p_fc = load_p
load_q_fc = load_q
prod_p_fc = prod_p
prod_v_fc = prod_v

# Now create the env with custom timeseries data

env = grid2op.make(
    "l2rpn_case14_sandbox",
    backend=LightSimBackend(),
    chronics_class=FromNPY,
    data_feeding_kwargs={
        "i_start": 0,
        "i_end": 48,
        "load_p": load_p,
        "load_q": load_q,
        "prod_p": prod_p,
        "prod_v": prod_v,
        # "gen_p_for_handler": PerfectForecastHandler("prod_p_forecasted"),
        # "gen_v_for_handler": PerfectForecastHandler("prod_v_forecasted"),
        # "load_p_for_handler": PerfectForecastHandler("load_p_forecasted"),
        # "load_q_for_handler": PerfectForecastHandler("load_q_forecasted"),
        "load_p_forecast": load_p_fc,
        "load_q_forecast": load_q_fc,
        "prod_p_forecast": prod_p_fc,
        "prod_v_forecast": prod_v_fc,
        "h_forecast": [i * 30 for i in range(1, 48)],
        "time_interval": timedelta(minutes=30),
    },
)

It is also straightforward to update chronics dynamically:

env.chronics_handler.change_chronics(
    new_load_p=timeserie.load_p,
    new_load_q=timeserie.load_q,
    new_prod_p=timeserie.gen_p,
    new_prod_v=timeserie.gen_v,
)


env.chronics_handler.change_forecasts(
    new_load_p=load_p_fc,
    new_load_q=load_q_fc,
    new_prod_p=gen_p_fc,
    new_prod_v=gen_v_fc,
)

Limitation

Unfortunately as specified around line 460 in Chronics/FromNPY.py, the current implementation only supports 1-step ahead forecasts.

I think internally, forecasts are expected to match the shape: (n_timesteps, n_elements)

For my application I would like to use n-th step ahead forecasts.

Possible implementation :

To allow FromNPY to support multi-step forecasts :

• Change the default forecasts with shape (n_timesteps, n_horizons, n_elements) or (n_horizons, n_timesteps, , n_elements)
• Keep backwards compatibility with n=1 by converting originals (n_timesteps, n_elements) to (n_timesteps, 1, n_elements) automatically

Near the end of the time series, higher-horizon forecasts will not have valid targets. Possible approaches:

• Allow zero-padding (simplest)
• Allow NaNs

Given the presence of i_end, zero-padding is likely acceptable and consistent with current behavior.

[RalfCortes]

Here is an example of implementation :

from typing import Optional

import numpy as np
from grid2op.Chronics import FromNPY
from grid2op.Exceptions import ChronicsError


class FromNPYMultiHorizon(FromNPY):
    """
    Extension of FromNPY supporting multi-horizon forecasts.

    Forecast arrays can be:
        (T, N)         -> automatically converted to (T, 1, N)
        (T, H, N)      -> used directly

    Where:
        T = number of timesteps
        H = number of forecast horizons
        N = number of loads / generators
    """

    MULTI_CHRONICS = False

    def __init__(
        self,
        *args,
        load_p_forecast: Optional[np.ndarray] = None,
        load_q_forecast: Optional[np.ndarray] = None,
        prod_p_forecast: Optional[np.ndarray] = None,
        prod_v_forecast: Optional[np.ndarray] = None,
        **kwargs,
    ):
        # Prevent parent from building recursive forecast object
        super().__init__(
            *args,
            load_p_forecast=None,
            load_q_forecast=None,
            prod_p_forecast=None,
            prod_v_forecast=None,
            **kwargs,
        )

        # Format forecast arrays
        self._load_p_forecast = self._format_forecast_array(load_p_forecast)
        self._load_q_forecast = self._format_forecast_array(load_q_forecast)
        self._prod_p_forecast = self._format_forecast_array(prod_p_forecast)
        self._prod_v_forecast = self._format_forecast_array(prod_v_forecast)

        # Determine number of horizons
        if self._load_p_forecast is not None:
            self.n_forecast_horizons = self._load_p_forecast.shape[1]
        else:
            self.n_forecast_horizons = 0

    @staticmethod
    def _format_forecast_array(arr: Optional[np.ndarray]) -> Optional[np.ndarray]:
        """
        Convert forecast array to shape (T, H, N)

        Accepts:
            (T, N)
            (T, H, N)
        """
        if arr is None:
            return None

        if arr.ndim == 2:
            # (T, N) → (T, 1, N)
            return arr[:, np.newaxis, :]

        if arr.ndim == 3:
            return arr

        raise ChronicsError("Forecast arrays must have shape (T, N) or (T, H, N)")

    def forecasts(self):
        """
        Return list of (datetime, dict) for each horizon.
        """
        if self._load_p_forecast is None:
            return []

        t = self.current_index

        if t >= self._load_p_forecast.shape[0]:
            return []

        results = []

        for h in range(self.n_forecast_horizons):
            dict_ = {}

            if self._load_p_forecast is not None:
                dict_["load_p"] = 1.0 * self._load_p_forecast[t, h, :]

            if self._load_q_forecast is not None:
                dict_["load_q"] = 1.0 * self._load_q_forecast[t, h, :]

            if self._prod_p_forecast is not None:
                dict_["prod_p"] = 1.0 * self._prod_p_forecast[t, h, :]

            if self._prod_v_forecast is not None:
                dict_["prod_v"] = 1.0 * self._prod_v_forecast[t, h, :]

            forecast_datetime = self.current_datetime + (h + 1) * self.time_interval

            results.append((forecast_datetime, {"injection": dict_}))

        return results

    def change_forecasts(
        self,
        new_load_p: Optional[np.ndarray] = None,
        new_load_q: Optional[np.ndarray] = None,
        new_prod_p: Optional[np.ndarray] = None,
        new_prod_v: Optional[np.ndarray] = None,
    ):
        """
        Update forecast arrays (effective after env.reset()).
        """
        if new_load_p is not None:
            self._load_p_forecast = self._format_forecast_array(new_load_p)

        if new_load_q is not None:
            self._load_q_forecast = self._format_forecast_array(new_load_q)

        if new_prod_p is not None:
            self._prod_p_forecast = self._format_forecast_array(new_prod_p)

        if new_prod_v is not None:
            self._prod_v_forecast = self._format_forecast_array(new_prod_v)

        if self._load_p_forecast is not None:
            self.n_forecast_horizons = self._load_p_forecast.shape[1]
        else:
            self.n_forecast_horizons = 0

    def check_validity(self, backend=None):
        super().check_validity(backend)

        if self._load_p_forecast is not None:
            assert self._load_p_forecast.shape[0] == self._load_p.shape[0]
            assert self._load_p_forecast.shape[2] == self.n_load

        if self._prod_p_forecast is not None:
            assert self._prod_p_forecast.shape[0] == self._prod_p.shape[0]
            assert self._prod_p_forecast.shape[2] == self.n_gen

But we also need to add a check in grid2op core :

L 614 MakeFromPath.py

  ((chronics_class_used != ChangeNothing) and 
         (chronics_class_used != FromNPY) and 
         (chronics_class_used != FromChronix2grid) and
         (chronics_class_used != FromHandlers)
         )
    ) and exc_chronics is not None:
        raise EnvError(
            f"Impossible to find the chronics for your environment. Please make sure to provide "
            f'a folder "{NAME_CHRONICS_FOLDER}" within your environment folder.'
        )

[RalfCortes]

Unfortunately while coding this Handler I encountered another problem : It seems that obs.simulate can have side effects and 'break' the initial observation we are working on.

Example :

Setup :

# Load initial env (no forecasts)

env = grid2op.make(
    "l2rpn_case14_sandbox",
    backend=LightSimBackend(),
)

init_obs = env.reset()

# Artificially create a timeseries :
load_p, load_q, prod_p, prod_v = (
    init_obs.load_p,
    init_obs.load_q,
    init_obs.prod_p,
    init_obs.prod_v,
)

load_p = np.concatenate([[load_p] * 48], axis=0)
load_q = np.concatenate([[load_q] * 48], axis=0)
prod_p = np.concatenate([[prod_p] * 48], axis=0)
prod_v = np.concatenate([[prod_v] * 48], axis=0)

load_p_fc = load_p
load_q_fc = load_q
prod_p_fc = prod_p
prod_v_fc = prod_v

# create a fake moving timeserie : 

load_p_fake = (load_p + (np.linspace(0, 1, len(load_p)) * 100).reshape(-1, 1)).round(1)

T, N = load_p_fake.shape
H = 48

indices = np.arange(T)[:, None] + np.arange(1, H + 1)[None, :]
indices = np.clip(indices, 0, T - 1)

load_p_fc = load_p_fake[indices]
load_q_fc = load_q[indices]
prod_p_fc = prod_p[indices]
prod_v_fc = prod_v[indices]

load_p_fc[1:, :, :] = 0  # Forecast available at t =0. We only need the first row for each forecast

env = grid2op.make(
    "l2rpn_case14_sandbox",
    backend=LightSimBackend(),
    allow_detachment=True,
    n_busbar=6,
    chronics_class=FromNPYMultiHorizon,
    data_feeding_kwargs={
        "i_start": 0,
        "i_end": 48,
        "load_p": load_p_fake,
        "load_q": load_q,
        "prod_p": prod_p,
        "prod_v": prod_v,
        "load_p_forecast": load_p_fc,
        "load_q_forecast": load_q_fc,
        "prod_p_forecast": prod_p_fc,
        "prod_v_forecast": prod_v_fc,
        "h_forecast": [i * 30 for i in range(1, 48)],
        "time_interval": timedelta(minutes=30),
    },
)


env.chronics_handler.change_chronics(
    new_load_p=load_p_fake,
    new_load_q=load_q,
    new_prod_p=prod_p,
    new_prod_v=prod_v,
)
env.chronics_handler.change_forecasts(
    new_load_p=load_p_fc,
    new_load_q=load_q_fc,
    new_prod_p=prod_p_fc,
    new_prod_v=prod_v_fc,
)

Issue :

This works fine !

obs_init = env.reset()
print(obs_init.get_time_stamp())
print(obs_init.load_p)

2019-01-01 00:30:00
array([21.7, 85.6, 43.6,  6.9, 11.6, 27.9,  8.6,  3.4,  5.3, 12.3, 14.6],
      dtype=float32)

obs_step, *_ = env.step(env.action_space())
print(obs_step.get_time_stamp())
print(obs_step.load_p)

2019-01-01 01:00:00
array([23.8, 87.7, 45.7,  9. , 13.7, 30. , 10.7,  5.5,  7.4, 14.4, 16.7],
      dtype=float32)

This also !

obs_simulate, *_ = obs_init.simulate(env.action_space(), time_step=1) # t = 0 ... 8 is Okay !

print(obs_simulate.get_time_stamp())
print(obs_simulate.load_p)
print(obs_simulate.rho)

2019-01-01 01:00:00
[23.8 87.7 45.7  9.  13.7 30.  10.7  5.5  7.4 14.4 16.7]

array([0.42673743, 0.35696393, 0.32683387, 0.2515636 , 0.72660244,
       0.19574797, 0.43717062, 0.71456206, 0.59599006, 0.8611846 ,
       0.06517944, 0.34333247, 0.36912158, 0.46527022, 0.49368683,
       0.49882945, 0.4822724 , 0.59464705, 0.46799862, 0.43246743],
      dtype=float32)

But here we have an issue with t = 9 :

obs_simulate, *_ = obs_init.simulate(env.action_space(), time_step=9)
# this causes a DC powerflow failure (env is done)

2019-01-01 05:00:00
[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
       0., 0., 0.], dtype=float32)

Now this doesn't "fail" but give always 0

obs_simulate, *_ = obs_init.simulate(env.action_space(), time_step=1) # t = 0 ... 8 is NOT Okay !

2019-01-01 01:00:00
[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
       0., 0., 0.], dtype=float32)

After t=9 was computed once it seems that we can not 'go back in time' and we messed up the original observation.

Explored workaround :

• I tried using obs_init.copy().simulate(...) but It did no change anything

Conclusion

I am not sure if this is a desired property. Indeed I am going through an 'impossible' state (= Loadflow divergence) but I feel like this shouldn't 'break' the initial observation ....

I have yet to find a functional workaround. It seems to me that the 'companion' obs_init._obs_env is broken and fails to reset properly when called from the simulate

[BDonnot]

Hello,

Can you create a sub issue for the last part ? (The fact that simulate can break subsequent calls). Also, does this happens with pandapower backend ?

And also for the "bug" before that. It's probably better to use "isinstance" rather than "==" in the make function when checking for time series folder.

I'll have a look when I'll be back to fix all that.

[RalfCortes]

Hello,

Indeed the issue does not arise when using Pandapower.
I will create a sub issue.

Thanks!

[BDonnot]

Then the first bug would rather be in lightsim2grid, thanks.

[RalfCortes]

Link to the issue : Grid2op/lightsim2grid#128