Impact of yttrium hydride formation on multi-isotopic hydrogen retention by a getter trap for the DONES lithium loop

Compliance with imposed hydrogen concentration limits in the lithium loop of the DEMO-Oriented Neutron Source (DONES) requires the installation of an yttrium-based hydrogen trap. To determine an appropriate H-trap design, it is essential to have access to a numerical tool capable of simulating hydrogen transport in the DONES lithium loop connected to an yttrium pebble-bed. In the past, a simplified model was created that allows such calculations when hydrogen concentrations in the lithium are low. However, in certain DONES operating phases, the concentration in the lithium is high and in a range where yttrium dihydride (YH2) formation is likely. Due to the anticipated great impact of YH2 formation on the H-trap performance a new model is developed that includes the mechanism of hydride formation. It is based on a mathematical reproduction of complete pressure-composition isotherms of the Li–H and Y–H systems. Thus, the conditions that trigger YH2 formation are determined and the variation of hydrogen solubility in different yttrium hydride phases is deduced. An approximate concentration-dependent relationship of hydrogen diffusivity in yttrium is derived and incorporated into the model. Simulations are performed to analyze the dynamics of the concentration decrease during purification of the lithium circuit prior to the experimental DONES phase by varying design parameters of the trap. It is found that hydride formation greatly increases the hydrogen gettering capacity of the H-trap and limits the maximum concentration in the lithium. Indeed, YH2 formation may be purposefully triggered to exploit its beneficial properties for DONES. Simulations of the hydrogen purification process during the experimental phase of DONES show that the H-trap must be replaced at least every 28 days to meet tritium limits. This work sets the conditions for the required pebble-bed mass of the H-trap at a given temperature to comply with the DONES safety requirements. Finally, the model is validated by numerical reproduction of experimental results.

Compliance with imposed hydrogen concentration limits in the lithium loop of the DEMO-Oriented Neutron Source (DONES) requires the installation of an yttrium-based hydrogen trap. To determine an appropriate H-trap design, it is essential to have access to a numerical tool capable of simulating hydrogen transport in the DONES lithium loop connected to an yttrium pebble-bed. In the past, a simplified model was created that allows such calculations when hydrogen concentrations in the lithium are low. However, in certain DONES operating phases, the concentration in the lithium is high and in a range where yttrium dihydride (YH 2 ) formation is likely. Due to the anticipated great impact of YH 2 formation on the H-trap performance a new model is developed that includes the mechanism of hydride formation. It is based on a mathematical reproduction of complete pressure-composition isotherms of the Li-H and Y-H systems. Thus, the conditions that trigger YH 2 formation are determined and the variation of hydrogen solubility in different yttrium hydride phases is deduced. An approximate concentration-dependent relationship of hydrogen diffusivity in yttrium is derived and incorporated into the model. Simulations are performed to analyze the dynamics of the concentration decrease during purification of the lithium circuit prior to the experimental DONES phase by varying design parameters of the trap. It is found that hydride formation greatly increases the hydrogen gettering capacity of the H-trap and limits the maximum concentration in the lithium. Indeed, YH 2 formation may be purposefully triggered to exploit its beneficial properties for DONES. Simulations of the hydrogen purification process during the experimental phase of DONES show that the H-trap must be replaced at least every 28 days to meet tritium limits. This work sets the conditions for the required pebble-bed mass of the H-trap at a given temperature to comply with the DONES safety requirements. Finally, the model is validated by numerical reproduction of experimental results.
Keywords: DONES, hydrogen trap, yttrium getter, liquid lithium, hydrogen diffusion, tritium transport modeling, metal hydrides (Some figures may appear in colour only in the online journal) * Author to whom any correspondence should be addressed.
Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

Introduction
An essential step on the way to safe and reliable fusion power plants is the execution of high-energy neutron irradiation experiments with fusion-relevant materials. The experimental facilities IFMIF (International Fusion Materials Irradiation Facility) and its scaled-down earlier version DONES (DEMO-Oriented Neutron Source) are designed for this purpose and will be constructed in the near future [1]. In these experiments, a 40 MeV deuteron beam colliding with a flowing liquid lithium (Li) target produces a neutron ray by undergoing nuclear stripping reactions. The reactions cause a gradual accumulation of protium ( 1 H), deuterium ( 2 H), tritium ( 3 H) and 7 Be in the liquid Li which flows through a closed loop system [2,3]. A particular threat to the safety of the plant is the production of radioactive 3 H and 7 Be. Both isotopes could be liberated in case of an accident or escape into the environment by permeation or leakages [4][5][6]. In addition, hydrogen isotopes and other non-metallic impurities such as C, N and O can lead to increased corrosion of the structural loop material and to an erosion of precipitating solid compounds [7]. For this reason, a DONES safety analysis decided to establish stringent contamination limits to be met by installing an impurity control and an impurity monitoring side loop system [8]. While N is absorbed by a titanium getter in the dump tank [9], the impurities O, C and part of the accumulating 7 Be are removed by a cold trap [2,10]. In order to meet the safety requirements for the hydrogen isotope concentrations in the loop, it was decided to install an yttrium pebble-bed as a separate hydrogen getter trap [2]. Yttrium metal (Y) has a much higher solubility of hydrogen isotopes than Li. As a result, a significant fraction of the dissolved isotopes in Li is absorbed by the Y pebbles [2]. The hydrogen retention rate of the pebble-bed is greatly reduced in N-or O-contaminated Li, which provides another reason why purification of the Li from these impurities is essential [11]. In addition to the corrosion-related maximum allowable total hydrogen isotope concentration of 10 wppm in the Li [2], radiation safety requires that the tritium inventory be kept below 0.3 g, in both the Y and Li [6].
The determination of appropriate trap designs and operating conditions to meet these safety limits requires the development of a simulation tool capable of realistically modeling the transport of hydrogen isotopes from flowing liquid Li into an arbitrarily sized Y getter bed. This was the purpose of a previous numerical study [12] using the simulation software EcosimPro © [13]. It allows performing dynamic simulations of the hydrogen isotope concentrations and fluxes in the trap and in the surrounding loop system by continuously solving a system of coupled differential transport equations and algebraic boundary conditions. The previous model was built on the assumption that the yttrium pebbles in the trap always remain in the low concentration α-Y solid solution phase regardless of the concentration in the surrounding Li. As long as the total concentration in the Li is sufficiently low, this approximation is valid [14]. However, in DONES the concentration in the Li is likely to exceed critical values and YH 2 formation is triggered. Yttrium dihydride formation strongly determines the thermodynamics of metal-hydrogen systems and thus the hydrogen isotope absorption dynamics of a Y getter trap [15,16]. Therefore, to ensure correct simulation results, the development of a new improved model that implements the physical mechanisms of hydride formation is essential and forms the objective of this work.
The new numerical tool comprises a model of the DONES Li loop with a parallel-connected purification side loop containing the hydrogen hot trap. To include the processes of hydride formation into the model, complete pressurecomposition isotherms of the Li-H and Y-H binary systems are mathematically reconstructed by using recent data from the literature. Thus, the solubility of hydrogen isotopes in different yttrium hydride phases as a function of the concentration in the Li is obtained and the conditions for YH 2 formation are derived. Moreover, the new model considers a theoretically determined concentration dependency of the hydrogen isotope diffusivity in different yttrium hydride phases. The integration of the processes of metal hydride formation into a dynamic hydrogen transport model is a novelty and is demonstrated in this work for the first time on this scale.
A detailed description of the developed numerical model is given in section 2. Section 3 presents simulation results of the time-evolving hydrogen isotope concentrations in the DONES Li loop encountered during a proposed initial purification run and throughout the DONES experimental phase. The simulations demonstrate the influence of YH 2 formation on the hydrogen absorption behavior of the trap. In addition, the simulation results and the underlying mathematical model are used to derive conditions for the minimum required pebblebed mass and appropriate operating modes to meet the safety requirements during the different stages of DONES operation. Finally, the numerical model is tested against experimental values of a previous deuterium retention experiment which was performed in the past [17].

Numerical model
To simulate the process of hydrogen isotope retention by a Y getter trap in the liquid Li loop of DONES, a numerical model of the circuit is developed from scratch using the software EcosimPro © [13]. In this study, absorption and permeation processes of hydrogen isotopes in loop components other than the getter trap are neglected. Therefore, it is sufficient to model hydrogen transport in the DONES Li loop with a simplified layout compared to the full piping and instrumentation diagram shown in [18].
A component flow chart of the created model is presented in figure 1. It consists of a main loop with a volume flow rate F main and a purification side loop connected in parallel. In this model, the purification side loop merely contains a Y getter trap which is exposed to a Li flow rate of F trap . The total volume of liquid Li in the loop system is labeled with V Li . The model of the Li circuit comprises a set of algebraic and coupled differential transport equations that determine the temporal evolution of the concentrations and particle fluxes at various spatial positions in the main loop and the trap.

