Navigation : Index # CC3 : SmartGrid - Ordonnancement de la production energetique sous incertitude

Étude de cas interdisciplinaire combinant programmation par contraintes (OR-Tools CP-SAT), inference probabiliste (modèle bayesien), optimisation multi-objectif et un jumeau numérique de reseau electrique.

Contexte metier

Un opérateur de reseau electrique doit, pour chaque heure, decider quelles centrales activer et a quel niveau de production (le unit commitment / dispatch) afin de satisfaire la demande tout en minimisant le cout economique et les emissions de CO2. Ce problème est rendu difficile par :

  • l’incertitude sur la production renouvelable (eolien, solaire) et sur la demande ;
  • les contraintes dures : capacites min/max des centrales, temps de montee/descente, reserve tournante (securite n-1) ;
  • les objectifs multiples conflictuels : cout, emissions, fiabilite.

Une seule technique ne suffit pas : un solveur déterministe ignore l’incertitude, un modèle probabiliste seul ne respecte pas les contraintes physiques. Cette étude de cas materialise leur composition ordonnee : filtrer les dispatches impossibles (contraintes), modeliser l’aleatoire (incertitude), puis optimiser (multi-objectif).

Architecture en 4 couches

Couche Rôle Technique Partie
Contraintes dures Filtrer les dispatches impossibles OR-Tools CP-SAT Partie 1
Incertitude Modeliser l’aleatoire renouvelable Modèle bayesien Partie 2
Optimisation Choisir le meilleur compromis Multi-objectif pondere Partie 3
Decision Synthese + comparaison Workflow + scoring Partie 4

Ce notebook est une fondation : le jumeau numérique (couche 0) est operationnel et les 4 couches de solveurs (CP-SAT, bayésien, multi-objectif, analyse comparative) sont implémentees et exécutées (Parties 1 à 4). Les 3 exercices d’extension (Exercices 1 à 3) sont volontairement laissés en stub pour les étudiants.

0. Installation et jumeau numérique

Le jumeau numérique represente un reseau electrique : un ensemble de centrales (avec cout, capacites, facteur d’emission) et une demande horaire a satisfaire, avec une production renouvelable stochastique.

# Installation des dependances (solveurs des Parties 1-4)
# %pip install ortools numpy scipy  # decommenter pour execution complete

from dataclasses import dataclass, field
from typing import List, Dict
import numpy as np


@dataclass
class PowerPlant:
    """Centrale de production electrique."""
    name: str
    pmin: float          # puissance minimale stable (MW)
    pmax: float          # puissance maximale (MW)
    cost_per_mw: float   # cout variable (EUR/MWh)
    co2_per_mw: float    # emission (kg/MWh)
    is_renewable: bool = False


@dataclass
class PowerNetwork:
    """Jumeau numerique d'un reseau electrique horaire."""
    plants: List[PowerPlant]
    demand_mw: List[float]               # demande par heure (MW)
    renewable_forecast_mean: List[float]  # production renouvelable moyenne attendue (MW)
    renewable_forecast_std: List[float]   # incertitude (ecart-type, MW)

    def total_capacity(self) -> float:
        return sum(p.pmax for p in self.plants)

    def net_demand(self, hour: int) -> float:
        """Demande a couvrir par les centrales pilotables apres renouvelable."""
        return max(0.0, self.demand_mw[hour] - self.renewable_forecast_mean[hour])


def demo_network() -> PowerNetwork:
    """Reseau de demonstration : 3 centrales pilotables + renouvelable sur 6 heures."""
    plants = [
        PowerPlant('Charbon', 50, 400, cost_per_mw=30, co2_per_mw=900),
        PowerPlant('Gaz',    20, 250, cost_per_mw=60, co2_per_mw=400),
        PowerPlant('Hydro',   0, 150, cost_per_mw=10, co2_per_mw=0,  is_renewable=False),
    ]
    return PowerNetwork(
        plants=plants,
        demand_mw=[500, 520, 480, 540, 600, 580],
        renewable_forecast_mean=[80, 60, 40, 50, 120, 150],
        renewable_forecast_std=[15, 12, 10, 12, 25, 30],
    )


