gpu_simulator.py

"""
Simulateur batch GPU — évalue N allocations en parallèle.

Principe
--------
Pour chaque allocation on maintient deux tenseurs d'état :
  vehicle_end_inter[N, V]  — fin d'intervention (dispo pour redirection)
  vehicle_at_base[N, V]    — retour caserne (dispo « en caserne »)

On itère séquentiellement sur les interventions (dépendance temporelle imposée),
mais chaque étape opère sur N simulations simultanément via des tenseurs PyTorch.

Approximations par rapport au simulateur séquentiel
----------------------------------------------------
- Redirection en retour : un véhicule qui a fini son intervention mais
  rentre à la caserne est disponible pour une nouvelle mission (comme le
  séquentiel), mais le temps de trajet est calculé depuis la caserne et non
  depuis sa position intermédiaire réelle.
- Contrainte plage_interdite_nuit pour PSE/VSAV non implémentée.
  (Effet marginal sur le score global ; peut être ajouté ultérieurement.)
- Indisponibilités : la matrice vehicle_unavail est pré-calculée une fois.
"""

from __future__ import annotations

from typing import Literal

import numpy as np
import torch

from data_loader import TYPE_ELIGIBLE, POMPE_TYPE, VSAV_TYPE   # constantes de type


MetricName = Literal["temps_moyen", "taux_couverture", "manque_local", "mixte", "couverture_0", "p95_vsav", "p99_vsav", "std_q95_pression", "composite"]
SEUIL_COUVERTURE = 8 * 60.0   # 8 minutes en secondes
WINDOW_8H = 24 * 3600.0       # fenêtre 8h→8h = 24h en secondes


