Simulation probabiliste de rupture de barrage et de rupture de digue : modélisation géotechnique dynamique de l'érosion et routage hydrodynamique des vagues

Les barrages, les réservoirs et les systèmes de digues de protection contre les inondations représentent certains des actifs d'infrastructure civile les plus critiques de la société moderne. Ils garantissent l’approvisionnement en eau pour l’irrigation agricole, génèrent de l’énergie hydroélectrique renouvelable et protègent des inondations des millions d’habitants vivant dans les plaines inondables historiques. Cependant, la rupture catastrophique d’un barrage en terre ou d’une digue libère de vastes volumes d’énergie potentielle en quelques minutes, générant des vagues de crues violentes et catastrophiques capables d’anéantir des communautés entières en aval, des ponts et des couloirs industriels.

Historiquement, les évaluations des risques de rupture de barrage reposaient sur des formules de rupture statiques et empiriques (telles que les équations classiques de Froehlich, MacDonald-Langridge-Monopolis ou Von Thun-Gillette) associées à un tracé hydraulique unidimensionnel (1D). Dans ces frameworks existants :

• L'ingénieur a supposé une géométrie de brèche fixe (forme de brèche trapézoïdale, largeur du fond et pentes latérales) et un temps de développement de brèche fixe (par exemple 1,5 heure) sur la base de diagrammes de dispersion de régression historiques.
• L'écoulement de la brèche a été acheminé vers l'aval à travers des coupes transversales 1D qui ignoraient complètement l'étendue de la vallée de la plaine inondable en 2D, la densité apparente du transport des sédiments et les réflexions des vagues sur les parois du canyon.

Cette approche empirique statique présente de dangereuses limites techniques. Les remblais en terre ne s'effondrent pas selon des programmes géométriques arbitraires ; ils échouent à cause de problèmes complexes et non linéairesmécanique de l'érosion géotechnique: érosion progressive par débordement vers l'arrière, érosion par cisaillement des conduits de tuyauterie internes et effondrements de cercles de glissement de stabilité de pente en rotation. De plus, les hypothèses statiques ignorent l’incertitude des paramètres géotechniques : cohésion du compactage du sol, distribution granulométrique (D_{50}), contrainte de cisaillement critique (tau_c) et les fluctuations du niveau d’eau des réservoirs.

Fin 2026, les normes de sécurité réglementaires (mandatées par la Federal Energy Regulatory Commission, FERC, et le US Army Corps of Engineers, USACE) exigentSimulations probabilistes de rupture de barrage et de rupture de digue basées sur des données physiques. En couplant directement la physique de l'érosion des brèches géotechniques transitoires avec des solveurs d'équations en eaux peu profondes bidimensionnelles (2D) accélérées par GPU et la quantification de l'incertitude d'ensemble de Monte Carlo, les ingénieurs en sécurité des barrages peuvent cartographier les risques d'inondation en aval avec une précision mathématique sans compromis.

Ce traité technique complet explore la physique géotechnique de l'érosion, l'hydrodynamique des vagues chargées de sédiments non newtoniens, l'échantillonnage d'ensemble de Monte Carlo et les protocoles de cartographie des risques d'évacuation d'urgence nécessaires pour exécuter des évaluations de rupture de barrage de pointe.

1. Mécanique géotechnique des brèches : débordement ou érosion interne des canalisations

Les barrages en remblai en terre échouent principalement à cause de l'un des deux mécanismes physiques suivants : le franchissement hydraulique ou la tuyauterie géotechnique interne.

L’équation d’érosion par contrainte de cisaillement excessive

Dans les modèles de brèches physiques (tels que USDA-ARS WinDAM ou DL-Breach), le détachement du sol le long de la limite de brèche est régi par leFormulation de contrainte de cisaillement excessive:

ε_r = k_d (τ_e - τ_c)^a

Où :

