Planification de la Production d'Électricité
Optimisation de la production électrique par programmation linéaire en nombres entiers mixtes (MILP)
Introduction
Ce projet analyse et optimise la planification journalière de la production d’électricité en utilisant la programmation linéaire en nombres entiers mixtes (MILP).
Ce travail a été réalisé dans le cadre du cours de Sophie DEMASSEY à l’École des Mines de Paris (PSL).
Le projet étudie différentes configurations de centrales électriques (thermiques et hydroélectriques) afin de minimiser les coûts de production tout en respectant les contraintes de demande et de réserve.
L’objectif principal est de déterminer le nombre optimal de centrales à activer et leurs puissances à chaque période de la journée, en tenant compte des coûts de fonctionnement, des coûts de démarrage et des contraintes opérationnelles.
Il faut tenir en compte que ce projet simplifie la nature réelle du problème, en supposant que la demande est connue et fixe avant la planification. En réalité, la demande est une variable aléatoire, et la planification doit tenir compte de la réserve pour assurer la stabilité de la réseau électrique.
Complément théorique sur l’optimisation linéaire
La programmation linéaire (PL) consiste à minimiser (ou maximiser) une fonction linéaire soumise à un ensemble de contraintes linéaires. Sous sa forme standard, un problème primal s’écrit :
\[\min_{x} ; c^T x \quad \text{sous} \quad A x = b, ; x \ge 0\]À ce problème est associé un problème dual, exprimé par :
\[\max_{y} ; b^T y \quad \text{sous} \quad A^T y \le c\]La faisabilité primale signifie qu’il existe au moins un vecteur \(x\) respectant les contraintes \(A x = b\) et \(x \ge 0\), tandis que la faisabilité duale impose l’existence d’un vecteur \(y\) satisfaisant \(A^T y \le c\). Lorsqu’un problème admet une solution réalisable pour les deux formulations (primal et dual), la dualité forte garantit que leurs valeurs optimales coïncident :
\[c^T x^* = b^T y^*\]Cette propriété fonde une grande partie de la théorie économique de l’optimisation : les multiplicateurs de Lagrange \(y_i\) associés au dual peuvent être interprétés comme des coûts marginaux ou prix d’ombre, représentant la valeur économique d’une unité supplémentaire de ressource \(b_i\).
Dans la pratique, les conditions dites de complémentarité relient les solutions primales et duales à l’optimum :
\[x_j (c_j - A_j^T y) = 0, \quad y_i (A_i x - b_i) = 0\]Ces relations démontrent que :
- une variable \(x_j\) strictement positive correspond à une contrainte duale saturée,
- et une contrainte primale saturée implique un multiplicateur dual non nul.
Ainsi, l’étude conjointe du primal et du dual permet non seulement d’identifier la solution optimale, mais aussi de comprendre la sensibilité du système aux variations de la demande, des coûts ou des capacités de production.
Données du Problème
Caractéristiques des Centrales Thermiques
| Catégorie | Nombre | Pmin (MW) | Pmax (MW) | Cmwh (€/MWh) | Cbase (€/h) | Cstart (€) |
|---|---|---|---|---|---|---|
| A | 12 | 850 | 2000 | 2.0 | 1000 | 2000 |
| B | 10 | 1250 | 1750 | 1.3 | 2600 | 1000 |
| C | 5 | 1500 | 4000 | 3.0 | 3000 | 500 |
Profil de Demande Journalière
| Période | Heures | Demande (MW) |
|---|---|---|
| 00h - 06h | 6 | 15000 |
| 06h - 09h | 3 | 30000 |
| 09h - 15h | 6 | 25000 |
| 15h - 18h | 3 | 40000 |
| 18h - 00h | 6 | 27000 |
Évolution des Modèles
Modèle 1.5 : Modèle de Base avec Puissance Minimale
Le premier modèle intègre la contrainte de puissance minimale pour chaque centrale activée :
\[\min \sum_{t=1}^{3} \sum_{p=1}^{5} P_{tp} \cdot C_{t} \cdot T_{p}\]Sous les contraintes :
\[\sum_{t=1}^{3} P_{tp} \geq D_{p}, \quad \forall p \in \{1, 2, 3, 4, 5\} \quad (C1)\] \[N_{tp} P^{\text{min}}_t \leq P_{tp} \leq N_{tp} P^{\text{max}}_t, \quad \forall t, p \quad (C2)\] \[N_{tp} \leq Nb_t, \quad \forall t, p \quad (C3)\]où :
- \(P_{tp}\) : puissance totale du type \(t\) à la période \(p\)
- \(N_{tp}\) : nombre de centrales actives du type \(t\) à la période \(p\)
- \(C_t\) : coût marginal (€/MWh)
- \(T_p\) : durée de la période \(p\)
- \(D_p\) : demande de la période \(p\)
Résultat : Coût optimal = 869 k€
En Gurobi, les variables sont déclarées comme suit :
# Variables
N = {} # Nombre de centrales actives (entières)
P = {} # Puissance produite (continues)
for t in USINES.keys():
for p in HOURS.keys():
N[f"{t},{p}"] = m.addVar(lb=0, name=f"N_{t},{p}", vtype=gp.GRB.INTEGER)
P[f"{t},{p}"] = m.addVar(lb=0, name=f"P_{t},{p}", vtype=gp.GRB.CONTINUOUS)
Modèle 2.1 : Intégration des Coûts Fixes (Cbase)
Le modèle évolue pour intégrer les coûts fixes de fonctionnement :
\[\min \sum_{t=1}^{3} \sum_{p=1}^{5} \left[(P_{tp} - N_{tp} \cdot P^{\text{min}}_t) \cdot C_{t}^{mwh} + N_{tp} \cdot C_{t}^{base}\right] \cdot T_{p}\]Cette formulation sépare :
- Le coût variable proportionnel à la production au-delà du minimum
- Le coût fixe proportionnel au nombre de centrales actives
Résultat : Coût optimal = 979 k€ (+12.7% par rapport au modèle 1.5)
Modèle 2.2 : Ajout des Coûts de Démarrage (Cstart)
L’introduction de la variable \(X_{tp}\) (nombre de centrales démarrées) nécessite une nouvelle contrainte de cohérence :
\[N_{tp} \leq X_{tp} + N_{t(p-1)} \quad (C4)\]La fonction objectif devient :
\[\min \sum_{t=1}^{3} \sum_{p=1}^{5} \left[(P_{tp} - N_{tp} \cdot P^{\text{min}}_t) \cdot C_{t}^{mwh} + N_{tp} \cdot C_{t}^{base}\right] \cdot T_{p} + \sum_{t=1}^{3} \sum_{p=1}^{5} X_{tp} \cdot C_{t}^{start}\]Résultat : Coût optimal = 1014 k€ (+3.6% par rapport au modèle 2.1)
En Gurobi, les variables sont déclarées comme suit :
# Contrainte de démarrage
for p in HOURS.keys():
for t in USINES.keys():
left_side = X[f"{t},{p}"] + Ntp(t, p, N, cyclique=False)
right_side = N[f"{t},{p}"]
m.addConstr(left_side >= right_side, name=f"allume_{t}_{p}")
Modèle 3.1 : Contrainte de Réserve (15%)
Pour garantir la sécurité du réseau, on impose une réserve de puissance de 15% :
\[\sum_{t=1}^{3} (N_{tp} \cdot P^{\text{max}}_t - P_{tp}) \geq 0.15 \cdot D_{p} \quad (C5)\]Résultat : Coût optimal = 1015 k€ (+0.1% par rapport au modèle 2.2)
L’impact marginal s’explique par le fait que la solution précédente disposait déjà d’une marge naturelle.
Modèle 4.1 : Planification Cyclique
Dans ce modèle, on impose que l’état final (période 5) soit cohérent avec l’état initial (période 1) :
\[N_{tp} \leq X_{tp} + N_{t(p-1)} \quad \text{avec } N_{t,0} = N_{t,5}\]Cette approche reflète une planification sur 24h répétable quotidiennement.
Résultat : Coût optimal = 989 k€ (-2.6% par rapport au modèle 3.1)
L’optimisation cyclique permet de réduire les coûts de démarrage en optimisant la transition jour-nuit.
En Gurobi, cette contrainte peut être implémentée comme suit :
def Ntp(t, p, N, cyclique=True):
if p == 1:
if cyclique:
Np_1 = N[f"{t},5"] # Boucle cyclique
else:
Np_1 = 0 # Démarrage à froid
else:
Np_1 = N[f"{t},{p-1}"]
return Np_1
Intégration des Centrales Hydroélectriques
Modèle 5.1 : Ajout de Deux Centrales Hydro
Deux centrales hydroélectriques sont disponibles :
| Type | Puissance (MW) | Coût base (€/h) | Coût démarrage (€) | Abaissement (m/h) |
|---|---|---|---|---|
| H1 | 900 | 90 | 1500 | 0.31 |
| H2 | 1400 | 150 | 1200 | 0.47 |
Ces centrales ne consomment pas de combustible (Cmwh = 0) mais abaissent le niveau du réservoir.
Résultat : Coût optimal = 890 k€ (-10.0% par rapport au modèle 4.1)
Les centrales hydro, plus économiques, remplacent partiellement les centrales thermiques coûteuses.
Modèle 5.2 : Contrainte de Pompage
Le niveau du réservoir doit être constant sur 24h. Le pompage nécessite 3000 MWh par mètre d’élévation :
\[\sum_{p=1}^{5}\sum_{t=1}^{5} P_{tp} \cdot T_p \geq \sum_{p=1}^{5} D_{p} \cdot T_p + \sum_{p=1}^{5}\sum_{t=4}^{5} 3000 \cdot M_t \cdot N_{tp} \cdot T_p \quad (C6)\]où \(M_t\) est la vitesse d’abaissement du réservoir (en m/h).
Résultat : Coût optimal = 987 k€ (+10.9% par rapport au modèle 5.1)
Le coût du pompage compense partiellement l’économie réalisée par l’hydro.
En Gurobi, cette contrainte peut être implémentée comme suit :
# Contrainte de pompage
left_side = sum(P[f"{t},{p}"]*HOURS[p] for t in USINES.keys() for p in HOURS.keys())
right_side = (sum(DEMAND[p]*HOURS[p] for p in HOURS.keys()) +
sum(3000*N[f"{t},{p}"]*USINES[t]['v_abai']*HOURS[p]
for t in USINES.keys() for p in HOURS.keys()))
m.addConstr(left_side >= right_side, name="Pompage")
Analyse Théorique : Sensibilité et Dualité
Merit Order et Solution Optimale
La méthode du Merit Order consiste à activer les centrales dans l’ordre croissant de leur coût marginal :
Pour \(p=1\) (demande = 15 GW) :
\[\mathbf{X_1} = \begin{bmatrix} P_{A1} = 0 \\ P_{B1} = 15 \\ P_{C1} = 0 \end{bmatrix}\]Cette solution est vérifiée numériquement par Gurobi.
Dualité et Coûts Marginaux
Le modèle primal en forme standard :
\[\begin{aligned} \min & \quad 9 P_{A1} + 8.28 P_{B1} + 16.5 P_{C1} \\ \text{s.t.} & \quad P_{A1} + \delta_A = 24 \\ & \quad P_{B1} + \delta_B = 17.5 \\ & \quad P_{C1} + \delta_C = 20 \\ & \quad P_{A1} + P_{B1} + P_{C1} = 15 \end{aligned}\]Son dual s’écrit :
\[\begin{aligned} \max & \quad 24 u_{A} + 17.5 u_{B} + 20 u_{C} + 15 v \\ \text{s.t.} & \quad u_A + v \leq 9 \\ & \quad u_B + v \leq 8.28 \\ & \quad u_C + v \leq 16.5 \end{aligned}\]Analyse de Sensibilité
Impact d’une augmentation de la demande (+1 GW)
La nouvelle solution est :
\[x_{\beta}^{'} = x_{\beta} + \Delta \cdot A_{\beta}^{-1} e_4\]Résultat : Coût additionnel = 8.28 k€
delta = 1
x_beta_ligne = x_beta + np.dot(delta*np.linalg.inv(A_beta1), e(4, 4))
cout_maintenant = np.dot(C_beta1, x_beta_ligne)[0][0]
print(f"Surcoût: {cout_maintenant - cout_avant:.2f} k€")
Impact d’une centrale supplémentaire de type B
- Période 1 : aucun impact (capacité excédentaire)
- Période 2 : économie de 630 €
Résultats Comparatifs
| Modèle | Variables | Contraintes | Nœuds explorés | Temps (s) | Coût (k€) |
|---|---|---|---|---|---|
| 1.5 | 30 | 50 | 1 | 0.02 | 869 |
| 2.1 | 30 | 50 | 1 | 0.02 | 979 |
| 2.2 | 45 | 65 | 1 | 0.06 | 1014 |
| 3.1 | 45 | 70 | 1 | 0.01 | 1015 |
| 4.1 | 45 | 70 | 1 | 0.01 | 989 |
| 5.1 | 75 | 110 | 1 | 0.04 | 890 |
| 5.2 | 75 | 111 | 16 | 0.15 | 987 |
Interprétation Économique
- Coûts fixes dominants : Le passage du modèle 1.5 au 2.1 (+12.7%) montre l’importance des coûts de base.
- Coûts de démarrage modérés : Le modèle 2.2 n’ajoute que 3.6% de coût supplémentaire.
- Planification cyclique avantageuse : Le modèle 4.1 réduit les coûts de 2.6% en optimisant les transitions.
- Hydroélectricité compétitive : Le modèle 5.1 réduit les coûts de 10% grâce à des coûts variables nuls.
- Pompage coûteux : Le modèle 5.2 montre que la contrainte de neutralité du réservoir ajoute 10.9% aux coûts.
Propriétés Mathématiques
Unimodularité Totale
La matrice des contraintes du modèle 1.5 possède la propriété d’unimodularité totale :
- Tous les coefficients sont dans \(\{-1, 0, 1\}\)
- Chaque colonne contient au plus deux éléments non nuls
Conséquence : La relaxation continue donne une solution entière, ce qui explique la résolution rapide (1 nœud exploré).
Écarts Complémentaires
À l’optimum, les conditions de Karush-Kuhn-Tucker (KKT) s’écrivent :
\[(c_j - u^T A_j) x_j = 0 \quad \text{et} \quad u_i (a_i x - b_i) = 0\]Ces conditions garantissent qu’à l’optimum :
- Les variables non basiques ont un coût réduit positif
- Les contraintes saturées ont un multiplicateur de Lagrange non nul
Conclusion
Ce projet démontre l’efficacité de la programmation linéaire en nombres entiers pour résoudre des problèmes complexes de planification énergétique. L’analyse progressive des modèles révèle :
- L’importance de modéliser tous les coûts (fixes, variables, démarrage)
- Les bénéfices de la planification cyclique pour réduire les coûts de transition
- L’intérêt économique de l’hydroélectricité, malgré les contraintes de pompage
- La stabilité numérique grâce à l’unimodularité de la structure du problème
Les résultats obtenus avec Gurobi montrent des temps de résolution très courts (< 0.2 s). Cependant, en réalité, ce type de problème contient des milliers de variables, ce qui augmente considérablement le temps de résolution.
Code Source
Jouer avec l’application Streamlit codé en Python utilisant Gurobi et Pandas sur:
👉 App