class GPUBatchSimulator:
    """
    Évalue un batch de N allocations en une seule passe.

    Parameters
    ----------
    interventions_np   : dict numpy retourné par data_loader.load_data()
    vehicles_np        : dict numpy retourné par data_loader.load_data()
    travel_matrix      : [N_cs, N_inter] float32 — CS→intervention
    return_travel_matrix : [N_cs, N_inter] float32 — x_final→CS
    vehicle_unavail    : [V, N_inter] bool
    window_start       : première intervention de la fenêtre (index)
    window_size        : nombre d'interventions à simuler (None = tout)
    metric             : métrique à minimiser
    device             : "cuda" | "cpu"
    """

    def __init__(
        self,
        interventions_np:        dict,
        vehicles_np:             dict,
        travel_matrix:           np.ndarray,
        return_travel_matrix:    np.ndarray,
        vehicle_unavail:         np.ndarray,
        window_start:            int  = 0,
        window_size:             int | None = None,
        metric:                  MetricName = "temps_moyen",
        device:                  str = "cuda",
        cstc_priority:           bool = True,
    ):
        self.metric  = metric
        self.cstc_priority = cstc_priority
        self.device  = torch.device(device if torch.cuda.is_available() else "cpu")
        dev = self.device

        # ── Données complètes (numpy) conservées pour set_window() ────────
        self._all = interventions_np          # référence légère, pas de copie
        self._N_total = len(interventions_np["timestamps"])

        # ── Matrices complètes sur GPU (fixes, ne bougent pas) ─────────────
        self.N_cs = travel_matrix.shape[0]
        self._travel_full  = torch.from_numpy(travel_matrix.astype(np.float32)).to(dev)        # [N_cs, N_total]
        self._ret_full     = torch.from_numpy(return_travel_matrix.astype(np.float32)).to(dev) # [N_cs, N_total]
        self._unavail_full = torch.from_numpy(vehicle_unavail).to(dev)      # [V,    N_total]

        # ── Véhicules (fixes) ─────────────────────────────────────────────
        self.V     = len(vehicles_np["ids"])
        self.vtype = vehicles_np["type"]
        self.vtype_t      = torch.from_numpy(self.vtype.astype(np.int64)).to(dev)
        self.type_eligible= torch.from_numpy(TYPE_ELIGIBLE).to(dev)  # [2, 3]

        vsav_pse = vehicles_np["vsav_pse_pair_idx"]
        pse_vsav = vehicles_np["pse_vsav_pair_idx"]
        self.vsav_has_pair = torch.from_numpy(vsav_pse >= 0).to(dev)
        self.pse_has_pair  = torch.from_numpy(pse_vsav >= 0).to(dev)
        self.vsav_pse_safe = torch.from_numpy(np.maximum(vsav_pse, 0).astype(np.int64)).to(dev)
        self.pse_vsav_safe = torch.from_numpy(np.maximum(pse_vsav, 0).astype(np.int64)).to(dev)

        # ── Fenêtre initiale ───────────────────────────────────────────────
        self.set_window(window_start, window_size)

    # ------------------------------------------------------------------
    def set_window(self, window_start: int, window_size: int | None = None) -> None:
        """
        Change la fenêtre d'interventions sans recharger les données.
        Appelable à chaque itération SA pour une fenêtre glissante aléatoire.

        Parameters
        ----------
        window_start : index de la première intervention
        window_size  : nombre d'interventions (None = jusqu'à la fin)
        """
        dev = self.device
        end = window_start + window_size if window_size else self._N_total
        end = min(end, self._N_total)
        sl  = slice(window_start, end)

        # Filtrer cstc inconnu
        valid = self._all["cstc_idx"][sl] >= 0
        global_indices = np.where(valid)[0] + window_start
        self.global_indices = global_indices
        self.N_inter        = int(valid.sum())

        def _f(arr):
            return arr[sl][valid]

        # Numpy (pour simulate_single_detailed et compute_all_metrics)
        self.timestamps_np       = _f(self._all["timestamps"])
        self.inter_type_np       = _f(self._all["type"])
        self.cstc_idx_np         = _f(self._all["cstc_idx"])
        self.duration_np         = _f(self._all["duration"])
        self.depart_np           = _f(self._all["depart"])
        self.proc_np             = _f(self._all["proc"])
        self.night_restricted_np = _f(self._all["night_restricted"])
        self.rouge_np            = _f(self._all["rouge"])

        # GPU — re-slice les matrices complètes (pas de transfer CPU→GPU)
        gi = torch.from_numpy(global_indices).to(dev)
        self.travel_t     = self._travel_full[:, gi]    # [N_cs, N_inter]
        self.ret_travel_t = self._ret_full[:, gi]        # [N_cs, N_inter]
        self.unavail_t    = self._unavail_full[:, gi]    # [V,    N_inter]

        self.timestamps_t = torch.from_numpy(self.timestamps_np.astype(np.float64)).to(dev)
        self.inter_type_t = torch.from_numpy(self.inter_type_np).to(dev)
        self.cstc_idx_t   = torch.from_numpy(self.cstc_idx_np).to(dev)
        self.duration_t   = torch.from_numpy(self.duration_np).to(dev)
        self.depart_t     = torch.from_numpy(self.depart_np).to(dev)
        self.night_t      = torch.from_numpy(self.night_restricted_np).to(dev)
        self.rouge_t      = torch.from_numpy(self.rouge_np).to(dev)

        self._t0 = float(self.timestamps_np[0]) - 1.0 if self.N_inter > 0 else 0.0

        # ── Fenêtres 8h→8h pour std_q95_pression ─────────────────────────
        if self.N_inter > 0:
            import pandas as pd
            dates = pd.to_datetime(self.timestamps_np, unit="s", utc=True)
            # Fenêtre 8h→8h : on décale de -8h pour que 8h→8h devienne 0h→0h
            shifted = self.timestamps_np - 8 * 3600.0
            first_day = np.floor(shifted[0] / 86400.0) * 86400.0
            self.inter_window_idx_np = np.floor((shifted - first_day) / 86400.0).astype(np.int32)
            self.N_windows = int(self.inter_window_idx_np.max()) + 1
            self.inter_window_idx_t = torch.from_numpy(self.inter_window_idx_np).to(dev)
        else:
            self.N_windows = 0
            self.inter_window_idx_t = torch.zeros(0, dtype=torch.long, device=dev)

    # ------------------------------------------------------------------
    def simulate_batch(
        self,
        alloc_cs_idx: torch.Tensor,   # [N_allocs, V] LongTensor
    ) -> torch.Tensor:
        """
        Simule N_allocs allocations en parallèle.

        Parameters
        ----------
        alloc_cs_idx : [N_allocs, V] LongTensor
            Pour chaque allocation n et chaque véhicule v, l'index du CS attribué.

        Returns
        -------
        scores : [N_allocs] float32 — score à minimiser (selon self.metric)
        """
        alloc = alloc_cs_idx.to(self.device)            # [N, V]
        N, V  = alloc.shape
        dev   = self.device
        INF   = float("inf")

        # État courant : deux timestamps par véhicule
        #   vehicle_end_inter : fin d'intervention (dispo pour redirection)
        #   vehicle_at_base   : retour caserne
        vehicle_end_inter = torch.full((N, V), self._t0, dtype=torch.float64, device=dev)
        vehicle_at_base   = torch.full((N, V), self._t0, dtype=torch.float64, device=dev)
        # Type de la dernière inter assignée à chaque véhicule (-1 = aucune)
        vehicle_last_inter_type = torch.full((N, V), -1, dtype=torch.long, device=dev)

        # Accumulation pour les métriques
        response_times = torch.full((N, self.N_inter), INF, dtype=torch.float32, device=dev)
        assigned_cs    = torch.full((N, self.N_inter), -1,  dtype=torch.long,    device=dev)

        # Accumulation pour std_q95_pression : temps engagé par CS par fenêtre 8h→8h
        engaged_per_cs_per_window = torch.zeros(
            (N, self.N_cs, max(self.N_windows, 1)), dtype=torch.float32, device=dev
        )

        row_idx = torch.arange(N, device=dev)  # [N] pour scatter

        # ── Couverture 0 : pré-calcul du nb de véhicules VSAV-éligibles par CS ──
        # Un véhicule est VSAV-éligible s'il peut répondre à une inter VSAV (VSAV ou PSE)
        vsav_eligible_mask = self.type_eligible[VSAV_TYPE][self.vtype_t]   # [V] bool
        # Compter par CS pour chaque allocation : nb_vsav_per_cs[n, c]
        nb_vsav_per_cs = torch.zeros((N, self.N_cs), dtype=torch.long, device=dev)
        for c in range(self.N_cs):
            nb_vsav_per_cs[:, c] = ((alloc == c) & vsav_eligible_mask.unsqueeze(0)).sum(dim=1)
        cov0_seconds = torch.zeros(N, dtype=torch.float64, device=dev)

        for i in range(self.N_inter):
            t         = self.timestamps_t[i].item()           # scalar float64
            inter_typ = self.inter_type_t[i].item()           # 0 ou 1
            dur       = self.duration_t[i].item()             # float32
            dep       = self.depart_t[i].item()               # float32 (temps préparation)

            # ── Disponibilité [N, V] ──────────────────────────────────────
            # Un véhicule est dispo dès qu'il a fini son inter (même en retour)
            avail = vehicle_end_inter <= t                     # [N, V] bool

            # ── Éligibilité par type [V] → broadcast [1, V] ──────────────
            type_ok = self.type_eligible[inter_typ][self.vtype_t]  # [V] bool
            type_ok = type_ok.unsqueeze(0)                          # [1, V]

            # ── Indisponibilité planifiée [V] → broadcast [1, V] ─────────
            unavail_i = self.unavail_t[:, i].unsqueeze(0)           # [1, V]

            # ── Masque éligible de base [N, V] ───────────────────────────
            eligible = avail & type_ok & ~unavail_i                 # [N, V]

            # ── Modularité PSE↔VSAV ───────────────────────────────────────
            is_night = self.night_t[i].item()
            is_rouge = self.rouge_t[i].item()

            # 1) VSAV couplé : interdit la nuit
            if is_night:
                eligible[:, self.vsav_has_pair] = False

            # 2) VSAV couplé : disponible si son PSE pair est en caserne
            #    OU si le PSE est absent mais parti sur une inter VSAV
            pse_base_of_vsav = vehicle_at_base[:, self.vsav_pse_safe]  # [N, V]
            pse_at_base      = pse_base_of_vsav <= t                   # [N, V]
            pse_last_type    = vehicle_last_inter_type[:, self.vsav_pse_safe]  # [N, V]
            pse_on_vsav      = (pse_last_type == VSAV_TYPE) & ~pse_at_base    # [N, V]
            vsav_pair_ok = ~self.vsav_has_pair | pse_at_base | pse_on_vsav    # [N, V]
            eligible = eligible & vsav_pair_ok

            # 3) PSE couplé + intervention POMPE : disponible seulement si VSAV pair en caserne
            if inter_typ == POMPE_TYPE:
                vsav_base_of_pse = vehicle_at_base[:, self.pse_vsav_safe]  # [N, V]
                vsav_at_base     = vsav_base_of_pse <= t                    # [N, V]
                pse_pair_ok = ~self.pse_has_pair | vsav_at_base          # [N, V] broadcast
                eligible = eligible & pse_pair_ok

            # 4) PSE couplé + intervention VSAV + nuit : seulement si rouge
            if inter_typ == VSAV_TYPE and is_night and not is_rouge:
                eligible[:, self.pse_has_pair] = False

            # ── Temps de trajet aller [N, V] ─────────────────────────────
            # alloc[n, v] = cs_idx → travel_t[cs_idx, i]
            travel_i   = self.travel_t[:, i]           # [N_cs]
            ret_travel_i = self.ret_travel_t[:, i]     # [N_cs]

            travel_nv    = travel_i[alloc]             # [N, V]   fancy indexing
            ret_travel_nv= ret_travel_i[alloc]         # [N, V]

            # ── Score = temps de trajet aller (minimiser distance) ────────
            scores = travel_nv.float().clone()
            scores[~eligible] = INF

            # ── Sélection du véhicule ─────────────────────────────────
            has_vehicle = eligible.any(dim=1)           # [N] bool
            best_v_global = scores.argmin(dim=1)        # [N]

            if self.cstc_priority:
                # Priorité CSTC : si un véhicule du CSTC est dispo, on le prend
                cstc_i = self.cstc_idx_t[i]                            # scalar int
                is_cstc = (alloc == cstc_i)                            # [N, V] bool
                cstc_eligible = eligible & is_cstc                     # [N, V]
                has_cstc = cstc_eligible.any(dim=1)                    # [N] bool

                scores_cstc = travel_nv.float().clone()
                scores_cstc[~cstc_eligible] = INF
                best_v_cstc = scores_cstc.argmin(dim=1)                # [N]
                best_v = torch.where(has_cstc, best_v_cstc, best_v_global)
            else:
                # Pas de priorité CSTC : toujours le véhicule le plus proche
                best_v = best_v_global

            # ── Temps de trajet du véhicule choisi ────────────────────────
            best_travel    = travel_nv[row_idx, best_v]       # [N] float32
            best_ret       = ret_travel_nv[row_idx, best_v]   # [N] float32
            best_cs        = alloc[row_idx, best_v]           # [N] int64

            # ── Mise à jour de l'état véhicule ───────────────────────────
            new_end_inter = (
                t
                + best_travel.to(torch.float64)
                + dur
            )                                                  # [N]
            new_at_base = (
                new_end_inter
                + best_ret.to(torch.float64)
            )                                                  # [N]

            valid_n = has_vehicle.nonzero(as_tuple=True)[0]   # indices des sims avec vhl
            if valid_n.numel() > 0:
                vehicle_end_inter[valid_n, best_v[valid_n]] = new_end_inter[valid_n]
                vehicle_at_base[valid_n, best_v[valid_n]]   = new_at_base[valid_n]
                vehicle_last_inter_type[valid_n, best_v[valid_n]] = inter_typ

                # Accumulation pression VSAV : seulement pour les véhicules VSAV-éligibles
                engaged = best_travel[valid_n] + dur + best_ret[valid_n]  # [valid] float32
                w = self.inter_window_idx_t[i].item()
                cs_of_valid = best_cs[valid_n]  # [valid] int64
                # Filtrer : ne garder que les véhicules VSAV ou PSE
                best_v_valid = best_v[valid_n]
                is_vsav_eligible = vsav_eligible_mask[best_v_valid]  # [valid] bool
                vsav_idx = is_vsav_eligible.nonzero(as_tuple=True)[0]
                if vsav_idx.numel() > 0:
                    engaged_per_cs_per_window[valid_n[vsav_idx], cs_of_valid[vsav_idx], w] += engaged[vsav_idx]

            # ── Couverture 0 : compter les occupés VSAV-éligibles par CS ──
            if i < self.N_inter - 1:
                t_next = self.timestamps_t[i + 1].item()
                dt = t_next - t
            else:
                dt = 0.0
            if dt > 0:
                # Un véhicule est "occupé" s'il n'a pas fini son inter
                busy = (vehicle_end_inter > t) & vsav_eligible_mask.unsqueeze(0)  # [N, V]
                # Compter les occupés par CS
                busy_per_cs = torch.zeros((N, self.N_cs), dtype=torch.long, device=dev)
                for c in range(self.N_cs):
                    busy_per_cs[:, c] = (busy & (alloc == c)).sum(dim=1)
                # Couverture 0 : tous les véhicules VSAV-éligibles occupés
                zero_mask = (busy_per_cs >= nb_vsav_per_cs) & (nb_vsav_per_cs > 0)  # [N, N_cs]
                cov0_seconds += zero_mask.sum(dim=1).to(torch.float64) * dt

            # ── Accumulation métriques ─────────────────────────────────────
            resp = dep + best_travel                           # [N] float32 (prép + trajet)
            resp = torch.where(has_vehicle, resp, torch.full_like(resp, INF))
            response_times[:, i] = resp
            assigned_cs[:, i]    = torch.where(has_vehicle, best_cs, torch.full_like(best_cs, -1))

        return self._compute_scores(response_times, assigned_cs, cov0_seconds,
                                    engaged_per_cs_per_window, alloc)

    # ------------------------------------------------------------------
    def _compute_scores(
        self,
        response_times:            torch.Tensor,   # [N, N_inter] float32
        assigned_cs:               torch.Tensor,   # [N, N_inter] int64
        cov0_seconds:              torch.Tensor,   # [N] float64
        engaged_per_cs_per_window: torch.Tensor,   # [N, N_cs, N_windows] float32
        alloc:                     torch.Tensor,   # [N, V] int64
    ) -> torch.Tensor:
        """Calcule le score scalaire à minimiser pour chaque allocation."""
        INF = float("inf")

        if self.metric == "couverture_0":
            return cov0_seconds.float() / 3600.0   # en heures, à minimiser

        if self.metric == "std_q95_pression":
            N = alloc.shape[0]
            dev = self.device
            # Nombre de véhicules VSAV-éligibles (VSAV+PSE) par CS [N, N_cs]
            vsav_elig = self.type_eligible[VSAV_TYPE][self.vtype_t]  # [V] bool
            nb_vsav_per_cs = torch.zeros((N, self.N_cs), dtype=torch.float32, device=dev)
            for c in range(self.N_cs):
                nb_vsav_per_cs[:, c] = ((alloc == c) & vsav_elig.unsqueeze(0)).sum(dim=1).float()
            # Ratio d'utilisation VSAV [N, N_cs, N_windows]
            denom = nb_vsav_per_cs.unsqueeze(2) * WINDOW_8H  # [N, N_cs, 1] × scalar
            denom = denom.clamp(min=1.0)
            ratio = engaged_per_cs_per_window / denom       # [N, N_cs, N_windows]
            # Mettre NaN pour les CS sans véhicules VSAV-éligibles
            no_vhl = (nb_vsav_per_cs == 0).unsqueeze(2)       # [N, N_cs, 1]
            ratio = ratio.masked_fill(no_vhl.expand_as(ratio), float("nan"))
            # Q95 par CS [N, N_cs]
            q95 = torch.nanquantile(ratio, 0.95, dim=2)     # [N, N_cs]
            # Std sur les CS (ignorer NaN)
            # torch.nanmean n'a pas de pendant nanstd, on calcule à la main
            q95_safe = q95.clone()
            valid_cs = ~torch.isnan(q95_safe)
            q95_safe[~valid_cs] = 0.0
            count = valid_cs.float().sum(dim=1).clamp(min=1)
            mean = q95_safe.sum(dim=1) / count
            diff2 = ((q95_safe - mean.unsqueeze(1)) ** 2) * valid_cs.float()
            std = (diff2.sum(dim=1) / count).sqrt()
            return std

        # P95/P99 VSAV : percentile des temps de réponse pour les inter VSAV
        if self.metric in ("p95_vsav", "p99_vsav"):
            q = 0.95 if self.metric == "p95_vsav" else 0.99
            vsav_mask = torch.from_numpy(self.inter_type_np == VSAV_TYPE).to(self.device)  # [N_inter]
            rt_vsav = response_times[:, vsav_mask]                  # [N, N_vsav]
            valid_vsav = rt_vsav < INF                               # [N, N_vsav]
            rt_for_q = rt_vsav.clone()
            rt_for_q[~valid_vsav] = float("nan")
            pq = torch.nanquantile(rt_for_q, q, dim=1)              # [N]
            return pq

        valid = response_times < INF                              # [N, N_inter] bool
        count = valid.float().sum(dim=1).clamp(min=1)            # [N]

        rt_safe = response_times.clone()
        rt_safe[~valid] = 0.0

        temps_moyen = rt_safe.sum(dim=1) / count                 # [N] float32

        if self.metric == "temps_moyen":
            return temps_moyen

        # Taux de couverture (% réponses < 8 min)
        covered      = (response_times <= SEUIL_COUVERTURE) & valid  # [N, N_inter]
        taux_couv    = covered.float().sum(dim=1) / count             # [N]

        if self.metric == "taux_couverture":
            return -taux_couv   # on minimise → négatif

        # Manque local (engin pas du CS compétent)
        cstc_expanded = torch.from_numpy(self.cstc_idx_np).to(self.device).unsqueeze(0)  # [1, N_inter]
        wrong_cs = (assigned_cs != cstc_expanded) & valid        # [N, N_inter]
        manque_local = wrong_cs.float().sum(dim=1) / count        # [N]

        if self.metric == "manque_local":
            return manque_local

        if self.metric == "composite":
            # ── Composite : toutes les métriques normalisées et pondérées ───
            # Valeurs de référence (allocation initiale année 2023)
            REF_TM     = 538.0    # temps moyen (s)
            REF_P95    = 879.0    # p95 temps de réponse VSAV (s)
            REF_COUV   = 0.422    # taux de couverture
            REF_ML     = 0.10     # manque_local ratio (~48992/478063)
            REF_STD_Q95= 0.09     # std q95 pression

            # Composantes normalisées (chacune ~1.0 pour l'allocation initiale)
            c_tm   = temps_moyen / REF_TM                           # [N]
            c_couv = (1.0 - taux_couv) / (1.0 - REF_COUV)          # [N] (à minimiser)
            c_ml   = manque_local / max(REF_ML, 1e-6)               # [N]

            # P95 VSAV
            vsav_mask = torch.from_numpy(self.inter_type_np == VSAV_TYPE).to(self.device)
            rt_vsav = response_times[:, vsav_mask]
            rt_for_q = rt_vsav.clone()
            rt_for_q[rt_vsav >= INF] = float("nan")
            p95_vsav = torch.nanquantile(rt_for_q, 0.95, dim=1)    # [N]
            c_p95  = p95_vsav / REF_P95                             # [N]

            # Std Q95 pression
            N_alloc = alloc.shape[0]
            vsav_elig = self.type_eligible[VSAV_TYPE][self.vtype_t]
            nb_vsav_cs = torch.zeros((N_alloc, self.N_cs), dtype=torch.float32, device=self.device)
            for c in range(self.N_cs):
                nb_vsav_cs[:, c] = ((alloc == c) & vsav_elig.unsqueeze(0)).sum(dim=1).float()
            denom = (nb_vsav_cs.unsqueeze(2) * WINDOW_8H).clamp(min=1.0)
            ratio = engaged_per_cs_per_window / denom
            no_vhl = (nb_vsav_cs == 0).unsqueeze(2)
            ratio = ratio.masked_fill(no_vhl.expand_as(ratio), float("nan"))
            q95_p = torch.nanquantile(ratio, 0.95, dim=2)
            q95_safe = q95_p.clone()
            valid_cs = ~torch.isnan(q95_safe)
            q95_safe[~valid_cs] = 0.0
            cnt = valid_cs.float().sum(dim=1).clamp(min=1)
            mean_q = q95_safe.sum(dim=1) / cnt
            diff2 = ((q95_safe - mean_q.unsqueeze(1)) ** 2) * valid_cs.float()
            std_q95 = (diff2.sum(dim=1) / cnt).sqrt()
            c_std  = std_q95 / max(REF_STD_Q95, 1e-6)              # [N]

            # Pondérations (somme = 1) — sans manque_local
            # temps_moyen: 45%, couverture: 35%, p95: 15%, std_q95: 5%
            score = 0.45 * c_tm + 0.35 * c_couv + 0.15 * c_p95 + 0.05 * c_std
            return score

        # Mixte : 0.5 × temps_moyen_normalisé + 0.3 × (1-couv) + 0.2 × manque
        tm_norm = temps_moyen / SEUIL_COUVERTURE   # normalisation approximative
        score   = 0.5 * tm_norm + 0.3 * (1.0 - taux_couv) + 0.2 * manque_local
        return score

    # ------------------------------------------------------------------
    def evaluate_single(self, alloc_1d: torch.Tensor) -> float:
        """Évalue une seule allocation (utile pour le score de référence)."""
        return self.simulate_batch(alloc_1d.unsqueeze(0))[0].item()

    # ------------------------------------------------------------------
    def simulate_single_detailed(
        self,
        alloc_1d:  torch.Tensor,   # [V] LongTensor — cs_idx par véhicule
        vehicles_np: dict,         # dict retourné par data_loader.load_data()
        cs_list:     list,         # liste des noms de CS
    ) -> "pd.DataFrame":
        """
        Simule une seule allocation et retourne le détail complet par intervention,
        équivalent GPU de la liste InterventionSimulee du simulateur séquentiel.

        Colonnes du DataFrame retourné
        -------------------------------
        inter_idx         int    index dans la fenêtre (0-based)
        timestamp         float  unix timestamp de l'intervention
        fem_mma           str    "POMPE" ou "VSAV"
        cstc              str    CS compétent
        proc              str    procédure ("R" = rouge)
        vehicle_idx       int    index interne du véhicule (0..V-1)
        vehicle_id               identifiant original du véhicule
        vehicle_type      str    "POMPE" / "VSAV" / "PSE"
        vehicle_cs        str    CS du véhicule dans cette allocation
        travel_time_s     float  temps de trajet aller (secondes)
        return_time_s     float  temps de trajet retour (secondes)
        depart_s          float  temps de préparation au départ (secondes)
        response_time_s   float  depart_s + travel_time_s
        covered           bool   response_time_s <= 480 s (8 min)
        manque_local      bool   vehicle_cs != cstc
        assigned          bool   False = aucun véhicule trouvé
        """
        import pandas as pd

        alloc = alloc_1d.to(self.device)   # [V]
        V     = alloc.shape[0]
        dev   = self.device
        INF   = float("inf")

        TYPE_NAMES = {0: "POMPE", 1: "VSAV", 2: "PSE"}
        INTER_TYPE_NAMES = {0: "POMPE", 1: "VSAV"}

        vehicle_end_inter = torch.full((V,), self._t0, dtype=torch.float64, device=dev)
        vehicle_at_base   = torch.full((V,), self._t0, dtype=torch.float64, device=dev)
        vehicle_last_inter_type = torch.full((V,), -1, dtype=torch.long, device=dev)

        rows = []

        for i in range(self.N_inter):
            t         = self.timestamps_t[i].item()
            inter_typ = self.inter_type_t[i].item()
            dur       = self.duration_t[i].item()
            dep       = self.depart_t[i].item()

            avail     = vehicle_end_inter <= t                               # [V]
            type_ok   = self.type_eligible[inter_typ][self.vtype_t]         # [V]
            not_unavail = ~self.unavail_t[:, i]                             # [V]
            eligible  = avail & type_ok & not_unavail                       # [V]

            travel_i     = self.travel_t[:, i]       # [N_cs]
            ret_travel_i = self.ret_travel_t[:, i]   # [N_cs]

            travel_v     = travel_i[alloc]            # [V]
            ret_travel_v = ret_travel_i[alloc]        # [V]

            scores = travel_v.float().clone()
            scores[~eligible] = INF

            has_vehicle = eligible.any().item()

            if has_vehicle:
                if self.cstc_priority:
                    # Priorité CSTC
                    cstc_i = self.cstc_idx_t[i].item()
                    is_cstc = (alloc == cstc_i)
                    cstc_eligible = eligible & is_cstc
                    has_cstc = cstc_eligible.any().item()
                    if has_cstc:
                        scores_cstc = travel_v.float().clone()
                        scores_cstc[~cstc_eligible] = INF
                        best_v = scores_cstc.argmin().item()
                    else:
                        best_v = scores.argmin().item()
                else:
                    best_v = scores.argmin().item()
                best_travel  = travel_v[best_v].item()
                best_ret     = ret_travel_v[best_v].item()
                best_cs_idx  = alloc[best_v].item()

                vehicle_end_inter[best_v] = t + best_travel + dur
                vehicle_at_base[best_v]   = t + best_travel + dur + best_ret
                vehicle_last_inter_type[best_v] = inter_typ

                vehicle_id   = vehicles_np["ids"][best_v]
                vehicle_type_str = TYPE_NAMES[int(vehicles_np["type"][best_v])]
                vehicle_cs_str   = cs_list[best_cs_idx]
            else:
                best_v = best_travel = best_ret = best_cs_idx = None
                vehicle_id = vehicle_type_str = vehicle_cs_str = None

            cstc_idx = self.cstc_idx_np[i]
            cstc_str = cs_list[cstc_idx] if cstc_idx >= 0 else "?"
            resp     = dep + best_travel if has_vehicle else INF

            rows.append(dict(
                inter_idx       = i,
                timestamp       = t,
                fem_mma         = INTER_TYPE_NAMES[inter_typ],
                cstc            = cstc_str,
                proc            = str(self.proc_np[i]),
                assigned        = has_vehicle,
                vehicle_idx     = best_v,
                vehicle_id      = vehicle_id,
                vehicle_type    = vehicle_type_str,
                vehicle_cs      = vehicle_cs_str,
                travel_time_s   = best_travel if has_vehicle else INF,
                return_time_s   = best_ret    if has_vehicle else INF,
                depart_s        = dep,
                response_time_s = resp,
                covered         = resp <= SEUIL_COUVERTURE if has_vehicle else False,
                manque_local    = (vehicle_cs_str != cstc_str) if has_vehicle else None,
            ))

        return pd.DataFrame(rows)

    # ------------------------------------------------------------------
    @staticmethod
    def compute_all_metrics(df: "pd.DataFrame") -> dict:
        """
        Calcule toutes les métriques à partir du DataFrame retourné par
        simulate_single_detailed() — même structure que eval_window() du
        simulateur séquentiel.

        Returns
        -------
        dict avec clés :
          "global" : métriques agrégées sur toutes les interventions
          "cs"     : dict {cs_name → métriques} par CS compétent (cstc)
        """
        import numpy as np
        import pandas as pd

        SEUIL = 480.0   # 8 min en secondes

        df_ok = df[df["assigned"]].copy()

        def _metrics_for(sub: "pd.DataFrame") -> dict:
            if sub.empty:
                return dict(
                    n=0,
                    temps_moyen=float("nan"),
                    temps_moyen_pompe=float("nan"),
                    temps_moyen_vsav=float("nan"),
                    temps_moyen_rouge=float("nan"),
                    taux_couverture=float("nan"),
                    p95=float("nan"),
                    p95_vsav=float("nan"),
                    p99_vsav=float("nan"),
                    manque_local=0,
                    manque_local_pompe=0,
                    manque_local_vsav=0,
                    manque_local_rouge=0,
                    n_non_assignees=0,
                )

            rt = sub["response_time_s"].values
            pompe  = sub[sub["fem_mma"] == "POMPE"]
            vsav   = sub[sub["fem_mma"] == "VSAV"]
            rouge  = sub[sub["proc"] == "R"]
            manque = sub["manque_local"].fillna(False)

            vsav_rt = vsav["response_time_s"].values if len(vsav) else None
            return dict(
                n                  = len(sub),
                temps_moyen        = rt.mean(),
                temps_moyen_pompe  = pompe["response_time_s"].mean() if len(pompe) else float("nan"),
                temps_moyen_vsav   = vsav["response_time_s"].mean()  if len(vsav)  else float("nan"),
                temps_moyen_rouge  = rouge["response_time_s"].mean() if len(rouge) else float("nan"),
                taux_couverture    = (rt <= SEUIL).mean(),
                p95                = float(np.percentile(rt, 95)),
                p95_vsav           = float(np.percentile(vsav_rt, 95)) if vsav_rt is not None else float("nan"),
                p99_vsav           = float(np.percentile(vsav_rt, 99)) if vsav_rt is not None else float("nan"),
                manque_local       = int(manque.sum()),
                manque_local_pompe = int(pompe["manque_local"].fillna(False).sum()),
                manque_local_vsav  = int(vsav["manque_local"].fillna(False).sum()),
                manque_local_rouge = int(rouge["manque_local"].fillna(False).sum()),
                n_non_assignees    = int((~df["assigned"]).sum()),   # sur tout df
            )

        global_metrics = _metrics_for(df_ok)
        global_metrics["n_non_assignees"] = int((~df["assigned"]).sum())

        cs_metrics = {}
        for cs, grp in df_ok.groupby("cstc"):
            cs_metrics[cs] = _metrics_for(grp)

        return {"global": global_metrics, "cs": cs_metrics}