Hydrogen transport in the lithium loop
The molar concentrations c pipe h,Li (z, t) of hydrogen isotopes h ∈ { 1 H, 2 H, 3 H} in liquid Li flowing with constant velocity v Li and flow rate F pipe along the longitudinal z-axis of a pipe of length l pipe and volume V pipe is described by the general mass continuity equation: Parameter σ h is a particle sink or source term with the unit (mol m −3 s −1 ). The second term describes the advection of hydrogen isotopes in z-direction. The hydrogen hot trap is considered as a cylindrical pipe of diameter d trap . It is densely filled with N peb round Y pebbles of radius r peb and individual surface area A peb = 4π r 2 peb . Consequently, the total Y surface area yields A Y = N peb A peb . The length l trap and inner volume V trap of the trap container are determined through the relations l trap = 4V trap /(π d 2 trap ) and Here, ϵ is the void fraction of the pebblebed. The volume and mass of the Y pebble-bed is given by V Y = 4/3 · N peb π r 3 peb and m Y = ρ Y V Y , respectively. Parameter ρ Y = 4469 kg m −3 is the density of pure Y. The volume of liquid Li filling the interstitial sites of the pebble-bed is determined by Equation (1) is used to calculate the concentrations c trap h,Li (z, t) in the Li flow through the interstitial sites of the pebble-bed. This is done by considering the retention of hydrogen isotopes into the Y pebbles as a particle sink. It is quantified by an isotope specific retention flux J ret,h (z, t) which is negative for fluxes orientated towards the pebble centers and varies along the z-direction of the trap. For its numerical treatment, the trap container is discretized into M segments z j = z 1 , . . . , z M , with ∆z = l trap /M. A sketch of the applied discretization is presented in figure 2. It allows equation (1) to be written in the following finite-difference form: label the inlet and outlet concentration of the trap. The isotope specific retention rateṅ h,Y (t) of the Y pebblebed is defined by: Consequently, the total retention rate is determined byṅ Also the temporal and spatial evolutions of the hydrogen isotope concentrations in the main loop are described by equation (1). Therefore, it is discretized along its longitudinal axis z' axis into U spatial segments z ′ k = z ′ 1 , . . . , z ′ U . The generation of hydrogen isotopes originating from the impelling deuterons into the Li target is quantified by a generation rateṅ gen,h with the unit (mol s −1 ) and considered to happen in the first segment z ′ 1 of the main loop pipe. Here, it is treated as a particle source with σ h =ṅ gen,h · U/V Li,main . This allows adapting equation (1) to the situation occurring in the Li target which is expressed by the following finite-difference form: where c main h,Li,in (t) is the inlet concentration of the main loop. Parameter V Li,main = V Li − V Li,trap is the volume of Li flowing through the main loop. Neither particle sinks nor sources are considered to occur in the remaining main loop sections z ′ k = z ′ 2 , . . . , z ′ U . This reduces the finite difference expression of the mass continuity equation in these segments to:  In this model, the Li density is calculated by ρ Li (T) = 562 − 0.1 · T(K) [19]. Both the temperature and the Li density are assumed constant everywhere in the loop. Consequently, mass continuity requires F main = F BP + F trap , where F BP designates the volume flow through the bypass parallel to the purification side loop (see figure 1). The value of the concentration at the main loop inlet results from the mixing of the Li flows through the trap and the bypass, according to: The average concentration c h,Li (t) in the entire Li system is calculated with:

Hydrogen transport into the yttrium getter bed
In both the lithium and the yttrium, hydrogen isotopes are dissolved in atomic form. In thermodynamic non-equilibrium, diffusion processes cause a net transport of the isotopes through the metal until thermodynamic equilibrium is attained. The occurring diffusion fluxes highly depend on the gradients of the concentrations at each location in Li and in Y. According to Fick's first law, a flux of hydrogen isotopes from the Li into the pebbles requires the existence of a concentration gradient in the Li which is oriented towards the pebble centers [20].
The model is based on the assumption that within a thin boundary layer surrounding each pebble, the concentration in the Li decreases from its core concentration value c trap h,Li (z j , t) down to a minimum value of c h,Li,I (z j , t) at the Li-Y interface [21]. The thickness of the boundary layer d BL depends on the hydrodynamic characteristics of the Li flow through the pebble-bed and can be approximated by a specific mass trans- [22]. Parameter Sh is the Sherwood number of the Li flow through the trap. The expression for the mass transfer coefficient in pebble-beds [23] which is used in this model is equivalent to the one derived in the equations (20)-(24) of the article [12]. Assuming a linear concentration decrease in the boundary layer allows writing an approximated expression for the retention flux occurring on the Li side of the Li-Y interface [24]: Parameter D Li is the diffusivity of hydrogen isotopes in Li. The model considers the diffusivity relation of deuterium in Li for all three hydrogen isotopes. It has the unit (m 2 s −1 ) and is determined from the tritium diffusivity measured in [25] by considering the isotope effect of hydrogen diffusion in metals [12,26,27]: Here, M2 H and M3 H label the molar masses of deuterium and tritium [26,27].
In contrast to the previous model presented in [12], this work considers a more accurate numerical description of hydrogen isotope transport inside of each Y pebble. From now on, a pebble is regarded as a sphere consisting of N + 2 spherical shells i. Figure 3 illustrates the considered discretized shell structure. The pebble center forms the first shell and is described as a sphere of radius r 1 = δ. The last shell is considered as a spherical layer of thickness δ with an inner radius of r N+2 = r peb − δ and an outer radius of r N+3 = r peb . The defined radial positions of the remaining shells r i = (i − 3/2) · ∆r with ∆r = r peb /N ≫ δ are situated in the middle between their inner and outer spherical surface planes. As a consequence, the volumes V i of the spherical shells i can be expressed by: (10) Using Fick's first law, the number of isotopes h that cross a spherical plane at radius r of a Y pebble per time yields: Figure 3. Discretized shell structure considered for the finite difference description of hydrogen transport inside of each Y pebble.
Integrating this equation on both sides from one to the next radial discretization node and assuming the same isotropic diffusion coefficient D Y (r i , z j , t) within the same shell i allows defining the following two finite difference relations [28]: and They approximately describe the time-evolving fluxes and concentrations of each hydrogen isotope h at the nodes i in the pebbles. Consequently, the retention flux in Y at the The diffusion relation for hydrogen isotopes in Y is concentrationdependent. Obtaining an expression for D Y (r, z, t) requires a previous theoretical analysis of the thermodynamics of the Li-H and Y-H system which is done in section 2.3. The final expression for D Y (r, z, t) is derived in section 2.5.

