Grain Boundary Segregation in Pd-Cu-Ag Alloys for High Permeability Hydrogen Separation Membranes

Dense metal membranes that are based on palladium (Pd) are promising for hydrogen separation and production due to their high selectivity and permeability. Optimization of alloy composition has normally focused on bulk properties, but there is growing evidence that grain boundaries (GBs) play a crucial role in the overall performance of membranes. The present study provides parameters and analyses of GBs in the ternary Pd-Ag-Cu system, based on first-principles electronic structure calculations. The segregation tendency of Cu, Ag, and vacancies towards 12 different coherent ∑ GBs in Pd was quantified using three different procedures for relaxation of supercell lattice constants, representing the outer bounds of infinitely elastic and stiff lattice around the GBs. This demonstrated a clear linear correlation between the excess volume and the GB energy when volume relaxation was allowed for. The point defects were attracted by most of the GBs that were investigated. Realistic atomic-scale models of binary Pd-Cu and ternary Pd-Cu-Ag alloys were created for the ∑5(210) boundary, in which the strong GB segregation tendency was affirmed. This is a starting point for more targeted engineering of alloys and grain structure in dense metal membranes and related systems.


Introduction
Cost-effective production of ultra-pure hydrogen can facilitate the widespread implementation of fuel cells and is one of the remaining bottlenecks before hydrogen can be introduced as an energy carrier on a large scale [1]. Dense metal membranes with high hydrogen permeance and selectivity have been identified as a promising enabling technology for efficiency improvement and cost reduction of hydrogen production. In particular, Pd-based hydrogen separation membranes are known to have 100% selectivity and high permeability, and thus allow for direct production of high purity hydrogen for use in fuel cells [2][3][4][5]. Combining these membranes with appropriate catalysts in membrane reactors to produce hydrogen from different sources has been described in numerous studies [1,[6][7][8].
It appears that the potential of binary Pd-based membranes has been exhausted in the literature, and several groups have recently started working on ternary compounds as the next generation membrane material [9][10][11][12][13][14][15][16][17][18][19][20]. This has many possible benefits: the surface can be engineered to enhance the tolerance to impurity gases [13], the permeability can be optimized beyond what is possible with binary alloys [14], and the mechanical strength can be increased (e.g., if the self-diffusion is hindered or the morphology is changed) [2]. There is also potentially a cost reduction that is involved if expensive elements are replaced with cheaper ones [21]. The challenge with this approach is that ternary compounds are difficult to engineer when using plating, rolling, etc. as processing techniques [14]. One solution is to use a non-equilibrium process like magnetron sputtering to synthesize the active

Methodology
The calculations made use of the Vienna ab initio simulation package (VASP) [42,43], employing plane-wave basis functions with the projector augmented wave method [44] and the density functional theory (DFT) at the Perdew-Burke-Ernzerhof generalized gradient approximation (GGA) level [45]. The self-consistency requirement was changes in the total electronic energy less than 10 −5 eV. The force relaxation criterion was 0.01 eV/Å, using the RMM-DIIS quasi-Newton method. Choosing a plane-wave cut-off energy of 500 eV and a k-point density of at least four points per Å −1 gave a numerical precision of better than 1 meV for relative total electronic energies. The GB models have been implemented in the supercell scheme, and the various supercells used are shown in Figure 1.
Three different schemes were employed as volume relaxation techniques (VRT): no relaxation, full relaxation of all cell parameters (including both cell size and shape), and only the relaxation of the cell length perpendicular to the GB plane (x). The atomic positions were relaxed in all the cases. Even if the present models do not physically resemble grains in a real material (the "grains" are sheets with infinite extension perpendicular to x), we can still learn about real materials by assessing how the different VRTs correspond to limiting cases of large or small domains in different directions. With this approach, the situation when only x is allowed to relax corresponds to the kind of relaxation that would find place for infinitely large grains, i.e., one single GB. With only one semi-infinite GB, the x direction will be fully relaxed, while the infinitely large lattice perpendicular to the GB is equivalent to an infinitely stiff lattice in the other directions. The relaxation along very large grains should thus be well represented by that found in our models when only x is allowed to relax. This option is not normally available in VASP but has been facilitated in a local version of the code. The relaxation procedures when either all or none lattice parameters are allowed to relax do not correspond directly to any physical distribution of GB in a real material. However, the two methods represent the outer limits of the kind of relaxation that would take place in a material with very small grains. In this case, relaxation perpendicular to any GB (corresponding to x in our models) will be countered by nearby GBs with other orientations. Similarly, relaxation along GBs will be partially allowed due to the finite size of the GB and the short distance to neighbor grains in all directions. Thus, the relaxation taking place in a material with very small grains should be in between the relaxation that is found when keeping the lattice fixed and when it is fully relaxed. We will in the following present results for all VRTs, hence describing relaxation effects that are likely to be seen in materials with large and small grain size.

Structure and Stability of Pure Pd GB Models
Initially, the structure and stability of pure Pd GB models are investigated. The 12 different GBs investigated in this work are listed in Table 1. They are based on the coincidence site lattice model [46] and constitute all coherent tilt GBs with ∑ up to 13, thus exhibiting a quite wide range of deviations from the pure crystal. They have been represented by periodic atomistic models, as shown in Figure 1, where pure Pd models without point defects are displayed. These models range between 44 and 100 atoms (see Table 1), and their size in the x direction is between 13 Å (∑5(310)) and 40 Å (∑13(510)). Each model has two cancelling GBs: One at the unit cell boundary and one in the middle, identified by the dotted, red lines. The relatively large variation in size is partially due to symmetry (the models are often based on the smallest possible supercells with two such cancelling GBs) and the differing interaction range between various GBs-complex GBs typically display strain fields with larger range than simple ones. Some inherent and calculated properties of the GB models are listed in Table 1. The excess volume per interface area V X is defined as the difference between the DFT relaxed unit cell volume and the corresponding volume of Pd atoms in the bulk, divided by the cross-section area A of the GB plane (perpendicular to x). The unit of V X is length and it can naively be interpreted as the accumulated extension of the lattice in the x direction due to the GB (when divided by 2, since there are two GBs per unit cell). The GB energy is defined as where the total energy E of the GB and of the bulk is that calculated by DFT, N GB is the number of atoms in the GB model, and the factor 2 is due to the presence of two GBs in each unit cell. The bond range ∆r b is the difference between the longest and shortest relaxed bond within the first coordination sphere of the atoms in each model. It is centered around the DFT equilibrium Pd-Pd distance of 2.80 Å with the shortest recorded relaxed bond distance being 2.43 Å and the longest one 3.46 Å; this gives ∆r b up to 1.04 Å, which is seen for the ∑5(210) model. The distinction between the coordination spheres gets unclear for some of the models, and we have used a cut-off of 3.5 Å on the far side of this definition, which typically is a minimum between the first and second coordination spheres in these systems. We first note that V X depends strongly on the VRT. This is not surprising, since the calculations with no volume relaxation are at the mercy of the initial model. When the unit cell is allowed to relax partially (only x) or fully, the spread in V X is smaller and the values of each model is consistent between the two procedures allowing for relaxation. V X is consistently smaller in the case of full relaxation, since a part of the volume expansion then is obtained within the GB plane. Keeping in mind that the only x relaxation corresponds to large grains while the relaxation of samples with small grains should be in between that of full relaxation and no relaxation, it appears that V X should be quite similar in magnitude, regardless of the typical grain size; the only x value is usually between the two others. The values vary between that of the ∑3(111) model (0.01-0.06 Å) and the ∑ 7(415) model (0.70-0.90 Å). These numbers may not be fully converged due to the limited size of the unit cells but are clear indicators of the deviation from perfect stacking. The former model can be viewed merely as a stacking fault, while the latter has quite large deviations from the perfect structure along the GB (see Figure 1). The different periodic grain boundary (GB) models included in this study. The dotted, red lines distinguish the GB planes. Unit cells are outlined with black lines, and the x direction is marked by the arrow. xgb is defined as the horizontal distance from the GB planes. The supercell size perpendicular to the figure plane (the z direction) is listed in Table 1. The atomic positions plotted in Figure 3 have been shown as green dots connected by a solid line for the Σ3(112 ) model.
The same pattern is found when studying the GB energy ; there is a clear correlation between VX and for the volume relaxed models, while the same is not the case when the unit cell was not relaxed. This correspondence is plotted in Figure 2a, showing that a very clear linear relation between VX and is obtained. The fit gives the empirical relation = 1.55 and = 1.19 for full and only x relaxation, correspondingly. The value is 0.97 and 0.93 for the two fits. This indicates that the volume misfit is an excellent predictor for the GB energy in the case of coherent defect-free GBs, which can potentially be used in experimental studies where VX may be easier to measure than . A similar correlation has been found previously for other materials, e.g., Ni [47].
The bond range ∆ is also following the excess volume and GB energy quite closely, if not as clearly as is the case between and (Figure 2b). Nevertheless, there is a clear correlation between and ∆ -not surprising, since both parameters are a measure of the deviation from the perfect bulk crystal structure. It is perhaps more surprising that ∆ is lower when no volume relaxation is allowed than in the relaxed case, as is seen e.g., for the Σ3(112 ) model. This can be explained from the larger freedom of the volume relaxed models to accommodate strain by changing the local coordination within some of the models. Despite the spread in bond range and accompanying variation of coordination number (number of nearest neighbors) in the near vicinity of the GBs, all of the models x Figure 1. The different periodic grain boundary (GB) models included in this study. The dotted, red lines distinguish the GB planes. Unit cells are outlined with black lines, and the x direction is marked by the arrow. x gb is defined as the horizontal distance from the GB planes. The supercell size perpendicular to the figure plane (the z direction) is listed in Table 1. The atomic positions plotted in Figure 3 have been shown as green dots connected by a solid line for the ∑ 3(112) model. The same pattern is found when studying the GB energy γ; there is a clear correlation between V X and γ for the volume relaxed models, while the same is not the case when the unit cell was not relaxed. This correspondence is plotted in Figure 2a, showing that a very clear linear relation between V X and γ is obtained. The fit gives the empirical relation γ = 1.55 V X and γ = 1.19 V X for full and only x relaxation, correspondingly. The R 2 value is 0.97 and 0.93 for the two fits. This indicates that the volume misfit is an excellent predictor for the GB energy in the case of coherent defect-free GBs, which can potentially be used in experimental studies where V X may be easier to measure than γ. A similar correlation has been found previously for other materials, e.g., Ni [47]. Table 1. The various grain boundary (GB) models investigated in this study, along with the number of atoms in their unit cells, the supercell size in the z direction l z (perpendicular to the figure plane in Figure 1), their excess volume divided by area V X , the calculated GB energy E GB , and the range of bond distances ∆r b within the first coordination shell. All results are for pure Pd GBs. Three volume relaxation techniques (VRT) were used: no relaxation of the unit cell, only relaxation of the x lattice constant (see Figure 1 for a definition of the x direction), and full relaxation of all cell parameters.

Segregation of Single Point Defects
The distribution of defects in the vicinity of GBs is determined by their relative energy at different positions and kinetics. We have only focused on energy in this work, since this is the most relevant in systems being allowed to equilibrate (which is the case in most membrane systems.) We have investigated three different defects with particular relevance for Pd membranes: substitutional Ag (AgPd) and Cu (CuPd), as well as vacancies (VPd). This was done with the three different VRTs described above. The energy of each vacancy was calculated at various sites with increasing distance from the GB; an example of this is shown with the green curve for the Σ3(112 ) model in Figure 1. All the sites with unique x coordinates, from the GB to the midpoint between neighbor GBs were  The bond range ∆r b is also following the excess volume and GB energy quite closely, if not as clearly as is the case between γ and V X (Figure 2b). Nevertheless, there is a clear correlation between V X and ∆r b -not surprising, since both parameters are a measure of the deviation from the perfect bulk crystal structure. It is perhaps more surprising that ∆r b is lower when no volume relaxation is allowed than in the relaxed case, as is seen e.g., for the ∑ 3(112) model. This can be explained from the larger freedom of the volume relaxed models to accommodate strain by changing the local coordination within some of the models. Despite the spread in bond range and accompanying variation of coordination number (number of nearest neighbors) in the near vicinity of the GBs, all of the models display local atomic structures very close to fcc in between the GB regions. This ensures that close to bulk behavior can be found furthest away from the GBs in the models.

Segregation of Single Point Defects
The distribution of defects in the vicinity of GBs is determined by their relative energy at different positions and kinetics. We have only focused on energy in this work, since this is the most relevant in systems being allowed to equilibrate (which is the case in most membrane systems.) We have investigated three different defects with particular relevance for Pd membranes: substitutional Ag (Ag Pd ) and Cu (Cu Pd ), as well as vacancies (V Pd ). This was done with the three different VRTs described above. The energy of each vacancy was calculated at various sites with increasing distance from the GB; an example of this is shown with the green curve for the ∑ 3(112) model in Figure 1. All the sites with unique x coordinates, from the GB to the midpoint between neighbor GBs were included (recall that cancelling GBs are present both at the unit size boundaries and at the red, dashed line.) The energy at the GB (x = 0) was used as the reference state, and the relative energy E gb between that of the impurity located at x = x gb and x = 0 was calculated for all 12 models. This is plotted in Figure 3 in the case of Cu Pd defects with relaxation only along the x direction. Since the number of unique sites in the x direction differs between the various models, this is also the case with the number of plot points in Figure 3. The same applies to the distance between neighbor GBs, which is why the extension of the various curves varies in Figure 3. Similar plots to the one shown in Figure 3 have been generated for the three different VRT using three different defects (Cu impurities, Ag impurities, and vacancies); in total, nine plots. For simplicity the only relax x with Cu impurities is the only series of plots shown here.
In order to quantify the tendency to segregate towards or away from the GB the segregation energy of a specific defect has been defined as follows. The lowest total energy (as calculated by DFT) of the defect among the three sites nearest the GB is taken as "the" energy at the GB, E GB (def). This is compared to the average value of the three energies at the farthest distances away from the GB, defined as the "bulk" value E bulk (def). The segregation energy of the defect is then defined as An example is shown in Figure 3: the three values with highest x gb (fitted with a black, dotted line) are used to calculate E bulk (Cu Pd ) for ∑13(510), and the point at x gb ≈ 1 Å (defining the lower black, dotted line) gives E GB (Cu Pd ). The resulting segregation energy of this example (marked with black arrows in the figure) is E segr (Cu Pd ) = −0.20 eV.
The model size should ideally be large enough for the energy to converge for large values of x gb . This can be seen for some of the models in Figure 3, but not all. Computational cost restricted the use of larger unit cells to achieve better convergence with respect to unit cell size. Nevertheless, some of the models exhibit excellent convergence as the x gb increases. The ∑13(510) model is one example, where the energy does not change more than 0.03 eV when x gb increases from 5 to 10 Å at 16 different sites, and the difference between the three sites with largest x gb is less than 1 meV. Other models fluctuate more, but most of them exhibit a clear trend when x gb increases. In some cases, there is no well-defined "bulk" energy; as an example, the energy of the ∑3(112) model varies significantly. This may be due to problems that are connected with the relatively small unit cells employed, and we have therefore in the following disregarded models where the average deviation from E bulk is larger than 0.03 eV. Larger deviations than this are designating models with severe relaxation effects that are deemed as unphysical.
well-defined "bulk" energy; as an example, the energy of the ∑3(11 _ 2) model varies significantly. This may be due to problems that are connected with the relatively small unit cells employed, and we have therefore in the following disregarded models where the average deviation from Ebulk is larger than 0.03 eV. Larger deviations than this are designating models with severe relaxation effects that are deemed as unphysical. Lines are drawn as guide to the eye. The segregation energy Esegr is taken as the difference between the average energy of the three sites furthest away from the GB and the lowest energy among the three nearest sites. The extension of a curve along the xGB axis corresponds to half the size of the unit cell in the x direction. See the text for details.
The segregation energy was used to characterize the behavior of the three defects CuPd, VPd, and AgPd for all three VRTs and all 12 GB models; this has been compiled in Figure 4. A number of models displayed unphysical relaxations, and their Esegr values have not been included in the figure. This leaves part of the figure without data, but enough results were generated to draw some general conclusions. The most striking feature of Figure 4 is that almost all the segregation energies are negative, which indicates a tendency for all three defects to segregate towards the GB. The only small Lines are drawn as guide to the eye. The segregation energy E segr is taken as the difference between the average energy of the three sites furthest away from the GB and the lowest energy among the three nearest sites. The extension of a curve along the x GB axis corresponds to half the size of the unit cell in the x direction. See the text for details.
The segregation energy was used to characterize the behavior of the three defects Cu Pd , V Pd , and Ag Pd for all three VRTs and all 12 GB models; this has been compiled in Figure 4. A number of models displayed unphysical relaxations, and their E segr values have not been included in the figure. This leaves part of the figure without data, but enough results were generated to draw some general conclusions. The most striking feature of Figure 4 is that almost all the segregation energies are negative, which indicates a tendency for all three defects to segregate towards the GB. The only small exceptions are for models where no volume relaxation was allowed, indicating that such relaxation is necessary in order to accommodate the defects at these specific GBs. The ∑3(111) twin boundary is special; since the GB is merely a stacking fault, there is nothing to gain from moving a point defect towards or away from the GB. E segr is thus very close to zero for all defects and VRTs for ∑3(111). All other GBs exhibit significant values of E segr for some or all defects and VRTs. The Cu Pd defect gives the smallest absolute values of E segr for most systems, indicating that the segregation tendency of Cu towards GBs in Pd is, in general, lower than that of Ag or vacancies.
Based on the computed segregation energies the equilibrium concentration at the GB of a point defect c def at a given temperature T can be found from the following equation [48]: where c 0 def is the overall ("bulk") equilibrium defect concentration, and k B is Boltzmann's number. This formula assumes negligible interaction between solutes and thermodynamic equilibrium; we shall see later in this paper that the higher concentration of solutes may actually increase the anticipated equilibrium concentration in many cases. The resulting concentration has been shown in Figure 5 for the V Pd , Cu Pd , and Ag Pd defects in Pd at T = 600 K, which is a relevant temperature for hydrogen separation membranes. The behavior is quite similar at higher temperatures (not shown), only with defect concentrations slightly closer to the bulk one (selected here to be 0.02). The negative segregation energies are reflected in defect concentrations significantly higher than c 0 def in most of the cases. In the example in Figure 3 (∑3(510)), this means that the site at layer 2 (directly next to the GB) has approximately 50% occupancy by the Cu solutes.
exceptions are for models where no volume relaxation was allowed, indicating that such relaxation is necessary in order to accommodate the defects at these specific GBs. The ∑3(111) twin boundary is special; since the GB is merely a stacking fault, there is nothing to gain from moving a point defect towards or away from the GB. Esegr is thus very close to zero for all defects and VRTs for ∑3(111). All other GBs exhibit significant values of Esegr for some or all defects and VRTs. The CuPd defect gives the smallest absolute values of Esegr for most systems, indicating that the segregation tendency of Cu towards GBs in Pd is, in general, lower than that of Ag or vacancies. Based on the computed segregation energies the equilibrium concentration at the GB of a point defect cdef at a given temperature T can be found from the following equation [48]: where def is the overall ("bulk") equilibrium defect concentration, and kB is Boltzmann's number. This formula assumes negligible interaction between solutes and thermodynamic equilibrium; we shall see later in this paper that the higher concentration of solutes may actually increase the anticipated equilibrium concentration in many cases. The resulting concentration has been shown in Figure 5 for the VPd, CuPd, and AgPd defects in Pd at T = 600 K, which is a relevant temperature for hydrogen separation membranes. The behavior is quite similar at higher temperatures (not shown), only with defect concentrations slightly closer to the bulk one (selected here to be 0.02). The negative segregation energies are reflected in defect concentrations significantly higher than def in most of the cases. In the example in Figure 3 (∑3(510)), this means that the site at layer 2 (directly next to the GB) has approximately 50% occupancy by the Cu solutes.    Table 1 and depicted in Figure 1 at T = 600 K. It is shown for the defects CuPd (red circles), AgPd (blue squares), and VPd (black triangles), using three VRTs: full volume relaxation (filled symbols), relaxation along the x direction only (halffilled symbols), and no volume relaxation (open symbols). The overall (bulk) defect concentration (chosen to be 0.02) is shown by the black, dotted line.
The variation of the defect concentration with temperature has been depicted for the ∑11(113) GB in Figure 6. The concentrations approach the equilibrium bulk concentration of 0.02 (black, longdashed line) as the temperature increases. The Cu impurity does not show a strong segregation tendency. In the case of only x relaxation (which corresponds to large grains) there is a weak segregation of Cu towards the GB. When no relaxation is allowed, however, there is a similarly weak segregation of Cu away from the GB. Since the situation with small grains corresponds to an   The variation of the defect concentration with temperature has been depicted for the ∑11(113) GB in Figure 6. The concentrations approach the equilibrium bulk concentration of 0.02 (black, long-dashed line) as the temperature increases. The Cu impurity does not show a strong segregation tendency. In the case of only x relaxation (which corresponds to large grains) there is a weak segregation of Cu towards the GB. When no relaxation is allowed, however, there is a similarly weak segregation of Cu away from the GB. Since the situation with small grains corresponds to an interpolation between no relaxation and full relaxation, we conclude that Cu may be weakly segregated away from the GB in this case. Ag is, on the other hand, quite strongly segregated towards the GB and there is virtually no difference between large (only x) and small (between no relaxation and full relaxation) grains. Vacancies are also strongly attracted to the GB, and the difference between materials with large and small grain size should be quite small-interpolated values between the curves based on no relaxation and full relaxation are likely to be very similar to the values of the curve based on only x relaxed. We can deduce the temperature dependence of the other GBs in Figure 5 from the behavior of the curves in Figure 6, since the segregation energy is the only parameter that determines the defect concentration in Equation (3).

Binary Systems with More Than One Impurity Atom
Even if there is a clear tendency of segregation of impurities towards the GB in many GB models, it is unclear from the above whether more than one atom can be attracted to the GBs simultaneously. Repulsive interactions between point defects may lead to lower concentrations than that predicted in Figure 5. This was therefore investigated in more detail for Cu in Pd within a ∑5(210) GB model, as shown in Figure 7. The ∑5(210) GB was selected as a typical representative of the coherent GBs investigated in this study; it has a characteristic GB energy (~1 eV), range of bond lengths ∆ ~1 Å, as well as segregation energies of Ag (~−0.3 eV) and Cu (~−0.15 eV). This leads to a clear segregation behavior of solitary solutes with at least 21% (Cu) or 57% (Ag) equilibrium occupancy near the GB at 600 K when the overall concentration is 2%.
Since Cu and Ag segregate towards different sites, one can expect the situation of both solutes segregating simultaneously towards the same GB. To investigate this possibility a periodic GB model with stoichiometry Pd3Cu was generated as a starting point; a fully relaxed model is shown in Figure 7b. The bond lengths of the relaxed model are color coded in this figure. The relatively large variation of bond lengths with both elongated and contracted bonds explain why both large and small atoms may be attracted to the GB from a geometric point of view; this creates sites where atoms of various radii could be fitted geometrically. The bond lengths of this particular model vary from 2.31 to 3.29 Å,

Binary Systems with More Than One Impurity Atom
Even if there is a clear tendency of segregation of impurities towards the GB in many GB models, it is unclear from the above whether more than one atom can be attracted to the GBs simultaneously. Repulsive interactions between point defects may lead to lower concentrations than that predicted in Figure 5. This was therefore investigated in more detail for Cu in Pd within a ∑5(210) GB model, as shown in Figure 7. The ∑5(210) GB was selected as a typical representative of the coherent GBs investigated in this study; it has a characteristic GB energy (~1 eV), range of bond lengths ∆r b~1 Å, as well as segregation energies of Ag (~−0.3 eV) and Cu (~−0.15 eV). This leads to a clear segregation behavior of solitary solutes with at least 21% (Cu) or 57% (Ag) equilibrium occupancy near the GB at 600 K when the overall concentration is 2%. results, it can be expected that all "small" sites (with short interatomic distances to the neighbor sites) close to this boundary are occupied by Cu. This can be quantified by ∆ , which indirectly indicates the size of the smallest sites. From Table 1, it is evident that small sites exist in all of the models of this study (except the ∑3(111) twin boundary), and in many cases to a larger degree than in the ∑5(210) model. We can thus expect that the segregation trend of Cu is global, and that these results can be transferred to almost all GBs.

Segregation in Ternary Pd-Cu-Ag Alloys
The results of the above studies clearly show that there is a strong segregation of various point defects towards most GBs, and that this is valid even for a large density of solutes. However, it is not clear how different solutes interact with each other. From the size of defects as compared to the available sites near the GB (the "size effect"), one could expect that many GBs display a combined segregation Since Cu and Ag segregate towards different sites, one can expect the situation of both solutes segregating simultaneously towards the same GB. To investigate this possibility a periodic GB model with stoichiometry Pd 3 Cu was generated as a starting point; a fully relaxed model is shown in Figure 7b. The bond lengths of the relaxed model are color coded in this figure. The relatively large variation of bond lengths with both elongated and contracted bonds explain why both large and small atoms may be attracted to the GB from a geometric point of view; this creates sites where atoms of various radii could be fitted geometrically. The bond lengths of this particular model vary from 2.31 to 3.29 Å, corresponding to ∆r b = 0.98 Å. This is significantly larger than ∆r b of the pure Pd ∑5(210) model, which only displayed bonds between 2.51 and 2.99 Å, thus ∆r b = 0.48 Å (Table 1). This means that the addition of Cu not only decreases the smallest bond length (which can be expected when a smaller atom is added), but it also increases the largest one. This reflects a higher flexibility of the lattice when atoms with more than one size are present.
Models with the Cu content reduced below the starting point of 25% were created by replacing Cu by Pd in the model in Figure 7b. The stability of these models was assessed by their formation energy defined as where E tot is the total electronic energy as calculated by DFT, N is the number of atoms in the GB model (listed in Table 1), and n is the number of Pd atoms being substituted by Cu. The reference energy of Pd and Cu is that of their standard state, fcc bulk. Due to the difference in standard state energy between Pd and Cu, increasing the Pd content typically increases the formation energy (it appears less stable). The Cu content is reduced to 23.75% when one Cu atom is replaced, and the most stable configuration of this model is with extra Pd placed at position 1 or 11 (Figure 7a, dashed blue curve with diamonds), i.e., at the very center of the GB. The models with less Cu also exhibit the most stable configuration when excess Pd is placed at position 1 or 11. As an example, the most stable configuration of the 22.5% Cu model is with extra Pd placed at position 1 and 11 (dashed red curve with empty squares). The latter configuration is used as the starting point for the 20% model, which displays the most stable configuration with an extra Pd placed at position 3 or 9. A number of models with 20% Cu and Pd placed according to this insight were then constructed. The most stable of those is shown in Figure 7c, where all of the Cu atoms were moved to positions 2 and 10, corresponding to E form = 1.9 eV. This is marked by the arrow at position 3 in (a). In conclusion, there is a very strong thermodynamic driving force for segregation of Cu towards the coherent ∑5(210) GB. From these results, it can be expected that all "small" sites (with short interatomic distances to the neighbor sites) close to this boundary are occupied by Cu. This can be quantified by ∆r b , which indirectly indicates the size of the smallest sites. From Table 1, it is evident that small sites exist in all of the models of this study (except the ∑3(111) twin boundary), and in many cases to a larger degree than in the ∑5(210) model. We can thus expect that the segregation trend of Cu is global, and that these results can be transferred to almost all GBs.

Segregation in Ternary Pd-Cu-Ag Alloys
The results of the above studies clearly show that there is a strong segregation of various point defects towards most GBs, and that this is valid even for a large density of solutes. However, it is not clear how different solutes interact with each other. From the size of defects as compared to the available sites near the GB (the "size effect"), one could expect that many GBs display a combined segregation of "small" and "large" defects (smaller and larger than the host atoms, respectively). We used a selection of Pd-Cu-Ag alloys with the composition Pd 0.8−δ Cu 0.2 Ag δ to investigate this hypothesis, again using the ∑5(210) GB as a representative model GB. The most stable models from Figure 7 were used as starting point, substituting Pd with Ag at various sites and concentrations.
The formation energy E form of a Pd-Cu-Ag alloy is defined as similar to that of Pd-Cu in Equation (4): Here, m is the number of Pd atoms that have been substituted by Ag, and the reference energy of Ag is that of the standard state, fcc bulk. This definition means that a negative E form indicates a stable compound compared to the pure metals.
The formation energy is plotted for different models with Ag substituted for Pd or Cu as a function of the site (which corresponds to the distance from the GB) in Figure 8. We recognize the trend from the binary alloy; Ag is most stable at the GB (represented by the positions 1 and 11) in the model where Pd has replaced Cu at positions 1 and 3. It is also the most stable in the model where two Cu atoms in addition have been segregated to site 2 and 10, but to a smaller extent. However, in the most stable situation when all the Cu atoms are moved towards the GB and only populate sites 2 and 10 (corresponding to Figure 7c), the Ag site with lowest energy is not anymore at the GB (position 1 or 11), but rather at sites 3 and 9, which are just outside the Cu sites. This is due to relatively large local relaxations around the Cu atoms, which reduce the size of the sites at the GB (position 1 and 11). These sites are thus not significantly smaller than other sites anymore. Adding to this is the reduced affinity between Ag and Cu when compared to that between Ag and Pd. Since the latter model is the most stable one, we can expect to find the enrichment of Cu in position 2 and of Ag in position 3 close to ∑5(210) GBs in Pd-Ag-Cu compounds. This is supported by the plot in Figure 8b, demonstrating that the formation energy decreases as the local Ag concentration increases.
These sites are thus not significantly smaller than other sites anymore. Adding to this is the reduced affinity between Ag and Cu when compared to that between Ag and Pd. Since the latter model is the most stable one, we can expect to find the enrichment of Cu in position 2 and of Ag in position 3 close to ∑5(210) GBs in Pd-Ag-Cu compounds. This is supported by the plot in Figure 8b, demonstrating that the formation energy decreases as the local Ag concentration increases.  (5)) of three fully relaxed ∑5(210) Pd64Cu16 models with Ag substituted at the sites defined in Figure 7c. One model was created from the periodic Pd60Cu20 model by substituting Pd for Cu at the positions 1, 3, 9, and 11 (blue line with diamonds), one had in addition moved two Cu atoms to positions 2 and 10 (red line with squares), and one had all Cu atoms located at positions 2 and 10 (green line with circles). The latter corresponds to the model shown in Figure 7c. The two most stable models from (a) were used to plot the formation energy as a function of local Ag concentration n in (b); Ag was then placed at position 1 (red circles) and position 3 (blue circles). The most stable model with n = 10% is shown in Figure 9b.  Figure 7c. The two most stable models from (a) were used to plot the formation energy as a function of local Ag concentration n in (b); Ag was then placed at position 1 (red circles) and position 3 (blue circles). The most stable model with n = 10% is shown in Figure 9b. How does the addition and segregation of Ag influence the bond distances at and near the ∑5(210) GB? This may be relevant for hydrogen solubility and diffusivity since both depend strongly on the interatomic distances. Figure 9 presents how this is quantified by the excess volume divided by GB area in Figure 9a and by actual bond distances in the most stable model with the composition Pd70Cu20Ag10 in Figure 9b. It is evident that an increasing amount of Ag at the GB leads to significantly increased excess volumes. Since the volumes are divided by the GB area (that of the relaxed unit cell), this primarily signifies the elongation of the GB unit cell in the x direction. The Pd-Pd-bonds in the near-GB region are thus clearly larger than those in the bulk Pd, and significantly more so than in the bulk Pd-Cu alloys, where the overall lattice constants are reduced due to the smaller size of Cu atoms.