net = demo_network()
print(f'Reseau : {len(net.plants)} centrales, capacite totale {net.total_capacity()} MW')
for h in range(len(net.demand_mw)):
    print(f'  H{h}: demande {net.demand_mw[h]} MW, renouvelable {net.renewable_forecast_mean[h]} '
          f'+/- {net.renewable_forecast_std[h]}, demande nette {net.net_demand(h):.1f} MW')
Reseau : 3 centrales, capacite totale 800 MW
  H0: demande 500 MW, renouvelable 80 +/- 15, demande nette 420.0 MW
  H1: demande 520 MW, renouvelable 60 +/- 12, demande nette 460.0 MW
  H2: demande 480 MW, renouvelable 40 +/- 10, demande nette 440.0 MW
  H3: demande 540 MW, renouvelable 50 +/- 12, demande nette 490.0 MW
  H4: demande 600 MW, renouvelable 120 +/- 25, demande nette 480.0 MW
  H5: demande 580 MW, renouvelable 150 +/- 30, demande nette 430.0 MW

Partie 1 : Le Dispatcher sous Contraintes (OR-Tools CP-SAT)

Théorie : le unit commitment est un problème NP-difficile. On modelise, pour chaque heure et chaque centrale pilotable, une variable binaire (activee) et une variable continue (niveau de production), soumises aux contraintes :

  • somme des productions = demande nette (equilibre offre/demande) ;
  • pmin <= production <= pmax quand la centrale est activee ;
  • reserve tournante : capacite excédentaire disponible >= marge de securite.

L’objectif primaire : minimiser le cout total.

# Partie 1 : modele CP-SAT du dispatch (unit commitment)
from ortools.sat.python import cp_model


def solve_dispatch_cp(network, cost_weights=None):
    """Resout le unit commitment par programmation par contraintes (CP-SAT).

    Pour chaque heure et centrale pilotable : binaire on (activee) + entier
    p (production MW). Contraintes : equilibre offre/demande nette + bornes
    de capacite (pmin <= p <= pmax quand active). Objectif pondere : on peut
    minimiser le cout (couts variables) ou le CO2 (facteurs d'emission).

    Renvoie un dict {hour: {plant_name: production_mw}}, '_cost', '_status'.
    cost_weights = (w_cost, w_co2) : poids de chaque terme dans l'objectif.
    """
    if cost_weights is None:
        cost_weights = (1.0, 0.0)
    H = len(network.demand_mw)
    plants = network.plants
    model = cp_model.CpModel()

    on, p = {}, {}
    for h in range(H):
        for i, plant in enumerate(plants):
            on[(h, i)] = model.NewBoolVar(f'on_{h}_{i}')
            p[(h, i)] = model.NewIntVar(0, int(plant.pmax), f'p_{h}_{i}')

    for h in range(H):
        nd = int(round(network.net_demand(h)))
        model.Add(sum(p[(h, i)] for i in range(len(plants))) == nd)
        for i, plant in enumerate(plants):
            model.Add(p[(h, i)] >= int(plant.pmin) * on[(h, i)])
            model.Add(p[(h, i)] <= int(plant.pmax) * on[(h, i)])

    w_cost, w_co2 = cost_weights
    # Objectif pondere cout+CO2. NB : on plie les poids DANS chaque terme et on somme
    # via cp_model.LinearExpr.Sum (jamais `poids * SumArray` ni le builtin sum()) : quand
    # un poids vaut 0, `0 * IntVar` produit un IntConstant(0) (objet) qu'on ne peut pas
    # additionner a une SumArray (__radd__ attend un float). Sum() accepte une liste
    # melangeant IntAffine et IntConstant sans probleme.
    wc, wco2 = int(round(w_cost)), int(round(w_co2))
    objective_terms = []
    for h in range(H):
        for i, plant in enumerate(plants):
            objective_terms.append(wc * int(round(plant.cost_per_mw)) * p[(h, i)])
            objective_terms.append(wco2 * int(round(plant.co2_per_mw)) * p[(h, i)])
    model.Minimize(cp_model.LinearExpr.Sum(objective_terms))

    solver = cp_model.CpSolver()
    solver.parameters.max_time_in_seconds = 10.0
    status = solver.Solve(model)
    if status not in (cp_model.OPTIMAL, cp_model.FEASIBLE):
        return None
    result = {}
    for h in range(H):
        result[h] = {plants[i].name: solver.Value(p[(h, i)]) for i in range(len(plants))}
    result['_cost'] = round(solver.ObjectiveValue(), 1)
    result['_status'] = solver.StatusName(status)
    return result