Determination of the Li-Y boundary condition
To complete the set of transport equations, two boundary conditions are required. The first condition yields from the principle of mass conservation across an interface: The second boundary condition is the relationship between the hydrogen isotope concentrations c h,Li,I and c h,Y,I at the Li-Y interface. It will be seen that this boundary condition mainly determines the hydrogen isotope absorption kinetics of the Y pebble-bed and describes the process of hydride formation in the getter material.
From now on, the term hydrogen (H) is considered as a unified expression for all three hydrogen isotopes, which leads to the definitions c Me ≡ h c h,Me and J Me ≡ h J h,Me . In general, when in this paper a quantity is not assigned to a specific hydrogen isotope by the index h, it refers to the sum of all hydrogen isotopes. The symbol Me stands for either Y or Li. Furthermore, in the following, it will be omitted to relate each quantity to the corresponding longitudinal position in the trap, since the boundary conditions apply to every trap segment.
Hydrogen is absorbed by the Y pebbles until the concentration and temperature-dependent chemical potential of hydrogen dissolved in lithium µ Li-H (c Li , T) equals the chemical potential of hydrogen dissolved in yttrium µ Y-H (c Y , T). When this equilibrium is established, all concentration gradients and net diffusion fluxes vanish. Thermodynamic equilibrium between the two metal-hydrogen systems Li-H and Y-H can be mathematically treated by considering a scenario in which both systems are in equilibrium with a third system, a diatomic hydrogen gas reservoir. Such a configuration satisfies the equilibrium condition: The equilibrium pressure p eq H2 of the diatomic gas would coexist with two distinct equilibrium concentrations c Li and c Y in the two metals. An increasing chemical potential implies a bigger equilibrium pressure and thus a higher equilibrium concentration. The higher the hydrogen solubility of the respective metal, the higher the adjusting equilibrium concentration at a given equilibrium pressure. Assuming an ideal gas, the chemical potential of the gas phase can be expressed by: where µ • H2 is the chemical potential of the gas phase at standard state and p • = 101 325 Pa is the standard state pressure. Consequently, relation (15) transforms to [15,29]: Here, R is the ideal gas constant. Hence, in thermodynamic equilibrium the equality of the chemical potentials between the Li-H and the Y-H systems may be represented by an equality of their temperature and concentration-dependent equilibrium pressures p eq Li (c Li , T) and p eq Y (c Y , T). For the model, it is assumed that within an infinitesimal layer crossing the Li-Y interface, the chemical potentials in both materials are always in equilibrium. According to equation (17), this allows defining the second boundary condition: where the variables c Li,I and c Y,I refer to the equilibrium concentrations in the corresponding metal at the Li-Y interface. It should be mentioned that in this model boundary condition (18) The equilibrium pressure plotted against the equilibrium concentration at constant temperature is represented by socalled pressure-composition (p-c) isotherms. They are characteristic for each metal-hydrogen solution. At very low concentrations the p-c isotherms of the Y-H and Li-H systems can be approximately described by the Sieverts' law: where K s,Me is the temperature-dependent but concentrationindependent Sieverts' constant of hydrogen isotopes in the observed metal. Consequently, if the concentrations in both metals are low, boundary condition (18) reduces to [14]: The model presented in [12] made use of this boundary condition. Thus, it enables hydrogen transport simulations in the low concentration regime. However, at higher concentrations, lithium and yttrium hydride phases form and the Sieverts' law as well as boundary condition (21) are no longer valid [30]. The formation of metal hydrides dictates the courses of the p-c isotherms beyond the Sieverts' regime. For this reason, hydride formation has a strong influence on the absorption characteristics of a Y getter trap exposed to flowing liquid Li of elevated H content. This work is dedicated to implementing a generalized expression of boundary condition (18) in the numerical model. Therefore, the relations which approximately describe the courses of the p-c isotherms of the Li-H and the Y-H systems at higher concentrations are derived from literature.
In the past, p-c isotherms covering the complete concentration range of the Li-H binary system could be measured for temperatures 700 • C < T < 900 • C [31,32]. However, experimentally determining their courses in the temperature range relevant for DONES (200 • C < T < 400 • C) is not possible since the equilibrium pressures are too small to be measured.
In the low-concentration liquid α-Li phase H is dissolved as a solute. As mentioned above, in this regime the p-c isotherms are described by relation (20). This model considers the same Sieverts' constant for all three hydrogen isotopes. Namely, the Sieverts' constant of deuterium in Li which was determined by Smith et al [32]: Hendricks et al [12] presents a collection of Sieverts' constants of the Li-H system which were measured in the experimental campaigns [30][31][32][33][34]. They differ within a narrower zone than the range, which is limited by the Sieverts' relations of 1 H and 3 H measured by Smith et al [32]. Hence, the model considers: as the upper and lower error limits of K s,Li . The maximum concentration c α-end Li of the pure solid solution α-Li phase is called terminal solubility of H in Li [15]. This critical temperature-dependent concentration could be measured in the experimental campaigns [35][36][37] for 1 H and 2 H between 250 • C and 500 • C. Yakimovich et al evaluated their data for 2 H dissolved in Li and derived the expression: with (mol H /(mol H + mol Li )) as the unit of c α-end Li [38]. All other equations in this work consider (mol m −3 ) as the unit of the concentration. In equation (25), T M ≈ 963 K is the monotectic temperature of the Li-H system [38]. Whenever a temperature T appears in a formula of this article, it is considered with the unit (K). The reported mean square error of the coefficients in equation (25) is 1.9%. It defines an error range of c α-end Li . The position of the phase boundary for protium dissolved in Li given in [38] is located within this error range. Figure 4 presents the calculated phase boundary c α-end Li plotted in the T-c binary phase diagram of the Li-H system.
Most of the graphs in this article are created using Matplotlib [39]. The corresponding equilibrium pressure p α-end Li at the end of the pure α-Li phase is estimated by inserting relation (25) into the Sieverts' law (20): This pressure is known as the α → β decomposition pressure of the Li-H system. It increases with temperature. At the decomposition pressure, the p-c isotherms reach a flat plateau area. Here, a very small increase in equilibrium pressure results in the precipitation of solid Li hydride compounds (β-LiH). The decomposition regime where the two hydride phases α-Li (liq.) + β-LiH (sol.) coexist is called the α → β phase transition. In fact, the concentration c β-start Li within the precipitating solid β-LiH compounds is significantly higher than the terminal solubility concentration c α-end Li that occurs in the coexisting liquid α-Li phase [35,40]. Nevertheless, in this model a concentration c Me of hydrogen in metals refers to the mean concentration averaged over the distinct concentrations occurring in coexisting metal-hydrogen phases. Consequently, the higher the proportion of a precipitating metalhydride phase, the higher the average concentration c Me in the metal.
Once the entire Li-H solution has transformed to solid β-LiH, the pure β-LiH phase is reached. The final saturation concentration in the pure β-LiH phase is here approximated by c β-start Li ≈ ρ Li /M Li . In this expression, M Li is the molar mass of Li. The calculated phase boundary of purely saturated β-LiH is plotted in the phase diagram of figure 4. For the model, it is assumed that when passing through the α → β phase transition the equilibrium pressure increases by about 1%, such that p β-start Li , with f αβ Li ≈ 1.01. The assumed slight slope of the plateau region agrees with previously measured pc isotherms of the Li-H system [31,32]. In the α → β phase transition regime as well as in the β-LiH phase, the p-c isotherms are here described by linear straight lines. To reproduce an almost infinite slope of the equilibrium pressure in the β-LiH phase the model assumes the hypothetical values p β-end Li = 10 10 Pa and c β-end Finally, the complete theoretically approximated relation describing p-c isotherms of the Li-H system can be expressed by: where the slopes of the straight lines are defined by:  In the relevant temperature range 200 • C < T < 400 • C, p-c isotherms of the Y-H system could only be measured for equilibrium pressures p eq Y > 100 Pa which led to very high equilibrium concentrations [41,42]. In contrast, the T-c phase diagram of the Y-H system could be experimentally explored in a much greater concentration and temperature range [42][43][44][45][46]. It reveals the ability of Y to absorb up to three H atoms per Y atom. Three hydride phases (α-Y, δ-YH 2±x and ϵ-YH 3−y ) and their transition regimes exist for the Y-H system.
In the low-concentration α-Y phase, hydrogen is interstitially dissolved in Y where it occupies tetrahedral sites without changing the hexagonal closed-packed (hcp) lattice structure of the pure metal [40]. To describe the p-c isotherms of the Y-H system in this regime, the model considers the same Sieverts' constant of deuterium in Y for each hydrogen isotope. Since it is estimated to be the most accurate for the relevant temperature range [12], it is made use of the Sieverts' relation measured by Begun et al [45]: A list of Sieverts' relations of the Y-H system found in literature [41,45,47] is given in [12]. As for the Li-H system, it is found that the difference between the Sieverts' relations measured for different hydrogen isotopes by Begun et al is greater than the disagreement among the Sieverts' relations measured in different experimental campaigns. Therefore, the upper and lower error limits of the Sieverts' constant are estimated by: Peng et al performed a thermodynamic study in which the entire T-c phase diagram of the Y-H system was theoretically calculated [48]. Their calculation is based on an improved thermodynamic description compared to the previously determined phase diagram by Fu et al [42]. For this reason, it is chosen to use Peng's theoretical results as the basis for the numerical reconstruction of the p-c isotherms in this model. The temperature-dependent equilibrium concentration c α-end Y which marks the end of the α-Y phase is determined by fitting the fifth order polynomial: (33) against the theoretically calculated phase boundary. Table 1 presents the obtained fitting coefficients. The phase boundary c α-end Y which is calculated in the fitting procedure is plotted in the Y-H phase diagram in figure 6. By comparing the obtained expression for c α-end Y with a relation derived by Beaudry et al, it is estimated that in the range 200 • C < T < 400 • C the uncertainty of the determined expression for c α-end Y is approximately [44]: The corresponding decomposition plateau pressure at the end of the α-Y phase is determined by inserting the obtained relation c α-end Y into the Sieverts' law (20), such that: In the α → δ phase transition regime, yttrium dihydride (YH 2 ) precipitates within the lattice. This goes along with a change of the hcp crystal structure into a face-centered cubic (fcc) structure in which all tetrahedral sites are occupied by hydrogen isotopes [40]. When most of the lattice has transformed to YH 2 , the solution enters the δ-YH 2±x phase. The increase in equilibrium pressure during the α → δ phase transition is estimated from measured p-c isotherms by Fu et al to be roughly 10%, such that: Using equation (35), the error of p δ-start Y yields: Fitting the polynomial (33) against the lower phase boundary of the δ-YH 2±x phase, which was theoretically calculated in [48], enables defining a relation for c δ-start Y .
It is plotted in figure 6 using the corresponding fitting parameters listed in table 1.
In the δ-YH 2±x phase, further hydrogen isotopes occupy octahedral sites of the fcc structure, which demands a relatively large amount of energy. Here, significantly rising the equilibrium concentration requires increasing the equilibrium pressure by several orders of magnitude. At some point, a second plateau area is reached, which is characterized by the precipitation of hcp yttrium trihydride (YH 3 ). A relation for the height of the δ → ϵ decomposition pressure p δ-end Y is estimated by applying a linear regression against the corresponding Van't Hoff plot which was measured by Peng et al in the relevant temperature range [48]: By comparing relation (38) with the Van't Hoff relation which yields from a fit to the equilibrium pressures p δ-end determined. Based on experimental data by Fu et al [42], the equilibrium pressure assigned to the end of the δ → ϵ transition regime and beginning of the δ-YH 2±x phase is here approximated by p δ-end In this model, the courses of the p-c isotherms in the α → δ and δ → ϵ phase transition regimes as well as in the δ-YH 2±x and ϵ-YH 3−y phases are described by linear straight lines linking the defined characteristic points p eq Y (c Y ) which mark the start and end points of the different hydride phases. This leads to the following definition of the p-c isotherms of the Y-H system: where the slopes are defined by: The maximum theoretically achievable concentration in the Y-H system is considered to be c ϵ-end To reproduce an almost infinite pressure increase in the ϵ-YH 3−y phase, the model assumes a hypothetical equilibrium pressure of   Figure 7 shows the fully assembled p-c isotherms of the Y-H system for T = 250 • C, T = 300 • C and T = 350 • C as mathematically defined by relation (42).
Finally, a generalized expression of boundary condition (18) is numerically determined after inserting the relations (27) and (42) into equation (18). The result is plotted in figure 8. It reveals the value of c Y,I which occurs at the pebble surface at a given value of c Li,I . In addition, the visualized boundary condition indicates the homogeneous equilibrium concentration that would occur inside the Y pebbles as soon as a certain constant concentration in the Li is established at thermodynamic equilibrium.  (21):

