Statistical analysis of the interaction between irradiation-induced defects and triple junctions

By using a generalized, spatially resolved rate theory, we systematically studied the irradiation-induced diffusion and segregation of point defects near triple junctions. Our model captured not only the formation, growth, and recombination of point defects but also the interaction of these defects with pre-existing defects. We coupled the stress field of the triple junction with defect diffusion via a modified chemical potential. The residual stress fields of grain boundaries and triple junctions are modeled via disclination mechanics theory. By assessing the behavior of 144 triple junctions with vacancy and interstitial defects, we correlated defect-sink efficiencies with key characteristics of triple junctions. For vacancies, the geometric configuration of triple junctions dominated sink efficiency, suggesting that equiaxed grains would resist the accumulation of vacancies more than elongated grains. For interstitials, the sink density of the grain boundaries composing the triple junctions dominated sink efficiency. Hence, the interstitial concentration may be managed by adjusting the structure of the grain boundaries. Overall, we illustrated the complex coupling between pre-existing defects and radiation-induced defects through interaction of their stress fields. This theoretical framework provides an efficient tool to rapidly assess defect management in microstructures.


Introduction
The ability of a material to manage and eliminate point defects, i.e., interstitials and vacancies, plays a significant role in determining how its mechanical properties will be affected by irradiation over extended periods of time. Immediately after the initial damage event, recombination eliminates most point defects. However, some residual defects remain, leading to the growth of point-defect clusters. These point defects and their associated clusters eventually diffuse toward different types of pre-existing microstructural sinks, which include dislocation loops, grain boundaries (GBs), and triple junctions (TJs).
An extensive amount of work has been conducted to understand and characterize the relationship between the density of these sinks in a given material (i.e., effect of grain morphology and grain size) and the associated ability of that material to resist radiation damage (see review by Beyerlein et al. [1]). For example, Harkness and Li [2] studied void formation in 304 stainless steel and suggested that increasing the presence of various microstructural sink types would increase resistance to damage. Following up on this hypothesis, Singh [3] further studied the relationship between grain size and void swelling, finding that smaller grains did indeed reduce radiation damage, likely due to the increased density of GB sinks. More recent work has focused on smaller grain sizes, specifically nanograined materials with grains much smaller than 1 µm. Yu et al. [4] have observed that reducing grain size lowers the density of He bubbles in Fe irradiated with He ions. Likewise, Han et al. [5] observed that the density and size of cavities were reduced in nanograined Cu as compared to coarse-grained Cu when subjected to He irradiation. Similarly, Alsabbagh et al. [6] observed a reduction in radiation-induced hardening and embrittlement in ultra-fine-grained steel as compared to coarse-grained steel, while Sun et al. [7] quantified the reduction in swelling and swelling rate in 304 L stainless steel when using an ultra-fine-grained structure.
While these examples show that nanograined materials exhibit improved resistance to radiation damage largely due to the fact that reduced grain size increases the volume fraction of GBs, the high density of GB networks (i.e., the high density of triple junctions) within these materials also affects their radiation tolerance. From the perspective of the microstructural defect management, GB networks, in terms of the character of triple junctions, represent unique features that serve as paths for diffusion/transport and regulate damage accumulation [8][9][10]. For example, Adlakha and Solanki [11] used atomistic simulations to illustrate the distinct role of TJs structures on defect binding and migration behavior in α-Fe for a couple of TJ configurations, demonstrating that interstitials preferably segregate at/near TJs over vacancy defects.
Theoretical and computational models characterizing point defects or solutes diffusion and interaction with microstructures are continuously developed [12][13][14][15][16][17][18][19]. However, due to the large number of parameters required to describe defect-defect interactions within a system, even in the case of a single element material such as Fe, most of the existing models typically consider only one type of point defect at a time, limiting studies to a small handful of cluster types. Furthermore, the degree to which most of the mean field rate theory formulations describe the coupling between sinks and defects is mostly uni-directional, i.e., only a limited number of models consider the chemo-mechanical coupling affecting point defect diffusion and segregation, in turn limiting the applicability of formulations describing sink strength of pre-existing microstructural sinks. In the case of GBs and TJs, the sink strength is not only impacted by its intrinsic residual stress field but also by the discrete structure of the sink considered. For example, we previously demonstrated that the density of sinks composing an infinite GB embedded into a bulk material can be directly correlated to the ability of that boundary to absorb defects [16]. However, sink density is only one of many parameters that influences the defect absorption. Identifying the full set of parameters and how they impact irradiation resistance is challenging. Some studies have taken a statistical approach to isolate the relevant variables affecting defect absorption [20][21][22].
In this article, we conduct such a statistical approach by using a spatially resolved rate theory to systematically study the interaction between irradiation-induced point defects and triple junctions. As detailed in the "Methods" section, our rate theory formulation  [47] builds upon a model coupling the intrinsic stress of pre-existing defects with point defect diffusion [16]. It captures not only the formation, growth, and recombination of point defects, but also the interaction of those defects with pre-existing microstructural sinks, allowing for a wide range of defect-defect interactions and defect types. New additions to this model comprise more complex clustering and dissociation reactions for defects, along with nonlinear elastic and long range effects of GBs and TJs. We performed a statistical analysis on 144 different triple junctions. The spatial configuration of the TJ constructs is illustrated in Fig. 1. A given TJ is described as a local geometrical and structural triplet of idealized GBs. Dihedral angles and sink densities of the triplet of GB composing the TJ are used as geometrical and structural (structure and energy of GBs composing the TJ) descriptors respectively, and vacancy and interstitial sink efficiencies are used as TJ sink descriptors. We used principal component analysis (PCA) on these descriptors to correlate the characteristics of TJ networks with their associated defect-sink efficiencies.
In what follows, we are only showing results for the first two principal components (PCs) of this analysis. The first two principle components were found to account for 40.2% and 25.5%, respectively, of the total variance in the behavior of the 144 TJs studies.