dispatch = solve_dispatch_cp(net)
if dispatch:
    print('Dispatch CP-SAT (min cout) - statut:', dispatch['_status'])
    for h in range(len(net.demand_mw)):
        prod = ', '.join(f'{k}={v}MW' for k, v in dispatch[h].items() if v > 0)
        print(f'  H{h}: {prod}')
    print(f'  Cout objectif: {dispatch["_cost"]}')
else:
    print('Aucune solution faisable (capacite < demande)')
Dispatch CP-SAT (min cout) - statut: OPTIMAL
  H0: Charbon=270MW, Hydro=150MW
  H1: Charbon=310MW, Hydro=150MW
  H2: Charbon=290MW, Hydro=150MW
  H3: Charbon=340MW, Hydro=150MW
  H4: Charbon=330MW, Hydro=150MW
  H5: Charbon=280MW, Hydro=150MW
  Cout objectif: 63600.0

Interprétation : le dispatch min-coût et ses limites

Le statut OPTIMAL confirme que CP-SAT a trouvé le dispatch de coût minimal satisfaisant toutes les contraintes (équilibre offre/demande, capacités min/max). Lecture du résultat :

  • Le charbon porte l’essentiel de la charge (270–340 MW) et s’ajuste heure par heure pour suivre la demande ; l’hydroélectricité reste plate à 150 MW, vraisemblablement son plafond de capacité (pmax).
  • Le coût optimal de 63 600 € est le plancher économique : c’est la stratégie « Économique » que l’on retrouvera dans la comparaison finale (Partie 4).
  • Limite de cette stratégie : minimiser le coût seul pousse vers le charbon (le moins cher par MWh), mais ce sont aussi les émissions de CO2 les plus élevées. C’est précisément ce déséquilibre que les Parties 2 (risque) et 3 (multi-objectif) viennent corriger en intégrant d’autres axes que le coût.

Exercice 1 : Reserve tournante (securite n-1)

Ajoutez au modèle CP-SAT une contrainte de reserve tournante : a chaque heure, la somme des capacites disponibles (non utilisees) des centrales activees doit couvrir la perte de la plus grande centrale en service (critere n-1).

Indice : pour chaque heure, la reserve = sum(pmax_i - prod_i pour i active) doit etre >= max(pmax_i pour i active).

# Exercice 1 : contrainte de reserve tournante n-1
def solve_dispatch_with_spinning_reserve(network: PowerNetwork):
    """Dispatch CP-SAT avec reserve tournante n-1."""
    result = None  # TODO etudiant : etendre solve_dispatch_cp avec la reserve n-1
    return result

Partie 2 : Le Previsionniste Probabiliste (modèle bayesien)

Théorie : la production renouvelable est aleatoire. On modelise l’ecart demande - renouvelable comme une variable aleatoire et on infere sa distribution, pour quantifier le risque de defaillance (probabilite que la demande nette depasse la capacite pilotable disponible).

# Partie 2 : modele bayesien de l'incertitude renouvelable
from math import sqrt, erfc


def failure_risk(network, hour):
    """Probabilite que la demande nette reelle depasse la capacite pilotable.

    Modele : la production renouvelable reelle ~ N(mean, std). La demande nette
    reelle = demande - renouvelable_reelle ~ N(net_demand, renewable_std).
    Risque = P(demande nette reelle > capacite pilotable totale), calcule via
    la survivor function d'une gaussienne : SF(z) = 0.5 * erfc(z / sqrt(2)).
    """
    mean_net = network.net_demand(hour)
    std = network.renewable_forecast_std[hour]
    capacity = network.total_capacity()
    if std <= 0:
        return 1.0 if mean_net > capacity else 0.0
    z = (capacity - mean_net) / std
    return 0.5 * erfc(z / sqrt(2))


