Module REM (Règles pour les calculs d'Evolution avec Mcnp)

C'est le module qui contrôle l'ensemble du calcul d'évolution du réacteur, il dialogue avec les autres modules.

Construction de l'arbre d'évolution

REM connaît la composition de départ de chacune des cellules. A partir des données nucléaires - décroissances radioactives (obtenues de DPH), interactions possibles (obtenues de ACE) - REM peut déterminer, de proche en proche, quels noyaux vont intervenir au cours d'une simulation. Il construit ainsi un « arbre d'évolution » qui contient tous les noyaux présents initialement dans le réacteur, et ceux qui sont produits plus tard. Cet arbre peut contenir plusieurs milliers de noyaux, ce qui est beaucoup. On cherche à élaguer l'arbre afin de limiter la taille du système à intégrer (IME).

On peut éliminer les noyaux dont la période est inférieure à une valeur donnée. Cela revient à dire que le noyau se désintègre instantanément en son descendant, il ne figure dans l'arbre que par son fils.
On peut également en éliminer tous les noyaux dont la production suppose un nombre de réactions (profondeur de l'arbre) supérieur à une valeur donnée - la probabilité que de tels noyaux soient produits est considérée comme si faible que les éliminer du calcul n'a pas de conséquence significative.

Le seuil de période et la profondeur maximum de l'arbre sont des paramètres qui sont fournis à REM à son lancement. Au fur et à mesure que l'arbre se construit, ACE est invoqué, et les fichiers de section efficace d'interaction des différents noyaux sont complétés.

Création de xsdir

Une fois que l'arbre est complet, REM dit à ACE de construire le fichier xsdir qui contient la liste des noyaux et le chemin d'accès au fichier binaire des sections efficaces d'interaction du dit noyau. A ce point, tout est prêt pour lancer la boucle de calcul.

Organisation du calcul

REM fournit à CGM, pour chaque cellule du volume, la composition (le vecteur $\overrightarrow{N}$ courant de la cellule), la température, la densité dans la cellule, les tallys, et il lui demande de construire le fichier d'appel à MCNP. Puis REM déclenche le calcul MCNP. Il invoque ensuite LTM pour récupérer les tallys que MCNP a calculés. Il dispose ainsi des < $\sigma_{{i}}^{}$ $\varphi$ > et des < Yj $\scriptstyle \rightarrow$ i $\sigma^{{cap}}_{{j}}$ $\varphi$ > dans toutes les cellules du volume. Il va pouvoir calculer les éléments de matrice d'évolution pour les passer à IME, en s'appuyant aussi sur les $\lambda_{{i}}^{}$ et les Y'j $\scriptstyle \rightarrow$ i $\lambda_{{j}}^{}$ qui viennent de DPH. REM doit passer à IME :

La matrice d'évolution A pour chaque cellule.
Le vecteur composition $\overrightarrow{N}$ pour chaque cellule.
Le vecteur d'alimentation (dans le cas d'un RSF) $\overrightarrow{B}$ pour chaque cellule.
Les consignes à respecter.
Le temps d'ajustement initial, et la procédure à invoquer chaque fois que le temps d'ajustement est atteint (le temps d'ajustement est variable, contrôlé par IME en fonction de la précision demandée).
L'instant de synchronisation et la procédure à invoquer quand il est atteint.
Le temps total d'intégration pour IME.

C'est REM qui déclenche un calcul IME.

Chaque fois qu'un temps d'ajustement est atteint, la procédure d'ajustement est invoquée. Elle mémorise les résultats (la puissance P, la réactivité k, la composition $\overrightarrow{N}$ de chaque cellule et, dans le cas des réacteurs à sels fondus, l'alimentation $\overrightarrow{B}$ de chaque cellule). La procédure peut invoquer une procédure du programme appelant (le programme principal) qui, lui-même, peut modifier les compositions, les alimentations. Finalement, la procédure rend la main à IME qui reprend son calcul, là où il en était, avec les nouveaux paramètres.

Lorsque la procédure de synchronisation est invoquée (le temps de synchronisation est atteint dans IME), les données sont mémorisées comme lors des ajustements. Des éléments de la matrice d'intégration, les compositions, les alimentations peuvent être ajustés par une procédure du programme appelant qui peut être invoquée comme dans le cas des ajustements. Un nouveau temps de synchronisation, et une nouvelle procédure de synchronisation sont passés à IME avant de lui redonner la main pour qu'il reprenne, là encore, le calcul là où il en était.

Il s'agit, avec ces « temps de synchronisation » de donner au programme principal la possibilité de faire des modifications des compositions, de la matrice d'intégration, du vecteur d'alimentation à des instants prédéterminés, simulant ainsi une intervention sur le réacteur et son mode de fonctionnement.

Il s'agit aussi, avec les « temps d'ajustement » et les « temps de synchronisation », de mémoriser régulièrement les paramètres caractérisant le fonctionnement du réacteur.

Quand le « temps total d'intégration » est atteint, IME rend la main à REM.

Contrôle du « pas MCNP »

Afin de limiter le nombre de calculs MCNP, (très longs par rapport aux calculs IME), il est tentant de faire une interpolation sur les < $\sigma$ $\varphi$ > fournis à IME. Nous avons choisi d'utiliser une spline pour interpoler les < $\sigma$ $\varphi$ > entre deux calculs MCNP, ce qui nous permet d'obtenir une dérivée continue autour d'un calcul MCNP.

Exemple de spline utilisée pour l'interpolation entre deux calculs MCNP

Ainsi, on donne à IME non pas la valeur des < $\sigma$ $\varphi$ >, mais leur évolution, c'est-à-dire les paramètres (a,b,c,d) de la fonction spline at3 + bt2 + ct + d qui approxime l'évolution des < $\sigma$ $\varphi$ >. Encore faut-il connaître les deux valeurs de < $\sigma$ $\varphi$ > entre lesquelles on cherche à faire l'interpolation. Voici comment cela se passe :

A un instant t1, on a fait un calcul MCNP d'où résultent des < $\sigma$ $\varphi$ > 1.
REM lance un calcul IME à < $\sigma$ $\varphi$ > 1 constant, (seul le paramètre d est non nul). Il demande l'intégration jusqu'à un intant t2 (pas MCNP) pour lequel IME trouve une composition $\overline{{N_{2}}}$ (la mémorisation des résultats intermédiaires - ajustements, synchronisations s'effectue pendant cet intervalle).
REM lance un calcul MCNP pour l'instant t2 avec la composition $\overline{{N_{2}}}$ , ce qui donne les < $\sigma$ $\varphi$ > 2 recherchés.
REM calcule les paramètres du spline d'interpolation entre < $\sigma$ $\varphi$ > 1 et < $\sigma$ $\varphi$ > 2 et relance le calcul IME entre les instants t1 et t2 à < $\sigma$ $\varphi$ > variable (il passe les paramètres a, b, c, d du fit du spline dans la matrice d'intégration). Il en résulte une nouvelle composition pour l'instant t2, soit $\overline{{N_{2}}}{^\prime}$ .
REM vérifie si l'écart quadratique absolu entre $\overline{{N_{2}}}$ et $\overline{{N_{2}}}{^\prime}$ [ $\chi^{{2}}_{}$ = ( $\overline{{N_{2}}}$ - $\overline{{N_{2}}}{^\prime}$ )2] est inférieur à une consigne fournie par le programme principal.
Si oui, la composition $\overline{{N_{2}}}{^\prime}$ est conservée. Le calcul IME reprend à partir de t2 jusqu'à un temps t3 tel que ( t3 - t2) = (t2 - t1), avec la composition $\overline{{N_{2}}}{^\prime}$ et < $\sigma$ $\varphi$ > 2 constant, afin de trouver < $\sigma$ $\varphi$ > 3, interpoler, etc.
Dans le cas où le $\chi^{{2}}_{}$ est très bon (inférieur à une deuxième consigne), le temps total d'intégration fixé à IME (le pas MCNP) peut être allongé. REM spécifie alors t3 tel que ( t3 - t2) > (t2 - t1).
Si, par contre, l'écart entre $\overline{{N_{2}}}$ et $\overline{{N_{2}}}{^\prime}$ est supérieur à la consigne, un nouveau calcul MCNP est lancé pour l'instant t2, avec la composition $\overline{{N_{2}}}{^\prime}$ , ce qui donne un < $\sigma$ $\varphi$ > 2'.
REM calcule les paramètres du spline d'interpolation entre < $\sigma$ $\varphi$ > 1 et < $\sigma$ $\varphi$ > 2' et relance le calcul IME entre les instants t1 et t2 à < $\sigma$ $\varphi$ > variable . Il en résulte une nouvelle composition pour l'instant t2, soit $\overline{{N_{2}}}{^\prime}{^\prime}$ .
Si l'écart quadratique entre $\overline{{N_{2}}}{^\prime}$ et $\overline{{N_{2}}}{^\prime}{^\prime}$ est inférieur à la consigne, on garde $\overline{{N_{2}}}{^\prime}{^\prime}$ . Le calcul IME reprend à partir de t2 jusqu'à un temps t3 tel que ( t3 - t2) = (t2 - t1), avec la composition $\overline{{N_{2}}}{^\prime}{^\prime}$ et < $\sigma$ $\varphi$ > 2' constant, afin de trouver < $\sigma$ $\varphi$ > 3, interpoler, etc.
Si l'écart quadratique entre $\overline{{N_{2}}}{^\prime}$ et $\overline{{N_{2}}}{^\prime}{^\prime}$ est supérieur à la consigne, on réduit t2 : on diminue l'intervalle t2 - t1 et on recommence l'intégration IME à partir de t1, en procédant comme ci-dessus pour les vérifications de précision, et en allongeant le pas MCNP dès que possible, c'est-à-dire quand on trouve un très bon écart quadratique entre $\overline{{N_{k}}}$ et $\overline{{N_{k}}}{^\prime}$ .

Avec ce système de contrôle de précision et d'allongement du pas MCNP dès que c'est possible, la fréquence des calculs MCNP est adaptée au comportement du réacteur. Si le réacteur est dans une phase transitoire, les pas MCNP sont courts, comme le sont les pas d'intégration d'IME, comme le sont aussi les temps d'ajustement. Lorsque le réacteur tend vers l'équilibre, le temps entre deux calculs MCNP augmente, les durées d'intégration IME sont de plus en plus longues. Cependant, grâce aux ajustements et synchronisations imposées à IME, on dispose de mesures intermédiaires qui rendent compte de l'évolution des compositions, de la puissance, de la réactivité, on dispose de la possibilité de modifier des paramètres du réacteur.

La façon dont s'intercalent les différentes étapes de calcul IME et MCNP est illustrée dans la figure suivante.

Quelques étapes du calcul d'évolution d'un réacteur montrant que les pas d'ajustement, les pas MCNP et les durées d'intégration RK4 changent au cours du temps, et pas forcément de façon régulière.

Simplification de l'intégration

Certains noyaux du vecteur composition ont une durée de vie suffisamment courte pour qu'on puisse considérer qu'ils sont en équilibre avec leur père quand le réacteur a fonctionné assez longtemps pour que leur nombre ait eu le temps de se stabiliser avec la population du père ; on peut accélérer le calcul en les mettant en équilibre avec le père, c'est-à-dire qu'ils restent dans le vecteur $\overrightarrow{N}$ mais leur population, au lieu d'être calculée en fonction des éléments de la matrice d'intégration, est fixée en fonction de la population du père.

Prenons comme exemple le 135Xe dans un réacteur à sels fondus. C'est un gaz, il sort rapidement du coeur du réacteur grâce à un bullage d'hélium. Sa période dans le coeur est de 30 secondes. Cependant, il a une section efficace de capture grande, si bien que son taux de capture, n < $\sigma$ $\varphi$ > ne peut être négligé : les captures sur le Xenon font perdre des neutrons proportionnellement à la quantité de Xenon présent, ne pas en tenir compte changerait la réactivité du réacteur. On est donc obligé de suivre le Xenon. Mais, si on le suivait comme les autres noyaux, les bruits de calcul introduiraient des fluctuations du vecteur $\overrightarrow{N}$ telles qu'on ne parviendrait jamais à allonger le pas d'intégration dans IME. D'où ce compromis, qui consiste à estimer la population du noyau en fonction de la population de son père, quand c'est possible, d'autant que la population peut être mieux estimée par l'hypothèse de l'équilibre avec le père que par le calcul. Bien entendu, le Xenon n'est pas le seul noyau qui peut être traité ainsi.

C'est IME qui assure la mise à l'équilibre d'un noyau avec son père, sous le contrôle de REM. Celui-ci donne à IME un temps, le seuil de demie vie effective (calculée en combinant le taux de réaction et la décroissance radioactive du noyau) en dessous duquel IME met le noyau en équilibre avec son père. Cependant, il faut être sûr que le réacteur a fonctionné assez longtemps pour que cette mise à l'équilibre soit valable. C'est pourquoi REM donne à IME un temps minimal de fonctionnement, à partir duquel la mise à l'équilibre est autorisée.

La boucle de calcul

En résumé, la boucle de calcul est une suite de calculs MCNP, suivis de calculs , le tout étant piloté par la précision, comme illustré dans la figure. Dans REM, c'est le « pas MCNP » (qui est aussi la durée totale d'intégration IME) qui est ajusté en fonction de la variabilité des vecteurs $\overrightarrow{N}$ .