Vacancy-sink efficiency of triple junctions
Let's first examine the behavior of vacancies in the vicinity of TJs. Figure 2 presents the correlations and variance between the configurational and structural variables and defectsink efficiencies around any TJ. Explanations on how to interpret correlation circles are provided in the "Methods" section. By visually inspecting Fig. 2, we can infer that the vacancy-sink efficiency shows a positive correlation with the smallest dihedral angle β 3 and a negative correlation with the largest dihedral angle β 2 . This suggests that having a −0. 6 Fig. 2 a Comparisons of the first two principal components (PCs) for the geometric/structural and sink descriptors for the triple junction (TJ). The variable correlations within the first two principal components space are tabulated in Table 1, while the correlation coefficients between the TJ geometrical and structural descriptors are listed in Table 2. b PC scores for all 144 TJ configurations with respect to the first two principal components. Open symbols corresponds to three illustrative TJ regions with low, medium and high sink efficiencies. The dashed elliptical region corresponds to seven TJ with the highest PC score. Their sink efficiencies are listed in Table 3   Table 1 Variable loadings for the first and second principal components  smaller range in the dihedral angles increases a TJ vacancy-sink efficiency. In other words, the geometric configuration of the TJ play a dominant role for determining vacancy-sink efficiency, suggesting that an equiaxed grain structure is preferable to an elongated grain structure for mitigating vacancies (Tables 1, 2). To further illustrate this point, we identified a group (circled in Fig. 2b) of TJ configurations that have a high vacancy-sink efficiency. Their vacancy-sink efficiencies as a function of the dihedral angles are listed in Table 3. Upon inspection, the trend is consistent with the PCA results: the TJs with the highest vacancy-sink efficiency tend to have a smaller range of dihedral angles. This is also evident by looking at the vacancy concentration profiles for three different TJ configurations as shown in Fig. 3. The two TJ configurations with obtuse GB configurations (Figs. 3a-b and data marked with either a diamond symbol or square symbol in all sink efficiency plots) have relatively few vacancies in the vicinity of the TJ. Conversely, the equiangular configuration (see Fig. 3c and data marked with a circle symbol in all sink efficiency plots), shows a large depleted zone around the TJ.  The nonlinear relationship between the geometrical dihedral descriptors and the vacancy-sink efficiency is plotted in Fig. 4a. There is a sharp transition in vacancy-sink efficiency when the range of dihedral angles is greater that 70 • . The TJs with the highest vacancy-sink efficiency all lie in the "equiaxed" region (i.e., β 2 − β 3 60 • ). Conversely, the TJ configurations with high dihedral angle range (i.e., β 2 − β 3 > 140 • ) have relatively low vacancy-sink efficiencies. However, the vacancy-sink efficiency profile of the 144 TJs tested is relatively "flat" when plotted as a function of the structural descriptor (see Fig. 4b), and thus no special correlation can be drawn on the vacancy-sink efficiency with respect to sink density. Overall, the role of equiaxiality in the TJ configuration is to be expected, since if GBs form a more or less symmetrical configuration, all vacancies (vacancies being the slowest species as compared to interstitials) have a chance to reach all sinks in their vicinity. If the configuration is highly asymmetric, vacancies have a lower chance to find sinks at GBs in the grain with an acute angle at the triple junction as compared to interstitials. The sinks on the GB of this grain near triple junctions do not participate as much in the point-defect removal and the overall efficiency drops.

