- Email: [email protected]
Phase-field simulation of irradiated metals
Computational Materials Science 50 (2011) 960–970
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
Phase-field simulation of irradiated metals Part II: Gas bubble kinetics Paul C. Millett a,⇑, Anter El-Azab b, Dieter Wolf c a
Idaho National Laboratory, Idaho Falls, ID 83415, USA Florida State University, Tallahassee, FL 32310, USA c Argonne National Laboratory, Argonne, IL 60439, USA b
a r t i c l e
i n f o
Article history: Received 14 September 2010 Received in revised form 14 October 2010 Accepted 27 October 2010 Available online 4 December 2010 Keywords: Phase-field simulation Irradiation damage Radiation damage
a b s t r a c t The phase-field model developed in Part I of this work is expanded to include fission gas generation, diffusion, and segregation within bubbles nucleated both homogeneously and heterogeneously along grain boundaries. Illustrative results are presented that characterize bubble growth and shrinkage, as well as the bubble density, size and nucleation rate as a function of the irradiation conditions. Finally, intergranular bubble characteristics such as shape, pinning energy and bubble density are investigated. Ó 2010 Elsevier B.V. All rights reserved.
1. Introduction The accumulation of gaseous species in nuclear fuel and structural materials leads to deleterious microstructural changes that are of utmost concern to both fission and fusion reactor applications [1]. Inert gas elements – such as helium, xenon and krypton – are generated due to fission reactions in nuclear fuel [2,3]. In addition, transmutation via (n,a) reactions in the metallic cladding leads to helium production [4]. In general, the precipitation of gaseous species into gas-filled bubbles leads to volumetric swelling in both structural materials and fuels, as well as high-temperature embrittlement and overall loss of mechanical strength in metallic components [1]. The processes of nucleation, growth, migration, and coalescence of gas-filled bubbles are all governed by (i) the dynamic production of point defects (i.e., vacancies, self-interstitials, and gas atoms) by irradiation and (ii) the pre-existing microstructural features of the material (i.e., grain boundaries (GBs), dislocation networks, and precipitates). Such features can themselves evolve as a result of interaction with point defects. A precise understanding of the inter-dependency of the point and extended defects, and their kinetics, as they relate to gas behavior still remains a very important subject in the field of nuclear reactor materials. According to conventional nucleation theory, in the early stages of irradiation, gas atoms diffuse throughout the crystal and begin
⇑ Corresponding author. E-mail address: [email protected] (P.C. Millett). 0927-0256/$ - see front matter Ó 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.commatsci.2010.10.032
to form clusters due to their insolubility in the solid [5]. This clustering can initiate either by chance encounters of the wandering gas atoms (i.e., homogeneous nucleation), or by gas-atom trapping at microstructural segregation sites such as GBs (i.e., heterogeneous nucleation). The subsequent growth or shrinkage of a gas bubble is governed by the fluxes of vacancies, interstitials, and gas atoms to or away from the bubble surface [6]. Basically, the time rate of change of the volume of a bubble is equal to the difference in the rates at which vacancies and interstitials are absorbed at the bubble surface times the (atomic) volume carried by each of these point defects [5]. In addition, gas-atom re-solution can limit bubble growth whereby fission fragments directly or indirectly cause gas atoms to be lost from the bubble [7,8]. If, however, gas bubbles are located on sinks (such as GBs) that recollect the lost gas atoms, they may grow indefinitely to the point of interconnection, thus leading to rapid fission gas release to the exterior of the material in the case of fuel [9]. Conversely, the existence of intergranular gas bubbles can, in turn, influence GB migration due to Zener pinning [10], thus ultimately affecting the microstructural evolution of the material. In the current paper, we expand the phase-field model developed in Part I to include gas atom evolution, thus allowing investigations of bubble nucleation and growth in irradiated metals. This phase-field model is based on the dynamics of generation, diffusion and reactions of vacancies, gas atoms, and self-interstitials (to be referred to merely as interstitials from here on) in a single component, polycrystalline metal. A polycrystal phase-field model [11] is incorporated, thus capturing lenticular-shaped intergranular bubble formation as well as bubble-induced Zener pinning of
961
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
GB migration. This model differs from that of Hu et al. [12], in that the order parameter g (introduced in Part I) differentiates the bubbles from the solid, and allows different thermodynamic and kinetic descriptions to be applied to each phase throughout space.
vacancies, interstitials, and gas atoms in the matrix; it is written in the form
f solid ¼ Efv cv þ Efi ci þ Efg cg þ kB T½cv lnðcv Þ þ ci lnðci Þ þ cg lnðcg Þ þ ð1 cv ci cg Þ lnð1 cv ci cg Þ;
2. Simulation approach In this model, we consider a polycrystalline metal in which bubbles are considered to be clusters of vacancies that contain a finite amount of gas atoms. (The current model is an extension of our previous work [13], which also considered single component metals.) The phase-field approach is based on expressing the total free energy functional F of the heterogeneous material in terms of the free energy of each of its phases and interfaces. Consistent with the Cahn–Hilliard definition of free energy for a non-uniform system [14], we describe the free energy functional of a polycrystalline solid consisting of both matrix and bubble phases, in terms of the vacancy concentration field cv, the interstitial concentration field ci, gas atom concentration field cg, an order parameter field g that distinguishes solid and void phases, and a set of order parameters /1?P that represent discrete crystallographic (grain) orientations
Z
hðgÞf solid ðca Þ þ jðgÞf bubble ðca Þ þ f pc ðg; /1!P Þ þ f ðca ; g; /1!P Þ dV:
F¼N
V grad
ð1Þ
In this expression, ca represents the set of concentration fields (cv, ci, and cg), and N is the number of lattice sites per unit volume of the crystal and thus all terms under the integral sign represent free energy per lattice site. The first three terms in the righthand-side correspond to the bulk system energy. The interfacial energy contributions, due to the bubble free surfaces and GBs, are accounted for by the gradient term. The excess free energy of the solid due to point defects fsolid(ca) is derived in terms of the enthalpy and entropic contributions of
ð2Þ
where Efv , Efi , and Efg are the vacancy, interstitial, and gas formation energies, respectively, kB is Boltzmann’s constant, and T is the absolute temperature. The entropic term in the brackets in Eq. (2) represents the configurational entropy of a defect-disordered crystal, shown schematically in Fig. 1, consisting of sites that are occupied either by a host atom, gas atom, vacancy, or a dumbbell self-interstitial defect. Note that this equation assumes that the gas atoms are substitutional defects that occupy regular lattice sites, which is typical of large fission gases. The free energy of the bubble phase is defined as
h i f bubble ¼ ðcv 1Þ2 þ c2i h i þ log cg þ cg kB T ln cg þ cg kB T lnðkB TÞ :
ð3Þ
The term in the first set of brackets was presented in Part I, in which an energy minimum is located precisely at cv = 1 and ci = 0. At this minimum, the free energy is exactly equal to zero. The term in the second set of brackets is the free energy of the gas phase within a bubble according to the ideal gas law. Its derivation is given in Appendix A. Here, the free energy is described by a gas concentration rather than the conventional gas pressure to be consistent with the use of a concentration variable in the matrix. The use of the ideal gas law is perhaps not appropriate for bubbles with a high density of gas, in which case the Van der Waals equation is more accurate [15], nevertheless, we make use of this derivation for the current work. To incorporate GBs and their mutual interactions with gas bubbles into the model, we employ a conventional polycrystal phase-field model, in which a set of order parameters /1?P are used to define the set of crystallographic grain orientations. The free energy of the polycrystalline material is defined as
f
polycrystal
ðg; /1;...;P Þ ¼ 3A
P X k¼1
!2 /2k
2
þg
4A
P X
! /3k
3
þg
:
ð4Þ
k¼1
For g = 0 (i.e., in the solid phase), this equation creates P energy wells located at {uk=l = 1, uk–l = 0}. For g = 1 (i.e., in the bubble phase), this equation has one energy well located at /1?P = 0, thus all grain order parameters go to zero inside a bubble. We note that this polycrystal model is qualitatively similar (the functional form differs) to that of Moelans et al. [16] that captures the pinning of GBs by particles; however, the current model allows for a dynamically-evolving bubble phase. To our knowledge, this is the first phase-field model explicitly taking into account dynamically-evolving gas bubbles interacting with planar defects, such grain boundaries. The gradient energy contribution, contributing primarily to gas bubble free surface energies and GB energies, are accounted for by the term
f grad ¼
Fig. 1. Schematic plot illustrating the quasi-lattice representation of vacancy, selfinterstitial, and gas atom defects. The configurational entropy of this system is included in Eq. (2), in which each lattice site can be occupied by either an atom, a vacancy, a gas atom, or a dumbbell self-interstitial.
X jw 2
jrwj2 ;
ð5Þ
where w represents the set of concentration fields (cv, ci, cg) and order parameter fields (g,/1?P). Following the standard procedure in the phase-field approach, the kinetic equations for the spatial and temporal evolution of the gas concentration field (in addition to the vacancy and interstitial concentrations derived in Part I) has the form
962
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
@cg 1 dF þ ng ðr; tÞ þ Pg ðr; tÞ: ¼ r Mg r N dcg @t
ð6Þ
The gas atom mobility is defined as
Mg ¼
Dg c g ; kB T
ð7Þ
where Dg is the gas atom diffusivity in the matrix. We note that, for a substitutional gas atom (as assumed in Eq. (2)), the gas mobility should depend on cv in addition to cg. This paper does not account for this point. Because these concentrations remain fairly dilute in the solid, this approximation is valid to a first order analysis. In the bubble phase, the gas can be assumed to instantly equilibrate when considering solid-state diffusion time scales. To account for this, the gas concentration inside each bubble is evenly distributed to provide a flat profile within each bubble. The gas concentration within the bubbles will vary from bubble to bubble. Bubble re-solution due to fission-fragment bombardment [5] is not included in the current model. This process can conceivably be including by adding a term that re-distributes fission gas from the bubbles to the surrounding solid regions. The exclusion of fission gas re-solution will result in higher bubble nucleation rates and growth rates predicted by the simulations. The order parameter fields that define the distinct grain orientations /i also evolve according to a set of Allen–Cahn equations
@/i ðr; tÞ dF ¼ L/ ; @t d/i
ð8Þ
where L/ is the grain boundary mobility. In the present work, Eq. (6) can be expressed in non-dimensional units following Part I leading to
@cg ~ 1 dF þ ~fð~r; ~tÞ þ P ~ M ~ g ð~r; ~tÞ ~ gr ¼r ~ dcg @~t N
ð9Þ
~ ¼ l r is the dimensionless gradient operator. Energy is where r kept in units of eV. In all of the simulations below, the grid spacing is Dx = k, the time step is Dt = 0.001 s, and the gradient terms are jv = ji = jg = jg = j/ = 1.0 eV k2. In Sections 4 and 5 below, we vary the Frenkel pair production rate and the gas atom production rate independently. In nuclear reactor conditions, these two are coupled to some degree; however, we disregard this coupling which is more applicable to dual-beam ion irradiation. Varying these production rates independently, as will be shown below, can lead to interesting results. 3. Gas bubble growth and shrinkage The size evolution of a gas bubble existing in a crystalline solid is governed by the flux rate of vacancies and interstitials to its surface. A net absorption of vacancies by the bubble will result in bubble growth, while conversely a net absorption of interstitials results in bubble shrinkage. Furthermore, the density of gas atoms within the bubble will fluctuate in time based on (i) the flux of gas atoms to the bubble and (ii) the change in the bubble volume itself. In order to validate that the current model can appropriately capture bubble growth and shrinkage based on the vacancy and selfinterstitial concentrations in the surrounding solid, as well as the direct dependence of this on the bubble’s gas concentration, we have simulated three benchmark cases, which are described below. For the first benchmark case, we have initialized the system to consist of a single void with a diameter of 20 k located at the center of the domain with periodic boundary conditions, as shown in Fig. 2a. Initially, the void itself contains no gas atoms, however, the surrounding solid does contain a finite concentration of gas atoms. The vacancy and interstitial concentrations in the solid are at their equilibrium values, respectively, thus we do not expect
Fig. 2. Increase of gas concentration inside an initially-empty void surrounded by a solid with a supersaturated gas concentration. The top figures show the evolution of (a) vacancy and (b) gas atom concentration fields. (c) Plots of the conserved concentration fields, cv, ci, cg, along a cross-sectional slice through the centerline of the bubble throughout time. The absorption of gas atoms into the void creates a gradient in the gas concentration in the adjacent solid region.
to observe any change in bubble size. For the simulations shown in Figs. 2–4, the energies of formation of the vacancies and interstitials in the crystal are taken as Efv ¼ Efi = 0.8 eV corresponding to eq equilibrium intrinsic defect concentrations of ceq v ¼ ci = 6.9 4 10 (the thermal energy in all simulations within this paper is taken as kBT = 0.11 eV, corresponding to a temperature of T = 1276 K). The gas atom defect energy in the solid is also taken as Efg = 0.8 eV. The reference chemical potential of the gas atoms in the gas phase is log = 0.4 eV throughout this paper (this value was arbitrarily chosen, with the condition that it be less than Efg ). The bulk recombination term is turned off (i.e., Rbulk = 0) so that cv in the bulk remains constant, however surface recombination is activated by assigning a value of Rsurf = 5 s1. The supersaturation of the gas in the solid is Sg = 120 (this supersaturation is calcuf lated with the supposed ceq g based on Eg , although there is no actual equilibrium concentration of impurity gas atoms in a crystal). Fig. 2a and b are contour plots that illustrate how the vacancy concentration (3a) and the gas atom concentration (3b) evolve in time. As can be seen, the gas concentration within the void increases in time due to the diffusional flux of gas atoms through the solid towards the bubble surface. The bubble volume, however, does not change in time due to the fact that there is not an excess quantity of vacancies or interstitials in the solid. Therefore, as more gas atoms enter the bubble, which has a constant size, the gas density in the bubble increases. This is shown quantitatively in Fig. 2c
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
963
which shows plots of the cv, ci, and cg concentration fields along a cross-sectional slice through the centerline of the bubble throughout time. In addition to the increase in gas concentration inside the bubble, in the nearby solid, a gradient in gas concentration develops due to the depletion of gas atoms from the solid to the bubble. This gradient drives an on-going flux of gas atoms to the bubble surface, which continues until the chemical potential of gas atoms in the bubble equals that in the solid. The last snapshot in Fig. 2 does not correspond to the fully equilibrated system, i.e., more time is needed before the gas concentration is done evolving (this is also true of Figs. 3 and 4 below). The second benchmark case, shown in Fig. 3, consists of a gas bubble that, initially, has a smaller diameter of 10 k. The initial vacancy concentration in the solid is supersaturated by Sv = 120. In addition, the gas concentration in the bubble is initialized at a high value corresponding to the case of an over-pressurized bubble. As shown in Fig. 3a and b, the bubble grows in size due to the absorption of vacancies from the solid. Furthermore, as the bubble’s area increases, the concentration of gas inside the bubble decreases. The total concentration of gas is conserved, thus this decrease in bubble gas concentration is solely due to the increase in bubble size, which is in turn governed by the diffusional flux of vacancies to the bubble surface. Contrary to the above case of void growth in a supersaturated vacancy field, a system that is supersaturated with interstitials will
Fig. 4. Shrinkage of a bubble due to a supersaturated self-interstitial concentration in the solid. The gas concentration in the bubble, which is initially under-saturated, increases due to the decrease in bubble volume. The top figures show the evolution of (a) vacancy and (b) gas atom concentration fields. (c) Plots of the conserved concentration fields, cv, ci, cg, along a cross-sectional slice through the centerline of the bubble throughout time.
Fig. 3. Growth of a bubble due to a supersaturated vacancy concentration in the solid. The gas concentration in the bubble, which is initially over-saturated, decreases due to the increase in bubble volume. The top figures show the evolution of (a) vacancy and (b) gas atom concentration fields. (c) Plots of the conserved concentration fields, cv, ci, cg, along a cross-sectional slice through the centerline of the bubble throughout time.
drive a void to shrink in size due to the continual occupation of empty lattice sites on the void surface by interstitials. This is essentially identical to vacancy-interstitial recombination in the bulk, albeit the rate may be higher due to the unique structure of the free surface. As mentioned in Part I, we consider the overall recombination term to have contributions from both the bulk and the free surfaces. Fig. 4 shows the shrinkage of a bubble surrounded by a solid containing a supersaturation of interstitials. The vacancy concentration in the bulk is initialized to equilibrium, where it roughly remains throughout time (again, bulk recombination is turned off in this simulation). As can be seen in Fig. 4b, as time progresses, the interstitial field becomes depleted adjacent to the bubble due to surface recombination, which in turn shrinks the bubble. The resultant gradient in ci drives more interstitials to the bubble surface, where they continue to annihilate the vacancies that constitute the bubble. Furthermore, the gas density in the bubble increases in time due to the decreasing bubble size. The plots in Fig. 4c illustrates the evolution in both cv, ci, and cg across the centerline of the bubble as it shrinks in time, indicating the gradient in ci, the concurrent bubble shrinkage denoted by the vacancy field, as well as the increasing gas density in the bubble. The use of an order parameter g to define the bubble phase allows us to distinguish the bubble and solid phases, and allows for
964
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
different conditions to apply in each region. For example, in the solid (where Eq. (2) applies), we assume gas atoms are substitutional, and the total vacancy, interstitial, gas, and atom concentrations add up to one. However, in the bubble phase (where Eq. (3) applies), we assume gas atoms occupy off-lattice sites, hence the concentrations can add up to larger than one.
4. Gas bubble evolution during irradiation The on-going production of point defects in irradiated materials eventually leads to the clustering of individual defect species, i.e., vacancies cluster to form voids and interstitials cluster to form dislocation loops. We have verified in Part I that the current phasefield model captures accurately the growth (shrinkage) rates of pre-existing voids due to the net fluxes of vacancies (interstitials) to their free surfaces, without irradiation. In the previous section of this paper, we have verified that the gas concentration within a bubble is allowed to increase or decrease according to the shrinkage or growth of the bubble, as well as the absorption of excess gas
atoms from the solid. In the current section, we will investigate the development and evolution of bubble microstructures during the on-going production of defects during irradiation. In these simulations, the initial domain is entirely solid material, with the equilibrium vacancy and interstitial concentrations (there is zero gas concentration initially). The defect formation energies are assigned values of Efv = 0.8 eV, Efi = 1.5 eV, and Efg = 1.2 eV; and the defect dif~v ¼ D ~ g ¼ 1, and D ~ i ¼ 4; in most fusivities are assigned values of D metals, Efi > Efv and Di > Dv, and in this section we assign these values accordingly, although the values do not correlate with a specific material. Furthermore, the bulk recombination term is turned on, Rbulk = 1 s1, in addition to the surface recombination term that was solely activated in the simulations of Section 3. Fig. 5 shows the evolution of g, cv, ci, and cg at progressive instances in time during an irradiation simulation. These vacancies and interstitials are produced stochastically in both space and time according to Part I. The gas production rate is assumed to by spatially uniform within the solid (g < 0.8) and Pgas has units of cg per time. As seen in Fig. 5, the ongoing production of vacancies and interstitials eventually leads to the nucleation of gas bubbles.
Fig. 5. Snapshots showing the nucleation and growth of bubbles in the presence of on-going cascades and gas production. The top-to-bottom rows represent the (a) g, (b) cv, (c) ci, and (d) cg fields. During the simulation, an interstitial production bias of 0.9 is assumed (i.e., there are 10% more vacancies introduced into the system), the net effect of which is similar to a dislocation bias. Within the voids, the vacancy concentration cv = 1, the interstitial concentration ci = 0, and the gas concentration varies.
965
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
Fig. 6. The evolution of (a), (d), (g) porosity, (b), (e), (h) bubble density, and (c), (f) average gas concentration inside the bubbles throughout time. The top row ((a), (b), (c)) represents a gas production rate of 1.0 104 cg/s, while the middle row ((d), (e), (f)) represents a gas production rate of 2.0 105 cg/s. In (g) and (h), the gas production rate is zero. As shown in (i), for increasing gas production rate, the bubble density at the end of stage II increases while the average bubble size decreases. Also, as shown in (c) and (f), the gas concentration within bubbles decreases during stage II (nucleation and growth) due to the relatively rapid bubble growth while the bubbles are still small.
The top-to-bottom rows in Fig. 5 correspond to the g, cv, ci, and cg contour plots at progressive instances in time (the contour intervals are indicated on the left side of the figure). As the voids nucleate and grow, the interstitial concentration inside the voids is reduced to zero, while the vacancy concentration evolves to a value of one. The gas concentration in the solid and within the bubbles varies from bubble to bubble and throughout the solid based on the local history of cascades and the spatial proximity of nearby bubbles. The production bias mentioned above ensures that a higher net flux of vacancies compared to interstitials arrives at the bubble surfaces, thus allowing bubble growth. During the intermediate stages, the processes of bubble nucleation and growth are occurring simultaneously while the bubble density is still low. In the latter stages, when the bubble density has increased, nucleation ceases and the evolution is controlled by growth and coarsening. As we show below, the gas concentration within the bubbles also evolves with time according to the production rate of gas atoms as well as the average bubble growth rate. In order to quantify the bubble microstructure throughout the simulation, we have calculated the porosity, void number density,
and average gas concentration within the bubbles (cg ) at regular intervals in time. This data is shown in Fig. 6. We vary two key parameters in this set of simulations: (i) the cascade-induced displacement rate (taking on values of either K = 2.3 103 dpa/s or K = 1.0 103 dpa/s), and (ii) the gas production rate (taking on values of Pgas = 0.0, 2.0 105, and 1.0 104 cg/s). The top row of Fig. 6 depicts the overall porosity (Fig. 6a), the bubble density (Fig. 6b), and the average gas concentration within Table 1 Values of k from Eq. (10) from fitting to the data shown in Fig. 6 and calculated from Eq. (11) for varying dpa rates and gas production rates. Also shown are the gas bubble nucleation rates. K (dpa/s) 3
2.3 103
1.0 10
Pgas (cg/s)
k (fit) (s3)
k (calculated) (s3)
J (k2 s1)
0.0 2.0 105 1.0 104
7
1.55 10 2.97 107 2.40 107
7
1.28 10 1.83 107 3.38 107
1.22 105 1.75 105 3.23 105
0.0 2.0 105 1.0 104
4.97 107 4.60 107 8.46 107
3.21 107 3.48 107 3.65 107
3.07 105 3.32 105 3.49 105
966
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
Fig. 7. (a) The nucleation rate of gas bubbles increases with increasing gas production rate for both displacement rates. This effect is more dramatic for the lower displacement rate (K = 1.0 103 dpa/s). (b) The average gas concentration within the bubbles decreases with increasing bubble diameter.
the bubbles (Fig. 6c) for a gas production rate of 1.0 104 cg/s. The porosity versus time plot (Fig. 6a) reveals that the bubble kinetics can be subdivided into three distinct regimes: stage I (incubation), stage II (nucleation and growth), and stage III (coarsening). This evolution is similar to our previous simulations in Part I. During stage II, the porosity follows an S-shaped curve characteristic of the Johnson–Mehl–Avarami equation [17]
h i 3 p ¼ po 1 exp kt ;
ð10Þ
where po is the porosity at the end of stage II, and k is a constant characteristic of the nucleation kinetics in the system, which for a 2D domain, is defined as
k¼
p _2 3
R J:
ð11Þ
In this equation, R_ is the void radius growth rate and J is the nucleation rate. The dashed lines in Fig. 6a, d, and g within the stage II regimes represent fits of the data to Eq. (10). For K = 2.3 103 dpa/s, these fits yield values of k = 4.97 107, 4.60 107, and 8.46 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. For K = 1.0 103 dpa/s, these fits yield values of k = 1.55 107, 2.97 107, and 2.40 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. This data is displayed in Table 1. The bubble nucleation rates, J, are obtained from Fig. 6b, e, and h, in which the dashed lines are linear fits representing the quasisteady-state nucleation regime. They yield nucleation rates of J = 1.22 105 k2s1 (Pgas = 0.0), 1.75 105 k 2s1 (Pgas = 2.0 105), and 3.23 105 k2s1 (Pgas = 1.0 104) for K = 1.0 103 dpa/s; and J = 3.07 105 k2s1 (Pgas = 0.0), 3.32 105 k2s1 (Pgas = 2.0 105), and 3.49 105 k2s1 (Pgas = 1.0 104) for K = 2.3 103 dpa/s. Please see Table 1. Taking into consideration the bubble growth rate of 0.1 k/s (obtained from measurements on the smallest voids in the simulations), the analytical values of k are obtained from Eq. (11). For K = 2.3 103 dpa/s, they are k = 3.21 107, 3.48 107, and 3.65 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. For K = 1.0 103 dpa/s, they are k = 1.28 107, 1.83 107, and 3.38 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. All of these values are within a factor of two of the fitted values obtained from Fig. 6a, d, and g. Fig. 6c shows the average gas concentration inside the bubbles (cg ) as a function of time (this value is averaged over all the bubbles in the domain). Interestingly, during the stage II regime, in which the bubble growth rate is highest, the gas density decreases from the initial value at the very onset of nucleation. This decrease occurs due to a faster bubble absorption rate of vacancies compared to gas atoms. Thus, the gas density decreases. During stage
III, however, cg begins to increase gradually and eventually levels off. This trend occurs for both displacement rates, however the magnitude of cg is lower for the higher displacement rate (K = 2.3 103 dpa/s) due to the fact that nucleation occurs earlier when the overall gas concentration is lower. The middle row of Fig. 6 displays similar results, although for a lower gas production rate, 2.0 105 cg/s. The results are qualitatively similar to the higher gas production rate, with the exception that the nucleation rates have decreased for both displacement rates. In addition, the cg values are generally lower, which is to be expected for a lower gas production rate. On the bottom row, Fig. 6g and h shows the porosity and bubble density evolution for the case of zero gas production. We see again that the bubble nucleation rates have decreased relative to the cases in which gas is produced. In Fig. 6i, the bubble density and the average bubble diameter are plotted as a function of gas production rate. These values are taken at the end of stage II, when nucleation ceases. We see a general trend of increasing bubble density and decreasing bubble diameter with increasing gas production rate. This is in agreement with experimental observations [1]. The generally-accepted explanation for this is that gas atoms act as trapping sites for vacancies, and therefore they represent nucleation sites for bubbles. In our simulations, the presence of gas atoms increases the free energy of the solid phase further away from equilibrium. Therefore, a smaller thermal fluctuation, or cascade event, is required to create a nucleation event. In addition to an increase in bubble density with increasing gas production rate, the bubble nucleation rates (or the slopes of Fig. 6b, e, and h) also increase with increasing gas production rate. This is shown in Fig. 7a for both displacement rates. Interestingly, the dependency of J (as well as the bubble density and bubble size) on the gas production rate is higher for the lower displacement rate. This is evident by the larger absolute values of the slopes in Figs. 6i and 7a for K = 1.0 103 dpa/s, relative to K = 2.3 103 dpa/s. Conceivably, this is perhaps due to the fact that bubble nucleation initiates early for the higher displacement rate, when the gas concentration is smaller. In Fig. 7b, the average gas concentration within the bubbles, cg , is plotted as a function of the bubble diameter. This data is for a gas production rate of 2.0 105 cg/s and a displacement rate of K = 2.3 103 dpa/s. For these conditions, cg decreases monotonically with increasing bubble size. This result is in accord with the data shown in Fig. 6f, in which cg decreases with time. As mentioned earlier, this decrease in cg with increasing bubble size is explained by a higher flux of vacancies than gas atoms to the growing bubbles. Had the gas production rate been considerably higher than the rate of increase in vacancy supersaturation (taking into consideration vacancy-interstitial recombination), this trend may
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
967
actually differ; although in this work we did not simulate such conditions.
5. Intergranular gas bubbles Extended microstructural defects such as GBs and dislocations strongly affect bubble structure evolution in irradiated materials due to their role as point defect sinks, as well as the inherent fact that impurity atom formation energies are lower at these features compared to the bulk. It is commonly believed that GBs act as heterogeneous nucleation sites for fission-gas bubbles; this being widely supported by many experimental observations of intergranular bubbles that are many times larger than their intragranular counterparts [2,3,18]. As intergranular bubbles grow to micronsizes and larger, their migration and coalescence on the GB faces, edges, and corners leads ultimately to an interconnected porosity that allows rapid fission gas release to the exterior of the material [2,3]. In addition, due to the triple junctions between bubble free surfaces and the GB, the inherent bubble shape is lenticular rather than spherical in order to minimize the interfacial energies of GB and free surfaces. Furthermore, the presence of intergranular bubbles effectively pins GB migration, according to the Zener drag theory [10]. As we show below, the current model is able to capture these important distinctions of intergranular bubble physics. Fig. 8 shows the interaction of a single GB with a stationary bubble that does not evolve (i.e., change shape) in time. (Here, the concentration fields and the g-field are not updated in time; only the /-fields representing grains 1 and 2 are solved for in time). In Fig. 8a, the equilibrium GB shapes are shown for different relative positions of the bubble. We see that, for this stationary bubble, the GB always intersects the bubble surface at a 90° angle. (Note that in these figures, the GB is not made to migrate across the domain, rather the initial position of the GB is changed from left to right). In Fig. 8b, the pinning energy that the bubble exerts on the GB is plotted as a function of bubble radius. This pinning energy is the difference between the system shown in the center of Fig. 8a and that shown on the far left of Fig. 8a. According to Zener [10,16], it is equal to the amount of GB that is eliminated by the presence of the bubble multiplied by the GB energy. For 2D, this is equal to 2rrgb. The data points in Fig. 8b correspond to the re-
Fig. 8. (a) GB pinning due to a stationary bubble (that does not evolve in time). (b) The pinning energy is equal to the GB energy, gb, multiplied by the length (for 2D) of the missing GB segment corresponding to the bubble diameter. The open circles represent values calculated from the model and the dashed line is the theory according to Zener [10].
Fig. 9. Evolution of the shape of a gas bubble on a GB from circular to lenticular (equilibrium) due to the GB energy. The (a) /, (b) cv, and (c) cg fields are shown throughout the evolution. (d) The dihedral angle, h, evolves exponentially from 90° to approximately 120°.
sults obtained from simulation, and the dashed line corresponds to the theory, showing a close agreement. This result is entirely due to the interaction in the free energy term of Eq. (4) between the order parameter representing the bubbles, g, and the order parameters representing the grains, /. With verification that the current model accurately captures the pinning effect of a gas bubble on a GB, we now focus on the influence of a GB on the shape of a gas bubble. As mentioned above, an intergranular bubble will assume a lenticular shape due to the Herring relation [19] governing the dihedral angles existing at the triple junction between the free surfaces and the GB. To verify that we correctly capture this, we have performed a simulation of a bicrystal with a single, initially-circular gas bubble residing on the GB. As shown in Fig. 9a, the bubble evolves from a circular shape to a lenticular shape as time progresses. Fig. 9b and c illustrates that the vacancy and gas atom concentrations, respectively, evolve accordingly. The evolution of the triple junction angle, h, is shown versus time in Fig. 9b, illustrating an exponential approach from h = 90° to h 120°. With the above analysis of the interaction between pre-existing gas bubbles and GBs, we now focus on bubble nucleation and growth in polycrystalline microstructures. Here, we re-activate the irradiation terms. Fig. 10 illustrates the nucleation and growth of intergranular gas bubbles due to the ongoing production of Frenkel pair defects and gas atoms. The simulation domain is again periodic in both directions (lx = ly = 6 lm) and contains 40 grains (i.e., P = 40 in Eq. (4), with isotropic GB properties, including energy and mobility) with an average diameter of d = 1.1 lm. The gas atom formation energy is reduced in the GBs by the following equation
f ¼ Ef U; E g g
ð12Þ
968
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
where the factor U is defined in Part I and is equal to 1.0 in the grain interiors and as low as 0.5 in the GBs and triple junctions. This results in preferential segregation of the gas atoms to the GB regions. The formation energies are Efv = 0.8 eV, Efi = 1.5 eV, an-
dEfg = 1.2 eV, and the recombination rates are Rbulk = 1.0 s1 and Rsurf = 5.0 s1. Fig. 10 shows sequential snapshots of the preferential nucleation of bubbles on GB and TJ regions, resulting from the enhanced
Fig. 10. (a)–(d) Snapshots throughout time of the nucleation and growth of intergranular gas bubbles in a polycrystalline grain structure, while (e) is a close-up view. Bubbles existing on GBs are lenticular shaped, whereas bubbles on triple junctions have a curved triangular shape. (f) A microscopy image of intergranular bubbles in UO2 taken from [2].
Fig. 11. Characteristics of intergranular gas bubbles for different grain sizes. The (a) porosity, (b) bubble density, and (c) gas concentration within the bubbles does not vary substantially for the different grain sizes. However, as the grain size increases, the (d) average bubble diameter increases, the (e) average GB bubble spacing decreases, and the (f) GB bubble coverage increases. (The data for (d), (e), and (f) are taken at the end of the simulations, t = 800 s).
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
gas atom concentration at these extended defects. As can be seen in the close-up view in Fig. 10e, bubbles existing on GBs are lenticular shaped, whereas bubbles on TJs have a curved triangular shape. In Fig. 10f, a micrograph image of intergranular bubbles in UO2 (taken from Ref. [2]) illustrates the qualitative similarities between the bubble structures in terms of the relative size of the bubbles with respect to the grain size, as well as the bubble shapes. Also note that, in the simulations, bubbles within the grain interiors have also begun to nucleate, although their shapes remain circular, as would be expected. In Fig. 11, we analyze the intergranular bubble characteristics for various grain sizes. The top row, Fig. 11a–c shows the porosity, bubble density, and cg throughout time for an average grain size between 37 k < d < 111 k. Each of these simulations are run with identical irradiation conditions of K = 1.0 103 dpa/s and Pgas = 8.2 106 cg/s. We observe that the porosity, bubble density, and cg do not vary substantially with the changing grain size (the exception is the smallest grain size, d = 37 k, in which the increase in porosity has a slightly slower rate, although bubble nucleation rate is similar to the other grain sizes). The average gas concentration inside the bubbles (Fig. 11c) grows most rapid at the onset of nucleation, and saturates at the later stages, in a very similar manner for each grain size. Furthermore, the bubble nucleation rates, obtained from Fig. 11b, do not vary substantially with grain size. However, when we plot the average GB bubble diameter (Fig. 11d), the average spacing on the GB lines (Fig. 11e), and the GB bubble coverage (percentage of GB length occupied by bubbles) (Fig. 11f), we see clear trends in each with respect to the grain size. With increasing grain size, the intergranular bubble diameter increases, the intergranular bubble spacing decreases, and the GB bubble coverage increases. Considering that the porosity is basically equal for each grain size at the end of the simulations, these trends can be rationalized by the fact that the total GB density is lower for the larger grain size. Therefore, for the same porosity and overall bubble density, and considering that virtually all bubbles reside on the GBs, bubbles should be closer spaced and larger for larger grain sizes. This, in turn, results in a higher GB coverage of bubbles. Interestingly, it appears that in the polycrystal specimens, the nucleation stage has yet to saturate, yielding only the first half of an S-shaped curve (Fig. 11a). We have observed that at about the time that nucleation ceases on the GBs, nucleation events begin to occur in the grain interiors. So, the system as a whole does not decelerate the porosity development. Eventually, nucleation in both the grain interiors and the GBs will stop, at which time the system will truly enter a coarsening stage.
6. Conclusions
turn, dictated by this change in bubble size, as well as by the flux of gas atoms to the bubble. In single crystal simulations, the porosity evolved in three distinct stages during irradiation: stage I (incubation), stage II (nucleation and growth), and stage III (coarsening). The bubble nucleation rate is shown to increase with increasing displacement rate as well as gas production rate. Increasing gas production rate is also shown to increase the bubble density and decrease the bubble size. Furthermore, at the early stages of bubble growth, the gas density within the bubbles (cg ) decreases due to a higher flux of vacancies, compared to gas atoms, to the bubble surfaces. In polycrystal simulations, the gas atom formation energy is decreased within the GB and triple junction regions, leading to preferential segregation and heterogeneous nucleation at these features. During intergranular bubble growth, the bubbles retain their equilibrium lenticular shapes on GB planes, and curved triangular shapes on triple junctions. As the grain size is increased, the porosity and nucleation rates remain constant; however, the average bubble size and GB bubble coverage increase, while the average GB bubble spacing decreases. Acknowledgements This work was supported through the INL Laboratory Directed Research and Development program under DOE Idaho Operations Office contract no. DE-AC07-051D14517V. We also acknowledge support from two DOE programs: the Energy Frontier Research Center (EFRC) for Materials Science of Nuclear Fuel and the Nuclear Energy Modeling and Simulation (NEAMS) program (Cetin Unal as Fuels IPSC program lead). In addition, this manuscript benefited from technical discussions with Michael Tonks (INL) and Srujan Rokkam (FSU). Appendix A. Free energy of ideal gas The free energy of the gas phase with respect to the gas concentration is derived below. We start with the total Gibb’s free energy [20] assuming an ideal gas law
Gðp; TÞ ¼ ngo þ nRT ln p;
ðA1Þ
where p is the pressure, n is the number moles, go is a reference molar chemical potential, and R is the universal gas constant. This equation can be re-written to obtain a volumetric energy density
gðp; TÞ ¼
n V
Av log þ bub
kB T ln p; bub
nAv
V
ðA2Þ
where Av is Avogadro’s number, Vbub is the volume of the bubble, and log is a reference atomic chemical potential of the gas. In terms of gas concentration, the equation is now written as
gðp; TÞ ¼ cg log þ cg kB T ln p: The phase-field model presented here extends our previous model in Part I for pure metals by incorporating fission gas atoms in a manner that distinguishes their point defect energetics and kinetics in the solid and gas bubble phases. This is accomplished by the use of a long-range order parameter, g(r,t), that spatially differentiates the two phases throughout the simulation domain. This is a key feature that distinguishes the current model with a previous phase-field model [12] for fission gas behavior in irradiated materials, in which a non-conserved order parameter is omitted. In addition, the g-field in the current model allows for a proper interaction between gas bubbles and GBs in terms of the GB pinning energies and the equilibrium bubble shapes. The rate of bubble growth or shrinkage, in the bulk or on GBs, is controlled solely by the flux of vacancies and interstitials to the bubble surfaces. The gas concentration within the bubbles is, in
969
ðA3Þ
Using the ideal gas law, the gas pressure is equal to
p¼
nAv
V
bub
kB T ¼ cg kB T;
ðA4Þ
and by substituting (A4) into (A3), we arrive at
f b ðcg Þ ¼ log cg þ cg kB T ln cg þ cg kB T lnðkB TÞ:
ðA5Þ
References [1] H. Trinkaus, B.N. Singh, J. Nucl. Mater. 323 (2003) 229. [2] I. Zacharie, S. Lansiart, P. Combette, M. Trotabas, M. Coster, M. Groos, J. Nucl. Mater. 255 (1998) 85. [3] S.B. Fisher et al., J. Nucl. Mater. 306 (2002) 153. [4] H. Ullmaier, Nucl. Fusion 24 (1984) 1039.
970
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
[5] D.R. Olander, Fundamental Aspects of Nuclear Reactor Fuel Elements (US-DOE, 1985). [6] K. Morishita, R. Sugano, J. Nucl. Mater. 353 (2006) 52. [7] J.A. Turnbull, J. Nucl. Mater. 38 (1971) 203. [8] M.H. Wood, J. Nucl. Mater. 82 (1979) 264. [9] J. Rest, J. Nucl. Mater. 321 (2003) 305. [10] C. Smith, Trans. Metall. Soc. AIME 175 (1948) 15. [11] L.Q. Chen, Y. Yang, Phys. Rev. B 50 (1994) 15752. [12] S. Hu, C.H. Henager, H.L. Heinisch, M. Stan, M.I. Baskes, S.M. Valone, J. Nucl. Mater. 392 (2009) 292. [13] P.C. Millett, S. Rokkam, A. El-Azab, D. Wolf, Mod. Sim. Mater. Sci. Eng. 17 (2009) 064003.
[14] J.W. Cahn, J.E. Hilliard, J. Chem. Phys. 28 (1958) 258. [15] L.D. Landau, E.M. Lifshitz, Statistical Physics, third ed., Elsevier, Burlington, MA, 1980. [16] N. Moelans, B. Blanpain, P. Wollants, Acta Mater. 53 (2005) 1771. [17] R.W. Balluffi, S.M. Allen, W.C. Carter, Kinetics of Materials, Wiley, Hoboken, 2005. [18] S. Fréchard, M. Walls, M. Kociak, J.P. Chevalier, J. Henry, D. Gorse, J. Nucl. Mater. 393 (2009) 102. [19] C. Herring, in: The Physics of Powder Metallurgy, McGraw-Hill, New York, 1951, pp. 143–179. [20] D.R. Gaskell, DR Introduction to Metallurgical Thermodynamics, second ed., Hemisphere, Publishing Corporation, 1981.
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
Phase-field simulation of irradiated metals Part II: Gas bubble kinetics Paul C. Millett a,⇑, Anter El-Azab b, Dieter Wolf c a
Idaho National Laboratory, Idaho Falls, ID 83415, USA Florida State University, Tallahassee, FL 32310, USA c Argonne National Laboratory, Argonne, IL 60439, USA b
a r t i c l e
i n f o
Article history: Received 14 September 2010 Received in revised form 14 October 2010 Accepted 27 October 2010 Available online 4 December 2010 Keywords: Phase-field simulation Irradiation damage Radiation damage
a b s t r a c t The phase-field model developed in Part I of this work is expanded to include fission gas generation, diffusion, and segregation within bubbles nucleated both homogeneously and heterogeneously along grain boundaries. Illustrative results are presented that characterize bubble growth and shrinkage, as well as the bubble density, size and nucleation rate as a function of the irradiation conditions. Finally, intergranular bubble characteristics such as shape, pinning energy and bubble density are investigated. Ó 2010 Elsevier B.V. All rights reserved.
1. Introduction The accumulation of gaseous species in nuclear fuel and structural materials leads to deleterious microstructural changes that are of utmost concern to both fission and fusion reactor applications [1]. Inert gas elements – such as helium, xenon and krypton – are generated due to fission reactions in nuclear fuel [2,3]. In addition, transmutation via (n,a) reactions in the metallic cladding leads to helium production [4]. In general, the precipitation of gaseous species into gas-filled bubbles leads to volumetric swelling in both structural materials and fuels, as well as high-temperature embrittlement and overall loss of mechanical strength in metallic components [1]. The processes of nucleation, growth, migration, and coalescence of gas-filled bubbles are all governed by (i) the dynamic production of point defects (i.e., vacancies, self-interstitials, and gas atoms) by irradiation and (ii) the pre-existing microstructural features of the material (i.e., grain boundaries (GBs), dislocation networks, and precipitates). Such features can themselves evolve as a result of interaction with point defects. A precise understanding of the inter-dependency of the point and extended defects, and their kinetics, as they relate to gas behavior still remains a very important subject in the field of nuclear reactor materials. According to conventional nucleation theory, in the early stages of irradiation, gas atoms diffuse throughout the crystal and begin
⇑ Corresponding author. E-mail address: [email protected] (P.C. Millett). 0927-0256/$ - see front matter Ó 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.commatsci.2010.10.032
to form clusters due to their insolubility in the solid [5]. This clustering can initiate either by chance encounters of the wandering gas atoms (i.e., homogeneous nucleation), or by gas-atom trapping at microstructural segregation sites such as GBs (i.e., heterogeneous nucleation). The subsequent growth or shrinkage of a gas bubble is governed by the fluxes of vacancies, interstitials, and gas atoms to or away from the bubble surface [6]. Basically, the time rate of change of the volume of a bubble is equal to the difference in the rates at which vacancies and interstitials are absorbed at the bubble surface times the (atomic) volume carried by each of these point defects [5]. In addition, gas-atom re-solution can limit bubble growth whereby fission fragments directly or indirectly cause gas atoms to be lost from the bubble [7,8]. If, however, gas bubbles are located on sinks (such as GBs) that recollect the lost gas atoms, they may grow indefinitely to the point of interconnection, thus leading to rapid fission gas release to the exterior of the material in the case of fuel [9]. Conversely, the existence of intergranular gas bubbles can, in turn, influence GB migration due to Zener pinning [10], thus ultimately affecting the microstructural evolution of the material. In the current paper, we expand the phase-field model developed in Part I to include gas atom evolution, thus allowing investigations of bubble nucleation and growth in irradiated metals. This phase-field model is based on the dynamics of generation, diffusion and reactions of vacancies, gas atoms, and self-interstitials (to be referred to merely as interstitials from here on) in a single component, polycrystalline metal. A polycrystal phase-field model [11] is incorporated, thus capturing lenticular-shaped intergranular bubble formation as well as bubble-induced Zener pinning of
961
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
GB migration. This model differs from that of Hu et al. [12], in that the order parameter g (introduced in Part I) differentiates the bubbles from the solid, and allows different thermodynamic and kinetic descriptions to be applied to each phase throughout space.
vacancies, interstitials, and gas atoms in the matrix; it is written in the form
f solid ¼ Efv cv þ Efi ci þ Efg cg þ kB T½cv lnðcv Þ þ ci lnðci Þ þ cg lnðcg Þ þ ð1 cv ci cg Þ lnð1 cv ci cg Þ;
2. Simulation approach In this model, we consider a polycrystalline metal in which bubbles are considered to be clusters of vacancies that contain a finite amount of gas atoms. (The current model is an extension of our previous work [13], which also considered single component metals.) The phase-field approach is based on expressing the total free energy functional F of the heterogeneous material in terms of the free energy of each of its phases and interfaces. Consistent with the Cahn–Hilliard definition of free energy for a non-uniform system [14], we describe the free energy functional of a polycrystalline solid consisting of both matrix and bubble phases, in terms of the vacancy concentration field cv, the interstitial concentration field ci, gas atom concentration field cg, an order parameter field g that distinguishes solid and void phases, and a set of order parameters /1?P that represent discrete crystallographic (grain) orientations
Z
hðgÞf solid ðca Þ þ jðgÞf bubble ðca Þ þ f pc ðg; /1!P Þ þ f ðca ; g; /1!P Þ dV:
F¼N
V grad
ð1Þ
In this expression, ca represents the set of concentration fields (cv, ci, and cg), and N is the number of lattice sites per unit volume of the crystal and thus all terms under the integral sign represent free energy per lattice site. The first three terms in the righthand-side correspond to the bulk system energy. The interfacial energy contributions, due to the bubble free surfaces and GBs, are accounted for by the gradient term. The excess free energy of the solid due to point defects fsolid(ca) is derived in terms of the enthalpy and entropic contributions of
ð2Þ
where Efv , Efi , and Efg are the vacancy, interstitial, and gas formation energies, respectively, kB is Boltzmann’s constant, and T is the absolute temperature. The entropic term in the brackets in Eq. (2) represents the configurational entropy of a defect-disordered crystal, shown schematically in Fig. 1, consisting of sites that are occupied either by a host atom, gas atom, vacancy, or a dumbbell self-interstitial defect. Note that this equation assumes that the gas atoms are substitutional defects that occupy regular lattice sites, which is typical of large fission gases. The free energy of the bubble phase is defined as
h i f bubble ¼ ðcv 1Þ2 þ c2i h i þ log cg þ cg kB T ln cg þ cg kB T lnðkB TÞ :
ð3Þ
The term in the first set of brackets was presented in Part I, in which an energy minimum is located precisely at cv = 1 and ci = 0. At this minimum, the free energy is exactly equal to zero. The term in the second set of brackets is the free energy of the gas phase within a bubble according to the ideal gas law. Its derivation is given in Appendix A. Here, the free energy is described by a gas concentration rather than the conventional gas pressure to be consistent with the use of a concentration variable in the matrix. The use of the ideal gas law is perhaps not appropriate for bubbles with a high density of gas, in which case the Van der Waals equation is more accurate [15], nevertheless, we make use of this derivation for the current work. To incorporate GBs and their mutual interactions with gas bubbles into the model, we employ a conventional polycrystal phase-field model, in which a set of order parameters /1?P are used to define the set of crystallographic grain orientations. The free energy of the polycrystalline material is defined as
f
polycrystal
ðg; /1;...;P Þ ¼ 3A
P X k¼1
!2 /2k
2
þg
4A
P X
! /3k
3
þg
:
ð4Þ
k¼1
For g = 0 (i.e., in the solid phase), this equation creates P energy wells located at {uk=l = 1, uk–l = 0}. For g = 1 (i.e., in the bubble phase), this equation has one energy well located at /1?P = 0, thus all grain order parameters go to zero inside a bubble. We note that this polycrystal model is qualitatively similar (the functional form differs) to that of Moelans et al. [16] that captures the pinning of GBs by particles; however, the current model allows for a dynamically-evolving bubble phase. To our knowledge, this is the first phase-field model explicitly taking into account dynamically-evolving gas bubbles interacting with planar defects, such grain boundaries. The gradient energy contribution, contributing primarily to gas bubble free surface energies and GB energies, are accounted for by the term
f grad ¼
Fig. 1. Schematic plot illustrating the quasi-lattice representation of vacancy, selfinterstitial, and gas atom defects. The configurational entropy of this system is included in Eq. (2), in which each lattice site can be occupied by either an atom, a vacancy, a gas atom, or a dumbbell self-interstitial.
X jw 2
jrwj2 ;
ð5Þ
where w represents the set of concentration fields (cv, ci, cg) and order parameter fields (g,/1?P). Following the standard procedure in the phase-field approach, the kinetic equations for the spatial and temporal evolution of the gas concentration field (in addition to the vacancy and interstitial concentrations derived in Part I) has the form
962
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
@cg 1 dF þ ng ðr; tÞ þ Pg ðr; tÞ: ¼ r Mg r N dcg @t
ð6Þ
The gas atom mobility is defined as
Mg ¼
Dg c g ; kB T
ð7Þ
where Dg is the gas atom diffusivity in the matrix. We note that, for a substitutional gas atom (as assumed in Eq. (2)), the gas mobility should depend on cv in addition to cg. This paper does not account for this point. Because these concentrations remain fairly dilute in the solid, this approximation is valid to a first order analysis. In the bubble phase, the gas can be assumed to instantly equilibrate when considering solid-state diffusion time scales. To account for this, the gas concentration inside each bubble is evenly distributed to provide a flat profile within each bubble. The gas concentration within the bubbles will vary from bubble to bubble. Bubble re-solution due to fission-fragment bombardment [5] is not included in the current model. This process can conceivably be including by adding a term that re-distributes fission gas from the bubbles to the surrounding solid regions. The exclusion of fission gas re-solution will result in higher bubble nucleation rates and growth rates predicted by the simulations. The order parameter fields that define the distinct grain orientations /i also evolve according to a set of Allen–Cahn equations
@/i ðr; tÞ dF ¼ L/ ; @t d/i
ð8Þ
where L/ is the grain boundary mobility. In the present work, Eq. (6) can be expressed in non-dimensional units following Part I leading to
@cg ~ 1 dF þ ~fð~r; ~tÞ þ P ~ M ~ g ð~r; ~tÞ ~ gr ¼r ~ dcg @~t N
ð9Þ
~ ¼ l r is the dimensionless gradient operator. Energy is where r kept in units of eV. In all of the simulations below, the grid spacing is Dx = k, the time step is Dt = 0.001 s, and the gradient terms are jv = ji = jg = jg = j/ = 1.0 eV k2. In Sections 4 and 5 below, we vary the Frenkel pair production rate and the gas atom production rate independently. In nuclear reactor conditions, these two are coupled to some degree; however, we disregard this coupling which is more applicable to dual-beam ion irradiation. Varying these production rates independently, as will be shown below, can lead to interesting results. 3. Gas bubble growth and shrinkage The size evolution of a gas bubble existing in a crystalline solid is governed by the flux rate of vacancies and interstitials to its surface. A net absorption of vacancies by the bubble will result in bubble growth, while conversely a net absorption of interstitials results in bubble shrinkage. Furthermore, the density of gas atoms within the bubble will fluctuate in time based on (i) the flux of gas atoms to the bubble and (ii) the change in the bubble volume itself. In order to validate that the current model can appropriately capture bubble growth and shrinkage based on the vacancy and selfinterstitial concentrations in the surrounding solid, as well as the direct dependence of this on the bubble’s gas concentration, we have simulated three benchmark cases, which are described below. For the first benchmark case, we have initialized the system to consist of a single void with a diameter of 20 k located at the center of the domain with periodic boundary conditions, as shown in Fig. 2a. Initially, the void itself contains no gas atoms, however, the surrounding solid does contain a finite concentration of gas atoms. The vacancy and interstitial concentrations in the solid are at their equilibrium values, respectively, thus we do not expect
Fig. 2. Increase of gas concentration inside an initially-empty void surrounded by a solid with a supersaturated gas concentration. The top figures show the evolution of (a) vacancy and (b) gas atom concentration fields. (c) Plots of the conserved concentration fields, cv, ci, cg, along a cross-sectional slice through the centerline of the bubble throughout time. The absorption of gas atoms into the void creates a gradient in the gas concentration in the adjacent solid region.
to observe any change in bubble size. For the simulations shown in Figs. 2–4, the energies of formation of the vacancies and interstitials in the crystal are taken as Efv ¼ Efi = 0.8 eV corresponding to eq equilibrium intrinsic defect concentrations of ceq v ¼ ci = 6.9 4 10 (the thermal energy in all simulations within this paper is taken as kBT = 0.11 eV, corresponding to a temperature of T = 1276 K). The gas atom defect energy in the solid is also taken as Efg = 0.8 eV. The reference chemical potential of the gas atoms in the gas phase is log = 0.4 eV throughout this paper (this value was arbitrarily chosen, with the condition that it be less than Efg ). The bulk recombination term is turned off (i.e., Rbulk = 0) so that cv in the bulk remains constant, however surface recombination is activated by assigning a value of Rsurf = 5 s1. The supersaturation of the gas in the solid is Sg = 120 (this supersaturation is calcuf lated with the supposed ceq g based on Eg , although there is no actual equilibrium concentration of impurity gas atoms in a crystal). Fig. 2a and b are contour plots that illustrate how the vacancy concentration (3a) and the gas atom concentration (3b) evolve in time. As can be seen, the gas concentration within the void increases in time due to the diffusional flux of gas atoms through the solid towards the bubble surface. The bubble volume, however, does not change in time due to the fact that there is not an excess quantity of vacancies or interstitials in the solid. Therefore, as more gas atoms enter the bubble, which has a constant size, the gas density in the bubble increases. This is shown quantitatively in Fig. 2c
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
963
which shows plots of the cv, ci, and cg concentration fields along a cross-sectional slice through the centerline of the bubble throughout time. In addition to the increase in gas concentration inside the bubble, in the nearby solid, a gradient in gas concentration develops due to the depletion of gas atoms from the solid to the bubble. This gradient drives an on-going flux of gas atoms to the bubble surface, which continues until the chemical potential of gas atoms in the bubble equals that in the solid. The last snapshot in Fig. 2 does not correspond to the fully equilibrated system, i.e., more time is needed before the gas concentration is done evolving (this is also true of Figs. 3 and 4 below). The second benchmark case, shown in Fig. 3, consists of a gas bubble that, initially, has a smaller diameter of 10 k. The initial vacancy concentration in the solid is supersaturated by Sv = 120. In addition, the gas concentration in the bubble is initialized at a high value corresponding to the case of an over-pressurized bubble. As shown in Fig. 3a and b, the bubble grows in size due to the absorption of vacancies from the solid. Furthermore, as the bubble’s area increases, the concentration of gas inside the bubble decreases. The total concentration of gas is conserved, thus this decrease in bubble gas concentration is solely due to the increase in bubble size, which is in turn governed by the diffusional flux of vacancies to the bubble surface. Contrary to the above case of void growth in a supersaturated vacancy field, a system that is supersaturated with interstitials will
Fig. 4. Shrinkage of a bubble due to a supersaturated self-interstitial concentration in the solid. The gas concentration in the bubble, which is initially under-saturated, increases due to the decrease in bubble volume. The top figures show the evolution of (a) vacancy and (b) gas atom concentration fields. (c) Plots of the conserved concentration fields, cv, ci, cg, along a cross-sectional slice through the centerline of the bubble throughout time.
Fig. 3. Growth of a bubble due to a supersaturated vacancy concentration in the solid. The gas concentration in the bubble, which is initially over-saturated, decreases due to the increase in bubble volume. The top figures show the evolution of (a) vacancy and (b) gas atom concentration fields. (c) Plots of the conserved concentration fields, cv, ci, cg, along a cross-sectional slice through the centerline of the bubble throughout time.
drive a void to shrink in size due to the continual occupation of empty lattice sites on the void surface by interstitials. This is essentially identical to vacancy-interstitial recombination in the bulk, albeit the rate may be higher due to the unique structure of the free surface. As mentioned in Part I, we consider the overall recombination term to have contributions from both the bulk and the free surfaces. Fig. 4 shows the shrinkage of a bubble surrounded by a solid containing a supersaturation of interstitials. The vacancy concentration in the bulk is initialized to equilibrium, where it roughly remains throughout time (again, bulk recombination is turned off in this simulation). As can be seen in Fig. 4b, as time progresses, the interstitial field becomes depleted adjacent to the bubble due to surface recombination, which in turn shrinks the bubble. The resultant gradient in ci drives more interstitials to the bubble surface, where they continue to annihilate the vacancies that constitute the bubble. Furthermore, the gas density in the bubble increases in time due to the decreasing bubble size. The plots in Fig. 4c illustrates the evolution in both cv, ci, and cg across the centerline of the bubble as it shrinks in time, indicating the gradient in ci, the concurrent bubble shrinkage denoted by the vacancy field, as well as the increasing gas density in the bubble. The use of an order parameter g to define the bubble phase allows us to distinguish the bubble and solid phases, and allows for
964
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
different conditions to apply in each region. For example, in the solid (where Eq. (2) applies), we assume gas atoms are substitutional, and the total vacancy, interstitial, gas, and atom concentrations add up to one. However, in the bubble phase (where Eq. (3) applies), we assume gas atoms occupy off-lattice sites, hence the concentrations can add up to larger than one.
4. Gas bubble evolution during irradiation The on-going production of point defects in irradiated materials eventually leads to the clustering of individual defect species, i.e., vacancies cluster to form voids and interstitials cluster to form dislocation loops. We have verified in Part I that the current phasefield model captures accurately the growth (shrinkage) rates of pre-existing voids due to the net fluxes of vacancies (interstitials) to their free surfaces, without irradiation. In the previous section of this paper, we have verified that the gas concentration within a bubble is allowed to increase or decrease according to the shrinkage or growth of the bubble, as well as the absorption of excess gas
atoms from the solid. In the current section, we will investigate the development and evolution of bubble microstructures during the on-going production of defects during irradiation. In these simulations, the initial domain is entirely solid material, with the equilibrium vacancy and interstitial concentrations (there is zero gas concentration initially). The defect formation energies are assigned values of Efv = 0.8 eV, Efi = 1.5 eV, and Efg = 1.2 eV; and the defect dif~v ¼ D ~ g ¼ 1, and D ~ i ¼ 4; in most fusivities are assigned values of D metals, Efi > Efv and Di > Dv, and in this section we assign these values accordingly, although the values do not correlate with a specific material. Furthermore, the bulk recombination term is turned on, Rbulk = 1 s1, in addition to the surface recombination term that was solely activated in the simulations of Section 3. Fig. 5 shows the evolution of g, cv, ci, and cg at progressive instances in time during an irradiation simulation. These vacancies and interstitials are produced stochastically in both space and time according to Part I. The gas production rate is assumed to by spatially uniform within the solid (g < 0.8) and Pgas has units of cg per time. As seen in Fig. 5, the ongoing production of vacancies and interstitials eventually leads to the nucleation of gas bubbles.
Fig. 5. Snapshots showing the nucleation and growth of bubbles in the presence of on-going cascades and gas production. The top-to-bottom rows represent the (a) g, (b) cv, (c) ci, and (d) cg fields. During the simulation, an interstitial production bias of 0.9 is assumed (i.e., there are 10% more vacancies introduced into the system), the net effect of which is similar to a dislocation bias. Within the voids, the vacancy concentration cv = 1, the interstitial concentration ci = 0, and the gas concentration varies.
965
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
Fig. 6. The evolution of (a), (d), (g) porosity, (b), (e), (h) bubble density, and (c), (f) average gas concentration inside the bubbles throughout time. The top row ((a), (b), (c)) represents a gas production rate of 1.0 104 cg/s, while the middle row ((d), (e), (f)) represents a gas production rate of 2.0 105 cg/s. In (g) and (h), the gas production rate is zero. As shown in (i), for increasing gas production rate, the bubble density at the end of stage II increases while the average bubble size decreases. Also, as shown in (c) and (f), the gas concentration within bubbles decreases during stage II (nucleation and growth) due to the relatively rapid bubble growth while the bubbles are still small.
The top-to-bottom rows in Fig. 5 correspond to the g, cv, ci, and cg contour plots at progressive instances in time (the contour intervals are indicated on the left side of the figure). As the voids nucleate and grow, the interstitial concentration inside the voids is reduced to zero, while the vacancy concentration evolves to a value of one. The gas concentration in the solid and within the bubbles varies from bubble to bubble and throughout the solid based on the local history of cascades and the spatial proximity of nearby bubbles. The production bias mentioned above ensures that a higher net flux of vacancies compared to interstitials arrives at the bubble surfaces, thus allowing bubble growth. During the intermediate stages, the processes of bubble nucleation and growth are occurring simultaneously while the bubble density is still low. In the latter stages, when the bubble density has increased, nucleation ceases and the evolution is controlled by growth and coarsening. As we show below, the gas concentration within the bubbles also evolves with time according to the production rate of gas atoms as well as the average bubble growth rate. In order to quantify the bubble microstructure throughout the simulation, we have calculated the porosity, void number density,
and average gas concentration within the bubbles (cg ) at regular intervals in time. This data is shown in Fig. 6. We vary two key parameters in this set of simulations: (i) the cascade-induced displacement rate (taking on values of either K = 2.3 103 dpa/s or K = 1.0 103 dpa/s), and (ii) the gas production rate (taking on values of Pgas = 0.0, 2.0 105, and 1.0 104 cg/s). The top row of Fig. 6 depicts the overall porosity (Fig. 6a), the bubble density (Fig. 6b), and the average gas concentration within Table 1 Values of k from Eq. (10) from fitting to the data shown in Fig. 6 and calculated from Eq. (11) for varying dpa rates and gas production rates. Also shown are the gas bubble nucleation rates. K (dpa/s) 3
2.3 103
1.0 10
Pgas (cg/s)
k (fit) (s3)
k (calculated) (s3)
J (k2 s1)
0.0 2.0 105 1.0 104
7
1.55 10 2.97 107 2.40 107
7
1.28 10 1.83 107 3.38 107
1.22 105 1.75 105 3.23 105
0.0 2.0 105 1.0 104
4.97 107 4.60 107 8.46 107
3.21 107 3.48 107 3.65 107
3.07 105 3.32 105 3.49 105
966
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
Fig. 7. (a) The nucleation rate of gas bubbles increases with increasing gas production rate for both displacement rates. This effect is more dramatic for the lower displacement rate (K = 1.0 103 dpa/s). (b) The average gas concentration within the bubbles decreases with increasing bubble diameter.
the bubbles (Fig. 6c) for a gas production rate of 1.0 104 cg/s. The porosity versus time plot (Fig. 6a) reveals that the bubble kinetics can be subdivided into three distinct regimes: stage I (incubation), stage II (nucleation and growth), and stage III (coarsening). This evolution is similar to our previous simulations in Part I. During stage II, the porosity follows an S-shaped curve characteristic of the Johnson–Mehl–Avarami equation [17]
h i 3 p ¼ po 1 exp kt ;
ð10Þ
where po is the porosity at the end of stage II, and k is a constant characteristic of the nucleation kinetics in the system, which for a 2D domain, is defined as
k¼
p _2 3
R J:
ð11Þ
In this equation, R_ is the void radius growth rate and J is the nucleation rate. The dashed lines in Fig. 6a, d, and g within the stage II regimes represent fits of the data to Eq. (10). For K = 2.3 103 dpa/s, these fits yield values of k = 4.97 107, 4.60 107, and 8.46 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. For K = 1.0 103 dpa/s, these fits yield values of k = 1.55 107, 2.97 107, and 2.40 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. This data is displayed in Table 1. The bubble nucleation rates, J, are obtained from Fig. 6b, e, and h, in which the dashed lines are linear fits representing the quasisteady-state nucleation regime. They yield nucleation rates of J = 1.22 105 k2s1 (Pgas = 0.0), 1.75 105 k 2s1 (Pgas = 2.0 105), and 3.23 105 k2s1 (Pgas = 1.0 104) for K = 1.0 103 dpa/s; and J = 3.07 105 k2s1 (Pgas = 0.0), 3.32 105 k2s1 (Pgas = 2.0 105), and 3.49 105 k2s1 (Pgas = 1.0 104) for K = 2.3 103 dpa/s. Please see Table 1. Taking into consideration the bubble growth rate of 0.1 k/s (obtained from measurements on the smallest voids in the simulations), the analytical values of k are obtained from Eq. (11). For K = 2.3 103 dpa/s, they are k = 3.21 107, 3.48 107, and 3.65 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. For K = 1.0 103 dpa/s, they are k = 1.28 107, 1.83 107, and 3.38 107 s3 for Pgas = 0.0, 2.0 105, and 1.0 104 cg/s, respectively. All of these values are within a factor of two of the fitted values obtained from Fig. 6a, d, and g. Fig. 6c shows the average gas concentration inside the bubbles (cg ) as a function of time (this value is averaged over all the bubbles in the domain). Interestingly, during the stage II regime, in which the bubble growth rate is highest, the gas density decreases from the initial value at the very onset of nucleation. This decrease occurs due to a faster bubble absorption rate of vacancies compared to gas atoms. Thus, the gas density decreases. During stage
III, however, cg begins to increase gradually and eventually levels off. This trend occurs for both displacement rates, however the magnitude of cg is lower for the higher displacement rate (K = 2.3 103 dpa/s) due to the fact that nucleation occurs earlier when the overall gas concentration is lower. The middle row of Fig. 6 displays similar results, although for a lower gas production rate, 2.0 105 cg/s. The results are qualitatively similar to the higher gas production rate, with the exception that the nucleation rates have decreased for both displacement rates. In addition, the cg values are generally lower, which is to be expected for a lower gas production rate. On the bottom row, Fig. 6g and h shows the porosity and bubble density evolution for the case of zero gas production. We see again that the bubble nucleation rates have decreased relative to the cases in which gas is produced. In Fig. 6i, the bubble density and the average bubble diameter are plotted as a function of gas production rate. These values are taken at the end of stage II, when nucleation ceases. We see a general trend of increasing bubble density and decreasing bubble diameter with increasing gas production rate. This is in agreement with experimental observations [1]. The generally-accepted explanation for this is that gas atoms act as trapping sites for vacancies, and therefore they represent nucleation sites for bubbles. In our simulations, the presence of gas atoms increases the free energy of the solid phase further away from equilibrium. Therefore, a smaller thermal fluctuation, or cascade event, is required to create a nucleation event. In addition to an increase in bubble density with increasing gas production rate, the bubble nucleation rates (or the slopes of Fig. 6b, e, and h) also increase with increasing gas production rate. This is shown in Fig. 7a for both displacement rates. Interestingly, the dependency of J (as well as the bubble density and bubble size) on the gas production rate is higher for the lower displacement rate. This is evident by the larger absolute values of the slopes in Figs. 6i and 7a for K = 1.0 103 dpa/s, relative to K = 2.3 103 dpa/s. Conceivably, this is perhaps due to the fact that bubble nucleation initiates early for the higher displacement rate, when the gas concentration is smaller. In Fig. 7b, the average gas concentration within the bubbles, cg , is plotted as a function of the bubble diameter. This data is for a gas production rate of 2.0 105 cg/s and a displacement rate of K = 2.3 103 dpa/s. For these conditions, cg decreases monotonically with increasing bubble size. This result is in accord with the data shown in Fig. 6f, in which cg decreases with time. As mentioned earlier, this decrease in cg with increasing bubble size is explained by a higher flux of vacancies than gas atoms to the growing bubbles. Had the gas production rate been considerably higher than the rate of increase in vacancy supersaturation (taking into consideration vacancy-interstitial recombination), this trend may
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
967
actually differ; although in this work we did not simulate such conditions.
5. Intergranular gas bubbles Extended microstructural defects such as GBs and dislocations strongly affect bubble structure evolution in irradiated materials due to their role as point defect sinks, as well as the inherent fact that impurity atom formation energies are lower at these features compared to the bulk. It is commonly believed that GBs act as heterogeneous nucleation sites for fission-gas bubbles; this being widely supported by many experimental observations of intergranular bubbles that are many times larger than their intragranular counterparts [2,3,18]. As intergranular bubbles grow to micronsizes and larger, their migration and coalescence on the GB faces, edges, and corners leads ultimately to an interconnected porosity that allows rapid fission gas release to the exterior of the material [2,3]. In addition, due to the triple junctions between bubble free surfaces and the GB, the inherent bubble shape is lenticular rather than spherical in order to minimize the interfacial energies of GB and free surfaces. Furthermore, the presence of intergranular bubbles effectively pins GB migration, according to the Zener drag theory [10]. As we show below, the current model is able to capture these important distinctions of intergranular bubble physics. Fig. 8 shows the interaction of a single GB with a stationary bubble that does not evolve (i.e., change shape) in time. (Here, the concentration fields and the g-field are not updated in time; only the /-fields representing grains 1 and 2 are solved for in time). In Fig. 8a, the equilibrium GB shapes are shown for different relative positions of the bubble. We see that, for this stationary bubble, the GB always intersects the bubble surface at a 90° angle. (Note that in these figures, the GB is not made to migrate across the domain, rather the initial position of the GB is changed from left to right). In Fig. 8b, the pinning energy that the bubble exerts on the GB is plotted as a function of bubble radius. This pinning energy is the difference between the system shown in the center of Fig. 8a and that shown on the far left of Fig. 8a. According to Zener [10,16], it is equal to the amount of GB that is eliminated by the presence of the bubble multiplied by the GB energy. For 2D, this is equal to 2rrgb. The data points in Fig. 8b correspond to the re-
Fig. 8. (a) GB pinning due to a stationary bubble (that does not evolve in time). (b) The pinning energy is equal to the GB energy, gb, multiplied by the length (for 2D) of the missing GB segment corresponding to the bubble diameter. The open circles represent values calculated from the model and the dashed line is the theory according to Zener [10].
Fig. 9. Evolution of the shape of a gas bubble on a GB from circular to lenticular (equilibrium) due to the GB energy. The (a) /, (b) cv, and (c) cg fields are shown throughout the evolution. (d) The dihedral angle, h, evolves exponentially from 90° to approximately 120°.
sults obtained from simulation, and the dashed line corresponds to the theory, showing a close agreement. This result is entirely due to the interaction in the free energy term of Eq. (4) between the order parameter representing the bubbles, g, and the order parameters representing the grains, /. With verification that the current model accurately captures the pinning effect of a gas bubble on a GB, we now focus on the influence of a GB on the shape of a gas bubble. As mentioned above, an intergranular bubble will assume a lenticular shape due to the Herring relation [19] governing the dihedral angles existing at the triple junction between the free surfaces and the GB. To verify that we correctly capture this, we have performed a simulation of a bicrystal with a single, initially-circular gas bubble residing on the GB. As shown in Fig. 9a, the bubble evolves from a circular shape to a lenticular shape as time progresses. Fig. 9b and c illustrates that the vacancy and gas atom concentrations, respectively, evolve accordingly. The evolution of the triple junction angle, h, is shown versus time in Fig. 9b, illustrating an exponential approach from h = 90° to h 120°. With the above analysis of the interaction between pre-existing gas bubbles and GBs, we now focus on bubble nucleation and growth in polycrystalline microstructures. Here, we re-activate the irradiation terms. Fig. 10 illustrates the nucleation and growth of intergranular gas bubbles due to the ongoing production of Frenkel pair defects and gas atoms. The simulation domain is again periodic in both directions (lx = ly = 6 lm) and contains 40 grains (i.e., P = 40 in Eq. (4), with isotropic GB properties, including energy and mobility) with an average diameter of d = 1.1 lm. The gas atom formation energy is reduced in the GBs by the following equation
f ¼ Ef U; E g g
ð12Þ
968
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
where the factor U is defined in Part I and is equal to 1.0 in the grain interiors and as low as 0.5 in the GBs and triple junctions. This results in preferential segregation of the gas atoms to the GB regions. The formation energies are Efv = 0.8 eV, Efi = 1.5 eV, an-
dEfg = 1.2 eV, and the recombination rates are Rbulk = 1.0 s1 and Rsurf = 5.0 s1. Fig. 10 shows sequential snapshots of the preferential nucleation of bubbles on GB and TJ regions, resulting from the enhanced
Fig. 10. (a)–(d) Snapshots throughout time of the nucleation and growth of intergranular gas bubbles in a polycrystalline grain structure, while (e) is a close-up view. Bubbles existing on GBs are lenticular shaped, whereas bubbles on triple junctions have a curved triangular shape. (f) A microscopy image of intergranular bubbles in UO2 taken from [2].
Fig. 11. Characteristics of intergranular gas bubbles for different grain sizes. The (a) porosity, (b) bubble density, and (c) gas concentration within the bubbles does not vary substantially for the different grain sizes. However, as the grain size increases, the (d) average bubble diameter increases, the (e) average GB bubble spacing decreases, and the (f) GB bubble coverage increases. (The data for (d), (e), and (f) are taken at the end of the simulations, t = 800 s).
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
gas atom concentration at these extended defects. As can be seen in the close-up view in Fig. 10e, bubbles existing on GBs are lenticular shaped, whereas bubbles on TJs have a curved triangular shape. In Fig. 10f, a micrograph image of intergranular bubbles in UO2 (taken from Ref. [2]) illustrates the qualitative similarities between the bubble structures in terms of the relative size of the bubbles with respect to the grain size, as well as the bubble shapes. Also note that, in the simulations, bubbles within the grain interiors have also begun to nucleate, although their shapes remain circular, as would be expected. In Fig. 11, we analyze the intergranular bubble characteristics for various grain sizes. The top row, Fig. 11a–c shows the porosity, bubble density, and cg throughout time for an average grain size between 37 k < d < 111 k. Each of these simulations are run with identical irradiation conditions of K = 1.0 103 dpa/s and Pgas = 8.2 106 cg/s. We observe that the porosity, bubble density, and cg do not vary substantially with the changing grain size (the exception is the smallest grain size, d = 37 k, in which the increase in porosity has a slightly slower rate, although bubble nucleation rate is similar to the other grain sizes). The average gas concentration inside the bubbles (Fig. 11c) grows most rapid at the onset of nucleation, and saturates at the later stages, in a very similar manner for each grain size. Furthermore, the bubble nucleation rates, obtained from Fig. 11b, do not vary substantially with grain size. However, when we plot the average GB bubble diameter (Fig. 11d), the average spacing on the GB lines (Fig. 11e), and the GB bubble coverage (percentage of GB length occupied by bubbles) (Fig. 11f), we see clear trends in each with respect to the grain size. With increasing grain size, the intergranular bubble diameter increases, the intergranular bubble spacing decreases, and the GB bubble coverage increases. Considering that the porosity is basically equal for each grain size at the end of the simulations, these trends can be rationalized by the fact that the total GB density is lower for the larger grain size. Therefore, for the same porosity and overall bubble density, and considering that virtually all bubbles reside on the GBs, bubbles should be closer spaced and larger for larger grain sizes. This, in turn, results in a higher GB coverage of bubbles. Interestingly, it appears that in the polycrystal specimens, the nucleation stage has yet to saturate, yielding only the first half of an S-shaped curve (Fig. 11a). We have observed that at about the time that nucleation ceases on the GBs, nucleation events begin to occur in the grain interiors. So, the system as a whole does not decelerate the porosity development. Eventually, nucleation in both the grain interiors and the GBs will stop, at which time the system will truly enter a coarsening stage.
6. Conclusions
turn, dictated by this change in bubble size, as well as by the flux of gas atoms to the bubble. In single crystal simulations, the porosity evolved in three distinct stages during irradiation: stage I (incubation), stage II (nucleation and growth), and stage III (coarsening). The bubble nucleation rate is shown to increase with increasing displacement rate as well as gas production rate. Increasing gas production rate is also shown to increase the bubble density and decrease the bubble size. Furthermore, at the early stages of bubble growth, the gas density within the bubbles (cg ) decreases due to a higher flux of vacancies, compared to gas atoms, to the bubble surfaces. In polycrystal simulations, the gas atom formation energy is decreased within the GB and triple junction regions, leading to preferential segregation and heterogeneous nucleation at these features. During intergranular bubble growth, the bubbles retain their equilibrium lenticular shapes on GB planes, and curved triangular shapes on triple junctions. As the grain size is increased, the porosity and nucleation rates remain constant; however, the average bubble size and GB bubble coverage increase, while the average GB bubble spacing decreases. Acknowledgements This work was supported through the INL Laboratory Directed Research and Development program under DOE Idaho Operations Office contract no. DE-AC07-051D14517V. We also acknowledge support from two DOE programs: the Energy Frontier Research Center (EFRC) for Materials Science of Nuclear Fuel and the Nuclear Energy Modeling and Simulation (NEAMS) program (Cetin Unal as Fuels IPSC program lead). In addition, this manuscript benefited from technical discussions with Michael Tonks (INL) and Srujan Rokkam (FSU). Appendix A. Free energy of ideal gas The free energy of the gas phase with respect to the gas concentration is derived below. We start with the total Gibb’s free energy [20] assuming an ideal gas law
Gðp; TÞ ¼ ngo þ nRT ln p;
ðA1Þ
where p is the pressure, n is the number moles, go is a reference molar chemical potential, and R is the universal gas constant. This equation can be re-written to obtain a volumetric energy density
gðp; TÞ ¼
n V
Av log þ bub
kB T ln p; bub
nAv
V
ðA2Þ
where Av is Avogadro’s number, Vbub is the volume of the bubble, and log is a reference atomic chemical potential of the gas. In terms of gas concentration, the equation is now written as
gðp; TÞ ¼ cg log þ cg kB T ln p: The phase-field model presented here extends our previous model in Part I for pure metals by incorporating fission gas atoms in a manner that distinguishes their point defect energetics and kinetics in the solid and gas bubble phases. This is accomplished by the use of a long-range order parameter, g(r,t), that spatially differentiates the two phases throughout the simulation domain. This is a key feature that distinguishes the current model with a previous phase-field model [12] for fission gas behavior in irradiated materials, in which a non-conserved order parameter is omitted. In addition, the g-field in the current model allows for a proper interaction between gas bubbles and GBs in terms of the GB pinning energies and the equilibrium bubble shapes. The rate of bubble growth or shrinkage, in the bulk or on GBs, is controlled solely by the flux of vacancies and interstitials to the bubble surfaces. The gas concentration within the bubbles is, in
969
ðA3Þ
Using the ideal gas law, the gas pressure is equal to
p¼
nAv
V
bub
kB T ¼ cg kB T;
ðA4Þ
and by substituting (A4) into (A3), we arrive at
f b ðcg Þ ¼ log cg þ cg kB T ln cg þ cg kB T lnðkB TÞ:
ðA5Þ
References [1] H. Trinkaus, B.N. Singh, J. Nucl. Mater. 323 (2003) 229. [2] I. Zacharie, S. Lansiart, P. Combette, M. Trotabas, M. Coster, M. Groos, J. Nucl. Mater. 255 (1998) 85. [3] S.B. Fisher et al., J. Nucl. Mater. 306 (2002) 153. [4] H. Ullmaier, Nucl. Fusion 24 (1984) 1039.
970
P.C. Millett et al. / Computational Materials Science 50 (2011) 960–970
[5] D.R. Olander, Fundamental Aspects of Nuclear Reactor Fuel Elements (US-DOE, 1985). [6] K. Morishita, R. Sugano, J. Nucl. Mater. 353 (2006) 52. [7] J.A. Turnbull, J. Nucl. Mater. 38 (1971) 203. [8] M.H. Wood, J. Nucl. Mater. 82 (1979) 264. [9] J. Rest, J. Nucl. Mater. 321 (2003) 305. [10] C. Smith, Trans. Metall. Soc. AIME 175 (1948) 15. [11] L.Q. Chen, Y. Yang, Phys. Rev. B 50 (1994) 15752. [12] S. Hu, C.H. Henager, H.L. Heinisch, M. Stan, M.I. Baskes, S.M. Valone, J. Nucl. Mater. 392 (2009) 292. [13] P.C. Millett, S. Rokkam, A. El-Azab, D. Wolf, Mod. Sim. Mater. Sci. Eng. 17 (2009) 064003.
[14] J.W. Cahn, J.E. Hilliard, J. Chem. Phys. 28 (1958) 258. [15] L.D. Landau, E.M. Lifshitz, Statistical Physics, third ed., Elsevier, Burlington, MA, 1980. [16] N. Moelans, B. Blanpain, P. Wollants, Acta Mater. 53 (2005) 1771. [17] R.W. Balluffi, S.M. Allen, W.C. Carter, Kinetics of Materials, Wiley, Hoboken, 2005. [18] S. Fréchard, M. Walls, M. Kociak, J.P. Chevalier, J. Henry, D. Gorse, J. Nucl. Mater. 393 (2009) 102. [19] C. Herring, in: The Physics of Powder Metallurgy, McGraw-Hill, New York, 1951, pp. 143–179. [20] D.R. Gaskell, DR Introduction to Metallurgical Thermodynamics, second ed., Hemisphere, Publishing Corporation, 1981.