REM reste dans sa boucle de calcul jusqu'à ce que le temps total de simulation du réacteur, fixé par le programme principal, soit atteint et là, il lui rend la main. Notons que le programme principal a pu suivre le processus grâce à l'activation de ses procédures en cours de calcul, lors des ajustements et/ou synchronisations).

Gestion des compositions et des cellules volumiques

Voyons comment sont gérés les éléments de volume tels qu'ils sont vus par REM, IME, CGM.

Conceptuellement, un réacteur se découpe en zones dont la fonction, donc le type de matériaux qu'elles contiennent, diffèrent. On est amené, par exemple, à caractériser le combustible, le modérateur, le caloporteur, les éléments de structure, les réflecteurs. Ces matériaux, que nous appellerons des « composites » sont à des températures différentes suivant leur emplacement dans le réacteur. Pour un composite donné, par exemple un réflecteur, la composition initiale est la même, quelle que soit la température de fonctionnement. Par contre, elle évolue de façon différente selon la zone de température.

Par ailleurs, les éléments de volume, les cellules au sens CGM, que nous appellerons les « cellules CGM », qui décrivent le découpage géométrique du réacteur demandé par le programme principal, peuvent fractionner une zone de température. La correspondance entre les « cellules CGM » et les zones de température doit être connue de REM, elle est fournie par le programme principal.