Interstitial-sink efficiency of triple junctions
While the sink efficiency for vacancies is dominated by the geometric configuration of the TJs, the sink efficiency of interstitials is for the most part uncorrelated with the dihedral angles (nearly 90 • angle between ζ I vectors and β vectors in Fig. 2a and as seen from scatter plots in Fig. 5a). Rather, as seen in Fig. 2a, interstitial sink efficiency positively correlates with the structural descriptors of the TJ network (i.e., the density of sinks for the TJ). The sink densities of the GBs comprising the triple junction account for most of the variance in interstitial sink efficiency. In particular, the GB with the highest density contributes most (in this case GB 31). Unlike the case with vacancies, when examining the PC scores for the 144 TJs tested in Fig. 2b, there is no distinct group of TJ configurations that possess a high interstitial sink efficiency.
The strong correlation between the interstitial sink efficiency and the structural descriptor of the TJ network shown in Figs. 5b, c reveals several insights on the origin of the defect-TJ interaction. First, TJs with the highest SIA sink efficiency are also those with the highest elastic energy (in the present case, GB "31", due to the geometric construct). This high efficiency could be attributed to the corresponding free volume associated with the GB and the high density of sinks along these boundaries (ρ (GB) S ∼ 1 nm −1 ). This remarks partially concurs with experimental observations: in general, high-angle grain boundaries preferably annihilate radiation-induced point defects and are usually biased to absorbed interstitials owing in part to an increased free volume [23,24]. The group of TJs that are most efficient at interstitial removal all have either a −15 • or −75 • misorientation along one of the GB (GB "31" in this specific configuration); however, these are not Interstitial sink efficiency versus a the difference between the largest and smallest dihedral angles of the TJ, b GB sink density (for GB "31"), and c GB energy. Open symbols correspond to TJ regions with low, medium and high sink efficiencies. The dashed elliptical region corresponds to TJs with the highest sink efficiency the same seven configurations that were highlighted as efficient sinks for vacancies. This observation has an important implication: in agreement with our previous work [16], for interstitials, the sink density of the TJ, and therefore the atomic structures of GBs composing a triple junction are key for governing the SIA sink efficiency. In this case, the sink region along the GBs is largely responsible for the elimination of defects at steady state, although defects are also removed via annihilation.

Vacancy vs. interstitial
The bias factor, used to compare interstitial and vacancy sink efficiencies (see Eq. (17)), show a similar (albeit weaker) correlation (see Fig. 6). with the structural descriptor of the TJ network rather than its geometrical descriptors. This difference between SIA and vacancy sink efficiency can be traced back to the difference in compositional expansion coefficient for the two defect types. The coefficient of compositional expansion (CCE) is the parameter that ultimately determines the strength of the stress-diffusion coupling for a particular defect type in the present framework, as seen in Eqs. (2) and (4). With a larger CCE, the stress coupling is stronger and thus would be more influential for SIAs. The transition from small to large CCE likely comes with a transition from negligible to significant stress coupling although this relationship cannot be characterized based on the present results. At low values of the CCE and very low defect concentrations, the stress coupling can likely be ignored; at high values and concentrations, doing so would lead to large inaccuracies.