• epsilon_rreprésente le taux d'érosion par décollement du sol (m/s ou cm/h).
• k_dest le coefficient d'érodabilité du sol (cm^3/N·s), déterminé via des tests d'érosion par jet in situ (JET) ou des tests de dispersion par sténopé immergés.
• tau_eest la contrainte de cisaillement hydrodynamique effective exercée par l'écoulement de l'eau sur la limite de brèche :

τ_e = γ_w R_h S_f

Oùgamma_west le poids unitaire de l'eau,R_hest le rayon hydraulique du canal de brèche, etS_fest la pente énergétique.

• tau_cest la contrainte de cisaillement critique du sol du remblai compacté (Pa), en dessous de laquelle aucun détachement ne se produit.
• $a$ est un exposant empirique, généralement pris comme $a = 1,0$.

Lorsque l'eau qui coule accélère à travers une brèche en développement,tau_epeut dépasser 500 Pa, dépassant largement la résistance critique au cisaillement du sol (tau_c ≈ 5 - 20 Pa) – entraînant des taux d’érosion volumétrique exponentiels.

2. Calcul dynamique de l'hydrogramme de débit sortant : système hydraulique de déversoir et d'orifice à large crête

À mesure que l'entaille de brèche s'élargit et s'approfondit, l'écoulement instantané s'écouleQ_{brèche}(t)bascule transitoirement entreécoulement à orifice immergé(au début de la tuyauterie interne) etécoulement à déversoir à large crête et surface libre(une fois la brèche du toit effondrée).

flowchart TD
    Pool["Reservoir Pool Level: H_res(t)"] --> Geom["Instantaneous Breach Notch Geometry<br/>• B_w(t), Z_s, H_b(t)"]
    Geom -->|Internal Piping Phase| Pipe["Pressurized Orifice Flow<br/>• Q = C_d · A_pipe(t) · √(2g · ΔH)"]
    Geom -->|Open Breach Crest Phase| Weir["Broad-Crested Weir Flow<br/>• Q = C_w · B_w(t) · (H_res - Z_b)^1.5 + C_s · Z_s · (H_res - Z_b)^2.5"]
    Pipe --> Wave["Downstream 2D Hydrodynamic Flood Wave Routing"]
    Weir --> Wave

La formulation combinée du déversoir trapézoïdal à large crête

Une fois qu'une tranchée ouverte est établie, le débit sortant est calculé via l'équation combinée du déversoir rectangulaire et triangulaire à large crête :

Q_{brèche}(t) = C_{w1} B_w(t) ≤ft( H_{res}(t) - Z_b(t) )^{3/2} + C_{w2} Z_s ≤ft( H_{res}(t) - Z_b(t) )^{5/2}

Où :