Discussion
The 12 GB models depicted in Figure 1 are by no means representing all the possible GBs in fcc Pd, even when being restricted to coherent ones without compensating dislocations or other defects. The relevance of the present results is thus restricted to a selection of such interfaces, which may not govern all important properties of these compounds. Nevertheless, we have seen that the selection of models gives a broad distribution of important parameters, like the interface energy, mismatch volume, range of bond lengths, etc. This indicates that the present results should be relevant at least for all coherent GBs in Pd alloys, even those with lower angles than the present ones. How does the addition and segregation of Ag influence the bond distances at and near the ∑5(210) GB? This may be relevant for hydrogen solubility and diffusivity since both depend strongly on the interatomic distances. Figure 9 presents how this is quantified by the excess volume divided by GB area in Figure 9a Figure 9b. It is evident that an increasing amount of Ag at the GB leads to significantly increased excess volumes. Since the volumes are divided by the GB area (that of the relaxed unit cell), this primarily signifies the elongation of the GB unit cell in the x direction. The Pd-Pd-bonds in the near-GB region are thus clearly larger than those in the bulk Pd, and significantly more so than in the bulk Pd-Cu alloys, where the overall lattice constants are reduced due to the smaller size of Cu atoms.

Discussion
The 12 GB models depicted in Figure 1 are by no means representing all the possible GBs in fcc Pd, even when being restricted to coherent ones without compensating dislocations or other defects. The relevance of the present results is thus restricted to a selection of such interfaces, which may not govern all important properties of these compounds. Nevertheless, we have seen that the selection of models gives a broad distribution of important parameters, like the interface energy, mismatch volume, range of bond lengths, etc. This indicates that the present results should be relevant at least for all coherent GBs in Pd alloys, even those with lower angles than the present ones.
The GBs in this study are all perfectly coherent, which of course is a simplification of the real situation. Many of the GBs found in real materials are quite complex, featuring all sorts of defects that to some extent compensate geometric mismatches that are inherent to the perfect GB. Nevertheless, many of these additional defects can be relatively far apart (e.g., in the case of dislocations), leaving near-perfect GBs over large areas, as demonstrated by several microscopy studies [49,50]. We therefore assume that our coherent models may be also relevant for a number of GBs where the lack of coherence is not too severe. We expect the correspondence to fail when going to truly amorphous GBs.
The results above suggest significant segregation of a variety of point defects towards virtually all GBs. However, the absolute numbers in Figure 5, summarizing the segregation tendency, should be applied with some care due to a number of reasons. The limited size of the models means that the strain originating from the GBs is not converged to zero at any place in the super cell. This challenge has been accommodated to some extent by using the three different volume relaxation methods-they represent the outer limits of how the unit cells should realistically be relaxed in the vicinity of GBs, and the real segregation tendencies should be somewhere between those limits. So even if there is no true bulk behavior between the GBs in the various models, the local relaxation effects and the resulting segregation should be correct within the boundaries that are described by the different VRTs.
Another potential source of error in these calculations is the lack of structurally compensating defects, most notably dislocations. They can accommodate significant parts of local strain and could as such counter some of the strongest segregation tendencies that were seen in this study. However, dislocations can attract point defects themselves, so this does not necessarily hinder segregation of defects, even if the nature of the segregation might be changed.
Another reason to take the absolute numbers of Figure 4 with some care is the rather large relaxation effects that were seen for some of the models. In some cases, the entire supercell was restructured, which led to a total energy being reduced by several eV in some of the cases. This can be understood as the starting point of a full relaxation to the lowest energy structure, which is the bulk without any GB. We have disregarded the points with largest restructuring effects, but it was difficult to distinguish between reasonable relaxations and unphysical effects due to the small size of the supercells; there was a continuous range from virtually no relaxation to almost complete reorganization of the GB model. We disregarded models where the average deviation from E bulk is larger than 0.03 eV, but some unphysical results due to limited unit cell size may still remain.
The calculations of the present study have all been performed without any explicit temperature being included. Temperature was included implicitly through the Arrhenius equation when obtaining equilibrium defect concentrations in Figure 5, but no other effect of temperature (thermal expansion, entropy, zero-point energy, phonon-based thermodynamics, etc.) was included. This has the potential to change the quantitative results significantly, but we expect that the qualitative trends remain unchanged.
Many of the same concerns apply when turning to the ternary compositions. We established that both larger (Ag) and smaller (Cu) atoms are attracted to a GB simultaneously, but the absolute value of the numbers and the actual sites of attraction may be different in reality than what is reported in this study. Furthermore, we did not consider simultaneous segregation of vacancies and solutes. With the knowledge from above, we expect that this is present in most boundaries, since vacancies and solutes have the ability to occupy different sites around a given GB.
The various sources of potential errors thus add up to a large and unknown uncertainty of the numbers presented in this study. The remaining main conclusion is unchanged, however: there is a clear and consistent tendency of many kinds of point defects (small and large substitutes, vacancies) to segregate towards GBs, due to the variety of local environments that are found there.
Despite the clear trends, it may be challenging to observe the segregation experimentally. It happens on the scale of single atomic layers, which makes high-resolution transmission electron microscopy the only viable way of directly probing such segregation. This relies on the ability to focus the electron beam along the GB plane, since this is where the change in concentration should be observed. It may also rely on the synthesis and heat treatment of the film; e.g., since sputtering is a non-equilibrium synthesis process, annealing (or operation under realistic conditions) may be required for the segregation to appear. Furthermore, Pd-membranes are typically aimed for hydrogen separation purposes, and the presence of hydrogen may influence the segregation behavior significantly, as has been seen in the case of segregation towards outer surfaces of similar systems [51].
Three different VRTs were compared for many of the calculations above, representing the boundaries of the likely volume relaxation regimes in real materials (from relaxation along x signifying low density of GBs to no relaxation representing very high density of GBs). However, some of the results indicate that not performing any relaxation of the volume is unreasonable. This is particularly clearly illustrated in Figure 2, where the linear trend between V X and γ can only be seen if the partial or full relaxation of the model is allowed. The conclusion from this observation is that some relaxation of the volume is required to move away from unreasonable situations arising from the somewhat arbitrary construction of the GB models. Most of the results in this study that are based on no volume relaxation should thus be neglected, apart from serving as a far-end borderline of the values.
Which of the VRTs allowing for volume relaxation to choose is less obvious. Both partial (only x) and full relaxation of the unit cell give results in reasonable correspondence with each other, indicating that relaxation of the coordinate perpendicular to the GB plane is most important. There are some differences between the two, depending on the particulars of the GB; most notably, the ∑ 11(233) GB exhibits large differences between only x and full relaxation, both for V X , γ, and ∆r b (Table 1). This reflects that this particular GB displays significant relaxation of the unit cell parameters parallel to the GB plane when allowed, in contrast to the other GB models. Such relaxation is most reasonable when the real GB resembles our simulation cell: infinitely long in the directions parallel to the GB plane but with a short distance between the GBs in the x direction. Partial relaxation (only x) describes situations where the average bulk lattice constant is maintained in these parallel directions by an infinitely stiff lattice, i.e., extending far away from the GB.