Ainsi, REM reçoit du programme principal une liste de matériaux - les composites - avec leur composition initiale, les zones de température (et leur température), la liste des cellules CGM correspondantes.

REM construit les arbres d'évolution pour chaque composite par zone de température, et non pas pour chaque cellule CGM, ce qui serait inutile : on obtiendrait une duplication d'arbres d'évolution dans le cas où plusieurs cellules CGM recouvrent une même zone de température.

Les calculs MCNP sont toujours faits individuellement sur les cellules CGM. Par contre, les calculs IME peuvent être faits sur un groupe de cellules qui recouvre tout ou partie d'une zone de température pour un composite donné. Un tel groupement de cellules CGM suppose qu'elles ont la même composition à tout instant (calcul IME), bien qu'elles ne soient pas forcément soumises exactement au même flux (calcul MCNP).

Une fois que les arbres d'évolution sont construits, REM peut déterminer quels sont les tallys qu'il voudra passer à CGM pour qu'il les demande à MCNP (ils figurent dans le fichier de commande MCNP). En ce qui concerne les sections efficaces, c'est simple, c'est les < $\sigma$ $\varphi$ > de toutes les réactions possibles, pour chacun des noyaux susceptible d'interagir (pour lesquels ACE a trouvé des sections efficaces) qui figurent dans l'arbre d'évolution. En effet, on démontre que MCNP peut calculer les sections efficaces moyennes d'interactions sur des noyaux qui ne sont qu'à l'état de trace (ou pas encore apparus dans la composition fournie à MCNP) du moment que les fichiers de sections efficaces lui sont fournis. Il suffit que MCNP ait évalué le flux (qui n'est, a priori, pas influencé par des noyaux qui ne sont dans le milieu qu'à l'état de trace). Ainsi, il est possible de suivre l'évolution des noyaux qui ne figurent pas dans les compositions initiales au fur et à mesure qu'ils apparaissent : leurs < $\sigma$ $\varphi$ > étant évalués, ils participent à la matrice d'évolution passée à IME.

