Optimisation stochastique d'un mécanisme de réaction d'oxyde d'uranium à l'aide de mesures de réacteur à flux de plasma

Aug 27, 2023

Rapports scientifiques volume 13, Numéro d'article : 9293 (2023) Citer cet article

Détails des métriques

Dans ce travail, une approche couplée d'algorithme génétique de Monte Carlo (MCGA) est utilisée pour optimiser un mécanisme de réaction d'oxyde d'uranium en phase gazeuse basé sur des mesures de réacteur à flux de plasma (PFR). Le PFR produit un plasma Ar stable contenant des espèces U, O, H et N avec des régions à haute température (3000–5000 K) pertinentes pour l'observation de la formation d'UO via la spectroscopie d'émission optique. Un traitement cinétique global est utilisé pour modéliser l'évolution chimique dans le PFR et produire des signaux d'émission synthétiques pour une comparaison directe avec les expériences. L'espace des paramètres d'un mécanisme de réaction d'oxyde d'uranium est ensuite exploré via un échantillonnage de Monte Carlo en utilisant des fonctions objectives pour quantifier l'accord modèle-expérience. Les résultats de Monte Carlo sont ensuite affinés à l'aide d'un algorithme génétique pour obtenir un ensemble de voies de réaction et de coefficients de vitesse corroborés expérimentalement. Sur les 12 canaux de réaction ciblés pour l'optimisation, quatre canaux sont bien contraints dans toutes les exécutions d'optimisation, tandis que trois autres canaux sont contraints dans certains cas. Les canaux optimisés mettent en évidence l'importance du radical OH dans l'oxydation de l'uranium dans le PFR. Cette étude comprend une première étape vers la production d'un mécanisme de réaction complet validé expérimentalement pour la formation d'espèces moléculaires d'uranium en phase gazeuse.

La cinétique de réaction des oxydes métalliques en phase gazeuse est d'une grande importance pour de nombreux domaines de recherche, notamment l'astrophysique, la science de la combustion, le génie nucléaire et la chimie des matériaux dans des environnements extrêmes. Ces dernières années, ce dernier domaine a produit de nombreux travaux expérimentaux et informatiques sur la chimie de la vapeur d'oxyde d'uranium (\({{\mathrm{UO_x}}}\))1. Les produits en phase gazeuse d'oxydes réfractaires, comme \({{\mathrm{UO_x}}}\), ont toujours été difficiles à produire en raison des températures de vaporisation élevées des oxydes parents. Plus récemment, les systèmes à plasma thermique ont fourni une voie pour produire facilement des métaux en phase gazeuse et étudier leur chimie dans des environnements réactifs. Cependant, les temps d'extinction rapides, la présence de radicaux de fond et la formation d'oxydes intermédiaires volatils dans de tels systèmes rendent difficile l'isolement de canaux de réaction spécifiques à étudier. Des problèmes similaires se posent dans d'autres systèmes réactifs à haute température tels que les combustibles de combustion de métaux. En conséquence, les mécanismes d'oxydation des métaux en phase gazeuse sont souvent basés sur des données expérimentales rares et des estimations théoriques de premier ordre, comme pour la formation d'oxyde d'aluminium2,3,4. De même, un mécanisme de réaction \({{\mathrm{UO_x}}}\) a été construit en utilisant une méthodologie comparable dans nos travaux précédents5. Bien que de tels mécanismes produisent des résultats qualitativement raisonnables qui peuvent s'aligner sur certains observables expérimentaux, une validation expérimentale détaillée est difficile à réaliser. Cette étape de validation est cruciale pour garantir que le mécanisme cinétique chimique peut être utilisé de manière prédictive pour informer les modèles ultérieurs. Ici, nous explorons une méthode pour déduire les coefficients de taux d'oxyde d'uranium (\({{\mathrm{UO_x}}}\)) basée sur des mesures expérimentales d'un système à plasma thermique.

En raison de la nature fortement couplée et non linéaire de la cinétique chimique dans les plasmas d'uranium, l'extraction d'informations sur la vitesse de réaction nécessite la résolution d'un problème d'optimisation. Dans ce problème, les paramètres sous-jacents du modèle (coefficients de taux) sont déterminés sur la base des sorties observées (c'est-à-dire des informations spectroscopiques). Résoudre un tel problème par des méthodes déterministes basées sur le gradient est difficile en raison de l'espace de paramètres potentiellement complexe avec de nombreux minima locaux. Dans ce cas, il faut plutôt utiliser une méthode d'optimisation capable d'explorer en continu tout l'espace des paramètres tout en localisant le minimum global. Une telle méthode précédemment utilisée pour les problèmes de cinétique chimique est l'algorithme génétique de Monte Carlo (MCGA)6. Cette technique est bien adaptée au problème actuel en raison de son efficacité à éviter la convergence sur les minima locaux et de sa facilité de mise en œuvre. Indépendamment de la méthodologie, la résolution d'un problème d'optimisation nécessite des évaluations répétées du modèle associé, souvent au nombre de milliers à des millions d'exécutions. Alors que des temps de calcul raisonnables sont atteints lors de la résolution de la cinétique chimique dans un système spatialement uniforme, le problème devient rapidement irréalisable lorsque la chimie est couplée à un transport de fluide complexe. Cette considération devient importante lors du choix d'un système expérimental pour informer le problème d'optimisation.

Les systèmes expérimentaux appropriés pour informer l'optimisation du mécanisme de réaction uranium-oxygène à pression atmosphérique comprennent les systèmes d'ablation laser7,8,9,10,11,12,13,14 et un réacteur à flux de plasma (PFR)15,16. Les systèmes d'ablation au laser utilisent un laser pulsé à haute intensité pour volatiliser les échantillons de métal, produisant des panaches de plasma réactif en expansion rapide. Si elle est effectuée par un laser suffisamment puissant dans des conditions atmosphériques, l'ablation s'accompagne d'une onde de choc à l'interface panache-ambiant rappelant une onde de souffle boule de feu. La composition chimique du panache en fonction du temps peut alors être mesurée par spectroscopie optique. Le réacteur à flux de plasma, quant à lui, produit un plasma d'uranium à l'aide d'une torche à plasma à couplage inductif (ICP) fixée à un tube de quartz. Alors que le plasma RF est généré par un flux d'argon, une solution aqueuse d'uranium est introduite dans la torche, produisant un plasma uranifère qui se refroidit en s'écoulant vers l'aval. La spectroscopie optique est utilisée pour mesurer l'évolution chimique d'espèces sélectionnées dans le plasma à divers points le long du tube.

Bien que les deux systèmes expérimentaux puissent théoriquement être utilisés pour optimiser un mécanisme de réaction, un avantage majeur de l'utilisation du réacteur à flux de plasma est que le temps de séjour des espèces peut être corrélé à la distance le long du réacteur. Autrement dit, si nous suivons une parcelle de fluide avec une composition chimique de départ connue à travers le réacteur, sa composition chimique à une position donnée peut être liée à son temps de séjour via le débit. Cette approche lagrangienne permet d'approximer l'évolution chimique dans un réacteur à flux de plasma à l'aide d'un modèle purement transitoire, tel qu'un modèle cinétique global15. En revanche, le comportement de transport complexe d'un panache d'ablation laser nécessite l'utilisation d'un modèle plus coûteux en termes de calcul, tel qu'un modèle de fluide compressible réactif17. Sur la base de notre expérience précédente avec les simulations ci-dessus, un modèle cinétique global du réacteur à flux se termine en quelques secondes, alors qu'un modèle d'ablation laser fluide peut prendre des heures ou plus. La réduction significative de l'effort de calcul requis pour modéliser le réacteur à flux de plasma est essentielle pour l'optimisation et motive l'utilisation du système dans ce travail. Ci-dessous, nous décrivons comment les expériences PFR, la modélisation et l'optimisation MCGA sont effectuées pour produire un mécanisme de réaction \({{\mathrm{UO_x}}}\) contraint expérimentalement.

Schéma du réacteur à flux de plasma. La sous-figure (a) montre la région de l'injecteur de la buse en amont, l'emplacement du plasma et le tube de quartz en aval. La sous-figure (b) met en évidence les trois canaux d'écoulement d'entrée annulaires concentriques, l'emplacement de la bobine de plasma à couplage inductif (ICP) et l'injecteur d'écoulement annulaire facultatif de la région d'entrée.