Comparison to experimental observations and relevance of the present findings
Numerous experimental studies have focused on understanding the role of GB character on sink efficiency (see for example [25]), yet there are far fewer studies that explicitly isolate the effect of TJs. In a small number of detailed studies, the effect of TJs can be observed in the data presented, but it is generally not the focus of the investigation. The effects of TJs can be indirectly seen for example in a study on 10 keV He irradiation of nanocrystalline Fe of various grain sizes, [26] where several examples are provided on the defect-denuded zones in the proximity of GBs/TJs as an indicator of sink efficiency. Despite the scatter in the data, the experiments confirmed the expected trend that higher misorientation GBs (and hence higher sink density), generally have more pronounced defect-denuded zones. In that study, as with others, there are not a sufficient number of observations to evaluate the relative sink efficiency of TJs as a function of dihedral angle or grain shape, and such a detailed experimental assessment could be productive for future study. In addition, the large scatter in the data presented in the majority of experimental studies correlating grain boundary character to radiation defect density suggest that other unidentified factors in the microstructure or radiation environment also contribute to the interaction between TJs and defects. These effects are presently not captured in our model.
While the present study focuses on modeling the behavior of bcc Fe, it is relevant to consider how these results might extend to other metals and other crystal structures. The behavior described herein is partially controlled by the ratios between the binding and migration energies (see Eq. (12)), the GB character, and the geometrical configuration fo the TJ, so other material systems with comparable characteristics would be expected to behave similarly. In Samaras et al. [27], vacancy cluster formation was compared in nanocrystalline grain boundary networks of fcc Ni and bcc Fe. While the behavior of TJs was not isolated, the study showed that both crystal structures result in GBs and TJs with similar periodic tensile/compressive hydrostatic stress fields, and in both cases glissile SIAs were attracted to regions of high hydrostatic tension in the GB/TJ network, facilitating annihiliation in the free volume of the GB/TJ. The rate of SIA migration and subsequent annihilation in fcc is generally slower than vacancies due to the close packed nature, although given sufficient time both systems eventually consume SIAs in the GBs/TJs. However, the resulting vacancy-dominated grain interiors for the two crystals (with distinct stacking fault energies) then evolved into distinct defects, with fcc forming stacking fault tetrahedra and bcc forming large vacancy clusters. Similarly, our conclusions and observations are in agreement with those from Adlakha and Solanki [11] on the role of the TJ strain energy for the aggregation of point defects.
The present study is built around a few constraining assumptions: the annihilation of defects does not change the GB structure, and in real materials, the GB forms a complex interacting network. Indeed, neighboring triple junctions and GBs may change the stress state of a given triple junction, especially in the nano-grained materials, in turn affecting the magnitude of the elastic fields. This may influence the point defect diffusion but it is not considered here. Our study is also constrained to idealized symmetric tilt grain boundaries (STGB), which are sessile in nature. Yet practical metallic systems contain a broader range of grain boundary types that are dynamic and responsive to the local stresses. For idealized STGBs, sinks are crystallographically periodic and depend on the disorientation between the two crystals. The triplet of GBs and the resulting TJs consist of the minimum set of disclination defects needed to achieve the requisite misorientation. However, GBs are rarely in this lowest energy configuration, and a number of other metastable configurations can occur with only minor elevation in interfacial energy but nontrivial distortion in the GB structure [28].
Likewise, we expect that TJs can exist in a multitude of nearly degenerate states beyond the idealized case presented here. Indeed, the evolutionary formation of grain boundary networks can be described by crystallization or defect-mediated grain boundary migration [29] with more complex crystallographic constraints at the TJ leading to well-established TJ drag [30]. And it is reasonable that both crystallization and migration can trap complex metastable defects and excess free volume in the TJ, such as illustrated in the work by Poletaev et al. [31,32]. This dynamics becomes even more complex during radiation events that are can alter the shape, location, and character of grain boundaries and triple junctions [33][34][35]. As such, the present reported sink efficiencies may represent lower-bounds on the efficiency of defected TJs. Nevertheless, the effects of realistic, excessively defected TJs on local stress fields and sink efficiency would be a valuable subject for future study, and the present statistical methodology can be useful for such a many-dimensional investigation.

Conclusions
We presented a theoretical framework pointing to the necessity of including the geometrical and structural "fingerprints" of pre-existing defects in order to appropriately consider defect management in microstructures. For vacancies, the geometric configuration of triple junction dominated sink efficiency. For interstitials, the sink density and structure attributes of the grain boundaries composing the triple junctions dominated sink efficiency. While the proposed framework simplifies microstructural configurations and the complexity of pre-existing defects (our model is restricted to 2D configurations), it is an efficient tool to rapidly assess defect managements in microstructures that correlates with experimental observations. In real materials, the spatial distribution of triple junction types is not random, leading to a highly correlated processes for defect-TJ interactions. Additional descriptors (excess free energy of the triple line, specific defect absorption efficiency, 3D neighboring grain boundaries configurations) of triple junctions could be used in order to better assess and more accurately describe TJ characters on managing irradiation defects.