Décroissances fictives

Lors de la simulation d'un réacteur, on peut être amené à simuler l'extraction de certains noyaux : dans un RSF, par exemple, on utilise un bullage d'helium pour extraire des produits de fission gazeux. Pour cela, nous avons introduit la notion de décroissance radioactive fictive. Soit un noyau i, bien réél. Il a ses sections efficaces d'interaction (fournies par ACE) et/ou ses décroissances radioactives (fournies par DPH). Lui ajouter une décroissance fictive, c'est lui ajouter une décroissance radioactive vers un noyau, disons k, « fictif », c'est-à-dire, en l'occurence, un noyau qui n'interagira pas avec les neutrons. En effet, un noyau extrait est, par définition hors flux, il ne peut plus interagir avec les neutrons. Par contre, on veut suivre son évolution et, s'il est radioactif, on veut pouvoir suivre ses descendants qui, eux aussi, doivent être hors flux.

L'artifice que nous avons utilisé est d'attribuer un numéro d'isomère négatif aux noyaux « fictifs », tel que leur numéro d'isomère d'origine puisse être retrouvé et, par conséquent, leur éventuelle décroissance radioactive suivie (qui se fera vers un noyau dont l'isomère est négatif, puisqu'il faut rester hors flux). ACE n'est pas interrogé pour les noyaux dont l'isomère est négatif, ils n'ont aucune section efficace d'interaction. Mais on veut suivre leur évolution avec IME : sachant quels noyaux il faut « extraire » dans chaque composite, et avec quelle constante de temps chacun « décroît » vers le noyau fictif associé, REM ajoute ces noyaux (et leurs descendants - fictifs également) aux arbres d'évolution. Ils sont ensuite traités de façon standard.