Le PFR est un système ICP disponible dans le commerce qui a été modifié18 pour étudier la cinétique chimique en phase gazeuse ainsi que la formation et la croissance des nanoparticules. Un schéma du PFR utilisé dans ce travail est illustré à la Fig. 1. La région d'entrée du PFR se compose de trois canaux d'écoulement annulaires concentriques, chacun avec un débit et une composition distincts. Une solution aqueuse de nitrate d'uranyle (\({{\textrm{UO}}}_2({{\mathrm{NO_3}}})_2\cdot {{\mathrm{6H_2O}}}\)) est nébulisée dans des gouttelettes liquides et introduite par le canal le plus interne (marqué en rouge) à l'aide d'un gaz porteur (argon). Pour un débit typique d'argon gazeux de 1 L/min, l'uranium est environ quatre ordres de grandeur moins abondant que l'argon dans le flux le plus interne (c'est-à-dire \(\sim\)100 ppm). Pour améliorer la cinétique d'oxydation, de l'oxygène gazeux peut être ajouté par ce canal, avec des débits typiques de 10 à 50 mL/min. De plus, le canal le plus à l'extérieur (marqué en bleu) fournit un débit de gaz d'argon supplémentaire de 12 à 14,4 L/min pour entretenir le plasma et refroidir la paroi de quartz externe. Le canal central (marqué en vert) n'est pas utilisé dans ces expériences. Sur la base des débits ci-dessus, les concentrations d'analyte sont de l'ordre de 10 à 100 ppm dans le flux aval, en fonction de l'étendue du mélange radial et de la diffusion. Les densités en nombre, les débits et la composition des composants fluides avant d'entrer dans le plasma sont répertoriés dans le tableau 1, le nitrate d'uranium étant divisé en ses molécules constitutives pour plus de commodité.

Un plasma RF de 40 MHz est généré en aval des canaux d'entrée à l'aide d'une bobine inductive entourant le tube de quartz externe. Étant donné que la majorité du flux d'entrée est de l'argon et que le plasma est généré à la pression atmosphérique, les propriétés thermodynamiques et de transport du plasma peuvent être rapprochées de celles d'un plasma d'argon LTE19,20. Les températures du plasma et du flux aval peuvent être modifiées en ajustant le débit d'argon le plus externe et la puissance fournie par l'alimentation. Enfin, un injecteur de flux annulaire optionnel peut être utilisé pour introduire un flux d'argon supplémentaire plus en aval de la bobine RF, bien que cette fonctionnalité n'ait pas été utilisée dans ce travail. Alternativement, une extension de tube de quartz de diamètre constant peut être connectée à la torche lorsque l'injecteur à flux annulaire n'est pas nécessaire.

La spectroscopie d'émission optique (OES) est utilisée pour suivre l'évolution de l'U et de l'UO dans le PFR. La lumière émise par le plasma est acheminée vers un spectromètre à l'aide d'un câble à fibre optique positionné à divers emplacements axiaux le long du réacteur à flux. Un étage de translation linéaire motorisé est utilisé pour déplacer le câble à fibre optique le long de l'axe x indiqué sur la figure 1, en maintenant la fibre à une distance radiale fixe du centre du réacteur. L'extrémité de la bobine RF est utilisée comme emplacement axial de référence \(x=0\) pour toutes les mesures (comme illustré à la Fig. 1).

Les configurations d'écoulement annulaire et de diamètre constant décrites ci-dessus présentent toutes deux des régions optiquement opaques où l'émission d'écoulement est obscurcie. Ces régions couvrent respectivement 0 à 3 cm et 3 à 5 cm de la bobine RF pour l'extension à diamètre constant et l'injecteur à flux annulaire. De plus, comme la pointe de la fibre optique est conductrice et n'a pas été isolée, la distance axiale minimale de la bobine RF a été maintenue à 1 cm pour éviter la formation d'arcs. Étant donné que les caractéristiques d'écoulement dans la région de la torche doivent être identiques pour les deux configurations, les limitations d'observation sont surmontées en combinant les données en amont et en aval prises avec les configurations d'écoulement annulaire et de diamètre constant, respectivement.

L'étalonnage d'un mécanisme de réaction par rapport à des quantités mesurées expérimentalement est un exemple de problème inverse, c'est-à-dire un problème où les équations gouvernantes et la solution sont connues, mais les paramètres d'entrée ne le sont pas. Ce type de problème n'admet généralement pas de solution unique et se pose plutôt comme un problème d'optimisation où l'adéquation d'une solution est dictée par une fonction objectif. La fonction objectif quantifie l'écart statistique de la solution calculée par rapport à la vraie solution. Par exemple, une fonction objectif courante est la somme des carrés des résidus de la solution :

where \({\varvec{k}}\) is a vector containing the reaction rate coefficients and \(n^{exp}_i\) and \(n^{calc}_i({\varvec{k}})\) are the measured and calculated species number densities at time point i, respectively. The optimization problem is solved by employing an iterative procedure that finds an optimal parameter set \({\varvec{k}}\) that minimizes the objective function \(\phi\). In the context of the current problem, an optimized \({\varvec{k}}\) value would represent a set of rate coefficients that closely match the uranium oxide formation rates observed in the laser ablation or PFR experiments. Typically, deterministic nonlinear least squares methods, such as the Gauss-Newton or Levenberg-Marquadt methods3.0.CO;2-R (1998)." href="#ref-CR21" id="ref-link-section-d58307020e1787"> 21,22,23, sont employés pour de tels problèmes d'optimisation. Les techniques informatiques modernes, telles que les réseaux de neurones24, peuvent également être utilisées à cette fin.

En raison du grand espace de paramètres du problème d'optimisation du mécanisme de réaction \({{\mathrm{UO_x}}}\) , l'espace des solutions peut être complexe et peut contenir de nombreux minima locaux. Les méthodes d'optimisation déterministes conventionnelles ont du mal à localiser un minimum global pour un tel problème, convergeant plutôt vers des minima locaux adjacents au point d'initialisation. Les méthodes de recherche exhaustives sont également inefficaces en raison de la demande de calcul pour cartographier un grand espace de paramètres. Pour éviter ces problèmes, nous utilisons une approche d'algorithme génétique de Monte Carlo (MCGA)6 pour optimiser le mécanisme de réaction \({{\mathrm{UO_x}}}\). Cette approche combine les méthodes d'optimisation stochastique de Monte Carlo et de l'algorithme génétique pour obtenir une optimisation globale pour les problèmes avec de grands espaces de paramètres. La partie Monte Carlo de l'approche utilise un échantillonnage aléatoire des paramètres de vitesse de réaction pour localiser les régions de bonne fitness dans l'espace des solutions. L'algorithme génétique optimise ensuite ces régions pour trouver le minimum global en utilisant des processus évolutifs de migration, de sélection, d'accouplement et de mutation. La nature stochastique de ces processus maintient la diversité parmi les ensembles de paramètres optimisés, évitant ainsi la convergence vers des minima locaux. Ainsi, l'approche MCGA permet une optimisation globale en recherchant l'ensemble de l'espace des solutions à l'aide de l'échantillonnage de Monte Carlo couplé à la diversité paramétrique inhérente à l'optimisation des algorithmes génétiques. De plus, MCGA peut être facilement adapté à différents systèmes expérimentaux, puisque les équations gouvernantes du système modélisé n'ont pas besoin d'être reformulées comme un problème inverse. Enfin, MCGA est facile à paralléliser : la fonction objectif de chaque ensemble de paramètres peut être évaluée indépendamment, de sorte que le processus d'évaluation peut être librement réparti entre les processeurs. La robustesse, la flexibilité et la rapidité de l'approche MCGA en font un excellent outil pour produire un mécanisme de réaction calibré expérimentalement pour la formation d'oxyde d'uranium.

Un diagramme décrivant l'approche MCGA est illustré à la Fig. 2. Suivant le diagramme, le processus de Monte Carlo peut être séparé en plusieurs tâches clés : génération du mécanisme de réaction, modification du coefficient de vitesse, évaluation du modèle et évaluation de l'ajustement. Chacune de ces tâches, et la manière dont elles s'intègrent dans la partie Monte Carlo de l'algorithme, sont détaillées ci-dessous dans des sous-sections nommées en conséquence. Toutes ces tâches, sauf la première, sont également présentes dans l'algorithme génétique, comme indiqué dans la dernière sous-section ici.

Diagramme des parties Monte Carlo (à gauche de la ligne pointillée) et de l'algorithme génétique (à droite de la ligne pointillée) de l'algorithme génétique de Monte Carlo (MCGA). L'ensemble des mécanismes de réaction produits par le processus de Monte Carlo sert de population de départ pour l'algorithme génétique, où chaque compagnon est un mécanisme de réaction. Les opérations génétiques responsables de la production des générations suivantes sont détaillées plus loin dans la Fig. 4.

La première étape du processus MCGA consiste à spécifier l'ensemble des canaux de réaction (un mécanisme de réaction) qui seront utilisés pour évaluer le comportement chimique du système. Le mécanisme de réaction comprend deux parties : un ensemble de canaux de réaction de l'uranium qui sont la cible de l'optimisation et un ensemble de canaux de réaction de soutien responsables de la chimie de fond. En raison de la solution de nitrate d'uranyle utilisée dans le PFR, ce dernier mécanisme consiste en divers canaux de réaction en phase gazeuse25,26,27 et en phase plasma25,28,29,30 qui détaillent le comportement chimique d'un plasma OHN. Pour réduire le coût de calcul de l'évaluation de la chimie de fond, une étape de réduction du mécanisme de réaction est effectuée sur ce mécanisme OHN. La réduction est effectuée en excluant les molécules dont la formation est défavorable dans les conditions d'intérêt (2000–5000 K), ainsi qu'en éliminant les réseaux de réaction étendus qui suivent les états atomiques et moléculaires excités des espèces mineures. Par exemple, l'azote étant présent en faible quantité, les réactions le mettant en jeu sont considérablement réduites, affectant au minimum les résultats de calcul tout en réduisant sensiblement le temps de calcul. Chaque étape de la réduction a été vérifiée en exécutant une simulation de test 0D et en vérifiant que la formation d'espèces \({{\mathrm{UO_x}}}\) était peu affectée. De plus, le mécanisme réduit a été vérifié après chaque optimisation MCGA en testant les mécanismes UO résultants avec les mécanismes OHN réduits et complets, trouvant un bon accord. Le mécanisme OHN réduit final comprend 44 espèces et 166 canaux de réaction, contre 81 espèces et 796 canaux de réaction pour le mécanisme complet.

Les canaux de réaction d'uranium ciblés pour l'optimisation sont répertoriés dans le tableau 2. Seules les réactions qui sont contraintes par les données expérimentales disponibles et les conditions du système sont incluses pour éviter un éventuel surajustement à partir de réactions mal contraintes. Par exemple, étant donné que les espèces réactives sont diluées dans le flux, les réactions à trois corps avec un tiers corps non Ar seront peu fréquentes et peuvent être exclues. De plus, étant donné que la pression du système est maintenue constante, la dépendance de la pression de réaction n'est pas contrainte. Par conséquent, la liste des réactions possibles est limitée aux réactions bimoléculaires, ce qui réduit considérablement le nombre de voies potentielles. De plus, étant donné que les mesures d'émissions comprenant les ensembles de données sont limitées à U et UO, la chimie des oxydes d'uranium supérieurs n'est pas bien contrainte. Bien que les mesures fournissent certaines contraintes sur la formation de \(\hbox {UO}_2\) via le taux de consommation d'UO, elles ne contiennent aucune information concernant la consommation de \(\hbox {UO}_2\) et les taux de formation de \(\hbox {UO}_3\). Par conséquent, seules les réactions impliquant U ou UO dans la direction exothermique sont prises en compte pour l'optimisation. Notez que le tableau 2 comprend également deux réactions d'ionisation associatives de notre mécanisme de réaction précédemment construit5. Ces réactions ont un impact important sur la chimie du plasma d'uranium en raison de la vitesse de réaction presque sphérique dure pour le canal d'ionisation associatif \({{\mathrm{U + O}}}\)31. Cependant, à notre connaissance, ce comportement n'est pas bien validé. Par conséquent, nous incluons ces canaux dans l'optimisation pour déterminer l'importance des voies d'ionisation associatives pour la formation de \({{\mathrm{UO_x}}}\). Bien que \({{\mathrm{UO^+}}}\) et \({{\mathrm{UO^+_2}}}\) ne soient pas directement mesurés ici, ces réactions sont partiellement limitées par les données U et UO disponibles.

Les coefficients de vitesse de départ (\(k_{est}\)) pour chaque réaction du tableau 2 sont estimés à l'aide de diverses approximations du premier ordre32 exprimées sous une forme de type Arrhenius modifiée :

où A est la fréquence de collision, T est la température du gaz, n est une constante de puissance de température, \(E_A\) est l'énergie d'activation et R est la constante du gaz. Les méthodes de la théorie des collisions simples (SCT) et du modèle simplifié des triples collisions (SMTC) sont utilisées pour calculer respectivement les coefficients de taux binaire et à trois corps. Les sections efficaces de collision pour les molécules sont estimées à partir des longueurs de liaison et des volumes de Van der Waals combinés des atomes constitutifs. Ces estimations fournissent une borne supérieure sur la fréquence de collision A et une constante de puissance de température \(n=0,5\) due à une contribution de la vitesse thermique. Aucune estimation a priori n'est faite pour l'énergie d'activation \(E_A\); tous les canaux de réaction sont initialement supposés être sans barrière. Les canaux de réaction sont exprimés dans la direction exothermique pour éviter des taux de réaction inverse non physiquement élevés. Notez également que les réactions d'association bimoléculaire sont formulées dans la limite de haute pression pour fournir une estimation du taux limite supérieur comme point de départ.

En plus des voies de réaction \({{\mathrm{UO_x}}}\) soumises à optimisation par le MCGA, nous considérons également plusieurs voies de réaction supplémentaires de l'uranium qui ne sont pas ajustées par l'algorithme. Ces canaux de réaction sont présentés dans le tableau 3 et consistent principalement en des réactions chimiques plasmatiques (ionisation, recombinaison, échange de charge) ainsi qu'en des réactions entre \(\hbox {UO}_2\) et \(\hbox {UO}_3\). Les taux de réaction de ces canaux sont maintenus fixes en raison d'un manque de données expérimentales contraignantes. Néanmoins, ces réactions fournissent des voies pour la chimie du plasma d'uranium et la formation d'oxydes supérieurs dans le modèle.

Une fois que le mécanisme de réaction UO cible est généré, la boucle principale de Monte Carlo commence. Chaque itération de cette boucle est indépendante et consiste à évaluer et évaluer une version modifiée du mécanisme de réaction UO généré. Chaque mécanisme modifié est produit en ajustant les paramètres d'Arrhenius du mécanisme d'origine comme suit :

où f est un facteur compris entre \(10^{-4}\) et \(10^0\) échantillonné aléatoirement à partir d'une distribution uniforme logarithmique de base 10, m est un facteur compris entre \(-3\) et 0 échantillonné aléatoirement à partir d'une distribution uniforme, et e est un facteur compris entre \(10^0\) et \(10^{4.6}\) échantillonné aléatoirement à partir d'une distribution uniforme logarithmique de base 10. Les bornes de ces facteurs sont choisies de manière à ne pas dépasser la limite supérieure physique fournie par les estimations initiales du taux de sphère dure. Ainsi, les valeurs \(k_{est}\) et \(k_{mod,min}\) indiquées dans le tableau 2 représentent respectivement les limites supérieure et inférieure des taux modifiés. Le facteur f est destiné à compenser la surestimation de la fréquence de collision A fournie par l'estimation initiale du taux de sphère dure \(k_{est}\). Le facteur m représente un changement de la dépendance à la température n de la forme d'Arrhenius modifiée. Le facteur e représente un ajustement de l'énergie d'activation \(E_A/R\). La valeur limite supérieure de \(e=10^{4.6}\) est choisie de sorte que la partie exponentielle de l'Eq. (3) réduit la vitesse de réaction de quatre ordres de grandeur à la température maximale du plasma de \(\sim 4500\) K. Par conséquent, si la valeur maximale de l'énergie d'activation est utilisée, le canal de réaction est effectivement retiré du mécanisme de réaction, indiquant que l'énergie d'activation est trop élevée pour que la réaction se produise dans le système actuel. Inversement, une énergie d'activation de \(10^0\) indique que pratiquement aucune barrière d'activation n'est présente pour la réaction, et le terme exponentiel dans l'Eq. (3) aura peu d'effet sur la vitesse de réaction.

Le mécanisme de réaction modifié est évalué à l'aide d'un modèle cinétique global afin de calculer l'évolution chimique à l'intérieur du réacteur à flux de plasma. Le modèle résout l'équilibre chimique transitoire d'un système fermé (adiabatique) spatialement uniforme (0D). Le modèle consiste en un système d'équations différentielles ordinaires fortement couplées, chacune régissant la concentration d'une espèce chimique particulière :

où \(n_i\) est la densité numérique des espèces i et \({\dot{S}}_{ij}\) est un terme source décrivant la contribution d'une réaction élémentaire j à la production ou à la consommation des espèces i. \({\point{S}}_{ij}\) est donné par :

où \(\Delta c_{ij}\) est le coefficient stoechiométrique net pour l'espèce i dans la réaction j, \({\dot{R}}_j\) est la vitesse de réaction pour la réaction j, \(k_j\) est le coefficient de vitesse pour la réaction j, et c est le coefficient stoechiométrique du réactif s dans la réaction j. Si une réaction est réversible, les taux de réaction avant et arrière (\(k_+\) et \(k_-\), respectivement) sont déterminés en utilisant le principe de l'équilibre détaillé via un coefficient d'équilibre \(K_{eq} = k_+ /k_-\). Le coefficient d'équilibre est déterminé par les propriétés thermodynamiques des espèces participant à la réaction5.

En règle générale, le système d'équations d'équilibre des espèces est complété par une équation d'équilibre thermique qui tient compte du chauffage ohmique, de la libération d'énergie chimique et du refroidissement convectif/conducteur/radiatif dans le flux de plasma. Ici, nous interpolons plutôt la température du plasma à chaque pas de temps sur la base d'un profil de température calibré expérimentalement. Cela réduit à la fois la complexité de calcul du modèle et améliore la cohérence entre les conditions de température modélisées et expérimentales. Un modèle CFD développé précédemment18 est utilisé pour suivre la température dans le PFR le long d'une ligne de courant lagrangienne. La trace de colis de fluide lagrangien fournit à la fois un profil de température et une corrélation temps/distance pour exprimer les concentrations 0D dépendant du temps en termes de localisation axiale. La figure 3 compare ce profil de température lagrangien aux mesures de température expérimentales disponibles. Les températures expérimentales sont extraites des intensités relatives des raies atomiques de Fe (en utilisant un analyte de nitrate de fer) et des probabilités de transition de raie correspondantes en utilisant la méthode du diagramme de Boltzmann18. Notez que cette méthode de détermination de la température s'est avérée précédemment donner des valeurs cohérentes dans l'incertitude de mesure pour les nitrates de fer, d'aluminium et de cérium15,33, indiquant une insensibilité à l'analyte utilisé.

Historique de température lagrangien produit par un modèle CFD du PFR18 (ligne continue) par rapport aux valeurs de température expérimentales obtenues à partir de mesures précédentes de la ligne atomique Fe via la méthode du tracé de Boltzmann (points).

Une version modifiée du package ZDPlasKin34 est utilisée pour intégrer les ODE ci-dessus. Le système 0D suit une parcelle fluide de gaz parfait sous pression atmosphérique, où la loi des gaz parfaits est appliquée en ajustant la densité numérique totale en fonction du profil de température. La condition initiale correspond à un emplacement à l'intérieur du canal d'écoulement d'analyte en amont de la bobine ICP, où l'écoulement est à température ambiante et les molécules d'analyte ne sont pas dissociées. Les concentrations initiales d'espèces sont calculées à l'aide des débits moléculaires expérimentaux \({\dot{N}}\) donnés dans le tableau 1 comme suit :

où \(n_{0,i}\) est la densité numérique initiale de la molécule i, \(P_0=1\) atm et \(T_0=300\) K. Un débit de canal d'analyte Ar de 1 L/min est utilisé pour ce calcul. Nous supposons que la ligne de courant lagrangienne représentative subit un mélange limité avec le flux d'Ar provenant du canal externe, de sorte que le rapport de mélange Ar/analyte reste constant tout au long de la simulation.

Enfin, l'évolution chimique modélisée est comparée à l'évolution observée expérimentalement pour évaluer l'accord du mécanisme de réaction modifié. Cette comparaison nécessite la conversion des sorties de densité numérique de l'étape d'évaluation en signaux d'intensité d'émission correspondants. Ici, nous utilisons l'hypothèse simplificatrice selon laquelle le plasma PFR est optiquement mince et les effets d'auto-absorption peuvent être négligés. L'intensité synthétique (en unités de puissance par volume) due à une désexcitation électronique d'un état supérieur 2 à un état inférieur 1 est calculée comme suit :

où h est la constante de Planck, c est la vitesse de la lumière, \(\lambda\) est la longueur d'onde de la raie, \(n_{2}\) est la population de l'état excité et \(A_{21}\) est la probabilité de transition. Pour une transition atomique, la densité numérique des atomes pour le niveau électronique e est exprimée en fonction de la densité numérique totale des espèces n as35 :

où \(g_{el,e}\) et \(\Delta E_{e0}\) sont respectivement le poids statistique et l'énergie par rapport à l'état fondamental du niveau e, \(k_B\) est la constante de Boltzmann et \(q_{el}\) est la fonction de partition électronique :

qui est une somme pondérée des populations dans tous les états électroniques. La fonction de partition électronique de l'atome d'uranium est calculée à partir des niveaux d'énergie des bases de données de spectroscopie atomique36,37,38. L'uranium a une forte raie d'émission atomique à 591,5 nm, qui est utilisée pour des comparaisons expérimentales dans cette étude.

Pour une transition moléculaire diatomique hétéronucléaire, la densité numérique des particules à un niveau donné est calculée comme suit :

où le niveau excité est décrit par les nombres quantiques électroniques, rotationnels et vibrationnels e, J et v, respectivement. Ici, une expression approximative est utilisée pour la fonction de partition interne totale \(q_{int}\) dans ces états :

où les termes de rotor rigide et d'oscillateur harmonique \(q_{rot,rr}\) et \(q_{vib,ho}\) sont donnés par :

où \(B_e\) est la constante de rotation et \(\omega _e\) est la fréquence vibrationnelle harmonique, toutes deux exprimées en unités d'énergie. Le terme de correction approximatif du premier ordre \(q_{corr}\)35,39 tient compte du mouvement anharmonique non rigide et du couplage rovibrationnel comme :

où \(D_e\) est la constante de distorsion centrifuge, \(\alpha _e\) est la constante de couplage rotation-vibrationnelle et \(\chi _e\) est la constante d'anharmonicité. Dans l'expression ci-dessus, \(D_e\), \(\alpha _e\) et \(\omega _e\chi _e\) sont exprimés en unités d'énergie. En raison de la complexité du spectre d'émission UO encombré, seules les constantes spectroscopiques de l'état fondamental ont été précédemment estimées dans la littérature40,41,42,43,44. Les constantes spectroscopiques et les niveaux d'énergie estimés par Konings et al.44 seront utilisés ici pour estimer l'intensité d'émission de la bande UO à 593,55 nm. Cette bande est supposée dominée par la transition [16.845]5–X(1)4 (une transition 0–0) observée par Kaledin et al.43. Bien que de nombreuses raies rovibrationnelles étroitement espacées contribuent également à cette bande, elles ne sont pas traitées dans les travaux actuels en raison à la fois d'un manque de constantes spectroscopiques et de la résolution limitée du spectromètre utilisé.

En théorie, l'éq. (7) pourrait être utilisé pour relier les intensités d'émission mesurées aux densités numériques absolues en utilisant un signal de référence bien connu (c'est-à-dire une ligne Ar forte). Cependant, les probabilités de transition des espèces \({{\mathrm{UO_x}}}\) ont des incertitudes aussi élevées que 50 %45, ce qui interdit une détermination précise des densités numériques basées sur les émissions. Pour contourner cette limitation, nous normalisons à la place les profils d'intensité d'émission expérimentaux et modélisés de sorte que la fonction objective minimise la différence dans les formes, plutôt que dans les amplitudes, des deux profils. Le signal d'émission le plus fort (le plus en amont) est utilisé comme point de normalisation de sorte que les courbes d'émission expérimentales et modélisées sont mises à l'échelle pour avoir la même amplitude au point d'émission le plus fort. Des valeurs anormalement basses des signaux d'émission U synthétiques (c'est-à-dire trop basses pour produire un signal détectable) sont évitées en incluant un terme de pénalité dans la fonction objectif qui vérifie le rapport maximal des intensités d'émission synthétiques \(I_{U/UO}=I_U/I_{UO}\). Comme son nom l'indique, le terme de pénalité pénalise mathématiquement les solutions qui se situent en dehors de la plage souhaitée de valeurs \(I_{U/UO}\) et est donné plus loin. Sur la base des considérations ci-dessus, la fonction objectif d'erreur quadratique moyenne (RMSE) suivante (\(\phi\)) est formulée :

où \(\phi _p\) est la fonction de pénalité et \(R^2_{lin,c}\) et \(R^2_{log,c}\) sont les coefficients linéaires et logarithmiques pondérés de détermination pour la courbe d'émission c donnée par :

où \(w_i=W_i/\sum _i^{\textrm{N}}W_i\) est un poids statistique normalisé, \({\bar{I}}_{lin}^{exp}=\sum _i^{\textrm{N}}w_iI_i^{exp}\) et \({\bar{I}}_{log}^{exp}=\sum _i^{\textrm{N}}w_ i\log (I_i^{exp})\) sont les moyennes linéaires et logarithmiques pondérées des intensités d'émission normalisées mesurées, respectivement, et \(I^{exp}_i\) et \(I^{calc}_i({\varvec{k}})\) sont les intensités d'émission normalisées mesurées et calculées au temps i, respectivement. N représente le nombre de points de données expérimentales comprenant la courbe d'émission c, et C représente le nombre total de courbes d'émission utilisées. L'utilisation de coefficients linéaires et logarithmiques ici vise à garantir que les changements de grande amplitude (qui dominent l'ajustement linéaire) et les changements de petite amplitude (qui dominent l'ajustement logarithmique) dans les signaux d'émission sont bien ajustés. Le premier est limité par la chute rapide des émissions près de la bobine RF et le second par la décroissance plus progressive des émissions plus en aval. Le poids statistique \(W_i\) de chaque point de données est donné par :

où \(S_j^{bck}\) et \(S_l^{mol}\) sont les signaux de fond et d'émission moléculaire consistant respectivement en B et M points de données expérimentaux. Cette pondération quantifie l'intensité d'une raie d'émission donnée par rapport à l'intensité du fond. Par conséquent, les raies faiblement émettrices de l'ordre de l'arrière-plan seront pondérées plus légèrement lors de l'évaluation de la fonction objectif par rapport aux raies à forte émission.

En raison de la corrélation inverse entre \(I_{U/UO}\) et \(n_{UO}/n_U\), le terme de pénalité suivant est utilisé :

où \({\bar{I}}_{U/UO}=11\) et \(I^{min}_{U/UO}=1\). La fonction déconseille le plus fortement les solutions où \(n_{UO}/n_U \gg 10\) pour éviter de localiser des solutions avec une bonne fitness mais des densités U déraisonnablement faibles. Les valeurs de probabilité de transition utilisées pour calculer \(I_U\) et \(I_{UO}\) sont \(A_{U,21}=3,15 \times 10^6\) s\(^{-1}\) pour la raie U à 591,54 nm46 et \(A_{UO,21}=3,8 \times 10^9\) s\(^{-1}\) pour la tête 0–0 de la bande UO 593,55 nm47, respectivement. Une incertitude du pire scénario de \(\pm 50\%\)45 est supposée pour ces deux valeurs. À partir des valeurs \(I_{U/UO}\) maximales expérimentales comprises entre 3 et 7, une plage approximative de \(1< I_{U/UO} < 21\) pour le rapport modélisé est obtenue. En calculant les densités en nombre d'espèces correspondant à la borne inférieure \(I_{U/UO}=1\), on trouve \(n_{UO}/n_U\approx 120\). En conséquence, les valeurs moyennes et maximales de \(I_{U/UO}=11\) et \(I_{U/UO}=21\) donnent \(n_{UO}/n_U\environ 11\) et \(n_{UO}/n_U\environ 5,5\), respectivement.

La dernière partie de la boucle de Monte Carlo consiste à décider de conserver ou de rejeter le mécanisme de réaction candidat en fonction de l'adéquation de la solution. Cette étape construit l'ensemble des mécanismes de réaction qui seront transmis à l'algorithme génétique pour une optimisation plus poussée. Les critères de conservation des solutions candidates peuvent être ajustés selon les besoins, mais sont généralement basés sur le dépassement d'une valeur seuil pour une ou plusieurs des mesures de fitness calculées ci-dessus. Par exemple, un critère simple tel que \(R^2_{lin}>0.5\) pour toutes les courbes de solution peut être utilisé. Si une solution est retenue, alors les valeurs statistiques calculées dans le cadre de l'Eq. (14) et les ensembles correspondants de paramètres de taux modifiés de l'Eq. (3) sont stockés. Il est préférable d'avoir un critère de sélection laxiste à ce stade afin de constituer un nombre important de mécanismes candidats. Les mécanismes transmis à l'algorithme génétique sont ensuite sélectionnés parmi ces candidats en fonction de la métrique de fitness utilisée par l'algorithme.

Diagramme illustrant la sélection du tournoi, le croisement uniforme et les opérations de mutation de la partie algorithme génétique de l'algorithme génétique de Monte Carlo (MCGA). Dans cet exemple, chaque "individu" ou "chromosome" (mécanisme de réaction) est composé de 18 "gènes" (6 canaux de réaction avec 3 coefficients de vitesse chacun).

La mise en œuvre de l'algorithme génétique (GA) utilisée ici utilise un mélange d'opérations et d'idées issues de plusieurs approches publiées précédemment6,48,49. Comme c'est souvent le cas, le choix des paramètres GA est largement heuristique basé sur les taux observés de convergence, d'amélioration de la condition physique et de variabilité de la solution. Par conséquent, bien que la mise en œuvre GA employée ici ne soit pas nécessairement optimale, elle fonctionne néanmoins suffisamment bien pour améliorer de manière fiable la condition physique de la population tout en maintenant un niveau de base de diversité.

Dans l'AG, la population de chaque nouvelle génération est créée en effectuant plusieurs opérations génétiques sur la population actuelle, comme illustré à la Fig. 4. Ici, une population fait référence à un ensemble d'individus, dont chacun est un mécanisme de réaction avec des coefficients de vitesse modifiés via Eq. (3). Chaque mécanisme contient un nombre de paramètres égal au nombre de réactions multiplié par le nombre de coefficients d'Arrhenius (soit 3 paramètres par réaction). La nouvelle génération commence par transférer directement un certain nombre des mécanismes les plus aptes (élites) de la génération actuelle à la nouvelle. Ici, nous sélectionnons un nombre d'élites égal à 5% de la population totale. Ensuite, la sélection des partenaires est effectuée à l'aide d'un k-tournoi de taille \(k=2\) avec une probabilité de 80 % de sélectionner le mécanisme le plus apte dans le tournoi. Autrement dit, chaque tournoi choisit au hasard 2 mécanismes parmi la population actuelle, et le plus apte des deux est sélectionné avec une probabilité de 80 % (sinon, l'autre mécanisme est sélectionné). Toute la population participe au processus de sélection, y compris les élites. Chaque couple de partenaires distincts ainsi sélectionné a alors une chance de se reproduire et/ou de muter avant d'être ajouté à la nouvelle population. Les partenaires peuvent être convertis en progéniture par croisement uniforme, dans lequel chaque paramètre (coefficient de taux) est commuté entre les deux partenaires avec une probabilité égale, ce qui entraîne deux mécanismes de progéniture. Dans ce travail, la probabilité de reproduction croisée pour un couple de partenaires donné est de 65 %. Qu'il y ait ou non reproduction, les deux mécanismes peuvent alors également subir une mutation. Dans cette opération, chaque paramètre (coefficient de taux) d'un mécanisme est randomisé dans des limites préalablement spécifiées avec une probabilité de 0,8 %. La probabilité de mutation est maintenue faible car même une seule opération de re-randomisation peut produire un changement radical dans le comportement du mécanisme. Étant donné que les opérations de reproduction et de mutation peuvent ou non se produire pour une paire de partenaires donnée, il est possible qu'une paire donnée passe simplement à la génération suivante (c'est-à-dire survive). Cela ne les supprime pas du pool de sélection, de sorte que les paires de partenaires peuvent à la fois survivre et se reproduire/muter. Cependant, chaque mécanisme de la nouvelle population doit être unique, de sorte que les mécanismes en double sont vérifiés et supprimés si nécessaire. Les opérations ci-dessus de sélection, de reproduction et de mutation sont répétées jusqu'à ce qu'une nouvelle population de la même taille que l'original soit générée.

La condition physique des élites étant déjà connue, seuls 95% de chaque nouvelle génération doivent être évalués. Il existe différentes méthodes pour paralléliser ce processus50, dont certaines impliquent la création de sous-populations soit pour permettre une évaluation non synchronisée, soit pour permettre des opérations évolutives supplémentaires (par exemple la migration). Ici, nous utilisons une simple parallélisation contrôleur-travailleur dans laquelle seule la tâche d'évaluation (qui est la tâche la plus exigeante en termes de calcul) est répartie entre le nombre de processeurs disponibles. Notez que, pour plus de diversité génétique, la première génération est construite en complétant la population initiale (les mécanismes générés par Monte Carlo les plus adaptés) avec un nombre égal de mécanismes avec des coefficients de taux aléatoires (dans les limites autorisées)6.

Plusieurs ensembles de données d'émission optique d'uranium ont été collectés à partir du PFR pour éclairer l'optimisation. Les paramètres variaient entre les ensembles de données, notamment les régions d'observation, la température d'écoulement (variée par les débits et la puissance RF) et la concentration d'oxygène dans le flux d'analyte, comme résumé dans le tableau 4. Des exemples de spectres d'émission d'uranium et de fond collectés à 3 et 8 cm de la bobine RF sont illustrés à la Fig. Par conséquent, le fond mesuré comprend l'émission due à la désexcitation des espèces de fond, le rayonnement continu (thermique), le bruit inhérent de l'instrument et toute autre lumière de fond parasite. Les spectres de l'uranium montrent des signaux de fond significatifs même lorsqu'ils sont corrigés pour le fond mesuré, comme observé dans les précédentes études de spectroscopie de l'uranium1. Ce bruit de fond « excessif » se produit probablement en raison de la multitude de raies d'émission d'uranium étroitement espacées dans le spectre visible combinées à la résolution limitée des spectromètres typiques. En tant que tel, l'extraction de l'intensité d'émission d'une ligne ou d'une bande d'uranium donnée nécessite d'abord de soustraire ce signal de fond supplémentaire. Le faire avec précision nécessiterait soit d'utiliser un spectromètre à résolution beaucoup plus élevée, soit de tenter de déconvoluer les pics d'intérêt à l'aide d'un modèle spectral complet d'oxyde d'uranium, qui se situent tous deux en dehors de la portée de ce travail. Au lieu de cela, nous supposons que les pics de fond sont à la fois beaucoup plus faibles et plus nombreux que les bandes d'intérêt, de sorte qu'un simple décalage peut être utilisé pour séparer approximativement l'un de l'autre. Il est difficile d'évaluer l'incertitude introduite par cette hypothèse sans un modèle spectral complet d'émission d'oxyde d'uranium. Qualitativement, il devrait avoir un impact plus important sur la bande UO à 593,55 nm que sur la raie U atomique à 591,5 nm, car le premier signal est généralement plus faible que la seconde. De plus, la bande UO est constituée de plusieurs raies rovibrationnelles rapprochées qui nécessitent une résolution spectrale beaucoup plus élevée (de l'ordre de 0,004 nm14) pour se résoudre correctement. Enfin, le calcul de la fonction de partition UO est approximatif en raison des informations limitées sur les états internes du système44. Par conséquent, la plus grande incertitude dans la procédure d'optimisation actuelle réside probablement dans la mesure et le calcul du signal dû à la tête 0–0 de la bande UO 593,55 nm.

Graphiques montrant les spectres d'émission d'uranium, de fond (argon) et d'uranium soustrait de fond mesurés à 3 et 8 cm de la bobine RF. Les lignes verticales indiquent les emplacements de la ligne U à 591,5 nm et de la bande UO à 593,55 nm.

Le fond d'argon et d'eau mesuré est négligeable par rapport au fond d'uranium pour la plupart des emplacements, de sorte que sa correction n'a pas d'impact notable sur les résultats. Même pour les emplacements en aval où le fond d'argon et d'eau était du même ordre de grandeur que le fond d'uranium, la contribution du fond d'argon et d'eau dans les régions spectrales d'intérêt était minime. Par conséquent, les mesures de fond n'ont pas été effectuées pour les ensembles de données ultérieurs, et à la place, un décalage constant basé sur le fond d'uranium total a été appliqué. Ceci est accompli en sélectionnant une plage de longueurs d'onde (586 à 586,5 nm) comme région de fond de référence (en raison du manque de pics visibles à tous les emplacements) et en utilisant le signal moyen dans cette plage comme valeur de décalage. Après application de ce décalage, les intensités de la raie U 591,5 nm et de la bande UO 593,55 nm sont calculées par intégration sur les pics correspondants.

Exemple de spectres d'uranium (a) et intensités de raies soustraites de fond intégrées correspondantes (b) pour les ensembles de données 1, 2, 3 sans flux \(\hbox {O}_2\) ajouté.

La figure 6 montre des exemples de spectres d'uranium de chaque ensemble de données ainsi que les valeurs d'intensité soustraites de fond intégrées correspondantes pour U et UO sur tous les emplacements mesurés. Les tracés d'intensité intégrés montrent à la fois le signal U et UO diminuant de manière monotone sur toute la plage d'observation. Cela reste vrai même en tenant compte de la diminution de l'intensité des émissions due à la baisse de température sur cette distance. Étant donné que les températures dans la région en amont (1 cm) devraient être suffisamment élevées (\(\sim\)4500 à 5000 K) pour dissocier au moins partiellement l'UO, la concentration en UO et l'intensité des émissions devraient initialement augmenter en aval. La tendance UO observée, cependant, suggère que la formation d'UO se produit plutôt encore plus en amont du premier point d'observation ici (c'est-à-dire dans la région de la bobine). L'origine de ce comportement peut être expliquée en modélisant les profils d'espèces dans le flux, comme illustré à la Fig. 7. Au début, le flux d'analyte se compose uniquement des molécules réactives constitutives \(\hbox {UO}_2\), \(\hbox {H}_2\)O, \(\hbox {NO}_3\) (non représenté) et du gaz porteur Ar à température ambiante. En se déplaçant vers l'aval, le flux rencontre un gradient de température prononcé à \(-2\) cm qui diminue la densité du nombre de gaz (via la loi des gaz parfaits) et produit une dissociation, une excitation et une ionisation rapides des molécules d'analyte. Notez que l'ionisation de l'uranium dans ce cas est presque entièrement due au canal d'ionisation associatif \({\mathrm{U+O}}\) du mécanisme non optimisé. Alors que \(\hbox {H}_2\)O et \(\hbox {NO}_3\) sont effectivement entièrement dissociés en leurs composants atomiques, le flux traverse la région du plasma trop rapidement pour dissocier complètement les molécules d'oxyde d'uranium. En conséquence, on voit apparaître deux pics UO ; le premier en raison du gradient de température initial séparant \(\hbox {UO}_2\) dans l'analyte et le second en raison du refroidissement en aval permettant aux réactions \({\mathrm{U+O}}\) et à la recombinaison électronique de UO\(^+\) d'avoir lieu. Après ce point, un refroidissement progressif induit la formation d'oxydes d'uranium supérieurs (\(\hbox {UO}_2\) et \(\hbox {UO}_3\)) qui finissent par épuiser les molécules d'UO précédemment formées. La diminution monotone du signal UO observée dans les ensembles de données ci-dessus correspond à ce dernier régime. Notez également qu'en raison de la disponibilité d'oxygène libre à partir de la dissociation de \(\hbox {H}_2\)O et \(\hbox {NO}_3\), l'uranium se sature vers un oxyde plus élevé (\(\hbox {UO}_3\)) que sa forme d'analyte initiale (\(\hbox {UO}_2\)).

Résultats non optimisés du modèle 0D \({{\mathrm{UO_x}}}\)5 montrant des profils d'espèces sélectionnés dans le flux en fonction du profil de température illustré à la Fig. 3.

Les comparaisons des profils d'émission synthétiques et des données d'émission mesurées sont présentées à la Fig. 8. Le profil d'émission synthétique d'UO (2a) montre un accord étonnamment étroit avec les données expérimentales. Cependant, le modèle sous-estime légèrement les taux d'épuisement de l'UO en amont et surestime les taux en aval. Un accord plus faible est observé pour le profil U (1a), qui diminue plus rapidement dans le modèle dans les régions en amont et en aval. Une caractéristique intéressante des données d'émission d'U est observée dans les tracés semi-logarithmiques. À savoir, le signal d'émission U semble saturer vers une valeur minimale après environ 4 à 5 cm. Le rapport signal sur bruit de fond pour la raie U à 591,5 nm à ces emplacements reste constamment élevé, de sorte que ce comportement ne peut être attribué au bruit de l'instrument ou au bruit de fond de l'uranium. En effet, ce comportement n'est pas observé pour la bande UO, qui est généralement plus faible que la raie U et se rapproche du fond après 6 cm. Une explication possible est que ce signal U aval provient de la lumière diffusée émise dans la partie amont du PFR. Étant donné que la raie d'émission U à 591,54 nm est plus intense que la bande UO à 593,55 nm dans la région en amont, la lumière diffusée pourrait contribuer de manière disproportionnée à l'intensité de la raie U en aval. Étant donné que cette observation peut être sans rapport avec les réactions chimiques, elle mérite une attention particulière lors de la procédure d'optimisation, comme discuté plus loin.

Profils d'émission synthétiques de (1) U et (2) UO (avec axes linéaires et semi-logarithmiques) générés par le modèle 0D \({{\mathrm{UO_x}}}\) non optimisé5 par rapport à (a) les mesures combinées des ensembles de données 2 et 3 et (b) les mesures de l'ensemble de données 1. La transparence des points expérimentaux indique le rapport signal sur fond.

Comme discuté précédemment, le but de l'étape de Monte Carlo du MCGA est d'effectuer une étude préliminaire de l'espace des paramètres du problème et de générer une population de départ pour l'algorithme génétique. En tant que tel, un critère assoupli pour retenir les mécanismes candidats est utilisé, ne nécessitant qu'un coefficient de détermination d'échelle linéaire positif pour tous les profils d'espèces modélisés (c'est-à-dire \(R^2_{lin}>0\)). Ce critère a été satisfait par 8,61% des mécanismes générés à partir de 2,3 millions d'échantillons. Ces 200 000 mécanismes candidats pourraient être utilisés comme population de départ pour la partie GA de l'algorithme. Cependant, pour maintenir le temps d'exécution de l'AG raisonnable, seul un sous-ensemble des mécanismes candidats est utilisé. Cela se fait en triant d'abord les mécanismes selon deux versions de la fonction objectif \(\phi\) précédemment définie dans l'Eq. (14). Les deux fonctions objectif sont : \(\phi _1\), qui inclut tous les termes de la fonction objectif :

et \(\phi _2\), qui exclut les termes \(R^2_{log}\) pour l'émission atomique d'uranium :

où D fait référence aux ensembles de données (y compris diverses conditions de débit d'oxygène), s fait référence à l'espèce et N représente le nombre total de termes dans chaque cas. Cette dernière formulation est incluse ici en raison de la saturation anormale du signal d'émission d'uranium en aval vu sur la figure 8, qui domine la courbe U expérimentale de l'espace logarithmique. Si cet effet est dû à des processus chimiques, nous nous attendons à ce qu'il soit bien adapté par l'algorithme MCGA utilisant \(\phi _1\). Cependant, si le comportement de saturation n'est pas capturé lors de l'utilisation de \(\phi _1\), l'effet peut être de nature non chimique, auquel cas il ne doit pas être utilisé pour contraindre le mécanisme de réaction et \(\phi _2\) doit être utilisé à la place.

La figure 9 trace les mécanismes générés par MC triés selon la forme physique normalisés en utilisant soit \(\phi _1\) ou \(\phi _2\) comme fonction objectif. La forme physique est ici définie comme \(1/\phi\) car notre formulation de \(\phi\) est minimisée pour le meilleur ajustement. Notez que \(\phi\) et la fitness sont des quantités statistiques plutôt que physiques. De plus, les valeurs de fitness ne sont significatives que dans le contexte de la fonction objectif utilisée pour les calculer et ne peuvent pas être comparées entre les fonctions objectifs. Pour plus de commodité, les valeurs de fitness sont ici normalisées par rapport à la valeur de fitness maximale pour la fonction objectif correspondante. Pour les deux fonctions objectives, seuls environ 5 % du mécanisme généré ont des valeurs de fitness à moins de 67 % de la fitness maximale. De plus, seuls quelques centaines de mécanismes (sur 200 000) font partie des 20 % supérieurs de fitness. Ce sous-ensemble de mécanismes supérieurs sert de population initiale pour l'algorithme génétique, comme discuté ci-dessous.

MC a généré des mécanismes candidats triés par valeurs de fitness (\(1/\phi\)) évalués à l'aide de fonctions objectives (a) \(\phi _1\) (Eq. 19) et (b) \(\phi _2\) (Eq. 20), qui ne diffèrent que par l'exclusion de l'émission d'U dans \(\phi _2\), et normalisées par rapport à la valeur de fitness maximale pour chaque fonction objectif.

Après avoir choisi les fonctions objectives cibles, l'algorithme génétique peut maintenant être utilisé pour optimiser le mécanisme de réaction \({{\mathrm{UO_x}}}\) pour les conditions expérimentales utilisées dans cette étude. Pour tester la fiabilité de l'algorithme génétique, nous effectuons quatre optimisations GA distinctes. Ces optimisations représentent des combinaisons de deux paramètres d'entrée avec deux options chacun. Le premier paramètre d'entrée définit la fonction objectif utilisée, les deux options étant \(\phi _1\) et \(\phi _2\) comme décrit par les équations. (19) et (20), respectivement. Le deuxième paramètre d'entrée définit la population initiale utilisée (c'est-à-dire l'ensemble de départ des mécanismes de réaction candidats), les deux options étant les 200 meilleurs mécanismes générés par MC selon \(\phi _1\) ou \(\phi _2\). Par conséquent, deux optimisations utilisent les fonctions objectifs \(\phi _1\) et \(\phi _2\) initialisées avec les 200 mécanismes MC les plus adaptés correspondants selon cette fonction objectif. Nous les désignerons comme les parcours de population optimaux. Les deux autres optimisations échangent plutôt les populations initiales utilisées pour chaque fonction objectif, de sorte que la population GA de départ est sous-optimale dans chaque cas. Autrement dit, l'optimisation \(\phi _1\) est lancée avec la population la plus adaptée selon \(\phi _2\) et vice versa. Cela teste si l'AG peut obtenir de manière fiable la même forme physique optimale quelle que soit la population de départ. Le chevauchement entre les mécanismes les plus adaptés pour les deux fonctions objectives n'est que de 2 mécanismes sur 200. Enfin, pour les populations optimales et sous-optimales, 200 échantillons MC de fitness arbitraire (c'est-à-dire générés sans contraintes/critères de sélection) sont ajoutés à la population initiale pour fournir une diversité supplémentaire pour l'optimisation6. Par conséquent, la population totale de GA de départ pour chaque série se compose de 400 mécanismes (individus).

Résultats d'optimisation GA utilisant \(\phi _1\) (tous les termes \(R^2_{log}\)) initialisés avec (1) des mécanismes optimaux et (2) sous-optimaux générés par MC. Les graphiques (a) montrent la forme physique normalisée moyenne et maximale en fonction de la génération, tandis que les graphiques (b) montrent des distributions de forme physique normalisées triées pour certaines générations.

Les résultats des optimisations GA utilisant \(\phi _1\) (tous les termes \(R^2_{log}\)) avec des populations initiales optimales et sous-optimales sont présentés à la Fig. 10. Les sous-parcelles 1a et 1b illustrent les résultats pour le cas optimal tandis que les sous-parcelles 2a et 2b montrent les résultats correspondants pour le cas sous-optimal. Notez que le nombre total de générations tracées ici (400) n'est pas corrélé avec le nombre de mécanismes dans la population (également 400). Les populations de départ optimales et sous-optimales sont représentées par les courbes bleues (étiquetées 0) en 1b et 2b, respectivement. Le cas optimal (1b) montre clairement la différence de fitness entre les 200 meilleurs mécanismes MC et les 200 mécanismes restants de fitness arbitraire par la discontinuité au milieu. Cette différence n'est pas évidente dans le cas sous-optimal (2b), qui montre la faible fitness globale de la population de départ en dehors de quelques mécanismes de fitness élevés. En regardant l'évolution de l'aptitude maximale et moyenne de la population (1a et 2a), nous constatons que l'aptitude augmente le plus considérablement au cours des 40 à 60 premières générations. Pendant ce temps, la diversité de la population initiale est mise à profit par la reproduction croisée pour localiser rapidement les régions de haute fitness dans l'espace des paramètres. Cela améliore la condition physique globale de la population (comme en témoigne la courbe moyenne) et localise les individus les plus en forme pour déplacer les «élites» de départ (augmentant la condition physique maximale). Cela peut également être observé dans les distributions de fitness de la population (1b et 2b), où le passage d'une forme concave à une forme convexe reflète l'augmentation de la fitness moyenne. Après les 100 premières générations, une grande partie de la population est homogénéisée vers une valeur de fitness similaire, comme en témoigne la forme de plateau des distributions de fitness de la population (1b et 2b). Passé ce stade, les mutations maintiennent un niveau de base de diversité génétique, comme l'indique la persistance d'un sous-ensemble de population à faible fitness qui représente environ 10 à 25% de la population totale. Cela permet également de poursuivre l'exploration des paramètres, conduisant à une amélioration progressive de la condition physique de la population (1a et 2a) au cours des 300 générations restantes. Après 400 générations, les populations de départ optimales et sous-optimales saturent vers une valeur de fitness similaire. Notez que si le mécanisme de départ le plus adapté est partagé par les populations optimales et sous-optimales, les mécanismes finaux les plus adaptés diffèrent. Curieusement, le cas de départ sous-optimal possède une meilleure forme physique après 400 générations que le cas optimal, mais il s'agit probablement d'un événement stochastique.

Résultats d'optimisation GA utilisant \(\phi _2\) (pas de \(R^2_{log}\) pour U) initialisé avec (1) des mécanismes optimaux et (2) sous-optimaux générés par MC. Les graphiques (a) montrent la forme physique normalisée moyenne et maximale en fonction de la génération, tandis que les graphiques (b) montrent des distributions de forme physique normalisées triées pour certaines générations.

La figure 11 présente le même ensemble de tracés pour les optimisations GA utilisant \(\phi _2\) (\(R^2_{log}\) exclu pour U). L'évolution de cette population GA est très similaire à celle du cas \(\phi _1\) précédemment examiné. Malgré la population initialement moins en forme (à la fois en moyenne et en maximum), la course sous-optimale (2a et 2b) arrive à nouveau à une fitness similaire à la course optimale (1a et 1b). Dans ce cas, la fitness finale du cas optimal est plus élevée que celle du scénario sous-optimal, mais cet écart se comblerait probablement à mesure que l'évolution se poursuivrait. Quoi qu'il en soit, les quatre exécutions ci-dessus montrent que l'optimisation GA fonctionne de manière fiable et produit des valeurs de fitness similaires quelle que soit la population initiale utilisée. Cependant, les cas sous-optimaux et optimaux ne convergent pas vers une valeur de fitness identique pour l'une ou l'autre des fonctions objectives sur les 400 générations observées ici. Cela peut être dû au nombre limité de générations réalisées et au choix des propriétés GA dictant l'exploration des paramètres et les taux de convergence. Ce comportement peut également être inhérent au problème d'optimisation lui-même en raison de la gamme limitée de conditions fournies par les données contraignantes et de la non orthogonalité de certains canaux de réaction. Ces considérations sont illustrées en examinant les coefficients de taux prédits par les populations optimisées, comme indiqué ci-dessous.

Coefficients de vitesse de réaction moyens (\({\bar{k}}\)) des populations GA optimisées à l'aide de (1) \(\phi _1\) et (2) \(\phi _2\) initialisées avec (a) des mécanismes optimaux et (b) sous-optimaux générés par MC. Les moyennes sont calculées à partir des mécanismes se situant à moins de 0,1 % de la meilleure forme physique (mécanismes \(\sim 300\) dans chaque cas) pour les réactions qui satisfont \({\bar{k}}>\sigma _k\) (où \(\sigma _k\) est l'écart type de k).

La figure 12 trace les coefficients de vitesse moyens des principaux canaux de réaction à partir des quatre populations optimisées. Les réactions dominantes sont identifiées en examinant la variation statistique des coefficients de taux dans les populations optimisées. Seules les réactions pour lesquelles l'écart type du coefficient de vitesse ne dépasse pas la moyenne sont incluses (sur la base de mécanismes à moins de 0,1 % de la meilleure forme physique). Cette condition est satisfaite pour au plus sept des 12 canaux de réactions optimisés, ce qui signifie que les réactions restantes ne sont pas bien contraintes par les données actuelles. Cependant, quatre réactions (R1, R4, R6 et R11 dans le tableau 2) sont systématiquement contraintes, comme indiqué par leur apparition dans tous les cas d'optimisation. Les trois réactions restantes (R5, R9 et R10) n'apparaissent que dans certaines des populations optimisées et ne sont donc que partiellement contraintes. Cela est probablement dû à la non-orthogonalité des différentes voies de réaction. Autrement dit, des canaux comme \({\mathrm{U + O}}\), \({\mathrm{U + OH}}\) et \({\mathrm{U + H_2O}}\) ou \({\mathrm{UO + O}}\), \({\mathrm{UO + OH}}\) et \({\mathrm{UO + H_2O}}\) peuvent se compenser puisqu'ils effectuent la même opération sur le espèces contraintes (c'est-à-dire ajouter O à U ou UO). Bien que les différentes conditions de débit d'oxygène de l'ensemble de données 1 limitent quelque peu ce comportement, la plage de concentrations \(\hbox {O}_2\) est limitée et seuls les emplacements en aval (> 3 cm) sont couverts. Ces réactions pourraient être mieux limitées à l'avenir en effectuant des mesures en amont sur une gamme plus large de conditions \(\hbox {O}_2\) ou avec des concentrations réduites en \(\hbox {H}_2\)O (c'est-à-dire en utilisant un nébuliseur désolvatant). L'ensemble de données actuel est également limité dans la plage de températures et de taux de refroidissement couverts, ce qui empêche également la localisation d'un véritable optimum global. Ceci est lié à la dépendance à la température de l'expression de la vitesse d'Arrhenius, car des vitesses de réaction similaires peuvent être obtenues en utilisant différentes combinaisons de coefficients si la réaction est active sur une plage de température limitée dans le système. Sur la base des observations ci-dessus, le manque de convergence vers un optimum global singulier semble provenir davantage des données contraignantes limitées que d'une lacune de l'algorithme génétique. Même ainsi, quatre réactions (R1, R4, R6 et R11 dans le tableau 2) sont toujours bien contraintes par l'optimisation actuelle, démontrant la fiabilité de la méthode MCGA et mettant en évidence les canaux de réaction dominants pour la formation d'UO dans le PFR.

Profils d'émission synthétiques de (1) U et (2) UO (avec axes linéaires et semi-log) générés par un mécanisme optimisé GA utilisant \(\phi _1\) (tous les termes \(R^2_{log}\)) comparés à (a) les mesures combinées des ensembles de données 2 et 3 et (b) les mesures de l'ensemble de données 1.

Avant d'effectuer une analyse détaillée des résultats de l'optimisation, nous devons décider quelle fonction objectif utiliser en vérifiant si la saturation du signal U aval est capturée par le mécanisme optimisé \(\phi _1\). La figure 13 montre les profils d'émission synthétiques d'U et d'UO produits par ce mécanisme par rapport à l'ensemble de données expérimentales complet. Bien que le résultat \(\phi _1\) MCGA produise un meilleur ajustement des données en amont sur le mécanisme non optimisé (Fig. 8), le comportement en aval (> 3 cm) de U est mal capturé. Bien que les valeurs R\(^2\) semblent adéquates (en dehors du cas 1 de l'ensemble de données 0 \(\hbox {O}_2\), l'inspection visuelle révèle que ni la diminution de 3 à 5 cm ni le comportement de saturation de 5 à 8 cm ne correspondent bien à aucun ensemble de données. Cela suggère que la saturation peut être due à un effet non chimique qui n'est pas pris en compte dans le traitement 0D actuel du PFR. Comme discuté précédemment, ce comportement est probablement causé par la diffusion optique de la forte raie d'émission U en amont. Par conséquent, nous la considérerons comme une contrainte invalide pour le problème d'optimisation actuel et concentrerons notre analyse ultérieure sur le mécanisme optimisé \(\phi _2\). Comme indiqué dans la section suivante, le résultat \ (\ phi _2 \) capture la diminution de 3 à 5 cm de l'intensité U tout en ne tenant pas compte de la saturation du signal ultérieure due à l'exclusion du terme log U.

L'ensemble final de canaux de réaction \({\mathrm{UO_x}}\) optimisés et les coefficients de vitesse correspondants obtenus par le MCGA sont présentés dans le tableau 5. Le tableau comprend également les canaux de réaction pertinents du mécanisme \({\mathrm{UO_x}}\) précédemment construit5. Seuls les coefficients de débit bien contraints sont répertoriés pour le mécanisme optimisé MCGA (voir Fig. 12). Les coefficients de vitesse pour les deux mécanismes sont tracés à la Fig. 14 pour l'intervalle de température \(3000 \le T \le 4500\) K, qui représente la plage sur laquelle l'optimisation est effectuée (voir Fig. 3). Étant donné que le mécanisme de réaction non optimisé ne tenait pas compte des interactions avec les molécules \({\mathrm{H_xO_y}}\) , seuls quelques canaux de réaction sont présents dans les deux mécanismes. Ceci est exacerbé par l'élimination des canaux de réaction \({\mathrm{UO_x + O_2}}\) par l'optimisation MCGA en raison de l'abondance relative d'OH dans le flux PFR, comme indiqué précédemment sur la Fig. 7. Par conséquent, seules deux voies de réaction (R6 et R11) peuvent être directement comparées entre les deux mécanismes. D'après la figure 14, nous pouvons voir que le coefficient de vitesse optimisé MCGA pour la première réaction (R6) est au moins d'un ordre de grandeur inférieur à l'estimation de la littérature. Pour cette dernière réaction (R11), la différence est encore plus grande, autour de quatre ordres de grandeur. En général, les estimations non optimisées se situent à un ordre de grandeur ou plus au-dessus des coefficients de taux optimisés. Il faut s'y attendre étant donné que les estimations de la littérature consistent en grande partie en des estimations de taux de collision de sphères dures de premier ordre sans barrière, qui sont la limite supérieure théorique des taux de réaction.

Comparaison des coefficients de débit optimisés de la littérature5 (plein) et MCGA (pointillés) pour les canaux (a) \({\mathrm{U + H_xO_y}}\) et (b) \({\mathrm{UO + H_xO_y}}\) dans le tableau 5. Les coefficients de débit sont calculés sur la plage de température donnée à l'aide des paramètres d'Arrhenius du tableau et de l'équation. (2).

Le tableau 5 et la figure 14 fournissent un aperçu supplémentaire des similitudes et des différences entre les mécanismes de réaction optimisés et précédemment construits. Tout d'abord, nous notons que la voie \({\mathrm{U+O}}\) dans le mécanisme optimisé est dominée par le canal d'ionisation associatif (R11) sur la réaction d'association moléculaire (R1). Ce comportement est en accord partiel avec une étude précédente31, qui mesurait une section efficace beaucoup plus élevée pour la réaction d'ionisation associative (R11) par rapport à la réaction d'association moléculaire (R1). L'étude a suggéré que le canal d'association moléculaire est effectivement fermé en raison de la prédominance de la voie d'ionisation associative. Dans le mécanisme optimisé, cependant, les deux canaux ne diffèrent que d'environ un ordre de grandeur. De plus, le coefficient de taux R11 optimisé suggère une section efficace beaucoup plus faible que celle observée dans l'étude susmentionnée. La voie optimisée s'avère également être effectivement sans barrière, ce qui est cohérent avec l'étude.

Ensuite, nous notons une certaine similitude dans les canaux R4 optimisés et R3 non optimisés, notamment en ce qui concerne l'énergie d'activation51 (\({\mathrm{E_A/R}}\) dans le tableau 5). Cela peut indiquer que le mécanisme d'abstraction pour U entrant en collision avec OH ou \(\hbox {O}_2\) se déroule de la même manière. En comparant davantage les deux réactions, nous voyons que le taux de collision optimisé est inférieur d'environ deux ordres de grandeur à l'estimation de la sphère dure. Néanmoins, le coefficient de vitesse optimisé global est toujours élevé par rapport aux autres réactions optimisées, comme le montre la figure 14.

En passant à R6, nous observons une grande différence d'énergie d'activation entre le résultat optimisé et notre estimation précédente. L'énergie d'activation dans le coefficient de vitesse non optimisé provient d'une estimation d'Eyring ajustée selon un calcul semi-empirique analogue pour un mécanisme d'oxydation d'Al2,5. La valeur de barrière non optimisée est essentiellement un sous-produit de l'ajustement ci-dessus, car l'estimation d'Eyring elle-même est sans barrière. La forte dépendance à la température de ce canal pour tous les résultats MCGA (voir Fig. 12) suggère que la barrière d'activation importante du taux optimisé est physiquement significative. On observe que le taux de collision pour ce canal est proche de la limite de la sphère dure, ce qui compense l'effet de la grande barrière d'activation. Ainsi, le coefficient de vitesse global se situe dans un ordre de grandeur du canal non optimisé à 4500 K et avoisine la même valeur que R4 à cette température, comme le montre la Fig. 14.

Ensuite, la réaction d'abstraction R9 s'avère être effectivement sans barrière et a un faible taux de collision global comparable au canal d'ionisation associatif R11. L'absence de barrière d'activation pour le canal d'abstraction R8 non optimisé analogue est simplement due à un manque d'informations dans la littérature, de sorte qu'une comparaison définitive entre R8 et R9 n'est pas tentée ici.

Enfin, le canal de réaction R10 restant a un taux de collision similaire à R4, mais avec une énergie d'activation plus élevée. Malgré le coefficient de vitesse relativement élevé, cette voie de réaction ne devrait être importante que dans la partie aval du PFR en raison de sa dépendance à \(\hbox {H}_2\)O, qui n'atteint des concentrations comparables à OH qu'à \(> 7\) cm où la température descend en dessous de 3000 K.

Coefficients de sensibilité de premier ordre pour (a) U et (b) UO en utilisant le mécanisme optimisé MCGA pour les conditions des ensembles de données 2 et 3.

Pour analyser où les canaux de réaction dominants sont les plus actifs dans le flux PFR, nous avons calculé les coefficients de sensibilité du premier ordre pour le système d'équations cinétiques chimiques52. Les résultats pour U et UO pour les conditions de jeu de données 2 et 3 sont indiqués sur la figure 15. En plus des canaux optimisés, ces parcelles incluent les réactions de formation fixe \ (\ mathrm {UO_2 / UO_3}} \) R13 et R14 à partir de la table. En aval (> 3 cm). Les tracés de sensibilité UO montrent que l'évolution UO en amont est relativement insensible au mécanisme de réaction. Plus en aval, UO est le plus sensible à R6 et R13. Les canaux optimisés restants du tableau 5, R9 et R10, semblent apporter des ajustements plus fins à l'évolution de l'UO en aval. Dans l'ensemble, le comportement en amont du mécanisme est limité principalement par les données U, tandis que les mesures U et UO jouent un rôle dans la contrainte du comportement en aval.

Profils modélisés de densité numérique d'espèces sélectionnées prédits par le mécanisme \({{\mathrm{UO_x}}}\) optimisé MCGA.

Les densités de nombre d'espèces produites par le mécanisme de réaction optimisé sont illustrées à la Fig. 16. Par rapport aux résultats non optimisés de la Fig. 7, le mécanisme MCGA produit nettement moins d'ions uranium dans la région de la bobine en raison du taux d'ionisation associatif R11 plus faible. Bien que cela se traduise par une population plus élevée d'U neutre, il reste environ un ordre de grandeur moins abondant que l'UO. Un rapport d'intensité de pic en amont de \(I_{U/UO}\approx 6\) est obtenu en utilisant le mécanisme optimisé, ce qui correspond à \(n_{UO}/n_U\approx 20\). Une valeur similaire est obtenue pour toutes les autres optimisations MCGA. Bien que cela se situe dans la plage de valeurs \(I_{U/UO}\) autorisées par le terme de pénalité, un étalonnage approfondi des densités absolues serait nécessaire pour valider ce résultat.

Profils d'émission synthétiques de (1) U et (2) UO (avec axes linéaires et semi-log) générés par le mécanisme \({{\mathrm{UO_x}}}\) optimisé par GA par rapport à (a) les ensembles de données combinés 2 et 3 mesures et (b) les mesures de l'ensemble de données 1.

Enfin, pour examiner l'adéquation du mécanisme optimisé, les profils finaux d'émission synthétique d'uranium et les ensembles de données contraignants sont tracés à la Fig. 17. Comme prévu, le mécanisme MCGA produit un meilleur ajustement de la plupart des données expérimentales par rapport aux résultats non optimisés de la Fig. Cela découle de l'utilisation de la fonction objectif \(\phi _2\), qui exclut le comportement de saturation comme contrainte. Cependant, le mécanisme optimisé produit un excellent ajustement des données U en amont (\(<5\) cm). Le mécanisme produit également un bon ajustement des données UO en amont et en aval. Comme indiqué précédemment, l'ajustement partiel des données UO en amont fournies par l'ensemble de données 3 est limité par le profil de température représentatif approximatif utilisé ici.

Dans ce travail, nous avons étudié une approche basée sur MCGA pour calibrer un mécanisme de réaction \({{\mathrm{UO_x}}}\) en utilisant des mesures d'émission d'expériences PFR. La sélection des canaux de réaction cibles et de leurs coefficients de vitesse potentiels a été effectuée à l'aide d'un ensemble limité d'hypothèses a priori. La cohérence entre le modèle PFR 0D et les expériences a été atteinte en utilisant un profil de température représentatif en accord avec les mesures de température disponibles. Les étapes d'échantillonnage de Monte Carlo et d'algorithme génétique ont été utilisées pour explorer et optimiser l'espace des paramètres du problème, respectivement, permettant un raffinement vers un maximum de fitness. Le mécanisme de réaction \({{\mathrm{UO_x}}}\) optimisé résultant a été analysé, mettant en évidence quatre canaux de réaction dominants qui étaient systématiquement contraints dans toutes les optimisations tentées et trois canaux supplémentaires qui étaient au moins partiellement contraints. Parmi ces canaux, il convient de noter l'implication du radical OH, qui n'était pas auparavant pris en compte dans le mécanisme non optimisé. Le mécanisme optimisé a prédit une cinétique plus lente pour la formation d'U et d'UO (la cinétique d'U étant plus impactée que l'UO) par rapport au mécanisme non optimisé, produisant des taux qui étaient de 1 à 2 ordres de grandeur ou plus inférieurs aux estimations précédentes de la sphère dure. Une caractéristique notable du mécanisme optimisé est le rapport de ramification inférieur entre l'ionisation associative et les voies neutres du canal \({\mathrm{U + O}}\), ce qui suggère que la voie neutre n'est pas éliminée au profit de l'ionisation associative comme indiqué dans une étude précédente31.

Dans l'ensemble, cette étude démontre la viabilité de l'utilisation d'une approche MCGA pour optimiser les coefficients de vitesse cinétique chimique pour un modèle PFR 0D à l'aide de mesures d'émission optique. Cependant, en raison de diverses limitations de l'étude actuelle, nous mettons fortement en garde contre l'utilisation des coefficients de taux optimisés résultants tels quels jusqu'à ce qu'une optimisation et une validation supplémentaires soient effectuées. Un mécanisme de réaction bien validé pourrait être atteint via la méthode MCGA après avoir incorporé divers raffinements dans un travail futur, comme indiqué ci-dessous.

Premièrement, l'importance des molécules OH pour l'oxydation de l'uranium dans le PFR suggère que la fugacité de l'oxygène et de l'hydrogène devrait être modifiée dans les études futures. Ceci peut être réalisé en utilisant un nébuliseur de désolvatation, qui élimine la majeure partie de l'eau de la solution d'analyte avant l'introduction dans le plasma. Deuxièmement, des mesures de température supplémentaires dans la région de la bobine du PFR amélioreraient la cohérence entre le modèle et les expériences, ce qui est important en raison de la sensibilité de l'évolution chimique simulée à l'historique de la température. La cohérence peut être encore affinée en considérant également les effets du mélange et de la diffusion radiale en aval des entrées d'écoulement. Troisièmement, la bande UO à 593,55 nm pourrait être mieux résolue en utilisant un réseau de spectromètre plus élevé. Cela aiderait à éliminer une certaine incertitude dans l'étalonnage introduit en raison du traitement de la bande UO comme la tête 0–0 de la bande. Quatrièmement, l'ensemble de données expérimentales pourrait être élargi pour inclure des informations sur la formation d'oxyde d'uranium plus élevé (c'est-à-dire \({\mathrm{UO_2}}\) et \({\mathrm{UO_3}}\)), en utilisant la spectroscopie infrarouge à transformée de Fourier (FTIR), par exemple. Cependant, cela nécessiterait d'étendre les calculs d'émissions atomiques et diatomiques utilisés ici à des molécules plus grosses, ce qui est potentiellement difficile pour les espèces d'uranium. Cinquièmement, la modélisation pourrait être affinée en effectuant des calculs détaillés de transport de rayonnement pour quantifier l'impact de l'auto-absorption et de la diffusion sur l'intensité d'émission simulée dans tout le PFR. Enfin, un moyen relativement simple d'améliorer le mécanisme généré par MCGA consiste à utiliser un ensemble de données plus important pour l'optimisation. Cela comprend la mesure des émissions dans les régions en amont et en aval sur une plage plus large de débits, de températures et de concentrations d'analytes.

Les ensembles de données générés pendant et/ou analysés pendant l'étude en cours sont disponibles auprès de l'auteur correspondant sur demande raisonnable.

Kautz, EJ et al. Spectroscopie optique et modélisation de l'oxydation en phase gazeuse de l'uranium : avancées et perspectives. Spectrochim. Acta B 185, 106283. https://doi.org/10.1016/j.sab.2021.106283 (2021).

Article CAS Google Scholar

Catoire, L., Legendre, J.-F. & Giraud, M. Modèle cinétique pour la combustion d'un accélérateur à bélier sensibilisé à l'aluminium. J. Propul. Puissance 19, 196–202. https://doi.org/10.2514/2.6118 (2003).

Article CAS Google Scholar

Huang, Y., Risha, GA, Yang, V. & Yetter, RA Effet de la taille des particules sur la combustion de la poussière de particules d'aluminium dans l'air. Combustion. Flamme 156, 5–13. https://doi.org/10.1016/j.combustflame.2008.07.018 (2009).

Article CAS Google Scholar

Starik, AM, Kuleshov, PS, Sharipov, AS, Titova, NS & Tsai, CJ Analyse numérique de la combustion du nanoaluminium dans la vapeur. Combustion. Flamme 161, 1659–1667. https://doi.org/10.1016/j.combustflame.2013.12.007 (2014).

Article CAS Google Scholar

Finko, MS et al. Un modèle de formation précoce d'oxydes moléculaires d'uranium dans des plasmas ablatés par laser. J.Phys. Appl. Phys. 50, 485201. https://doi.org/10.1088/1361-6463/aa92f5 (2017).

Article CAS Google Scholar

Berkemeier, T. et al. Note technique : Algorithme génétique de Monte Carlo (MCGA) pour l'analyse de modèles de cinétique chimique multiphase afin de déterminer les coefficients de transport et de vitesse de réaction à l'aide de plusieurs ensembles de données expérimentales. Atmos. Chim. Phys. 17, 8021–8029. https://doi.org/10.5194/acp-17-8021-2017 (2017).

Article ADS CAS Google Scholar

Heaven, MC, Barker, BJ & Antonov, IO Spectroscopie et structure des liaisons actinides les plus simples. J.Phys. Chim. A 118, 10867–10881. https://doi.org/10.1021/jp507283n (2014).

Article CAS PubMed Google Scholar

Hartig, KC, Ghebregziabher, I. & Jovanovic, I. Détection à distance de l'uranium et de ses isotopes par spectrométrie isotopique moléculaire d'ablation laser à filament femtoseconde. Sci. Rép. 7, 1–9. https://doi.org/10.1038/srep43852 (2017).

Article Google Scholar

Hartig, KC, Harilal, SS, Phillips, MC, Brumfield, BE & Jovanovic, I. Évolution du monoxyde d'uranium dans les plasmas d'uranium induits par laser femtoseconde. Opter. Express 25, 11477–11490. https://doi.org/10.1364/OE.25.011477 (2017).

Article ADS CAS PubMed Google Scholar

Weisz, DG et al. Formation de \(^{238}\)U\(^{16}\)O et \(^{238}\)U\(^{18}\)O observée par spectroscopie d'émission résolue en temps après ablation laser. Appl. Phys. Lett. 111, 034101. https://doi.org/10.1063/1.4991824 (2017).

Article ADS CAS Google Scholar

Harilal, SS, Brumfield, BE, Glumac, N. & Phillips, MC Élucidation des caractéristiques spectrales du monoxyde d'uranium à partir d'un plasma produit par laser. Opter. Express 26, 20319–20330. https://doi.org/10.1364/OE.26.020319 (2018).

Article ADS CAS PubMed Google Scholar

Skrodzki, PJ et al. Suivi de la formation d'oxydes dans les plasmas d'uranium produits par laser. Opter. Lett. 43, 5118–5121. https://doi.org/10.1364/OL.43.005118 (2018).

Article ADS CAS PubMed Google Scholar

Weerakkody, EN et al. Formation résolue dans le temps d'oxydes d'uranium et de silicium suite à l'ablation laser de \(\text{ U}_3\)\(\text{ Si}_2\). Spectrochim. Acta B 170, 105925. https://doi.org/10.1016/j.sab.2020.105925 (2020).

Article CAS Google Scholar

Weerakkody, EN & Glumac, NG Spectroscopie d'absorption quantitative de plasmas produits par laser. J.Phys. Appl. Phys. 54, 125201. https://doi.org/10.1088/1361-6463/abd210 (2021).

Article ADS CAS Google Scholar

Koroglu, B. et al. Évolution chimique en phase gazeuse des oxydes d'uranium, d'aluminium et de fer. Sci. Rep. 8, 10451. https://doi.org/10.1038/s41598-018-28674-6 (2018).

Article ADS CAS PubMed PubMed Central Google Scholar

Koroglu, B. et al. Etude expérimentale de la volatilité de l'uranium lors de la condensation de vapeur. Anal. Chim. 92, 6437–6445. https://doi.org/10.1021/acs.analchem.9b05562 (2020).

Article CAS PubMed Google Scholar

Finko, MS & Curreli, D. Simulations de la dynamique du panache de plasma d'uranium dans l'oxygène atmosphérique produit par ablation laser femtoseconde. Phys. Plasmas 25, 083112. https://doi.org/10.1063/1.5034470 (2018).

Article ADS CAS Google Scholar

Koroglu, B. et al. Réacteur à flux plasma pour la surveillance en régime permanent de processus physiques et chimiques à haute température. Rev. Sci. Instrument. 88, 093506. https://doi.org/10.1063/1.5001346 (2017).

Article ADS CAS PubMed Google Scholar

Devoto, RS Coefficients de transport de l'argon partiellement ionisé. Phys. Fluides 10, 354. https://doi.org/10.1063/1.1762115 (1967).

Article ADS CAS Google Scholar

Devoto, RS Coefficients de transport de l'argon ionisé. Phys. Fluides 16, 616. https://doi.org/10.1063/1.1694396 (1973).

Article ADS CAS Google Scholar

Tadi, M. & Yetter, R. A. Evaluation of the rate constants in chemical reactions. Int. J. Chem. Kinet. 30, 151–159. 3.0.CO;2-R">https://doi.org/10.1002/(SICI)1097-4601(1998)30:2<151::AID-KIN7>3.0.CO;2-R (1998).

3.0.CO;2-R" data-track-action="article reference" href="https://doi.org/10.1002%2F%28SICI%291097-4601%281998%2930%3A2%3C151%3A%3AAID-KIN7%3E3.0.CO%3B2-R" aria-label="Article reference 21" data-doi="10.1002/(SICI)1097-4601(1998)30:23.0.CO;2-R">Article CAS Google Scholar

Lee, HS & Litzinger, TA Etude cinétique chimique de la décomposition HAN. Combustion. Flamme 135, 151–169. https://doi.org/10.1016/S0010-2180(03)00157-3 (2003).

Article CAS Google Scholar

Maraun, D. et al. Identification des constantes de vitesse et des spectres d'absorption non observables dans la dynamique des réactions biochimiques non linéaires. Int. J. Bifurc. Chaos 14, 2081–2092. https://doi.org/10.1142/S0218127404010473 (2004).

Article CAS MATH Google Scholar

Lemes, NHT, Borges, E. & Braga, JP Un algorithme général pour résoudre des problèmes inverses linéaires et non linéaires. J. Braz. Chim. Soc. 18, 1342–11347. https://doi.org/10.1590/S0103-50532007000700008 (2007).

Article CAS Google Scholar

Capitelli, M., Ferreira, CM, Gordiets, BF & Osipov, AI Plasma Kinetics in Atmospheric Gases, vol. 31 de la série Springer sur la physique atomique, optique et des plasmas (Springer, 2000).

Varga, T. et al. Optimisation d'un mécanisme de combustion d'hydrogène par mesures directes et indirectes. Proc. Combustion. Inst. 35, 589-596. https://doi.org/10.1016/j.proci.2014.06.071 (2015).

Article CAS Google Scholar

Zhang, Y., Mathieu, O., Petersen, EL, Bourque, G. & Curran, HJ Un modèle thermique pour l'ablation laser pulsée nanoseconde de l'aluminium. Combustion. Flamme 182, 122–141. https://doi.org/10.1063/1.4995972 (2017).

Article CAS Google Scholar

Starikovskaia, SM, Starikovskii, AY & Zatsepin, DV Oxydation de l'hydrogène dans un mélange stoechiométrique hydrogène-air dans l'onde d'ionisation rapide. Combustion. Modèle de théorie. 5, 97–129. https://doi.org/10.1088/1364-7830/5/1/306 (2001).

Article ADS CAS MATH Google Scholar

Murakami, T., Niemi, K., Gans, T., O'Connell, D. & Graham, WG Cinétique chimique et espèces réactives dans les plasmas hélium-oxygène à pression atmosphérique avec des impuretés d'air humide. Sources de plasma Sci. Technol. 22, 015003. https://doi.org/10.1088/0963-0252/22/1/015003 (2013).

Article ADS CAS Google Scholar

Schroter, S. et al. Cinétique chimique dans un plasma d'hélium à pression atmosphérique contenant de l'humidité. Phys. Chim. Chim. Phys. 20, 24263–24286. https://doi.org/10.1039/c8cp02473a (2018).

Article CAS PubMed PubMed Central Google Scholar

Fite, WL, Lo, HH & Irving, P. Ionisation associative dans les collisions U + O et U + \(\text{ O}_2\). J. Chem. Phys. 60, 1236-1250. https://doi.org/10.1063/1.1681188 (1974).

Article ADS CAS Google Scholar

Carr, RW Modélisation des réactions chimiques Vol. 42, 43–99 (Elsevier, 2007).

Google Scholar

Rodriguez, K. et al. Agrégation en phase vapeur de nanoparticules d'oxyde de cérium dans un plasma à refroidissement rapide. ACS Terre Espace Chim. 6, 1997–2008. https://doi.org/10.1021/acsearthspacechem.2c00073 (2022).

Article ADS CAS Google Scholar

Pancheshnyi, S., Eismann, B., Hagelaar, GJM et Pitchford, LC Code informatique ZDPlasKin. http://www.zdplaskin.laplace.univ-tlse.fr (2008).

McQuarrie, DA Mécanique statistique (Harper & Row, 1976).

MATH Google Scholar

Palmer, BA, Keller, RA et Engleman, R. Un atlas des intensités d'émission d'uranium dans une décharge à cathode creuse. Technologie. Rep. LA-8251-MS, Laboratoire scientifique de Los Alamos, Los Alamos, NM (1980).

McBride, BJ & Gordon, S. Programme informatique pour le calcul et l'ajustement des fonctions thermodynamiques. Technologie. Rep. NASA-RP-1271, Administration nationale de l'aéronautique et de l'espace, Washington, DC (1992).

Blaise, J. & Wyart, J.-F. Niveaux d'énergie de constantes sélectionnées et spectres atomiques des actinides. http://web2.lac.u-psud.fr/lac/Database/Contents.html.

Lewis, GN & Merle, R. Thermodynamique 2e éd. (McGraw-Hill, 1961).

Google Scholar

Huber, KP & Herzberg, G. Spectres moléculaires et structure moléculaire (Springer, 1979).

Réserver Google Scholar

Heaven, MC, Nicolai, J.-P., Riley, SJ & Parks, EK Spectres électroniques résolus en rotation pour le monoxyde d'uranium. Chim. Phys. Lett. 119, 229–233. https://doi.org/10.1016/0009-2614(85)80065-8 (1985).

Article ADS CAS Google Scholar

Kaledin, LA, McCord, JE & Heaven, MC Spectroscopie laser de l'UO : Caractérisation et affectation des états dans la gamme de 0 à 3 eV, avec une comparaison avec la structure électronique de ThO. J. Mol. Spectrosc. 164, 27–65. https://doi.org/10.1006/jmsp.1994.1055 (1994).

Article ADS CAS Google Scholar

Kaledin, LA & Heaven, MC Spectroscopie électronique d'UO. J. Mol. Spectrosc. 185, 1–7. https://doi.org/10.1006/jmsp.1997.7383 (1997).

Article ADS CAS PubMed Google Scholar

Konings, RJM et al. Les propriétés thermodynamiques des éléments f et de leurs composés. Partie 2. Les oxydes de lanthanides et d'actinides. J.Phys. Chim. Réf. Données 43, 013101. https://doi.org/10.1063/1.4825256 (2014).

Article ADS CAS Google Scholar

Doidge, PS Un recueil et un examen critique des forces de l'oscillateur de la ligne de résonance de l'atome neutre pour l'analyse de l'absorption atomique. Spectrochim. Acta B. 50, 209–263. https://doi.org/10.1016/0584-8547(94)01199-B (1995).

Annonces d'article Google Scholar

Mandal, PK et al. Durées de vie radiatives, rapports de branchement et probabilités absolues de transition de l'uranium atomique par mesures de photoionisation retardée. Appl. Phys. B 116, 407–413. https://doi.org/10.1007/s00340-013-5711-z (2014).

Article ADS CAS Google Scholar

Tyagi, R., Zhang, Z. & Pitzer, RM Spectre électronique des molécules UO et UO\(^+\). J.Phys. Chim. A. 118, 11758. https://doi.org/10.1021/jp505722y (2014).

Article CAS PubMed Google Scholar

Elliott, L. et al. Algorithmes génétiques pour l'optimisation des mécanismes réactionnels de la cinétique chimique. Programme. Combustion d'énergie. Sci. 30, 297–328. https://doi.org/10.1016/j.pecs.2004.02.002 (2004).

Article CAS Google Scholar

Sikalo, N., Hasemann, O., Schulz, C., Kempf, A. & Wlokas, I. Une méthode basée sur un algorithme génétique pour l'optimisation des mécanismes de cinétique réduite. Int. J. Chem. Kinet. 47, 695–723. https://doi.org/10.1002/kin.20942 (2015).

Article CAS Google Scholar

Cantú-Paz, E. Une enquête sur les algorithmes génétiques parallèles. Technologie. Rep. 97003, Université de l'Illinois à Urbana-Champaign, Champaign, IL (1997).

Zeng, XL, Huang, SQ & Ju, XH Etude ab initio sur la réaction de l'uranium avec l'oxygène. J. Radioanal. Nucl. Chim. 298, 481–484. https://doi.org/10.1007/s10967-013-2442-x (2013).

Article CAS Google Scholar

Turányi, T. Analyse de sensibilité des systèmes cinétiques complexes. Outils et applications. J. Math. Chim. 5, 203–248. https://doi.org/10.1007/BF01166355 (1990).

Article MathSciNetGoogle Scholar

Télécharger les références

Ce travail a été réalisé en partie sous les auspices du DoE américain par le Lawrence Livermore National Laboratory sous le contrat DE-AC5207NA27344. Le financement a été fourni par la subvention de recherche et développement dirigée par le laboratoire (LDRD) 20-SI-006 (K. Knight, PI). Le projet ou l'effort décrit a été ou est parrainé par le ministère de la Défense, Agence de réduction des menaces pour la défense dans le cadre de l'Alliance de recherche universitaire sur la science des matériaux dans les environnements extrêmes, HDTRA1-20-2-0001. Le contenu de l'information ne reflète pas nécessairement la position ou la politique du gouvernement fédéral, et aucune approbation officielle ne devrait en être déduite.

Laboratoire national Lawrence Livermore, Livermore, Californie, 94550, États-Unis

Mikhail Finko, Batikan Koroglu, Kate E. Rodriguez, Timothy P. Rose, Jonathan C. Crowhurst, Harry B. Radousky et Kim B. Knight

Département de génie nucléaire, plasma et radiologique, Université de l'Illinois Urbana-Champaign, Champaign, IL, 61820, États-Unis

David Curreli

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

Vous pouvez également rechercher cet auteur dans PubMed Google Scholar

MF a conçu le travail, effectué la modélisation et l'analyse des données, préparé les figures et rédigé le projet original. MF, BK et KR ont mené les expériences et acquis les données. MF, BK, TR et JC ont interprété les données et édité le manuscrit. DC, HR et KK ont aidé à la conception, obtenu des fonds et administré le travail. Tous les auteurs ont examiné le manuscrit.

Correspondance avec Mikhail Finko.

Les auteurs ne déclarent aucun intérêt concurrent.

Springer Nature reste neutre en ce qui concerne les revendications juridictionnelles dans les cartes publiées et les affiliations institutionnelles.

Libre accès Cet article est sous licence Creative Commons Attribution 4.0 International, qui permet l'utilisation, le partage, l'adaptation, la distribution et la reproduction sur n'importe quel support ou format, à condition que vous accordiez le crédit approprié à l'auteur ou aux auteurs originaux et à la source, fournir un lien vers la licence Creative Commons et indiquer si des modifications ont été apportées. Les images ou tout autre matériel de tiers dans cet article sont inclus dans la licence Creative Commons de l'article, sauf indication contraire dans une ligne de crédit au matériel. Si le matériel n'est pas inclus dans la licence Creative Commons de l'article et que votre utilisation prévue n'est pas autorisée par la réglementation légale ou dépasse l'utilisation autorisée, vous devrez obtenir l'autorisation directement du détenteur des droits d'auteur. Pour voir une copie de cette licence, visitez http://creativecommons.org/licenses/by/4.0/.

Réimpressions et autorisations

Finko, M., Koroglu, B., Rodriguez, KE et al. Optimisation stochastique d'un mécanisme de réaction d'oxyde d'uranium à l'aide de mesures de réacteur à flux de plasma. Sci Rep 13, 9293 (2023). https://doi.org/10.1038/s41598-023-35355-6

Télécharger la citation

Reçu : 18 octobre 2022

Accepté : 16 mai 2023

Publié: 07 juin 2023

DOI : https://doi.org/10.1038/s41598-023-35355-6

Toute personne avec qui vous partagez le lien suivant pourra lire ce contenu :

Désolé, aucun lien partageable n'est actuellement disponible pour cet article.

Fourni par l'initiative de partage de contenu Springer Nature SharedIt

En soumettant un commentaire, vous acceptez de respecter nos conditions d'utilisation et nos directives communautaires. Si vous trouvez quelque chose d'abusif ou qui ne respecte pas nos conditions ou directives, veuillez le signaler comme inapproprié.