Chemomechanical model of stress-induced diffusion
The chemomechanical diffusion model used in this work is based on the chemical potential proposed by Swaminathan et al. [36] and expanded on by Cui et al. [37]; The constant μ (0) α gives the chemical potential at standard state; R is the standard gas constant; T is absolute temperature; γ α is the activity coefficient; c α gives the mole fraction of point defects (either vacancy or SIA) per mole of atoms in the solid; and τ α is the mechanical contribution to the chemical potential for species "α". The subscript α indicates the defect type, e.g., a size n = 3 vacancy cluster. Species explicitly modeled are single vacancy and vacancy clusters of size up to 4, and single SIA and SIA clusters of size up to 3, resulting in the consideration of 7 different types of defect configurations. In the results presented here, c α is small enough that the activity coefficient γ α is taken to be approximately 1.
In our previous work [16], τ α was based on a reduced model that only considered the contribution from the hydrostatic stress in the chemical potential. The present work builds upon this method by considering not only the hydrostatic stress but also the effects of nonlinear elasticity. This is achieved by using the more general form of τ α derived by Cui et al. [37] to compute the stress coupling term; where V m is the molar volume of the pure material, C is the anisotropic elasticity tensor, and F c is the compositional deformation gradient attributed to the point defects. F in,e is the deformation gradient due to the intrinsic fields and the elastic response of the grain boundaries and triple junction, and E in,e = 1 2 F T in,e · F in,e − I is the corresponding Green strain tensor.
We used the finite element method to determine the combined displacementsû due to compositional effects and the elastic response [16], from which the corresponding deformation gradient tensorF is calculated. The total deformation is multiplicatively decomposed following the method used by Cui et al. [37] in deriving the stress-dependent chemical potential; The intrinsic deformation gradient is calculated from the known intrinsic fields, and the compositional deformation gradient is calculated as; where η V and η I are the CCE for vacancies and SIAs (respectively), and the values c V and c I indicate the total fractions of both point defect types. The tensor F in,e is then found by removing the compositional deformation gradient; The stress is calculated based on the intrinsic elastic fields describing the GB and TJ, and the elastic response to the eigenstrains introduced by point defects. By including eigenstrains in the calculations, the stress field is affected by defect concentrations just as the diffusion is influenced by the stress field. This mechanism is key to describing the two-way coupling between stress fields and point-defect concentrations. Once the chemical potential is known, the rate equation for defect concentration is derived via Fick's first and second laws: where D α is diffusivity of species α. The additional a αβ f β terms in Eq. (7) account for defect-defect interactions, specifically clustering and recombination. Each term f β describes an interaction rate between two types of defects, and the coefficient a αβ depends on the type of defects α and β that are interacting. For example, when two point defects of size n = 1 recombine, a αβ = −1 for both defect types since the reaction removes one defect from each population. These additional terms also include point-defect implantation rate, i.e., the rate at which point defects are added to the simulation. For all simulations included in this study, Frenkel pairs are implanted at a uniform, constant rate of f implentation = 0.001 dpa s −1 . Although the displacement rate reported in this study is in units of dpa s −1 , experimental measures of radiation dose rate are typically in units of particle flux, or particles cm −2 s −1 . In experimental studies using ion irradiation, the relationship between displacement rate and ion flux depends on many factors including the mass and energy of the incident ion as well as the thickness of the sample. In the present study, the implantation rate of 0.001 dpa s −1 would be equivalent to implanting 50 keV helium ions in a 100 nm Fe thin film with an ion flux of ∼ 10 14 He + cm −2 s −1 [38].

Defect-defect interactions
In this work, only small, mobile defect clusters are considered; large, immobile point defect clusters are excluded. Specifically, sizes n = 1, 2, 3, 4 are included for vacancy clusters, and sizes n = 1, 2, 3 are included for SIA clusters [39]. SIA clusters are treated as dislocation loops that diffuse in one-dimension (1D), and vacancy clusters are treated as spheres that migrate in three-dimensions (3D). The rates for defect-defect interaction f β depend on the type of defects that are interacting. The term f β describes annihilation or dissociation interactions of SIAs and vacancies via recombination depending on the nature of the defect cluster (3D-3D interaction vs. 3D-1D). The rate for 3D-3D or 1D-1D interactions (i.e., vacancy-vacancy or interstitial-interstitial) is given by [40]; The values n and m indicate the sizes of the interacting defects. Used as subscripts, they indicate a property related to that defect, e.g. D n refers to the diffusivity of the first defect, which as size n. The value Z int accounts for the tendency of SIA clusters to absorb other SIAs. For SIA-SIA interactions Z int = 1.15; otherwise, Z int = 1 [40]. The coefficient Z ann accounts for the increased recombination of point defects near a GB and is described below. For clustering rates, Z ann = 1. For 3D-1D interactions (i.e., interstitial-vacancy interactions), the rate is given by [40]; where b is the length of the Burgers vector. The geometric constants ω and ω 2D depend on the atomic volume of the solid; The annihilation coefficient Z ann accounts for increased mobility and recombination of point defects near GBs, and is computed as follows for a GB centered along x = 0 [16]; The parameter Z (0) ann defines the annihilation coefficient along the GB, and W determines the width of the region around each GB where annihilation is more frequent. Calculation of W is based on the elastic energy of the GB and follows the method used in our previous work [16].
The rate of dissociation of individual defects out of larger clusters is based on the binding energy of the cluster and the migration behavior of the resulting defects [40]; This chemomechanical method has been used to study point defect diffusion around dislocation loops and GBs [16], although it may be extended to any microstructural defect for which the intrinsic stress field may be computed. For the present work, the intrinsic stresses around triple junctions are considered, as described below.