Yttrium hydride formation in lithium
Its upper and lower error limits ∆ + c − Li,I and ∆ − c − Li,I are determined by error propagation: Its error range is given by: where ∆ ± p δ-start Y is determined through relation (37). In figure 9 the calculated critical concentrations c − Li,I and c + Li,I are plotted against temperature together with their error ranges. Figure 8 discloses that at c Li,I > c + Li,I , further increasing the concentration in the Li leads to an almost negligible rise in pebble surface concentration c Y,I . This can be attributed to the steep increase of the p-c isotherms in the δ-YH 2±x phase. By observing the figures 5 and 7, it is found that in the DONESrelevant temperature range the equilibrium pressure p β-start Li , at which a Li-H solution would fully transform to solid β-LiH, is not high enough to trigger a formation of YH 3 . For this reason, the formation of YH 3 surrounded by H loaded Li is not Hence, the concentration c δ-start Y can be considered as the approximate saturation concentration of H in Li-exposed Y at DONES-relevant temperatures.
The ratio between the concentrations that adjust in the Y pebbles and in the Li at thermodynamic equilibrium is called partitioning coefficient K Li−Y D of the Li-Y-H system. It is evident that the higher the partitioning coefficient, the higher the H gettering capacity of the Y getter at a certain equilibrium concentration in the Li. During thermodynamic non-equilibrium the partitioning coefficient is defined by the ratio of the interface concentrations, according to: Figure 10 presents the temperature-dependent partitioning coefficient of the Li-Y-H system as a function of c Li,I for different temperatures. It is calculated from the numerically obtained boundary condition (18) shown in figure 8. As the plot reveals, in the Sieverts' regime of the low concentration α-Y phase the partitioning coefficient takes a constant value K Li−Y D,0 . Inserting the low concentration boundary condition (21) into equation (51) Consequently, with the equations (22) and (30) it follows that: The error of this low concentration partitioning coefficient is calculated by error propagation: using the estimated errors (23)

Diffusion in different yttrium hydride phases
A physical quantity which, according to equation (10), strongly determines the retention flux of hydrogen isotopes into the Y pebbles, is the diffusion coefficient Several previous experimental studies have shown that the diffusion coefficient in Y is dependent on the hydride phase of the Y-H system [16,49]. This can be attributed to the fact that the driving force for chemical diffusion arises from the gradient of the chemical potential in the metal-hydrogen system. This is expressed by the generalized Fick's first law [15,50]: The quantity B Y is the H mobility in Y. Equation (54) allows writing a general expression for the temperature-and concentration-dependent diffusion coefficient of H in Y: The chemical potential µ Y-H of atomic H at a given position in the Y pebble where a certain concentration c Y occurs can be obtained by expressing µ Y-H in terms of the corresponding chemical potential of a diatomic H 2 gas that would adjust at thermodynamic equilibrium with the corresponding H-Y solution of concentration c Y . Therefore, equation (16) is inserted into equation (15) and the equilibrium pressure is replaced by the derived temperature-dependent p-c isotherm p eq Y (c Y , T), which is expressed by equation (42): Consequently, the derivative in equation (55) writes: At low concentrations, the diffusion coefficient is concentration-independent [15]. It is determined using equation (55) by inserting the Sieverts' law (20) into equation (57): For the value of D 0 Y , this model considers the diffusion relation of deuterium in Y for all three hydrogen isotopes. It is determined from the protium diffusivity reported in [27] by taking into account the isotope effect of diffusion [12,26,27]: Finally, substituting the equations (57) and (58) into equation (55) provides a relation for the concentrationdependent H diffusion coefficient in Y: In this description, the derivative of the equilibrium pressure p eq Y (c Y , T) which is described by relation (42) equals the phasespecific slopes s Y given by the equations (43)- (46). Figure 11 shows the calculated concentration-dependent diffusion coefficient in Y for three different temperatures. It can be seen that at concentrations corresponding to the phase transitions the diffusivity plunges. When entering a new phase, it abruptly rises by several orders of magnitude above the constant diffusivity D 0 Y of the α-Y phase. For the numerical description, the derived concentrationdependent diffusion coefficient is inserted into equation (12), z j , t)]. Here, c Y (r i , z j , t) refers to the summed concentration of all isotopes in the corresponding shell of the discretized pebble.
The previous model presented in [12] assumed the constant diffusivity D 0 Y for the entire concentration range. However, as It is important to mention that, in contrast to the new model described in this paper, the previous model presented in [12] used a different relation than equation (59) for the calculation of the diffusion coefficient D 0 Y . The relation used in [12] is reported by Talbot et al [51] and gives slightly smaller diffusivities than equation (59). Nevertheless, when in the following reference is made to a simulation with the old model, what is meant is a simulation in which the diffusion coefficient D 0 Y given by relation (59) and boundary condition (21) is considered to be valid over the entire concentration range. Figure 12 presents two simulated diffusion processes of H into a Y pebble of radius r peb = 5 × 10 −4 m assuming a constant interface concentration of c Y,I = 1.42 × 10 5 mol m −3 > c ϵ-start Y and a temperature of T = 300 • C. The upper plot shows the concentration profiles in the pebble after three different points in time that result if the old model is used for the simulation. Since the old model assumes a homogeneous diffusion coefficient, a continuous concentration profile is observed which rises smoothly over time. It follows the typical solution of Fick's second law for a diffusion process into a homogeneous sphere at constant interface concentration (see chapter 6 in [20]). The bottom graph displays a simulation of the same process, this time carried out with the new model. Areas shaded in color indicate the concentration ranges corresponding to the different hydride phases. The interrupted concentration scale in the lower plot illustrates the fact that the new model, unlike the old model, has taken into account a maximum H content of yttrium. It is found that a concentration-dependent diffusion coefficient results in a concentration profile which is no longer continuous. In fact, in the concentration ranges representing the α → δ and δ → ϵ transition regimes, the diffusing atoms pile up and form steep diffusion fronts. Due to the much slower diffusivity in the transition regimes compared to the pure hydride phases, the rate of diffusing atoms entering a transition range from an adjacent higher phase is bigger than the rate of atoms leaving the transition range into a lower phase. Moreover, it is observed that the α → δ front moves faster than the δ → ϵ front. These results are in good accordance with measurements made by den Broeder et al who experimentally visualized the same behavior of H diffusion in different yttrium hydride phases at lower temperatures [49].
It should be noted that the material-altering effects of radioactive decaying tritium on the H transport properties of the two metals are not taken into account in this model.

Simulation results
In this section, the developed numerical model is applied to simulate the retention of hydrogen isotopes by an H-trap installed in the DONES Li loop. Therefore, the volume flow rate through the main loop is set to F main = 1 × 10 −1 m 3 s −1 . The current design of the DONES Li loop envisages a total Li volume of approximately V Li = 14 m 3 which is considered as a reference value for the following simulations [52]. Unless otherwise noted, the simulations are executed assuming a flow rate through the trap of F trap ≈ 2 × 10 −3 m 3 s −1 , a trap diameter of d trap = 0.15 m and a pebble radius of r peb = 5 × 10 −4 m.
The pebble-bed void fraction of perfectly round pebbles packed in a cylindrical trap container with d trap ≫ r peb is usually determined to be approximately ϵ ≈ 0.4 [53]. However, in a real Y getter-bed the sizes and shapes of the pebbles are expected to strongly vary. Therefore, a packing density as high as ϵ ≈ 0.4 is likely not achievable. The mass transfer coefficient α f in equation (8) decreases with increasing void fraction [12]. Hence, to obtain conservative simulation results, the following calculations assume a higher void fraction of ϵ = 0.6.
The applied discretization density of the spatial subsystems considered in the model are set to N = 198, M = 6 and U = 3. To solve the system of equations introduced in the previous section, the finite difference equations (2), (5), (12) and (13) are manually solved in for-loops over the indexes i, j and k. The temporal parts of the differential equations (2), (4), (5) and (13) are automatically computed by EcosimPro © . It applies the DASSL algorithm [54] by making use of the Newton-Raphson method combined with the backward differentiation formula.
In order to cope with the imposed safety requirements for DONES, the installed H-trap needs to be capable of maintaining the average total hydrogen isotope concentration in the Li below c limit Li = 10 wppm. Further radiation safety requirements restrict the inventories of tritium in Y to m limit 3 H,Y = 3 × 10 −4 kg and in Li to m limit 3 H,Li = 3 × 10 −4 kg [6]. This is equivalent to c limit 3 H,Li = 7.1 × 10 −3 mol m −3 considering that V Li = 14 m 3 . The current design basis of the DONES impurity control loop envisages T = 350 • C as the design temperature of the H-trap [18].