• B_w(t)est la largeur inférieure instantanée de l'entaille de brèche en expansion (m).
• Z_b(t)est l'élévation instantanée du fond de l'entaille de brèche, s'érodant vers le bas jusqu'à ce qu'elle atteigne le substrat rocheux ou les strates de fondation non érodables.
• Z_sest la pente latérale de la brèche trapézoïdale (Z_s = (Deltax) / (Deltaz), généralement 0,5$ à 1,0$ pour les argiles cohésives, et jusqu'à 1,5$ pour les sables non cohésifs).
• C_{w1}, C_{w2}sont des coefficients de débit de déversoir calibrés pour les effets des eaux résiduaires de submersion.

Simultanément, l'élévation du bassin du réservoirH_{res}(t)chute dynamiquement selon l’équation de conservation stade-stockage du réservoir :

frac{dH_{res}}{dt} = frac{Q_{entrée}(t) - Q_{brèche}(t) - Q_{évacuateur}(t)}{A_{res}(H_{res})}

OùA_{res}(H_{res})est la superficie du bassin du réservoir en altitudeH_{res}.

3. Acheminement des vagues hydrodynamiques : vagues de déferlement de barrage non newtoniennes 2D

Une vague de rupture de barrage n’est pas une eau propre et claire ; c'est une boue violente et de haute densité composée d'eau, de gravier de remblai pulvérisé, de sédiments en suspension et d'arbres déracinés.

Dans les canyons de montagne escarpés ou les vallées fluviales étroites, la concentration volumique de sédiments (C_v) peut dépasser 20\% - 45\%$, transformant la vague de crue en unesuspension hyperconcentrée non newtonienne.

Hydrodynamique plastique de Bingham non newtonienne

>

Formulation de la contrainte de cisaillement apparente :

τ = τ_y + µ_B · (∂u / ∂z)

>

Paramètres :

• τ_y: Limite d'élasticité de la boue de sédiments (Pa) - doit être surmontée pour initier le mouvement du fluide

• µ_B: Viscosité plastique dynamique Bingham (Pa·s)

• Densité apparente: ρ_boue = ρ_w · (1 - C_v) + ρ_s · C_v

Le solveur hydrodynamique de capture de chocs 2D

Pour acheminer cette onde violente et à grande vitesse vers l’aval sans instabilité numérique ni diffusion non physique :

• Les modèles utilisent les équations 2D pour eaux peu profondes avecTermes sources de la résistance à l'écoulement du plastique Bingham ou de Herschel-Bulkley (tau_{bx}, tau_{par}).
• Les solveurs appliquent des schémas de capture de chocs TVD (Diminution de la variation totale) de second ordre pour résoudre les problèmes abrupts.forages à vagues de rupture de barrageavançant à des vitesses dépassant12 à 20 m/sen bas de la vallée en aval.
• Les équations de quantité de mouvement résolvent les forces centrifuges latérales lorsque la vague contourne des virages serrés de la rivière, prédisant avec précision la montée des vagues de super-élévation contre les parois du canyon.

4. Cadre d'ensemble probabiliste de Monte Carlo : quantification de l'incertitude géotechnique

Dans l’ingénierie médico-légale de la sécurité des barrages, les propriétés géotechniques exactes ne sont jamais connues avec une certitude à 100 %. La qualité du compactage, la teneur en humidité, la cohésion et l'homogénéité interne d'un remblai en terre varient sur sa longueur linéaire.

Plutôt que d'exécuter une seule simulation déterministe du « pire des cas », qui pourrait être irréaliste ou dangereusement optimiste, les ingénieurs modernes exécutentEnsembles probabilistes de Monte Carlo.

flowchart TD
  Params["Input Parameter Probability Distributions\n• Reservoir Head: Normal(μ = Full Pool, σ = 1.2m)\n• Soil Erodibility k_d: Log-Normal(μ = 10 cm³/N-s, σ = 0.4)\n• Critical Shear Stress τ_c: Weibull(λ = 8.5 Pa, k = 2.1)\n• Cohesion c': Normal(μ = 25 kPa, σ = 5 kPa)\n• Breach Width B_w: Beta distribution (Froehlich)"]
  LHS["Latin Hypercube Sampling (LHS)\nGenerates 500 Stratified Scenarios"]
  GPU["Automated GPU Solver Cluster\n500 Coupled 2D Simulations (4.2 GPU Hours)"]
  Hazard["Probabilistic Hazard Outputs\nP10, P50, P90 Inundation Maps & Arrival Time Bands"]

  Params --> LHS --> GPU --> Hazard

Mesures probabilistes des inondations

À partir de l’ensemble de 500 analyses, le système extrait des probabilités de dépassement spatialement distribuées :

• Carte des inondations P90: La limite d'inondation conservatrice avec un niveau de confiance de 90\%$, utilisée pour la conception de la protection des infrastructures civiles (centrales nucléaires, complexes hospitaliers, stockage de produits chimiques).
• Carte des inondations P50: Le scénario médian attendu, utilisé pour le zonage municipal, l’évaluation des taux d’assurance et l’itinéraire d’évacuation primaire.
• Incertitude de l’heure d’arrivée des vagues: Quantifie l'arrivée la plus précoce possible du front d'onde à grande vitesse dans des centres de population spécifiques en aval (T_{arrivée, P10}), offrant aux responsables des situations d'urgence des marges de sécurité défendables pour les fenêtres d'autorisation d'évacuation.

5. Matrice de comparaison complète des méthodologies de rupture de barrage

Le tableau suivant compare les méthodes empiriques existantes aux cadres de simulation de rupture de barrage 2D probabilistes modernes basés sur la physique :

6. Étude de cas d'ingénierie réelle : audit de sécurité des barrages en terre à haut risque

Pour valider ce cadre probabiliste sur un actif d'infrastructure civile critique, FloodSim a exécuté un audit de rupture de barrage dans le cadre d'un plan d'action d'urgence (PAE) pour un barrage en remblai en terre de 85 pieds de haut (26 mètres) retenant un réservoir d'eau potable et de contrôle des crues de 32 000 acres (39,5 millions de mètres cubes) situé à 6 milles en amont d'une vallée municipale urbaine de 95 000 habitants.

Découvertes critiques pour la sécurité des personnes :

1. La dangereuse illusion du « décalage »: L'ancien modèle 1D prévoyait que la vague de crue initiale nécessiterait 65 minutes pour parcourir 6 miles en aval jusqu'au centre-ville. En revanche, le solveur GPU 2D haute résolution, prenant en compte le rétrécissement de la vallée et la canalisation de l'écoulement supercritique le long du corridor fluvial principal, a révélé que le forage hydraulique à grande vitesse frapperait le noyau municipal dansseulement 38 minutes. Si les responsables des urgences s’étaient appuyés sur l’ancien modèle, les ordres d’évacuation seraient arrivés près de 30 minutes trop tard.
2. Super-affouillement du pont inter-États: La simulation 2D a cartographié des vitesses d'écoulement localisées dépassant6,2 m/ssur les piles d'un pont routier interétatique, calculant les contraintes de cisaillement inférieures qui affaibliraient les fondations des piles dans les 12 minutes suivant l'arrivée des vagues.
3. Réétalonnage du plan d’action d’urgence: Armée de courbes d'arrivée probabilistes exactes, l'agence municipale de gestion des urgences a déplacé quatre sirènes d'évacuation, établi des fermetures automatisées des portes d'autoroute et établi un protocole certifié d'autorisation d'évacuation de 22 minutes, garantissant ainsi la sécurité totale de la vie publique.

7. Classification des risques d'évacuation en aval : courbes australiennes de vulnérabilité aux risques

L'évaluation du risque pour la sécurité de la vie publique nécessite de traduire la profondeur brute des inondations ($h$) et la vitesse ($v$) en catégories de danger pour la stabilité humaine. Une faible profondeur d'eau de 0,4 mètre coulant à 3,0 m/s est bien plus mortelle qu'une flaque d'eau stagnante de 1,0 mètre, car les forces de traînée hydrodynamiques balayent les piétons adultes et les automobiles hors des routes.

Les plans d'action d'urgence modernes utilisentCourbes australiennes de vulnérabilité aux risques liés à la gestion des urgences (H1 à H6):

<!--ec:block {"type": "table", "hasHeaderRow": true, "rows": [{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"CATÉGORIES DE CLASSIFICATION DES RISQUES (H1 - H6)"}],"isHeader":true}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"Intensité hydrodynamique totale : D \\times V = h \\times \\sqrt{u^2 + v^2} (m²/s)"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"tapez": "span", "text": "}]}]}, {"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"[H1 : Sûr pour tous] : D \\times V \\le 0,4 \m²/s, h \\le 0,3m"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"Capable de marcher en toute sécurité ; les enfants et les personnes âgées en sécurité."}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"tapez": "span", "text": "}]}]}, {"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type": "span", "text": [H2 : dangereux pour les petits] : D \\times V \\le 0,6 \m²/s, h \\le 0,5m"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"Dangereux pour les jeunes enfants et les personnes âgées fragiles."}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"tapez": "span", "text": "}]}]}, {"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"[H3 : Dangereux pour les véhicules] : D \\times V \\le 0,8 \m²/s, h \\le 0,6m"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type": "span", "text": "Les berlines et les véhicules légers perdent de la traction et de la flottaison."}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"tapez": "span", "text": "}]}]}, {"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"[H4 : Dangereux pour tous] : D \\times V \\le 1,0 \m²/s, h \\le 1,0m"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type": "span", "text": "Adultes emportés par les pieds ; tous les véhicules flottent."}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"tapez": "span", "text": "}]}]}, {"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"[H5 : Dommages structurels] : D \\times V \\le 4,0 \m²/s, h \\le 4,0m"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"Les bâtiments en maçonnerie non renforcée subissent des dommages."}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"tapez": "span", "text": "}]}]}, {"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"type":"span","text":"[H6 : Effondrement catastrophique] : D \\times V > 4,0 \m²/s"}]}]},{"type": "tableRow", "cells": [{"type": "tableCell", "content": [{"_type": "span", "text": "Défaillance structurelle complète des bâtiments résidentiels"}]}]}]} -->