Representation of microstructural defects
Microstructural defects are represented using disclination dipoles. Such mathematical constructs have been used to represent equivalent edge dislocations and GBs [15,41]. Here, we use sets of disclination dipoles to construct GBs, and GBs are subsequently used to create TJs. We use disclination mechanics as opposed to edge dislocations since disclination dipoles can be used as a mathematical construct for a large class of GBs with arbitrary misorientations between 0 • and 90 • .
The intrinsic stresses and displacements of wedge disclinations and disclination dipoles are based on the linear-elastic theory of defects [16,42,43], and the stress fields predicted by the disclination models have been successfully validated against atomistic modeling [41].

Disclination structural unit model
GBs are represented in this work using the disclination structural unit model (DSUM) to achieve misorientations θ 12 , θ 23 , and θ 31 [15,[44][45][46]. Disclination dipoles are located periodically along each GB, and the overall misorientation is determined by the strength of the disclinations, the arm length in each dipole, and the spacing among dipoles.
When using DSUM to create infinite-length GBs, infinite series formulas are used to determine the intrinsic stresses [15]. In a TJ, however, each GB must terminate at the triple point, which disallows the use of infinite-series formulas. These GBs are instead constructed by using known formulas for disclination stresses and displacements along with a finite number of disclination dipoles. The contributions from each disclination are then superimposed to find the intrinsic, mechanical fields. These dipoles are placed such that the first disclination occurs at the triple point, and subsequent disclinations are spaced along each GB according to the DSUM construct.
For simpler DSUM structures, the GB may consist of a single wall of disclination dipoles spaced by period length H. For more complicated GBs, there may be two or more disclination dipoles in each period. This results in two or more "walls" of disclination dipoles,  Bias factor plotted against the sink density of GB "31". Open symbols correspond to TJ regions with low, medium and high sink efficiencies each constructed using the same Franks vector, dipole arm, and period, but each with a slightly different offset. As in our previous work, material and DSUM parameters were chosen to simulate BCC Fe [16].

Constructing triple junctions using disclinations
Triple junctions and their intrinsic stress fields are constructed in three main steps. The kinematics of the triple junction are first determined, i.e., the orientation of each GB and the misorientations between the three grains, therefore physically constraining the possible sets of realistic misorientations. The individual GBs are then represented using the disclination structural unit model (DSUM), in which a series of repeating disclination dipoles are used to achieve the desired misorientation along each GB as described in the previous section. In the last step, we introduce compensating disclinations in order to remove non-physical long-range stresses created by the use of finite-length (disclinationbased) GBs.

Triple-junction kinematics
The layout of each triple junction is determined by the kinematic equations from Upadhyay et al. [47]. With nine equations and seven variables, each TJ is determined by two degrees of freedom. Equation (13) The orientations of the three GBs β 12 , β 23 , and β 31 with respect to the reference x-axis are described as, The dihedral angles β 1 , β 2 , and β 3 between the GBs are then given by: For these simulations, the triple junctions are determined by varying the misorientation angles [θ 12 ] and [θ 23 ]. The angle between the orientation of the first grain and the  reference x-axis is assumed to be zero, i.e., the x-axis is aligned with the structure of the first grain. The remaining parameters are computed via Eqs. (13)- (15). The geometrical interpretation of these parameters is illustrated in Fig. 1.

Compensating disclinations
When constructing a finite-length GB using a repeating wall of disclination dipoles, it is necessary to include a compensating disclination at each end of the dipole wall [47,48] in order to suppress any non-physical, long-range stress field. The magnitude of the Franks vector associated with these compensating disclinations is equal to the misorientation associated with the wall of disclination dipoles. For GBs constructed using more than one wall of disclination dipoles, a set of compensating disclinations must be added to the ends of each wall. We illustrate in Fig. 7 the stress field around an isolated triple junction with and without these compensating disclinations. The associated point defect sinks arranged along the same triple junction is then shown in Fig. 8.
When modeling a larger network of grains, these compensating disclinations can account for instance, for non-perfect geometrical configurations where particular triple point where the sum of misorientations around any particular triple point is not a perfect multiple of 90 • . In this work, compensating disclinations allow a single triple point to be simulated using relatively short GBs without having to construct a surrounding network of compatible grains and GBs.

Simulation setup
Our spatially resolved rate theory framework combines the triple junction construct and defect-evolution governing equations described above. We used this framework in a series of 2D simulations of point-defect diffusion near triple junctions during irradiation. To simulate the addition of defects during irradiation, Frenkel pairs (i.e., size n = 1 vacancies and SIAs) were introduced throughout the control volume at a uniform rate of 0.001 dpa s −1 . Each simulation began with uniform defect concentrations and ran until steady state was reached. We use the sink efficiency ζ α (where the subscript α indicates the type of defect, e.g., a size n = 4 vacancy cluster) and bias factor B as metrics to evaluate the steady-state concentration profiles of point defects around each triple junction. These metrics measure the ability of the triple junction and its three associated GBs to act as a local "network" configuration (i.e., a triple junction configured as the assembly of three GBs as opposed to an isolated GB). Defect removal rate and average concentration are calculated within a 10-nm radius around the triple point in order to avoid influence from the simulation boundaries and size dependence from the GB lengths. The sink efficiency ζ α in the vicinity of the triple junction is defined as; wheref gives the average rate at which defects are removed from the computation domain via interaction with sinks, ρ s is the sink density is given by ρ s ,c is the average concentration, and c th is the concentration at thermal equilibrium. Because the concentration is significantly smaller at thermal equilibrium than during stead-state irradiation, we simplify Eq. (16) by assuming that c th α ≈ 0. The sink density is 1/ L x L y , indicating one triple junction per control volume with dimensions L x and L y . Sink density is equal across all simulations since we used the same control volume dimensions, and each simulation includes only a single triple junction. Similarly, the bias factor is defined as follows; Since this work includes multiple sizes of defect clusters, sink efficiency and bias factor are computed separately for each defect type and size. Size n = 1 vacancy clusters and size n = 3 interstitial clusters accounted for the majority of the point defect concentrations at steady state. The list of parameters used for this work can be found in Table 4. The range of triple junctions simulated in this work was determined by varying the two misorientation parameters θ 12 and θ 23 , which correspond to the misorientation angles of the first and second GBs. The angle θ 12 varied from 5 • to 85 • using 5 • increments. For each value of θ 12 , θ 23 was varied from 5 • up to the current value of θ 12 . Based on this range of values, the resulting misorientation θ 31 is always negative, and it is calculated along with the remaining parameters according to Eqs. (13)- (15). TJ configurations for which θ 12 + θ 23 = 90 • were excluded since this case is physically unrealistic. In total, we simulated 144 physically and geometrically constrained TJ configurations.
Constant-valued Dirichlet boundary conditions were used for all defect concentrations such that concentrations did not evolve along the simulation boundaries. The boundary values were determined by simulating defect implantation in a volume with no stress field and no defect sinks until equilibrium was reached between implantation and recombination. In this way, these boundaries assume a far-away bulk of defects at steady-state with Annihilation coefficient ann 20 -no sinks. Zero-traction boundaries were used for the mechanical equations, allowing the control volume as a whole to expand and contract without constraint. In determining the ideal size for our simulations, we repeated the same triple junction at both 50 × 50 nm 2 and 100 × 100 nm 2 . We found that doubling the simulation dimensions to 100 nm caused no change in sink efficiency for SIAs or vacancies, so we are confident that our use of 100×100 nm 2 for all triple junctions in this work is more than large enough to avoid boundary effects.

Numerical implementation
The mechanical equilibrium equations and the elastic response to the eigenstrains introduced by point defects are solved using 4-node linear finite elements, and a central finite difference scheme is used for the diffusion equations. Deformation gradients are computed numerically from the displacement fields using a central finite difference scheme. Each finite element used to solve the mechanical equations corresponds with a node in the diffusion grid; thus the finite element nodes are offset from the diffusion nodes. To calculate stresses at the diffusion nodes, the stresses at the Gauss points in each element are computed and averaged. The defect-defect interaction rates are calculated for each volume element using the concentration values from the corresponding finite difference node. The diffusion and mechanical equations are solved separately at each time step, with periodic boundary conditions used for both. Mechanical equilibrium is first solved in order to determine the stresses in each element, during which concentrations remain constant.
Diffusion rates are then calculated for the current time step, along with the defect-defect interaction rates. Updated concentration values are computed using an explicit integration scheme. The weak form, finite element formulation, periodic displacement boundaries, and finite difference scheme are described in further detail in our previous work [16].

Principal component analysis
The relatively large number of TJ configurations simulated, as well as the large number of variables, presented a challenge in analyzing the resulting data. To identify trends among the large number of variables characterizing the 144 TJs, we used Principal Component Analysis [49,50]. We consider each TJ as a local geometrical and structural network. As such, for each TJ, the TJ configuration was determined using two independent "configuration" input variables per Eqs. (13)- (15). Additionally, within the TJ, we characterized each individual GB by the density of sinks along its length, yielding three "structural" input variables. We considered the following "radiation tolerance" sink efficiency output variables, which quantify point defect absorption for each TJ: vacancy sink efficiency ζ V , SIA sink efficiency ζ I , and bias factor B. Since sink efficiency is calculated separately for each defect type and cluster size, there eight more variables: seven for sink efficiency and one for bias factor. This leads to at least thirteen variables, with 144 unique observations for each one.
Size n = 1 vacancies (V 1 ) and size n = 3 SIAs (I 3 ) were used when computing the output variables because these represented the majority of the point defects at steady state. Additionally, the sink efficiency values and bias factor at various sizes were found to correlate strongly, and including all sizes in the PCA would be redundant.
Though misorientation angles were used to create the TJ configurations, any two of the nine variables can be assigned independently when evaluating Eqs. (13)- (15). As such, the two independent geometrical variables can be freely chosen when analyzing the TJ sink efficiency data. For this work, the largest (β 2 ) and smallest (β 3 ) dihedral angles were used since they show the most compelling relationship between grain structure and sink efficiency.
As the dihedral angles represent the TJ structure, the sink density was chosen to represent the structure of the GBs themselves. The equivalent sink density was determined for each GB by considering the density of disclination dipoles from DSUM and then converting to the equivalent density of edge dislocations. This conversion is carried out by finding the number of edge dislocations with Burgers vector b that are equivalent to each disclination dipole [15]; Then, the sink density of a particular GB is; where N is the number of disclination dipoles per DSUM period, which corresponds to the number of minority structural units in one period. Including GB elastic energy proved to be redundant due to its close relation with the sink density, thus it was not included in the PCA. Before performing the PCA, the input and output variables were normalized such that each has a mean of zero and a standard deviation of one for the full set of 144 observations. Note that the goal of this analysis is to investigate correlations between input and output variables, rather than between input variables which are uncorrelated due to their random choice. Using the thirteen variables describing the defect behavior and defect accumulation near TJs, a n × m matrix X of data points is constructed (n = 144 and m = 13), where X ij represents the value of the input or output variable j obtained in dataset i. This matrix is then converted into a deviation matrix D, with its components D ij defined as, where,X j is the average value of the variable j over all datasets, and σ X j is the standard deviation of variable j over all datasets. The m × m correlation matrix R is constructed to provide the correlation between all input and output variables such that, The value for the correlation R ij between any two variables i and j is simply: using the Einstein convention to imply a sum over k.
The eigenvectors of the correlation matrix between all the variables are the principal components of the system. The first principal component PC 1 represents the direction in our thirteen dimensional space along which the value of PC 1j D ij has the maximum variance. The second principal component PC 2 is then defined as the direction orthogonal to PC 1 with maximum variance of PC 2j D ij , and so on. The variance of principal component PC i is given by the eigenvalue of that eigenvector.
To better understand the relationships between defect behaviors and damage accumulation near the TJs, rather than using a scatter plot or correlation matrix, we represented each variable as a function of the first two principal components [49,50] to visualize their relationships, also known as the correlation circle. The correlation PC plot shows vectors pointing away from the origin to represent the TJ descriptors. The angle between the vectors is an approximation of the correlation between the variables. A small angle indicates the variables are positively correlated, an angle of 90 • indicates the variables are not correlated, and an angle close to 180 • indicates the variables are negatively correlated. The length of the line indicates how well the variable is represented in the plot.