Conception d’un outil de simulation d’activation neutronique à vocation industrielle

Conception d’un outil de simulation d’activation
neutronique à vocation industrielle

Les sections efficaces

Dans la section 3.3.1, les probabilités utilisées afin de calculer les différentes valeurs nécessaires à la simulation ont été calculées grâce aux libres parcours moyens dans le matériau. Le libre parcours moyen est lui aussi issu d’autres valeurs qui sont les sections efficaces macroscopiques. L’équation 3.3 donne ainsi la relation entre le libre parcours moyen et la section efficace macroscopique. λlpm,j = 1 Σj (3.3) où · λlpm,j est le libre parcours moyen pour le processus j (en cm) · Σj est la section efficace macroscopique du processus j dans la matériau observé (en cm−1 ). La section efficace est donc la valeur de base permettant de calculer l’ensemble des probabilités nécessaires au transport de particules par méthode Monte-Carlo. D’après l’équation 3.3, il est possible de constater que le libre parcours moyen est inversement proportionnel à la section efficace macroscopique. Cela se traduit par le fait que plus les chances d’interaction (entre le neutron et la matière) sont importantes, plus la distance entre deux collisions est faible. Ainsi, plus un matériau sera dense et contiendra des isotopes ayant de fortes chances d’interagir avec des neutrons (i. e. ayant une section efficace importante), plus les neutrons feront de collisions.

Estimateurs de la méthode Monte-Carlo

La section 3.3.1 a permis de décrire le principe de fonctionnement du transport de particules par méthode Monte-Carlo. Cependant, sans résultat de sortie (appelé score), la méthode Monte-Carlo n’aurait pas d’intérêt. Il est donc nécessaire de comptabiliser une ou des grandeurs en sortie du calcul. Pour cela, l’espérance de la valeur souhaitée en sortie doit être calculée. En supposant que f(x) soit la fonction de densité de probabilité de contribution à la valeur de sortie par échantillonnage aléatoire du transport et ayant pour contribution x, alors l’espérance du résultat est définie par l’équation 3.4 E(x) = Z xf(x)dx (3.4) où · f(x) est la fonction densité de probabilité de contribuer de x au résultat en échantillonnant aléatoirement le transport · x est la contribution au score. La fonction f(x) n’étant généralement pas connue, la méthode Monte-Carlo échantillonne implicitement celle-ci en transportant aléatoirement les particules. Les scores sont alors calculés à l’aide d’un ou plusieurs estimateurs de l’espérance qui sont définis dans cette partie du manuscrit. De plus, afin de prendre en compte le fait que l’échantillon de particules simulées doit être représentatif du vrai nombre de particules issues des sources, on associe à chaque particule simulée et provenant de la source un poids numérique ws défini lors du lancement de la particule par l’équation 3.5 ws = Isource Isim (3.5) où · Isource est l’intensité des sources de particules (en s −1 ) · Isim est le nombre de particules issues de la source et simulées (sans unité). Le poids défini par l’équation 3.5 permettra donc de pondérer les résultats observés afin de prendre en compte le fait que seul un échantillon de particules est simulé. Le poids des particules pourra également être modifié judicieusement lors du transport des particules afin d’accélérer les calculs. Le poids de la particule avant d’effectuer la collision i (noté wi) sera donc égal à ws ∗ Qi−1 j=1 wb,j avec wb,j le poids de biaisage associé à la collision j. Si aucun biaisage n’est utilisé, alors wb,j = 1. Ce procédé sera décrit plus précisément dans la partie 3.3.4. Il existe différentes manières de comptabiliser un résultat issu d’un transport par méthode Monte-Carlo. Étant donné que le transport par méthode MonteCarlo est basé sur des interactions ponctuelles (en termes de position, d’énergie et de direction), les résultats sont généralement moyennés selon un découpage plus ou moins fins de l’espace des phases défini par l’utilisateur. Le but de cette section sera donc de présenter les principales façons de procéder afin de calculer ces résultats.

Estimateur collision