DONES initial purification phase
It is assumed that the Li which is going to be purchased to fill the DONES Li loop for its future experimental campaigns will be contaminated with a relatively high protium concentration. Furthermore, filling the DONES Li loop the first time will cause the stainless steel walls of the pipes to release an additional significant amount of protium into the Li. Part of the dissolved protium will be retained by a cold trap which is working at a temperature of T = 190 • C. The cold trap retains solidified Li hydride (β-LiH) which precipitates above the temperature-dependent terminal solubility c α-end Li . Equation (25) reveals that at T = 190 • C the terminal solubility of the solid solution α-Li phase is c α-end Li ≈ 28.7 mol m −3 . Hence, it is to expect that after the first purification run using only the cold trap, the liquid Li will have a maximum protium concentration of approximately c init Li, 1 H = 28.7 mol m −3 = 57 wppm. To meet the DONES safety requirements, prior to continuous DONES operation, the Li should undergo another purification run using the H-trap in which its concentration is further reduced to c Li, 1 H < 10 wppm = 5.1 mol m −3 . This initial purification process using the H-trap is simulated and analyzed in the following paragraphs. Since protium is assumed to be the only isotope present during the initial purification phase, its concentration c Li, 1 H equals the total concentration c Li . Figure 13 presents the simulated temporal evolution of the average concentration c Li in the loop during initial purification It is found that after a period of several minutes in which the empty Y pebbles absorb the protium atoms as a solute, the surface concentration c Y,I surpasses the maximum concentration c α-end Y of the pure α-Y phase and hydride formation sets in (see figures 14 and 15). The precipitation of high concentrated YH 2 at the surface is characterized by an abruptly accelerated increase in surface concentration (see figure 15). As figure 16 reveals, during hydride formation the protium concentration in Li at the Li-Y interface abides in the concentration range c − Li,I < c Li,I < c + Li,I . Figure 17 demonstrates how the rapidly growing concentration c Y,I at almost constant value of c Li,I causes the partitioning coefficient K Li−Y D to increase above its constant value K Li−Y D,0 in the pure α-Y phase (see figure 7). Approximately 1 h after the start of the purification process, the increase in interface concentration at the pebble surface comes to a halt at the saturation concentration c Y,I = c α-end Y (see figure 15) and the partitioning coefficient reaches its maximum value (see figure 17). By this time, the Y surface has entered the pure δ-YH 2±x phase. Once surface  figure 10, the sudden rise in the concentration c Li,I implies that the partitioning coefficient decreases. At about the same time, the quickly diffusing protium atoms in the α-Y phase have homogeneously distributed in the entire pebble. In contrast, the much slower elapsing diffusion process occurring in the α → δ phase transition regime closer to the surface leads to the formation of a relatively steep diffusion front (see red shaded area in figure 14). It slowly flattens out over time while the pebble surface remains at the saturation concentration c Y,I = c δ-start Y . Figure 16 reveals that after approximately 3 h, the H absorption flux induced by the concentration gradient that occurs at the pebble interface has lowered the average concentration in the loop so much that c Li,I reaches a maximum before decreasing alongside with the average concentration c Li . According to figure 10, this results in a slow but persistent rise in the partitioning coefficient back towards its maximum value. The transport processes come to an end once thermodynamic equilibrium establishes and the concentrations in both the Li-H and the Y-H systems have reached steady-state values. At this time the concentration profiles in the pebbles are flat. Figure 14 shows that at equilibrium the Y pebbles of a getter-bed with a mass of m Y = 15 kg have entirely transformed to pure δ-YH 2±x and occur at a homogeneous concentration of c Y = c δ-start Y .
The influence of the implemented mechanism of hydride formation on the retention dynamics of the trap becomes  figure 17). For this reason, the final equilibrium concentration c eq Li ≈ 19 wppm which yields from the new model is much lower than the adjusting equilibrium concentration c eq Li ≈ 30 wppm predicted by the old model. The great difference between the two simulated equilibrium concentrations in Li obtained with and without the consideration of hydride formation illustrates the significant magnitude by which YH 2 formation increases the gettering capacity of Y. Ignoring the processes of hydride formation as done in the old model might therefore lead to an overestimation of the required trap mass to meet the DONES safety requirements. This example demonstrates the importance of incorporating the thermodynamics of hydride formation into the model. Observing figure 13, it is found that regardless of the pebble-bed mass, with the examined parameter set (T = 350 • C, F trap = 2 × 10 −3 m 3 s −1 , r peb = 5 × 10 −4 m, d trap = 0.15 m, ϵ = 0.6) the purification process takes about 12 h to reach thermodynamic equilibrium. It is found that the higher the trap mass, the lower the average equilibrium concentration c eq Li in the loop that adjusts at thermodynamic equilibrium. Figure 15 reveals that different pebble-bed masses lead to specific equilibrium concentrations in the pebbles. In fact, three distinct scenarios are identified.
The first scenario describes a case like the one discussed above (m Y = 15 kg, T = 350 • C), where, at thermodynamic equilibrium, the entire pebble-bed converts to pure δ-YH 2±x . This means that the concentration everywhere in the pebbles is c Y = c δ-start Y . As the derived concentration boundary condition illustrated by figure 8 indicates, full saturation occurs only if the establishing equilibrium concentration in Li satisfies c eq Li > c + Li,I . In this case, the final concentration c eq Li can be calculated assuming that at thermodynamic equilibrium the mole numbers of hydrogen isotopes in the yttrium n eq Y and in the lithium n eq Li are distributed according to n init Li = n eq Y + n eq Li . Here, n init Li is the mole number of hydrogen isotopes present in the Li prior to purification. With n init Li = c init Li V Li , n eq Li = c eq Li V Li and n eq Y = c δ-start Y V Y the equilibrium concentration yields: According to figure 6, the concentration c δ−start Y decreases with increasing temperature. Therefore, in the first scenario, the gettering capacity slightly reduces with increasing temperature. Solving equation (61) for m Y yields the pebble-bed mass required to reduce the equilibrium concentration from a given initial concentration c init Li down to a certain equilibrium value: The error of the required pebble-bed mass is determined by the error of the phase boundary ∆c δ-start Y which was estimated in equation (39). Consequently, error propagation yields: The second scenario refers to a case where during the concentration decrease the value of c Li falls below the minimum concentration required for hydride formation c − Li,I , such that c eq Li < c − Li,I . As figure 13 shows, such a situation is observed for a trap with a mass of m Y = 100 kg operating at T = 350 • C. The pink dash-dotted lines in the figures 15 and 17 reveal that if a sufficiently heavy trap is installed the concentration c Y,I never exceeds c α-end Y and no hydride formation takes place. Hence, if c eq Li < c − Li,I the pebble-bed occurs exclusively in the α-Y phase and the concentration boundary condition is described by relation (21). With equation (51), this implies that at thermodynamic equilibrium the mole number of H atoms in the pebbles is given by In this case, the relation n init Li = n eq Y + n eq Li transforms to: Also in the second scenario, the gettering capacity decreases with increasing temperature. This is attributed to the temperature dependency of K Li−Y D,0 which is visualized in figure 10 and mathematically expresses in equation (52). According to equation (64), in the second scenario the required pebble-bed mass to lower the concentration in Li down to a certain equilibrium value is determined by: Its error originates from the error of K Li−Y D,0 given by equation (53), such that: A third scenario arises when, for the observed temperature and pebble-bed mass, the values of c eq Li resulting in the equations (62) and (65) do not satisfy the respective conditions for c eq Li . In this case, the Li-Y-H system levels out at an equilibrium concentration c − Li,I < c eq Li < c + Li,I . According to figure 8, in this regime, the pebbles have only partially transformed to YH 2 once thermodynamic equilibrium is reached. Figure 9 shows that both c − Li,I and c + Li,I strongly increase with The equations (62) and (65) allow calculating the minimum pebble-bed mass which is required to reduce the concentration in the DONES Li loop from c init Li = 28.7 mol m −3 = 57 wppm to the limit value c limit Li = 10 wppm = 5.1 mol m −3 . It is plotted in figure 18 against the temperature.
From figure 9 it can be derived that if the operating temperature is smaller than a critical value T < T − crit = ( 347 ± 60 45 ) • C the DONES limit concentration satisfies c limit Li > c + Li,I and the pebble bed follows the defined first scenario. In this case, the minimum required trap mass is significantly lower than at higher operating temperatures and is described by the equations (62) and (63). They yield an almost constant value of m Y ≈ (18 ± 1) kg throughout the temperature regime T < T − crit . According to figure 9, if T > T + crit = ( 350 ± 60 45 ) • C the limit concentration satisfies c limit Li = 5.1 mol m −3 < c − Li,I . Consequently, in this regime, any pebble-bed will behave as described in the defined second scenario where the minimum required trap mass and its error range are given by the equations (65) and (66). It can be seen that if T > T + crit a much larger minimum pebble-bed mass m Y > 78 kg is required to comply with the concentration limit than at lower operating temperatures. Such a large pebble-bed mass is necessary because, without hydride formation, the gettering capacity remains relatively low. Figure 18 reveals that the defined third scenario occurs if T − crit < T < T + crit . In this regime, the minimum required pebblebed mass is in the range of 18 kg < m Y < 78 kg.
The old model did not take hydride formation into account which is why it suggests a minimum required trap mass given by equation (65) regardless of the operating temperature (see orange dashed line). Using the old model as a tool to determine appropriate design parameters of the trap at temperatures T > (265 ± 5) • C could hence only lead to a non-problematic overestimation of the minimum required trap mass. However, at temperatures T < (265 ± 5) • C, using the old model would lead to a severe underestimation of the necessary pebble-bed mass. This is because below this temperature the orange dashed line falls below the curve which yields from the new model. Therefore, implementing hydride formation in the model is essential.
In figure 18, the green shaded area represents an appropriate combination of operating temperature and pebble-bed mass which ensures meeting the safety requirement c Li < c limit Li during initial purification of the DONES Li loop. Picking the design parameters from the red area must be avoided, as this would lead to equilibrium concentrations above the concentration limit. The relatively large error range which enters deep into the green shaded area arises from the error ranges of the critical temperatures T − crit and T + crit . An optimal H-trap design for the DONES purification phase should be based on the intention to exploit the gettering capacity enhancing effect of yttrium dihydride formation. Therefore, the optimal trap should have an operating temperature of T < T − crit . To ensure hydride formation to occur and the DONES limit to be met under all circumstances, the conservatively estimated error range in figure 18 must be taken into account when choosing the design parameters. Considering this, it is found that if T < T − crit , a trap with a pebble-bed mass of m Y = 20 kg would ensure meeting the requirements with the least amount of yttrium. Figure 19 displays the simulated protium concentration decrease during initial purification at three different temperatures that would result from a trap with a mass of m Y = 20 kg. As predicted by figure 18, the two temperatures that satisfy T < T − crit lead to a concentration decrease well below the safety limit. In contrast, T = 350 • C is above the critical temperature. At this temperature, a pebble-bed mass of m Y = 20 kg can be assigned to the red shaded area in figure 18 and therefore does not comply with the concentration safety limit (see figure 19). It is found that the higher the temperature, the less time it takes for the system to reach thermodynamic equilibrium. This can be attributed to the fact that the diffusion coefficients in Li and in Y increase with temperature (see figure 11).
According to the error range in figure 18, an operating temperature of T = 290 • C is the maximum temperature where hydride formation in a pebble-bed with a mass of m Y = 20 kg is guaranteed. Furthermore, with this combination of design parameters, the purification run would be completed in a reasonable period of time of between one and two days (see figure 19). Having determined the optimal trap mass and working temperature it remains to be investigated how varying the flow rate, pebble diameter and aspect ratio AR ≡ l trap /d trap affects the purification dynamics. In the following, these quantities are referred to as secondary design parameters. Figure 20 depicts the simulated temporal evolution of the average concentration in the Li yielding from the determined optimal trap mass and temperature for different sets of the secondary parameters. The reference parameter set which was assumed for all previous simulations is F trap = 2 × 10 −3 m 3 s −1 , r peb = 5 × 10 −4 m and d trap = 0.15 m. Given that m Y = 20 kg and ϵ = 0.6, the trap length yields l trap = 0.63 m and the aspect ratio is AR = 4.2. Figure 20 discloses that operating a trap at a lower Li flow rate F trap causes the purification process to finish later. This has the following reasons. As discussed in detail in [12], the mass transfer coefficient of the pebble bed appearing in equation (8) follows the proportionality peb for Reynolds numbers 1.6 × 10 −3 < Re < 55 and the relation α f ∝ v 0.69 Li r −0.31 peb if 55 < Re < 1500 [23]. Since a smaller flow rate F trap implies a lower Li flow velocity v Li it leads to a reduced mass transfer coefficient and thus to a lower retention flux into the pebbles. Moreover, at a lower flow rate the residence time of Li in the trap container is higher. As a result, the Li flowing past is exposed to the Y getter for a longer time leading to smaller concentrations c trap h,Li in the trap. Consequently, the concentration gradient in Y at the pebble surface and thus the retention flux is smaller, so that thermodynamic equilibrium is reached later.
Furthermore, the plot in figure 20 shows that using a trap with wider pebbles also leads to a delay in the purification time. As equation (3) discloses, this is due to the fact that a reduced total surface area of the Y pebbles at equal pebblebed mass implies a reduced H retention rate. Moreover, wider pebbles imply a smaller mass transfer coefficient.
Increasing the trap diameter and thus reducing the aspect ratio of the trap container does not increase the residence time. However, it does decrease the Li flow velocity and thus the mass transfer coefficient of the pebble bed. Therefore, one would expect the time to reach equilibrium to increase with trap diameter. However, figure 20 shows no noticeable change in the shape of the concentration decrease as the trap diameter is increased from d trap = 0.15 m to d trap = 0.3 m. This is due to a second overlapping effect that increases the retention flux when moving towards wider trap diameters. It arises from the effect a reduced mass transfer coefficient has on the hydride formation dynamics. During α → δ phase transition, the pebble surface concentration in Y increases abruptly and the concentration gradients on both sides of the Li-Y interface increase alongside. This is directly related to an increase in H retention flux. However, the increased retention flux lasts only until the pebble surface reaches the saturation concentration. A smaller mass transfer coefficient causes hydride formation to proceed slower and surface saturation to happen later. Therefore, in wider traps the increased retention flux occurring during hydride formation lasts over a longer period of time.
Finally, it is found that in contrast to the pebble-bed mass and operating temperature, neither the Li flow rate, nor the pebble diameter or aspect ratio has an influence on the establishing concentration c eq Li at thermodynamic equilibrium.