print('Risque de defaillance (capacite pilotable =', net.total_capacity(), 'MW) :')
for h in range(len(net.demand_mw)):
    r = failure_risk(net, h)
    print(f'  H{h}: demande nette moyenne {net.net_demand(h):.0f} MW, '
          f'incertitude +/-{net.renewable_forecast_std[h]} MW -> risque {r:.2e}')
Risque de defaillance (capacite pilotable = 800 MW) :
  H0: demande nette moyenne 420 MW, incertitude +/-15 MW -> risque 6.86e-142
  H1: demande nette moyenne 460 MW, incertitude +/-12 MW -> risque 6.72e-177
  H2: demande nette moyenne 440 MW, incertitude +/-10 MW -> risque 4.18e-284
  H3: demande nette moyenne 490 MW, incertitude +/-12 MW -> risque 1.87e-147
  H4: demande nette moyenne 480 MW, incertitude +/-25 MW -> risque 8.20e-38
  H5: demande nette moyenne 430 MW, incertitude +/-30 MW -> risque 3.00e-35

Interprétation : un risque de queue négligeable, piloté par l’incertitude

Les risques de défaillance sont astronomiquement faibles (de l’ordre de 10⁻³⁵ à 10⁻²⁸⁴) : avec 800 MW de capacité pilotable pour une demande nette de 420–490 MW, la marge est si large que la probabilité d’un dépassement (queue de distribution gaussienne via erfc) est numériquement négligeable.

  • Le risque ne suit pas le niveau de demande mais l’incertitude de la prévision renouvelable : H5 (demande 430 MW, incertitude ±30 MW) affiche le risque le plus élevé (≈ 3·10⁻³⁵), très au-dessus de H0 (demande 420 MW, ±15 MW, ≈ 7·10⁻¹⁴²) — alors que la demande de H0 est plus faible. C’est l’écart-type renewable_forecast_std, pas la demande, qui pilote la probabilité de queue.
  • Enseignement : tant que la capacité pilotable surpasse largement la demande nette, la fiabilité est quasi gratuite. Le risque ne devient pertinent qu’en cas de marge tendue ou d’incertitude élevée — c’est ce que l’Exercice 2 explore en amplifiant l’incertitude renouvelable.

Exercice 2 : Sensibilite au prior sur la variabilite eolienne

Analysez comment le risque de defaillance (Partie 2) evolue quand on modifie l’incertitude du renouvelable (multipliez renewable_forecast_std par 0.5, 1.0, 2.0). Quelle heure devient la plus risqueuse ? Concluez sur la robustesse du dispatch determine en Partie 1.

# Exercice 2 : analyse de sensibilite au prior d'incertitude
def sensitivity_to_renewable_std(network: PowerNetwork, factors=(0.5, 1.0, 2.0)):
    """Retourne le risque max par heure pour chaque facteur d'incertitude."""
    result = None  # TODO etudiant : boucler sur les facteurs, retourner {factor: max_risk}
    return result

Partie 3 : L’Optimiseur Multi-Objectif

Théorie : on veut minimiser simultanement le cout economique, les emissions CO2 et le risque de defaillance. Ces objectifs sont conflictuels (le charbon est bon marche mais très emissif). On construit un score pondere S = a*cout + b*co2 + c*risque et on cherche le dispatch minimisant S, en explorant le front de Pareto.

# Partie 3 : score multi-objectif (cout + CO2 + risque)
def dispatch_metrics(dispatch, network):
    """Calcule cout total, CO2 total et risque max d'un dispatch CP-SAT."""
    if dispatch is None:
        return {'cost': float('inf'), 'co2': float('inf'), 'max_risk': 1.0}
    cost = co2 = 0.0
    for h in range(len(network.demand_mw)):
        for plant in network.plants:
            prod = dispatch[h].get(plant.name, 0)
            cost += prod * plant.cost_per_mw
            co2 += prod * plant.co2_per_mw
    risks = [failure_risk(network, h) for h in range(len(network.demand_mw))]
    return {'cost': cost, 'co2': co2, 'max_risk': max(risks)}