En exportant des couches SIG dynamiques deH_1à traversH_6Grâce aux classifications par incréments de 5 minutes en aval du barrage rompu, les responsables des urgences peuvent identifier quelles routes d'évacuation restent praticables en toute sécurité et à quel moment précis elles se transforment en pièges mortels.

8. Protocole de conformité réglementaire pour les ingénieurs en sécurité des barrages

Les ingénieurs en sécurité des barrages, les consultants géotechniques et les directeurs des services publics municipaux qui mèneront des évaluations des risques de rupture de barrage en 2026 doivent respecter cette liste de contrôle de conformité en cinq étapes :

1. Effectuer des tests d'érosion par jet in situ (JET): Ne jamais estimer les coefficients d'érodibilité (k_d) à partir des seules cartes de classification des sols des manuels. Effectuez des tests standardisés ASTM D5852/D6665 pour calibrer les véritables paramètres de cisaillement excessif.
2. Modèle d’élargissement de brèche géotechnique dynamique: Éliminez les hypothèses trapézoïdales statiques. Utilisez des algorithmes d'érosion couplés basés sur la physique qui modélisent la migration vers l'arrière, l'effondrement rotationnel des pentes latérales et le rabattement des réservoirs.
3. Déployez des solveurs de capture de chocs 2D accélérés par GPU: Acheminez les ondes de sortie de brèche à travers des maillages 2D non structurés dérivés d'un LiDAR haute résolution, garantissant ainsi que la topographie critique de la vallée en aval, les remblais routiers et les ponceaux sont hydroconditionnés.
4. Quantifier l’incertitude des paramètres géotechniques: Exécutez un minimum de 250 exécutions d'ensemble de Monte Carlo sur des distributions de paramètres stratifiées de l'Hypercube latin pour générer des contours de danger défendables P10, P50 et P90.
5. Classer les routes d'évacuation en aval via des courbes profondeur-vitesse: Traduire les sorties hydrodynamiques en standardH_1 - H_6classifications des dangers pour les personnes et les véhicules afin de fournir des renseignements exploitables et vitaux aux directeurs de gestion des urgences.

Grâce à l'intégration disciplinée de la mécanique géotechnique des sols, de l'hydrodynamique non newtonienne des eaux peu profondes accélérée par GPU et de la quantification probabiliste de l'incertitude de Monte Carlo, la simulation moderne de rupture de barrage permet aux ingénieurs civils de protéger la société contre les scénarios de défaillance les plus catastrophiques imaginables.