DONES experimental phase
As soon as the DONES experimental phase is started, nuclear stripping reactions in the target component of the loop lead to an accumulation of hydrogen isotopes in the Li flow. Current estimations expect the generation rates [6,8,18]: The total hydrogen isotope generation rate is defined byṅ gen = hṅ gen,h . Due to the preliminary purification procedure which was discussed in the last section, the 1 H concentration at the start of the DONES experimental phase is expected to satisfy the condition c Li < c limit Li = 10 wppm. To provide conservative simulations, this work assumes an initial protium concentration of c1 H,Li = 10 wppm = 5.1 mol m −3 in the Li.
The generation rates yield that without the use of an H-trap, the 3 H inventory limit in Li (m limit 3 H,Li = 3 × 10 −4 kg) would be exceeded after the time: Without a trap, the total concentration in the Li would grow over time as indicated by the olive dashed line in figure 21.

Case A: non-replaced H-trap.
To demonstrate the effect of hydride formation on the behavior of a Y getter trap exposed to a liquid Li flow with continuous hydrogen isotope generation, a scenario is considered in which a single H-trap is used throughout the entire lifetime of DONES. Figure 21 showcases the simulated development of the average total hydrogen isotope concentration in the DONES Li loop over an operating time of 6 years at a fix temperature of T = 350 • C considering different pebble-bed masses of the H-trap. The zoom shown in figure 21 discloses that the initial 1 H concentration of c init Li, 1 H = 5.1 mol m −3 drops sharply at the beginning before reaching an initial thermodynamic equilibrium some hours later. The establishing initial equilibrium concentration can be calculated using the equations (61) and (64) or c − Li,I < c eq Li < c + Li,I , depending on whether the purification process belongs to the defined first, second or third scenario, respectively (see section 3.1). Whether c eq Li < c − Li,I or c eq Li > c + Li,I also determines the hydride phase in which the yttrium pebbles occur at initial equilibrium. In the example shown in figure 21, all three equilibrium concentrations satisfy c eq Li < c − Li,I (350 • C), which means that the yttrium pebbles must occur in the α-Y phase. Indeed, when examining the process over the scale of a year a concentration increase in the Li is observed which arises from the proceeding isotope generation in the target. The zoom demonstrates how the increase in concentration proceeds slowly compared to the time required to reach initial equilibrium, and therefore has a negligible effect on the value of the equilibrium concentration. This is because compared to the scale of a year, the time required for the hydrogen isotopes to diffuse from a pebble surface to the center is negligibly short. Therefore, it can be assumed that despite the continuous isotope generation, the Li-Y-H system remains at thermodynamic equilibrium and the concentration in the pebbles is approximately homogeneous. This means that the boundary condition visualized in figure 8 applies to both the Li system and the Y system as a whole, not just to their interfaces. Consequently, according to the equations (21) and (51), the isotope specific concentrations satisfy the condition n h,Y /V Y = K Li−Y D,0 n h,Li /V Li as long as the Y pebbles occur in the α-Y phase. It is obvious that the rates at which the generated isotopes are distributed between the Li and Y systems respect the relation:ṅ Therefore, the rates at which the molar numbers of the dissolved isotopes in the two metals increase while the pebbles occur in the α-Y phase are: Figure 22 displays the simulated time evolution of the approximately isotropic total hydrogen isotope concentration that Li,I to c + Li,I . Therefore, as long as hydride formation progresses the concentration in the Li occurs in a quasistationary state. By far most part of the generated isotopes is absorbed by the Y pebbles where they are directly converted to YH 2 . Saturation occurs once the nearly homogeneous concentration in the pebbles reaches c δ-start Y . Consequently, the mole numbers of added hydrogen isotopes to the Y and the Li during hydride formation are given by: The transformation in expression (74) is performed by making use of the equations (35), (36), (47) and (49). From the equations (74) and (75) it becomes apparent that n HF Y ≫ n HF Li ≈ 0. The time τ HF from the beginning of hydride formation to reaching saturation is obtained by considering that τ HFṅgen = n HF Li + n HF Y , such that: (76) Equation (76) and the figures 21 and 22 demonstrate that τ HF increases linearly with the pebble-bed mass. The rates at which the mole numbers of hydrogen isotopes in the Li and Y grow during hydride formation are given by: Assuming that during DONES operation the pebbles are approximately in equilibrium with the Li system, such that c Y,I ≈ c Y , the derivative of equation (20) with respect to time provides a general expression for the growth rate of isotopes of type h in Y: By considering that during hydride formationṅ HF Y ≈ṅ gen anḋ n HF Y ≈ 0, inserting relation (79) into equation (71) yields the following expression for the growth rates of the mole numbers of the isotopes h in Li during hydride formation: It can be seen that if c h,Li /c Li >ṅ gen,h /ṅ gen , thenṅ HF h,Li < 0, and if c h,Li /c Li <ṅ gen,h /ṅ gen , thenṅ HF h,Li > 0. Once the pebbles have completely converted to pure δ-YH 2±x , the pebble-bed is exhausted because the concentration in the pebbles cannot get any higher. As a consequence, the total retention rateṅ sat Y = 0, which implies that after saturation occursṅ sat Li =ṅ gen . Substituting the conditions of a saturated trapṅ sat Y = 0 and n sat Y = c δ-start V Y into the equations (72) and (79) provides the growth rates of the mole numbers of isotope species h in the Li after the pebble-bed is exhausted: Hence, if after trap exhaustion the generation of isotopes in the loop keeps changing the ratio c h,Li /c Li over time there is an ongoing isotope exchange between the Y and the Li. Both figures 21 and 22 contain a simulation executed with the old model developed in [12]. Since it does not incorporate the physics of hydride formation it predicts the curves to continue as if the Y-H system forever remained in the α-Y phase. Moreover, the old model does not consider the existence of a maximum saturation concentration of the Y pebbles in Li. Therefore, according to the old model, the pebble-bed can take up an infinite amount of hydrogen isotopes. However, the possibility of the trap being exhausted above the saturation concentration is a critical safety consideration and must be taken into account in the design of the H-trap for DONES. This is another reason why the inclusion of hydride formation in the model is so important. Figure 23 presents the simulated evolution of the average total hydrogen isotope concentration in the loop at different temperatures, considering a trap with a fix pebble-bed mass of m Y = 15 kg. It is found that at T = 200 • C and T = 250 • C, the initial concentration drop proceeds according to the third scenario defined in section 3.1, meaning that hydride formation starts right at the beginning. The figure shows that the height of the steady state plateau increases with temperature. This is in accordance with the temperature dependencies of the critical concentrations c − Li,I and c + Li,I depicted by figure 9. Moreover, it is observed that the time period hydride formation lasts decreases with increasing temperature. According to equation (76) and figure 6, this is due to the fact that the width of the concentration range assigned to the α → δ phase transition regime of the Y-H system narrows with increasing temperature. Figure 24 displays the simulated evolution of the tritium content dissolved in Li which results from a trap with a mass of m Y = 15 kg. Instead of abruptly entering into a flat stationary state, during hydride formation, the tritium inventory increase slows down rather smoothly until saturation occurs. This behavior is explained by equation (80) considering that due to the large amount of protium in the system c3 H,Li /c Li < n gen, 3 H /ṅ gen . It is found that at T = 200 • C the establishing stationary concentration is low enough to keep the 3 H inventory below the DONES safety limit m limit 3 H,Li = 3 × 10 −4 kg for more than 5 years.