Ce dont on dispose à la fin de REM

On a déjà vu que les compositions trouvées à chaque ajustement et synchronisation lors du calcul d'évolution sont mémorisées, pour toutes les cellules que connaît IME, qui peuvent correspondre à des groupes de cellules CGM. En outre, on mémorise tous les tallys demandés à MCNP. Afin d'en faciliter l'exploitation, ces informations sont organisées de manière à regrouper en un ensemble

les résultats de l'évolution durant tout un pas MCNP
les résultats du calcul MCNP de fin de ce pas

Nous appellerons un tel ensemble un inventaire.

Pour chaque inventaire, pour chaque cellule, on a donc mémorisé :

les compositions de chaque instant caractéristique (ajustement, synchronisation)
les tallys déterminés par MCNP
la réactivité telle que déterminée par MCNP (qui la calcule différemment de IME) et, pour chaque instant, des données globales du réacteur, soit :
le nombre de neutrons créés
le nombre de neutrons absorbés
la réactivité évaluée par IME
le nombre de neutrons source (le paramètre $\Phi$ sur lequel joue IME pour ajuster la puissance)

spline : courbe du 3 eme degré couramment utilisée pour l'interpolation entre des points expérimentaux. Elle assure la continuité des dérivées aux points de jonction entre deux splines (points expérimentaux).