Conclusions
Atomic-scale calculations based on density functional theory were used to investigate various properties of low-number coherent ∑ grain boundaries (GBs) in fcc Pd, Pd-Cu, and Pd-Cu-Ag alloys. Their excess volumes, grain boundary energies, and ranges of bond lengths were computed while using three different volume relaxation techniques: no relaxation of volume, full relaxation of all lattice parameters, and partial relaxation of volume only allowing for one lattice constant to change. A linear correlation between the excess volume and grain boundary energy was found in pure Pd when partial or full relaxation of the volume was allowed. The tendency to segregate towards the GBs was assessed for three different point defects: Cu Pd , Ag Pd , and V Pd . Virtually all GBs exhibited strong segregation tendencies for all defects, quantified by the segregation energy and the corresponding equilibrium concentration of the defects at or near the GBs. An amplified segregation tendency was observed when increasing the local solute concentration, indicating that many of these sites might be nearly fully occupied by substitutional defects. This also demonstrated that the initial study at the dilute limit is clearly relevant for alloys with higher concentrations of the alloying element. Ternary Pd-Cu-Ag compounds were finally investigated, and simultaneous segregation of both Cu and Ag towards the ∑5(210) GB was observed; i.e., the results also hold for more complex alloy systems. The most stable model furthermore displayed a pronounced increase in the excess volume, indicating significantly increased local lattice parameters. In summary, this study demonstrates how insight from first principles calculations can be used to understand the complex behavior of point defects at and around grain boundaries of metals and alloys. Funding: Financial support from the Research council of Norway through the CLIMIT program (Contract No. 215666/E20) is gratefully acknowledged. The computational part was executed with a grant from the Notur metacenter for supercomputing.