def multi_objective_score(dispatch, network, weights=(1.0, 1.0, 1.0), norms=None):
    """Score = w_c*cout_norm + w_co2*co2_norm + w_r*risque_norm, normalise par les normes de reference (1.0 = composante a sa norme ; 0 = meilleur)."""
    m = dispatch_metrics(dispatch, network)
    if norms is None:
        norms = {'cost': m['cost'] or 1.0, 'co2': m['co2'] or 1.0, 'max_risk': 1.0}
    wc, wco2, wr = weights
    return (wc * (m['cost'] / max(norms['cost'], 1e-9))
            + wco2 * (m['co2'] / max(norms['co2'], 1e-9))
            + wr * (m['max_risk'] / max(norms['max_risk'], 1e-9)))


score_min_cost = multi_objective_score(dispatch, net, weights=(1.0, 1.0, 1.0),
                                       norms={'cost': 200000, 'co2': 2000000, 'max_risk': 1.0})
print(f'Score multi-objectif du dispatch min-cout (composantes normalisees par reference, 0 = meilleur): {score_min_cost:.3f}')
Score multi-objectif du dispatch min-cout (composantes normalisees par reference, 0 = meilleur): 1.137

Interprétation : décomposition du score multi-objectif

Le score 1,137 du dispatch min-coût se décompose, avec les poids (1, 1, 1) et les normes {cost: 200 000, co2: 2 000 000, max_risk: 1,0} :

  • coût : 63 600 / 200 000 ≈ 0,32 (faible → bon marché) ;
  • CO2 : 1 638 000 / 2 000 000 ≈ 0,82 (élevé → polluant) ;
  • risque : ≈ 0 (la capacité excède largement la demande, cf. Partie 2).

Le terme CO2 domine (0,82) : optimiser le coût seul produit un dispatch bon marché mais sale. C’est la traduction chiffrée du conflit d’objectifs. Le score multi-objectif rend ce compromis comparable entre stratégies (Partie 4) : un score plus bas indique un meilleur équilibre global coût/CO2/risque.

Exercice 3 : Ajouter un troisieme objectif (equite territoriale)

Ajoutez un troisieme objectif : equite territoriale (eviter qu’une même region subisse toutes les centrales les plus polluantes). Supposons que chaque centrale a un attribut region. Mesurez la concentration regionale de la pollution et ajoutez-la au score.

# Exercice 3 : objectif d'equite territoriale
def territorial_equity_score(dispatch, network: PowerNetwork) -> float:
    """Mesure de concentration de la pollution par region (0 = equitable)."""
    score = None  # TODO etudiant : mesurer la concentration (ex: variance/entropy regionale des co2)
    return score

Partie 4 : Decision finale et analyse comparative

Synthese : on compare les stratégies (CP-SAT pur, multi-objectif weighted-sum, front de Pareto) sur les 3 axes cout/CO2/fiabilite + l’equite territoriale. Le tableau comparatif ci-dessous est produit a partir des Parties 1-3 deja executees.

# Partie 4 : analyse comparative des strategies de dispatch
def comparative_analysis(network):
    """Compare 3 strategies CP-SAT (economique, ecologique, compromis)."""
    eco = solve_dispatch_cp(network, cost_weights=(1.0, 0.0))    # min cout
    green = solve_dispatch_cp(network, cost_weights=(0.0, 1.0))  # min CO2
    comp = solve_dispatch_cp(network, cost_weights=(1.0, 3.0))   # compromis

    m_eco = dispatch_metrics(eco, network)
    m_green = dispatch_metrics(green, network)
    m_comp = dispatch_metrics(comp, network)
    norms = {
        'cost': max(m_eco['cost'], m_green['cost'], 1.0),
        'co2': max(m_eco['co2'], m_green['co2'], 1.0),
        'max_risk': max(m_eco['max_risk'], m_green['max_risk'], 1e-9),
    }
    s_eco = multi_objective_score(eco, network, weights=(1.0, 1.0, 1.0), norms=norms)
    s_green = multi_objective_score(green, network, weights=(1.0, 1.0, 1.0), norms=norms)
    s_comp = multi_objective_score(comp, network, weights=(1.0, 1.0, 1.0), norms=norms)

    table = {}
    for name, d, m, s in [('Economique (min cout)', eco, m_eco, s_eco),
                          ('Ecologique (min CO2)', green, m_green, s_green),
                          ('Compromis (pondere)', comp, m_comp, s_comp)]:
        table[name] = {'cout_EUR': round(m['cost']), 'co2_kg': round(m['co2']),
                       'risque_max': m['max_risk'], 'score': round(s, 3)}
    return table