Case B: replaced H-traps.
Indeed, since the amount of generated tritium in the loop follows the relationṅ gen, 3 H · t = n3 H,Li + n3 H,Y , using a single trap for the entire DONES lifetime to satisfy m3 H,Li < m limit 3 H,Li = 3 × 10 −4 kg would at some point violate the second imposed tritium safety requirement m3 H,Y < m limit 3 H,Y = 3 × 10 −4 kg. Nevertheless, the calculations in the previous section show that hydride formation in a sufficiently large trap prevents the concentration in Li from increasing above a certain temperature-dependent value for several years. This effect can be regarded a natural safety control system that is triggered in the event of an accidental, uncontrolled increase in concentration occurring in the DONES Li loop.
One proposed method to meet both limits at all times is to frequently replace the installed H-trap before the tritium contents in both the Li and in the Y exceed the limits [6].
To determine an appropriate operating temperature, pebblebed mass and trap replacement period τ trap it is considered that a non-exhausted H-trap operates in a state in which it is least efficient, i.e. when the Y pebbles occur in the α-Y phase where no hydride formation occurs.
As the DONES experimental phase starts, the tritium contents in both metals begin to grow. It is imagined that once an arbitrary tritium content m max 3 H,Li in the Li is reached, the trap is replaced by an empty trap. As a result, the tritium content in the Li decreases within a short period of time before it reaches an equilibrium, similar to the initial concentration decrease shown in figure 21. By hypothetically assuming that the concentration decrease occurs instantaneously, the value of the equilibrium tritium content in the Li can be calculated with equation (64): (82) While the tritium content in the Li plunges, the mass of absorbed tritium m3 H,Y in the new Y pebble-bed quickly rises until the concentration profiles in the pebbles are flat. Due to the reasons discussed in the previous section it can be assumed that, in spite of the proceeding hydrogen isotope generation, from this moment onwards the pebbles approximately remain in thermodynamic equilibrium with the Li system. Therefore, boundary condition (21) is approximately satisfied and yields: Moreover, the assumed thermodynamic equilibrium implies that the progressive production of hydrogen isotopes in the loop causes the tritium content in the two metals to increase approximately linearly, with slopes given by the equations (72) and (73). Therefore, the time evolution of the slowly increasing tritium content in the Li after the trap exchange is approximately described by: Inserting the equations (73) and (82) into equation (84) and replacing m3 H,Li with the maximum tritium content m max 3 H,Li that occurred in the Li just before the trap was replaced gives the time τ trap it will take for the tritium content in Li to rise back to the same maximum value m max 3 H,Li as prior to the trap replacement: Transforming equation (85) yields the maximum inventory of tritium in the Li m max 3 H,Li that establishes at a given value of τ trap : If τ trap = τ max trap , the adjusting maximum tritium content in the Li given by equation (86) is: The requirement for the pebble-bed mass to meet the other limit m3 H,Li < m limit 3 H,Li = 3 × 10 −4 kg is determined by inserting m limit 3 H,Li into equation (86): The minimum required pebble-bed mass is highest in case the maximum permissible replacement period of τ max trap is chosen. In this case, the requirement for the Y mass m Y is obtained by inserting equation (88) into equation (90): The relative error of m min Y is given by: Figure 25 presents a plot of the temperature-dependent minimum pebble-bed mass and its error range required to meet the imposed tritium limits of the DONES experimental phase considering that the trap is replaced every τ max trap = 28 days. The design parameters should be chosen such that they are located in the green shaded area, outside of the plotted error range. Traps with a mass and operating temperature in the red area would cause the tritium content in the Li to level out above the limit value m limit 3 H,Li . In contrast to the above thought experiment, in a real DONES operating scenario, it is not waited until the tritium contents in the Y and Li reach their maximum permissible values before the trap is exchanged. Instead, the trap replacement procedure is started from the beginning with a fix replacement period τ trap . As a result, the tritium contents in both Li and Y converge toward their maximum values determined by the equations (86) and (87) over the course of several replacement iterations.
From the equations (91) and (92) it follows that if τ trap = τ max trap , the minimum required Y mass at T = 350 • C is m Y = (16.4 ± 1.7) kg (see figure 25). The simulated temporal evolutions of the tritium contents in Y and in Li, that result from such a pebble-bed mass and operating temperature, are showcased in the figures 26 and 27 for different exchange periods, respectively. It can be seen that the shorter the trap exchange period, the faster the tritium contents in the two metals converge toward their maximum values. Figure 26 confirms the prediction from equation (91), that for an exchange period of τ max trap = 28 days the maximum tritium content in the pebble bed matches m limit 3 H,Y . Moreover, figure 27 verifies that if at T = 350 • C the Y mass is m Y = (16.4 ± 1.7) kg and τ trap = τ max trap the maximum adjusting concentration in the Li is indeed converging against the value m limit 3 H,Li . If the exchange period is longer, the tritium contents in the bed and the Li level out at higher values m3 H,Y > m limit 3 H,Y and m3 H,Li > m limit 3 H,Li . This must be prevented at all costs. A shorter replacement period allows keeping the tritium contents in both metals at values below the limit contents without changing the pebble-bed mass or temperature. The figures also demonstrate that if the pebble-bed mass has the value given by equation (91) and τ trap = τ max trap , the maximum tritium contents in the Li and the Y are equal. Figure 28 presents the simulated time evolution of the tritium content in Li at a fix temperature T = 350 • C and replacement period τ trap = τ max trap , varying only the pebble-bed mass. The plot validates equation (89) according to which the adjusting maximum tritium content in Li increases with decreasing pebble-bed mass. When choosing τ trap = τ max trap , the maximum tritium content in Y is m3 H,Y = m limit 3 H,Y , regardless of the observed temperature and Y mass (see equation (87)).
The variation of the adjusting maximum tritium content in Li that occurs when changing the temperature at a fix pebblebed mass of m Y = 7.5 kg and fix trap replacement period τ trap = τ max trap is depicted in figure 29. In this case, the value of m3 H,Li is determined by equation (89). It is found that m3 H,Li increases with temperature. The reason for this is that K Li−Y D,0 decreases with increasing temperature (see equation (52)).
The temporal evolutions of the protium, deuterium, tritium and total hydrogen isotope concentrations that occur when a trap with m Y = 16.4 kg operating at T = 350 • C is replaced with an exchange period of τ trap = τ max trap is shown in figure 30. The plots reveals, that during several trap replacements the protium concentration which starts at a value of c init Li, 1 H = 5.1 mol m −3 is gradually reduced before it keeps  oscillating at a constant level. In contrast, the deuterium and tritium concentrations approach their equilibrium values from below. Regardless of whether the initial concentration is zero or higher, the maximum concentration values below which the individual concentrations will oscillate after several  replacement iterations are determined by equation (86), in which the tritium generation rate and the molar mass must be replaced by the corresponding values of the observed isotope. It is evident that the values to which the concentrations of the different isotopes converge increase with increasing temperature and decreasing mass of the pebble-bed. Figure 30 shows that if the tritium limits for both the Li and the Y are satisfied, the total hydrogen isotope concentration in Li levels out considerably below a value where hydride formation would be triggered. However, at lower temperatures hydride formation might occur during the first trap replacement iterations when the total hydrogen isotope concentration is still relatively high. According to the findings of the last section, hydride formation would keep the concentration in Li stationary until the trap is exchanged.
The simulated effects a changing pebble diameter, Li flow rate through the pebble-bed or aspect ratio of the trap container have on the concentration dynamics can be seen in figure 31. It shows the evolution of the total hydrogen isotope concentration which is plotted in figure 30, magnified to day 336. First, it is found that a variation of the secondary parameters only affects the time to reach the intermediate equilibrium each time the trap is exchanged. It has no influence on the maximum establishing concentration in the Li. For the reasons discussed in section 3.1 the time that elapses before the intermediate equilibrium is established increases with increasing pebble diameter and decreasing flow rate. In contrast to figure 20, since no hydride formation takes place a reduced aspect ratio is found to slightly delay the concentration drop. Due to a smaller diffusion coefficient a lower temperature would delay the time to reach the intermediate equilibrium (see figure 19).
Whether the exhausted traps should be replaced by completely new Y pebble-beds or whether a thermal degassing process under vacuum conditions could regenerate the used pebble-beds for reinsertion remains to be investigated.

