Soutenance mémoire
Le code de calcul
La simulation numérique complète l’expérience physique lorsque celle-ci est difficile à mettre en œuvre et/ou trop coûteuse, voire impossible à effectuer. Elle repose sur la modélisation du système physique d’intérêt sous la forme d’équations mathématiques traduisant la réalité physique. Par exemple, les codes de calcul de mécanique des fluides s’appuient sur les équations de Navier-Stokes, tandis que les codes de calcul simulant les phénomènes électromagnétiques s’appuient sur les équations de Maxwell. La solution de ces équations est implémentée dans un langage de programmation informatique menant à la création d’un code de calcul 2, représenté dans ce mémoire par la fonction paramétrée yθ : (xc,xe) ∈ X c ×X e ∈ Rd ×Rq −→ yθ(xc,xe) ∈ R m, (2.2) où θ ∈ T ⊂ Rp est un vecteur de paramètres pouvant avoir une interprétation physique (voir chapitre 5). Notons que l’exécution du code nécessite d’attribuer une valeur à la variable environnementale xe en entrée, bien qu’elle soit incertaine dans la réalité. Nous verrons dans la section 2.1.4.1 comment propager son incertitude à la réponse du code. Nous nous intéressons dans ce mémoire à la validation des codes de calcul déterministes. Autrement dit, pour deux réalisations identiques des variables d’entrée et du paramètre (xc1,xe1,θ1) = (xc2,xe2,θ2), alors yθ1(xc1,xe1) = yθ2(xc2,xe2). De plus, nous utilisons le code comme une fonction boîte noire, ce qui signifie que les équations mathématiques reliant yθ(xc,xe) à xc et xe n’ont pas besoin d’être connues. Enfin, nos contributions sont destinées à la validation des codes de calcul pour lesquels le paramètre θ est incertain. Dans la section 2.1.3, nous détaillons les différentes sources d’incertitudes affectant à la fois le système physique et le code de calcul.
L’incertitude affectant la réponse du code de calcul
Pour propager l’incertitude aléatoire du vecteur Xe à la réponse du code de calcul, l’échantillonnage Monte-Carlo est la stratégie la plus naturelle à condition que l’on puisse disposer d’un nombre conséquent d’appels au code. Lorsque l’on veut propager ensemble, sans les confondre, l’incertitude épistémique sur θ et l’incertitude du vecteur Xe à la réponse du code, Roy et Oberkampf (2011) proposent de représenter l’incertitude sur θ par un ensemble de valeurs candidates et militent donc pour une représentation par intervalle de l’incertitude épistémique. Par ce biais, l’incertitude sur yθ(xc,Xe) est quantifiée par une boîte de probabilité, ou « p-box » dans la littérature anglo-saxonne (Ferson et al., 2003). La construction d’une boîte de probabilité repose sur la propagation de l’incertitude sur Xe à valeur de θ fixée, ce processus étant répété pour un grand nombre de valeurs plausibles pour θ. Ainsi, l’incertitude globale affectant yθ(xc,Xe) peut être représentée par une famille de fonctions de répartition. À valeur de x c fixée, la figure 2.1 illustre la boîte de probabilité obtenue en supposant que yθ(Xe) := yθ(xc,Xe) suit une loi gaussienne de moyenne θ ayant pour plage d’incertitude l’intervalle [3, 5]. La boîte de probabilité nous permet de distinguer la part d’incertitude provenant de l’aléa Xe , caractérisée par sa forme, de la part d’incertitude provenant du manque de connaissance sur la valeur de θ, caractérisée par sa largeur. On pourrait penser, à tort, qu’un résultat similaire serait obtenu en spécifiant d’un point de vue bayésien une loi uniforme a priori sur θ, puis en propageant à travers le code l’incertitude affectant Xe conditionnellement à θ, le tout nécessitant une double boucle Monte-Carlo. L’article de Ferson et Ginzburg (1996) illustre les différences entre ces deux approches
Impact du plan d’expériences numériques
Pour illustrer l’importance du plan d’expériences sur l’erreur commise lorsque l’on échantillonne la loi a posteriori approchée (2.32) à la place de la loi a posteriori cible (2.21), supposons que la fonction analytique yτ : x ∈ [0, 1] −→ (6x −2)2 ×sin (τx −4) (2.33) joue de nouveau le rôle du code de calcul. Soit Xf = (0.1, 0.3, 0.8)T. Nous simulons des mesures physiques zf suivant l’équation (2.20) avec θ = 12 et λ = 0.3. Dans la procédure d’estimation décrite dans la section précédente, nous supposons toujours que λ est connu (fixé ici à 0.3) si bien que seul le paramètre θ est à estimer. Dans un premier temps, nous calculons (2.21) à l’aide d’un algorithme de Métropolis-Hastings en choisissant la loi a priori uniforme π(θ) ∝ 1[5,15]. La figure 2.7 (à droite) montre une densité de la loi a posteriori cible à deux modes, l’un autour de la vraie valeur θ ≈ 12 à partir de laquelle ont été simulées les mesures zf et l’autre autour de θ ≈ 9.8. Maintenant, on échantillonne la loi a posteriori approchée (2.32) pour plusieurs tailles de plan d’expériences LHD maximin . La figure 2.8 illustre les densités de probabilités obtenues. On constate que pour un nombre restreint de simulations M = 45 et M = 60, elles sont très éloignées de la loi cible. Lorsque M augmente, dans notre exemple lorsque M = 75 et M = 90, alors les régions autour de θ = 12 et θ = 9.8 sont clairement identifiées. Ces résultats illustrent que plus la taille du plan d’expériences est élevée, plus la distribution a posteriori (2.32) ressemblera à (2.21).
Comparaison des méthodes et interprétation des résultats
Hypothèses de modélisation Deux caractéristiques rendent les approches de Higdon et al. (2004) et Bayarri et al. (2007b) plus cohérentes que celle de Kennedy et O’Hagan (2001a) pour le calage et la validation d’un code de calcul :
1. le paramètre ρ est fixé à 1,
2. l’erreur de code est modélisé par un processus gaussien de moyenne constante βb et les auteurs recommandent de poser βb = 0.
Ces deux hypothèses restrictives impliquent, d’une part que les réponses du code de calcul fournissent des prédictions à l’échelle de r (x), et d’autre part que le code de calcul peut prédire de façon non biaisée la moyenne du système physique sur X . On aimerait que ces deux hypothèses soient vérifiées en amont d’une procédure de calage et de validation. En effet, est-il bien raisonnable de calculer des prédictions d’un système physique à l’aide d’un code de calcul dont les réponses sont à une échelle différente et/ou translatées avec la réalité ? Nous croyons que si de tels problèmes subsistent, ils sont du ressort de la conception du code de calcul. Ensuite, si βb = 0 et qu’il existe une valeur τ = θ? telle que bθ?(x) = 0 pour tout x, alors le maximum de vraisemblance ˆt du modèle (2.65) converge vers θ? sous l’hypothèse que le code n’est pas coûteux (Loeppky et al., 2006). Autrement dit, l’estimation de θ dans le modèle (2.65) lorsque n est suffisamment grand coïncide avec l’estimation de θ dans (2.20) sous l’hypothèse que le code n’est pas remplacé par un émulateur processus gaussien. Techniques d’estimation utilisées La méthode de Higdon et al. (2004) propose d’estimer conjointement θ et les hyperparamètres des processus gaussiens mis en jeux, tandis que les techniques proposées dans Kennedy et O’Hagan (2001a) et Bayarri et al. (2007b) ont en commun l’estimation des hyperparamètres du code à partir de la vraisemblance partielle des simulations y(DM), appelée approche modulaire dans Liu et al. (2009). Cette démarche est assez intuitive dans le sens où seules les simulations y(DM) influent sur l’estimation des hyperparamètres (Ψy ,σ2y ) du processus gaussien modélisant le code. En ce qui concerne l’estimation des hyperparamètres (Ψb,σ2b,λ2 ), Kennedy et O’Hagan (2001a) maximisent la vraisemblance des mesures physiques zf conditionnellement aux simulations après avoir intégré par rapport à la loi a priori de θ (2.53). Bayarri et al. (2007b) optent pour une approche plus hybride dans laquelle Ψb est estimé en maximisant la vraisemblance des erreurs virtuelles (2.74) calculée à partir d’une valeur a priori θ?, puis les lois a posteriori de σ2b et λ2 sont échantillonnées avec celles de θ et des données latentes yθ(X f ) et b(X f) dans un algorithme de Gibbs basé sur les conditionnelles complètes (2.66) à (2.69). Cela permet la construction naturelle d’un prédicteur de l’erreur de code (2.71), ce qui ne figure pas dans la méthode de KOH. Notons enfin que chaque nouvelle évaluation de (2.48) demande d’inverser une matrice de variance-covariance de taille (n×M)×(n×M), ce qui peut être préjudiciable lorsque la taille du plan d’expériences M est importante. Combiné avec un paramètre θ de grande dimension, l’échantillonnage de (2.48) peut alors devenir très coûteux.
|
Table des matières
1 Introduction
2 État de l’art
2.1 La simulation numérique d’un système physique
2.1.1 Le système physique
2.1.2 Le code de calcul
2.1.3 Les sources d’incertitudes
2.1.4 Propagation des sources d’incertitudes
2.1.5 Systèmes physiques étudiés
2.2 Vérification et Validation d’un code de calcul
2.2.1 Vérification d’un code de calcul
2.2.2 Validation d’un code de calcul
2.3 Calage des codes de calcul
2.3.1 Modélisation statistique
2.3.2 Impact du plan d’expériences numériques
2.3.3 Incertitude de prédiction
2.4 Calage des codes de calcul avec erreur de code
2.4.1 Méthode de Kennedy et O’Hagan
2.4.2 Méthode de Bayarri
2.4.3 Méthode de Higdon
2.4.4 Comparaison des méthodes et interprétation des résultats
2.4.5 Modélisation alternative de l’erreur de code
2.4.6 Cas d’un code linéaire
2.5 Méthodes de calage et de validation pour la grande dimension
2.5.1 L’analyse bayésienne linéaire
2.5.2 Les méthodes d’assimilation de données
2.6 Conclusion
3 Bayesian model selection for the validation of computer codes
3.1 Introduction
3.2 Statistical modeling for code validation
3.3 Bayesian model selection
3.4 Application to code validation
3.5 Simulation study
3.6 Validation of an industrial code for power plant production control
3.7 Conclusion
4 Adaptive numerical designs for the calibration of computer codes
4.1 Introduction
4.2 Calibration framework
4.3 Adaptive designs for calibration
4.4 Simulation study
4.5 Conclusion
5 Validation of a computer code for the energy consumption of a building, with application to optimal electric bill pricing
5.1 Introduction
5.2 Overview of the study
5.3 Calibration of the dynamic thermal code
5.4 Validation of the dynamic thermal code
5.5 Optimal forecasts of consumption
5.6 Conclusion
6 Conclusions et perspectives
A Le processus aléatoire gaussien
A.1 Définition d’un processus aléatoire
A.2 Vecteur et processus gaussien
A.3 Modélisation par processus gaussien
A.4 Les plans d’expériences
A.4.1 Plans d’expériences remplissant l’espace des entrées
A.4.2 Plans d’expériences orientés vers un objectif spécifique
A.5 Inférence bayésienne du processus aléatoire gaussien
B L’analyse bayésienne linéaire
Télécharger le rapport complet