Une manière pouvant être utilisée afin de comptabiliser un résultat de calcul MonteCarlo est l’estimateur collision. À la différence des méthodes déterministes évoquées dans la section 3.2, les valeurs de sortie simulées ne sont pas directement calculées. En effet, dans le cas de la méthode Monte-Carlo, les résultats sont mis à jour à chaque fois qu’une particule atteint un détecteur. Dans le cas de l’estimateur collision, le flux scalaire pour un groupe énergétique g d’un détecteur de volume V est estimé par l’équation 3.6 φg = 1 V X i∈A wi Σt(Ei) (3.6) où · V est le volume du détecteur (en cm3 ) · A est l’ensemble de tous les évènements correspondant à une collision entre une particule et un noyau dans le détecteur et ayant une énergie comprise entre Einf,g et Esup,g · Einf,g la borne inférieure du groupe énergétique g (en MeV ) · Esup,g la borne supérieure du groupe énergétique g (en MeV ) · Σt(Ei) est la section efficace macroscopique totale du détecteur à l’énergie Ei (en cm−1 ) · wi est le poids de la particule avant d’effectuer la collision i (en s −1 ). L’estimateur collision permet donc de calculer le flux neutronique dans un détecteur à partir des collisions entre les neutrons et les noyaux. À noter cependant que dans le cas où la composition du détecteur n’autorise pas un nombre important de collisions, alors peu d’évènements se produiront et le résultat obtenu ne sera que peu judicieux. Afin de prendre en compte tous les évènements se produisant dans un détecteur, l’estimateur corde peut être utilisé.

Estimateur corde

L’estimateur corde permet de calculer le flux scalaire de particules lors d’un transport par méthode Monte-Carlo en utilisant la distance parcourue par les particules entre chaque interaction. L’estimateur corde peut ainsi être défini par l’équation 3.7 φg = 1 V X i∈B liwi (3.7) où · V est le volume du détecteur (en cm3 ) · B est l’ensemble de tous les évènements correspondant à un passage dans le détecteur d’une particule ayant une énergie comprise entre Einf,g et Esup,g · Einf,g la borne inférieure du groupe énergétique g (en MeV ) Esup,g la borne supérieure du groupe énergétique g (en MeV ) · li est la distance parcourue dans le détecteur par une particule entre la position de l’interaction i − 1 et la position de l’interaction i (en cm) · wi est le poids de la particule avant d’effectuer l’interaction i (en s −1 ). L’estimateur corde permet donc également de calculer le flux scalaire. Cependant, à la différence de l’estimateur collision, il tient compte de tous les évènements se produisant dans un détecteur comme par exemple le processus de transport (qui n’est pas une collision). Ainsi, avec un même historique de transport, un détecteur utilisant l’estimateur corde fournira un résultat plus convergé qu’un détecteur utilisant l’estimateur collision. Enfin, un autre estimateur pouvant être utilisé est l’estimateur surfacique.

Estimateur surfacique

L’estimateur surfacique permet de calculer le flux scalaire à partir des particules traversant la surface du détecteur. L’estimateur surfacique peut être défini par l’équation 3.8 φg = 1 S X i∈C wi ||Ω~ · ~n|| (3.8) où · S est la surface du détecteur (en cm2 ) · C est l’ensemble de tous les évènements correspondant à une traversée de la surface du détecteur par une particule entre l’interaction i − 1 et i et ayant une énergie comprise entre Einf,g et Esup,g · Einf,g la borne inférieure du groupe énergétique g (en MeV ) · Esup,g la borne supérieure du groupe énergétique g (en MeV ) · Ω~ est le vecteur direction de la particule avant l’interaction i · ~n est le vecteur normal à la surface du détecteur au point d’entrée de la particule orienté vers l’intérieur du détecteur · ||Ω~ · ~n|| est la norme du produit scalaire entre Ω~ et ~n · wi est le poids de la particule lorsqu’elle traverse la surface du détecteur avant l’interaction i (en s −1 ). L’estimateur surfacique permet ainsi de calculer le flux scalaire de particules à partir des évènements de traversée de la surface du détecteur. Ceci peut être utile en particulier dans le cas où le volume du détecteur ne peut pas être calculé ce qui peut arriver pour des formes complexes. Cependant, l’estimateur surfacique doit être utilisé avec précaution car comme cela peut être vu dans l’équation 3.8, un ratio de valeurs est utilisé dans la somme. Ainsi, dans le cas où la norme du produit scalaire tend vers 0, le flux divergera. Cela se produit donc quand la direction de la particule devient tangente à la surface du détecteur. En pratique, seules les particules ayant un angle d’incidence supérieur à une valeur définie par l’utilisateur sont prises en compte afin d’éviter la divergence. Ceci implique que le flux calculé par l’estimateur surfacique n’est pas forcément la vraie valeur du flux. C’est donc pour cette raison que les estimateurs collision et corde sont généralement les estimateurs utilisés pour les résultats du Monte-Carlo bien que dans certains cas l’estimateur surfacique puisse être plus judicieux.