analysis = comparative_analysis(net)
print('Analyse comparative des strategies de dispatch :')
print(f'{"Strategie":<25} {"Cout (EUR)":>12} {"CO2 (kg)":>12} {"Risque max":>12} {"Score":>8}')
print('-' * 75)
for name, m in analysis.items():
    print(f'{name:<25} {m["cout_EUR"]:>12} {m["co2_kg"]:>12} '
          f'{m["risque_max"]:>12.2e} {m["score"]:>8.3f}')
print()
best = min(analysis.items(), key=lambda kv: kv[1]['score'])
print(f'Strategie optimale (score multi-objectif minimal) : {best[0]}')
Analyse comparative des strategies de dispatch :
Strategie                   Cout (EUR)     CO2 (kg)   Risque max    Score
---------------------------------------------------------------------------
Economique (min cout)            63600      1638000     3.00e-35    1.596
Ecologique (min CO2)            106800       918000     3.00e-35    1.560
Compromis (pondere)             106800       918000     3.00e-35    1.560

Strategie optimale (score multi-objectif minimal) : Ecologique (min CO2)

Interprétation : le coût de la transition bas-carbone

Les trois stratégies sont ici scorées avec des poids égaux (1, 1, 1) et des normes tirées des extrêmes (coût et CO2 max). Décomposition des scores :

  • Économique (1,596) = coût 0,60 + CO2 1,00 + risque 0 : excellent sur le coût, maximal sur le CO2.
  • Écologique (1,560) = coût 1,00 + CO2 0,56 + risque 0 : inverse, bon sur le CO2, maximal sur le coût.

L’Écologique l’emporte (1,560 < 1,596) bien qu’elle coûte plus cher (106 800 € contre 63 600 €) : la réduction de CO2 (1 638 000 kg émis par l’Économique contre 918 000, soit 720 000 kg évités / −44 %) compense la perte économique. Compromis = Écologique : les deux stratégies convergent vers un dispatch identique (mêmes coût, CO2 et risque). Avec les poids de recherche (1, 3) donnant 3 fois plus de poids au CO2, la somme pondérée sélectionne en fait la solution min-CO2 — le compromis coïncide ici avec l’extrémité écologique du front.

Coût marginal d’évitement : passer du min-coût au min-CO2 coûte +43 200 € (+68 %) pour éviter 720 000 kg de CO2, soit environ 0,06 €/kg de CO2 évité. C’est le type d’arbitrage chiffré qu’un opérateur de réseau doit expliciter avant de décider.

Conclusion

Cette étude de cas illustre la composition ordonnee des paradigmes IA sur un problème metier reel (transition energetique) : on filtre les dispatches impossibles (CP-SAT) avant de modeliser l’incertitude (bayesien) avant d’optimiser (multi-objectif). Inverser l’ordre produit un système soit trop rigide (contraintes ignorent l’aleatoire), soit trop flou (decision sans contraintes physiques).

Etat de la fondation : jumeau numérique operationnel (couche 0) + 4 couches de solveurs implémentees (CP-SAT en Partie 1, bayésien en Partie 2, multi-objectif en Partie 3, analyse comparative en Partie 4 — toutes avec sorties réelles). Les 3 exercices d’extension (Exercices 1 à 3) restent en stub pour les étudiants.

Retour au sommet