Experimental validation of the model
Finally, the presented numerical model is tested against experimental results obtained by Yamasaki et al [17]. In their experiments, m Li = 0.6 kg of liquid Li at T = 300 • C is loaded with an initial deuterium concentration of c init 2 H,Li = 160 wppm = 40.1 mol m −3 before it is set in motion through an experimental liquid Li circuit connected in line with a Y bed. The trap container used for the experiments has a length of l trap = 0.3 m and contains m Y = 1 × 10 −2 kg of 2-3 mm wide yttrium chips. The Li flow rate through the trap is reported to be F trap = 2.5 × 10 −5 m 3 s −1 . As soon as the deuterium loaded Li passes the pebble-bed, its deuterium content is continuously reduced. Yamasaki et al measured the decrease in deuterium concentration that occurred after five time periods of different lengths using the chemical dissolution method [17]. Figure 32 shows the ratios between the measured concentration c2 H,Li of the extracted Li samples and the initial deuterium concentration in the Li c init 2 H,Li as a function of the Li volume that has passed through the trap [17].
The described experiment is numerically reproduced using the presented EcosimPro © model illustrated in figure 1. Therefore, it is assumed that F main = F trap = 2.5 × 10 −5 m 3 s −1 . Moreover, the simulation considers a pebble radius of r peb = 5 × 10 −4 m. The void fraction and diameter of the pebble-bed Figure 32. Ratios between the deuterium concentration of extracted Li samples that were exposed to a Y bed and the initial deuterium concentration in the Li prior to the contact with the Y, measured in deuterium retention experiments performed by Yamasaki et al [17]. The graph contains simulations that reproduce the experimental values with the old and the new model. is set to ϵ = 0.6 and d trap = 5 × 10 −3 m which implies that the simulated trap has a length of l trap ≈ 2.8 × 10 −1 m. The discretization densities of the trap and pipe segments are set to the same values like in the simulations discussed in the previous sections.
The simulation results obtained with the numerical reproduction of Yamasaki's experimental data is plotted with a solid blue line in figure 32. It can be seen that the simulation agrees very well with the first four data points. The last data point is located considerably below the simulated curve.
When analyzing the observed purification process, it can be certainly excluded that the pebbles in Yamasaki's experiment have completely transformed to pure δ-YH 2±x . There is just not enough deuterium in the system to load a Y mass of m Y = 1 × 10 −2 kg with a homogeneous saturation concentration of c δ-start Y . Therefore, the establishing equilibrium concentration in Yamasaki's experiment must have satisfied the condition c eq Li < c + Li,I (300 • C) (see figure 9). The two critical concentrations c − Li,I (300 • C) and c + Li,I (300 • C) are plotted together with their overlapped error range (blue shaded area) in figure 32. It is found that the simulated concentration levels out in the concentration range c − Li,I < c eq Li < c + Li,I . Hence, according to the simulation, at thermodynamic equilibrium, the Y chips in Yamasaki's experiment occur in the α → δ transition regime. Given that c eq Li < c + Li,I (300 • C), an overestimation of the critical concentrations c − Li,I and c + Li,I implies an overestimation of the simulated equilibrium concentration. In fact, the error range of the critical concentrations can be considered as the error range of the simulated equilibrium concentration. As the equations (37), (48) and (50) indicate, this error range is mainly determined by the errors of the Sieverts' constants of the Li-H and Y-H systems. Figure 32 shows that the error bar of the last experimental data point enters deep into the error range of the simulated equilibrium concentration (see figure 31). Therefore, it is reasonable to say that the simulation results and the experimental data match within their error ranges. The experimental data suggests that the model rather underestimates the H gettering capacity of Y in liquid Li. For this reason, designing the H-trap of DONES based on simulations performed with the presented model would probably lead to a slightly oversized trap. Nevertheless, from a safety point of view, this is beneficial. Figure 32 includes a former numerical reproduction of Yamasaki's experimental data which was presented in [12] (see dashed orange line). It was performed using the old model of H retention from flowing liquid lithium by an yttrium pebble-bed. A comparison of the curves produced by the old and the new model shows that the inclusion of the hydride formation process in the model has significantly improved the agreement of the simulation results with the experimental data.

Conclusion
An appropriate design of the hydrogen hot trap for the DONES Li loop depends on the development of an accurate numerical simulation tool that can predict the hydrogen absorption behavior of a Y pebble-bed exposed to hydrogen-loaded flowing liquid Li. This was accomplished by a previous work that focused merely on the low concentration range. However, in certain operating phases of DONES the hydrogen concentration in the Li is relatively high and likely leads to a formation of YH 2 in the pebble-bed.
To account for these phenomena, in this work a new numerical model of hydrogen transport in a Y getter trap connected to the liquid Li loop of DONES has been developed which incorporates the thermodynamic processes of hydride formation. Therefore, a generalized concentration boundary condition for the Li-Y interface was derived by approximating the courses of complete pressure-composition-isotherms of the Y-H and Li-H systems in the relevant temperature range of DONES. This allowed determining the solubility of hydrogen in Y as a function of the concentration in the surrounding liquid Li. It was found to highly increase when YH 2 forms. By estimating the profile of the chemical potential of hydrogen isotopes inside of the Y pebbles during thermodynamic non-equilibrium a relation for the concentration dependency of the diffusion coefficient in different yttrium hydride phases could be obtained.
By simulating the time-evolving hydrogen isotope concentrations in the Li and the Y bed during different DONES operating phases, the effects of YH 2 formation on the trap performance were analyzed in detail. It was found that when choosing an appropriate pebble-bed mass and operating temperature, hydride formation is triggered and may be exploited purposefully to increase the hydrogen gettering capacity of the trap. This is shown to be particularly important for the initial purification run of the DONES Li loop prior to the experimental phase, when the hydrogen concentration in the Li is highest. The simulations have shown that if during the DONES experimental phase the ongoing generation of hydrogen isotopes causes the concentration in the Li to exceed a certain value, the formation of YH 2 is triggered, maintaining the concentration in the Li at a constant level. In addition, a mode of operation relevant to the DONES experimental phase has been presented that involves replacing the H-trap after at least every 28 days so that the tritium contents in the Li and the Y remain below safety limits. Furthermore, this work established the required conditions for the mass of the Y pebble-bed and the operating temperature of the installed H-trap which allow meeting the safety limits for the hydrogen isotope concentrations in DONES. Algebraic formulas have been derived that enable calculating the minimum required trap mass at any given temperature for both the initial purification phase and the experimental phase of DONES, taking into account error ranges.
Finally, the developed model was tested against experimental results. Simulation and measurement data were found to be in good accordance within each other's error limits. Indeed, compared to the old model, the agreement between simulation and experimental data could be considerably improved.