Table des matières

Remerciements
Glossaire
Acronymes
Liste des symboles
Table des figures
Liste des tableaux
1 Introduction
1.1 Contexte
1.1.1 Contexte général
1.1.2 Modélisation de l’activation neutronique
1.1.3 Méthodes de calcul de l’activation neutronique
1.2 Objectifs
1.3 Organisation du manuscrit
I Théorie et méthodes numériques pour l’activation neutronique
2 Théorie de l’activation neutronique
2.1 Introduction
2.2 Les phénomènes mis en jeux
2.2.1 La transmutation
2.2.2 La fission
2.2.2.1 La fission spontanée
2.2.2.2 La fission induite
2.2.3 La radioactivité
2.2.3.1 La radioactivité α
2.2.3.2 La radioactivité β
2.2.3.3 L’émission spontanée
2.2.4 Le phénomène d’activation neutronique
2.3 Les grandeurs de dosimétrie
2.3.1 La dose absorbée
2.3.2 La dose équivalente
2.3.3 La dose efficace
2.3.4 Le débit d’équivalent de dose ambiant
2.4 Les équations à résoudre
2.4.1 Les grandeurs d’intérêt
2.4.2 Notion de section efficace
2.4.3 L’équation du transport stationnaire de Boltzmann
2.4.4 Les équations de Bateman généralisées
2.4.5 Le couplage des équations
2.5 Les données nucléaires
2.5.1 Les sections efficaces
2.5.1.1 Les sections efficaces de transport
2.5.1.2 Les sections efficaces d’activation
2.5.1.3 Le traitement des sections efficaces
2.5.2 Les données de décroissance radioactive
2.5.3 Les taux de fission
2.6 Conclusion
3 Méthodes numériques de résolution des équations liées à l’activation neutronique
3.1 Introduction
3.2 Résolution de l’équation de Boltzmann stationnaire par méthode déterministe
3.2.1 Méthode des caractéristiques
3.2.2 Méthode des probabilités de première collision
3.2.3 Méthode des harmoniques sphériques
3.2.4 Méthode des ordonnées discrètes
3.3 Transport de particules par méthode Monte-Carlo
3.3.1 Principe de la méthode Monte-Carlo appliquée au transport des particules
3.3.1.1 Activation neutronique et méthode Monte-Carlo
3.3.1.2 Description détaillée
3.3.1.3 Mode de calcul analogue
3.3.1.4 Mode de calcul non-analogue
3.3.2 Les sections efficaces
3.3.3 Estimateurs de la méthode Monte-Carlo
3.3.3.1 Estimateur collision
3.3.3.2 Estimateur corde
3.3.3.3 Estimateur surfacique
3.3.3.4 Incertitude statistique
3.3.4 Techniques de réduction de variance
3.3.4.1 Notion de biaisage
3.3.4.2 Le splitting
3.3.4.3 La roulette Russe
3.3.4.4 La capture implicite
3.3.4.5 Les fenêtres de poids
3.4 Résolution des équations de Bateman généralisées
3.4.1 Linéarisation des chaînes de réaction
3.4.2 Matrices exponentielles
3.4.3 Résolution des ODE
3.4.3.1 La méthode BDF
3.4.3.2 Schéma prédicteur/correcteur
3.4.3.3 La méthode Block BDF
3.5 Codes de calcul
3.5.1 Codes de calcul Monte-Carlo
3.5.2 Codes de calcul d’inventaire isotopique
3.6 Méthodes de modélisation de l’activation neutronique
3.6.1 Prise en compte des gradients de flux
3.6.2 Méthode Direct-1-Step
3.6.3 Calcul à la volée
3.6.4 Méthode Rigorous-Two-Step
3.7 Conclusion
II Mise en place d’une méthodologie d’analyse d’activation neutronique à vocation industrielle
4 Étude des problématiques géométriques liées à l’activation neutronique
4.1 Introduction
4.2 Prise en compte des gradients de flux grâce au maillage virtuel superposé
4.2.1 Les différents types de maillage
4.2.1.1 Motif de base et structure du maillage
4.2.1.2 Conformation au modèle
4.2.1.3 Système de coordonnées
4.2.2 La maillage cartésien structuré non conforme
4.2.2.1 Avantages
4.2.2.2 Limites
4.3 Découpage automatisé du modèle et calcul des inventaires isotopiques initiaux grâce à la CAO
4.3.1 Calcul Monte-Carlo des matériaux englobés
4.3.2 Modélisation CAO
4.3.3 Calcul CAO du volume des cellules englobées
4.3.3.1 Le procédé
4.3.3.2 Les limites
4.4 Le choix du maillage cartésien adapté automatiquement
4.4.1 Choix des paramètres du maillage pour le calcul des flux neutroniques
4.4.2 Comparaison au maillage tétraédrique
4.4.2.1 Mémoire vive nécessaire pour le maillage
4.4.2.2 Durée de mise en place des maillages
4.4.3 Performance et précision du maillage
4.4.4 Évaluation du gain de précision
4.4.4.1 Impact du maillage adapté automatiquement sur la précision du flux neutronique et des inventaires isotopiques
4.4.4.2 Exemple de sensibilité au volume de cellule englobée dans une maille
4.5 Conclusion
5 Traitement des aspects numériques de l’activation neutronique
5.1 Introduction
5.2 Calcul Monte-Carlo
5.2.1 Parallélisation et gestion de la RAM
5.2.1.1 Découpage énergétique des flux neutroniques
5.2.1.2 Limite à la parallélisation du transport Monte-Carlo
5.2.1.3 Mutualisation des détecteurs
5.2.2 Incertitudes statistiques
5.3 Aspects numériques du calcul d’inventaire isotopique
5.3.1 Le solveur d’EDO
5.3.1.1 La méthode BBDF
5.3.1.2 Recherche du zéro par itération de Newton-Raphson
5.3.1.3 Le solveur de système linéaire
5.3.2 Optimisations du solveur d’EDO
5.3.2.1 Contrôle de l’erreur et gestion du pas
5.3.2.2 Réduction de la taille du système
5.3.2.3 Parallélisation du calcul d’inventaire isotopique
5.3.3 Automatisation des calculs d’inventaire isotopique
5.3.3.1 Stabilité du schéma numérique
5.3.3.2 Configurations du solveur d’EDO
5.3.3.3 Linéarisation de chaînes
5.3.3.4 Algorithme de calcul des inventaires isotopiques
5.4 Échantillonnage des sources de décroissance
5.4.1 Échantillonnage de la position et de l’intensité des particules de décroissance
5.4.2 Énergies ponctuelles ou par histogramme
5.4.3 Biasing énergétique
5.5 Conclusion
III Résultats et développements futurs
6 Résultats et analyses
6.1 Introduction
6.2 Benchmark ITER port plug
6.2.1 Paramètres de simulation
6.2.1.1 Géométrie
6.2.1.2 Maillage et sources
6.2.1.3 Bases de données nucléaires
6.2.1.4 Stations de calcul
6.2.2 Autres codes d’analyse d’activation neutronique
6.2.2.1 Spécificités de l’outil R2Smesh
6.2.2.2 Spécificités de l’outil R2SUNED
6.2.2.3 Couplage MCNP/FISPACT
6.2.3 Résultats
6.3 Performances
6.3.1 Calcul des flux neutroniques
6.3.1.1 Limite à la parallélisation
6.3.1.2 Temps de calcul
6.3.2 Calcul des inventaires isotopiques
6.3.2.1 Calcul de l’inventaire initial
6.3.2.2 Génération des fichiers d’entrée FISPACT
6.3.2.3 Résolution des équations de Bateman
6.3.2.4 Temps de calcul complet
6.3.3 Calcul du débit d’équivalent de dose ambiante
6.3.3.1 Efficacité du transport des photons de décroissance
6.3.3.2 Biaisage de la source de décroissance radioactive .
6.4 Unification du processus de simulation de l’activation neutronique
6.4.1 Réduction de l’espace disque utilisé
6.5 Conclusion
7 Bilan et perspectives
7.1 Conclusions
7.1.1 Fonctionnalités de la méthodologie RayActive
7.1.2 Résultats de RayActive
7.2 Futurs développements
7.3 Bilan
IV Bibliographie et annexes
Bibliographie
A Découpage énergétique par défaut de RayActive
A.1 Découpage énergétique par défaut pour le flux neutronique dans RayActive
A.2 Découpage énergétique par défaut pour les sources de photons issues des décroissances radioactives dans RayActive
B Développement des coefficients du solveur d’EDO
C Composition des matériaux utilisés pour les calculs présentés
C.1 Composition massique des matériaux utilisés pour le calcul du débit d’équivalent de dose ambiante du modèle du dispositif expérimental
C.2 Composition massique des matériaux utilisés pour le benchmark du ITER port plug
Résumé