Alloy surface segregation in reactive environments: A first-principles atomistic thermodynamics study of Ag3Pd(111) in oxygen atmospheres

We present a first-principles atomistic thermodynamics framework to describe the structure, composition and segregation profile of an alloy surface in contact with a (reactive) environment. The method is illustrated with the application to a Ag3Pd(111) surface in an oxygen atmosphere, and we analyze trends in segregation, adsorption and surface free energies. We observe a wide range of oxygen adsorption energies on the various alloy surface configurations, including binding that is stronger than on a Pd(111) surface and weaker than that on a Ag(111) surface. This and the consideration of even small amounts of non-stoichiometries in the ordered bulk alloy are found to be crucial to accurately model the Pd surface segregation occurring in increasingly O-rich gas phases.


I. INTRODUCTION
The chemical properties of an alloy surface like its catalytic function or corrosion resistance depend intricately on the detailed surface composition. This compositional, and concomitant geometric and electronic structure is not only different from the corresponding surfaces of the parent metals, but differs often also substantially from the bulk alloy due to segregation of one species to the surface. The resulting segregation profile can depend on the surface orientation, the bulk composition, and temperature. A reactive environment, such as one where catalysis or corrosion occurs, can further modify the chemical composition of the surface, if one alloy component interacts more strongly with a gas phase species than the others, resulting in adsorbate-induced segregation. It is consequently usually not possible to predict the properties and functions of alloy surfaces in materials science applications solely on the basis of the formal bulk composition; the effects of the corresponding (possibly reactive) environment must explicitly be considered.
Previous theoretical attempts to address the segregation thermodynamics at an alloy surface from firstprinciples have mostly focused on the interdependence with the bulk alloy reservoir. This includes tabulated segregation energies of a transition metal impurity in another transition metal 1,2 , or approaches based on the coherent potential approximation 3 , screened generalized perturbation method 4 or surface cluster expansions [5][6][7][8] . The understanding emerging from these studies points at a crucial influence of the bulk reservoir structure and composition on the surface segregation profile, even if this concerns only slight excesses of one type of atom in samples with equal nominal bulk composition. This complexity has not yet been met in the scarce attempts to also explicitly account for the effect of a (reactive) environment [9][10][11] , but the results reveal nevertheless an at least equally sensitive dependence of the adsorbateinduced segregation on the detailed adsorption properties and gas phase conditions.
In view of these findings we set up a first-principles atomistic thermodynamics [12][13][14][15][16] framework based on density-functional theory that is able to address the structure and composition of an alloy surface while quantitatively accounting for the dependencies on both the bulk alloy and the gas phase reservoir. The method is applied to obtain the segregation profile of a Ag 3 Pd(111) surface in contact with an oxygen atmosphere, and we illustrate the intricate role of the bulk alloy reservoir by discussing the effect of already a small amount of nonstoichiometries (in form of anti-sites) in a nominally ordered Ag 3 Pd bulk alloy. Not surprisingly we find Agsegregation to occur under O-poor environmental conditions, e.g. typical for ultra-high vacuum experiments, while Pd-segregation is preferred with increasing oxygen content in the gas phase. Disentangling the various factors contributing to this segregation profile reveals, however, an intriguing complexity that needs to be captured in the modeling to obtain quantitatively accurate results. In the course of this analysis we can thus scrutinize various approximate treatments suggested in the literature, like the neglect of non-stoichiometries in the second or deeper substrate layers, the substitution of bulk segregation with that of dilute impurities in the parent metal hosts, or the interpolation of bonding properties at the alloy surface from those known at the parent metal surfaces.

A. Surface free energy
In the following we describe a first-principles atomistic thermodynamics framework to address alloy surface segregation in reactive environments. Our modeling aims at a system composed of a gas phase containing species j at partial pressures p j and a common temperature T , and a multi-component alloy A n A(bulk) B n B(bulk) . . . I n i(bulk) , where x i(bulk) = n i(bulk) / i n i(bulk) is the mole fraction of species i (i = A, B, ..., I) in the alloy bulk. To evaluate the surface free energy of a particular surface of this alloy in the given environment, we can concentrate on the part of the system that is affected by the surface, which contains (normalized to surface area A surf ) N j gas phase particles and a total of N solid particles at x i mole fractions. The most stable surface structure and composition is then the one, which minimizes the surface free energy, defined as where G(T, p j , N j , N, x i ) is the Gibbs free energy of the surface, and µ i and µ j are the chemical potentials of alloy species i and gas phase species j, respectively. For clarity, we will develop the methodological framework in this work for a binary alloy in contact with a single component gas phase, with obvious generalizations to more component systems. Specifically, we consider a Ag n Ag(bulk) Pd n P d(bulk) alloy and an oxygen atmosphere described by an oxygen pressure p O2 . In view of the total energy calculations described below, we will furthermore consider the surface affected part of the alloy to be described by an inversion-symmetric slab geometry exposing two identical surfaces. With A surf then denoting the area of the surface unit-cell, eq. (1) simplifies to where N slab refers now to the total number of metal atoms per surface unit-cell in the slab, N O similarly to the number of O atoms in the surface geometry, and in the second step we have exploited that the mole fractions sum to unity by construction (x Ag + x Pd = 1). Within our thermodynamic theory we assume the surface to be in equilibrium with the surrounding gas phase, and correspondingly the oxygen chemical potential is determined in the entire system by the gas phase reservoir, µ O = 1/2µ O2(gas) (T, p O2 ). The alloy surface is furthermore assumed to be in equilibrium with the underlying bulk alloy reservoir. This means that the chemical potentials µ Ag and µ Pd are not independent either, but are related by the bulk alloy chemical potential µ bulk , i.e. the Gibbs free energy per atom of the bulk alloy. For a bulk alloy of composition Ag n Ag(bulk) Pd n Pd(bulk) , the relationship between the chemical potentials is therefore where we have again used the constraint on mole fractions. With this equation and introducing the difference in chemical potentials of the two metal species, ∆µ Ag−Pd = µ Ag − µ Pd , we can rewrite eq. (2) to This equation has now an intuitive structure, with the first two terms giving the free energy difference between the surface structure and an equivalent number of bulk alloy atoms. The third term accounts for a possible difference in the surface and bulk stoichiometries and the fourth term for the additional number of O atoms present in the surface geometry. Finally, we note that throughout this text, we use the sign convention that a more negative Gibbs free energy will indicate a more stable state of the system. In the case of a gas phase chemical potential this translates to µ O2(gas) approaching −∞ in the limit of an infinitely dilute gas, since adding a particle will then yield an infinite gain in entropy. As a consequence, γ > 0 indicates the cost of creating the surface between the solid bulk phase and the homogeneous gas phase.

B. Dependence on gas phase and alloy bulk reservoir
The influence of the surrounding gas phase on the surface structure and composition enters eq. (3) through the term N O /2 µ O2(gas) (T, p O2 ). In increasingly O-rich environments, this term will favor surface structures containing an increasing number of oxygen atoms per surface area. In our approach we assess this effect by computing the surface free energy for a range of oxygen chemical potentials. At sufficiently low µ O2(gas) (O-poor conditions), only clean surface structures without oxygen atoms will be stable, which then forms a natural lower bound for the considered range of oxygen chemical potentials. As an upper bound, we use µ O2(gas) (O − rich) = E total O2(gas) , where E total O2(gas) is the total energy of an isolated O 2 molecule (including the zeropoint energy). This roughly corresponds to conditions where oxygen would start to condense at the sample at sufficiently low temperatures 16 .
At the accuracy level of interest to our study, the oxygen gas phase is well described by ideal gas laws, which then allows to readily relate the oxygen chemical potential to specific temperature and pressure conditions. For this we write 16 where k B is the Boltzmann constant, and ∆µ O contains now all temperature and pressure dependent free energy contributions of the O 2 molecule internal degrees of freedom. ∆µ O (T, p O2 ) can be obtained can be obtained within the ideal gas approximation, or alternatively, one can use tabulated enthalpy and entropy values at standard pressure p 0 = 1 atm to determine ∆µ O (T, p 0 ) 18 . For an oxygen gas phase, both approaches yield virtually identical results in the range of temperatures and pressures of interest to our study 16,17 . Apart from the dependence on the gas phase chemical potential, eq. (3) shows that the surface free energy depends also critically on the difference of the chemical potentials of the metal atoms ∆µ Ag−Pd . This variable reflects the change in bulk chemical potential due to changes in x Ag(bulk) , as can be seen by evaluating the derivative of eq. (2) x Ag(bulk) ∂µ Ag ∂x Ag(bulk) where the term in brackets is equal to zero by the Gibbs-Duhem equation 17 . The value of ∆µ Ag−Pd will therefore depend sensitively on the bulk reservoir structure and composition, e.g. whether the bulk reservoir is ordered or random, or whether there is a slight excess of one type of atom. This will vary from sample to sample, and can have a significant impact on the surface energetics, even for samples with equal nominal bulk composition 19 .
In our application below, we will demonstrate this crucial dependence on the detailed bulk reservoir structure and composition by considering the case of just a small amount of non-stoichiometries (in form of antisites) in the ordered Ag 3 Pd bulk alloy. Rather than focusing on specific data for a specific sample, the aim is thereby to discuss the general effects that such bulk antisites can have on adsorption-induced segregation. Instead of explicitly computing ∆µ Ag−Pd e.g. by means of cluster expansion techniques 6,8,20-23 , we therefore proceed by roughly estimating the range that ∆µ Ag−Pd can span, and discuss the results below for the obtained two bounds. For the estimate we use as limits the phase separation of Ag, when the bulk reservoir is rich in Ag (Ag-rich limit), and the phase separation of Pd, when the bulk reservoir is rich in Pd (Pd-rich limit). Phase separation of Ag will occur, when µ Ag > µ Ag,fcc , where µ Ag,fcc is the chemical potential of Ag atoms in the fcc bulk structure of the parent metal. Equivalently, phase separation of Pd will occur, when µ Pd > µ Pd,fcc , where µ Pd,fcc is the chemical potential of Pd atoms in the fcc bulk structure of the parent metal. Using eq. (2) and approximating µ bulk with its value for the perfectly stoichiometric alloy sample µ bulk(stoich) , we therefore arrive at and as a first estimate for the range of ∆µ Ag−Pd due to nonstoichiometries in the ordered bulk alloy reservoir.

C. Solid phase Gibbs free energies
In order to evaluate the surface free energy of a specific configuration with eqs. (3) and (6) we need to compute the solid phase Gibbs free energies, G slab , µ bulk(stoich) , µ Ag,fcc and µ Pd,fcc . For this it is useful to decompose them into several contributing terms, namely where E total is the total energy (excluding the zero point energy), F vib the vibrational free energy (including the zero point energy), and S conf. the configurational entropy. A crucial aspect that governs our analysis of all of these terms is that the quantity of interest to us, namely the surface free energy, does not depend on absolute Gibbs free energies. What enters into eq. (3) is a difference of Gibbs free energies of slab, bulk alloy and bulk metals. This can allow for quite some degree of cancellation, e.g. due to similar free energy contributions in the three systems or due to similar errors in the computed total energies. We start our analysis by noting that the pV term is completely negligible for the transition metal alloys discussed here 16 . As for the vibrational free energy, a proper evaluation would require a systematic computation of all vibrational frequencies at the surface and in the bulk systems 16 . Here, it can be particularly beneficial to realize that only the differences of the vibrational free energy contributions matter for eq. (3), i.e. the surface free energy is only affected by the changes of the vibrational modes. In the application discussed below, we will address the effect of on-surface O adsorption on the segregation profile of a Ag 3 Pd(111) surface. In this case, the vibrational modes of Ag and Pd atoms in the bulk alloy and bulk metal structures will be rather similar, as will be their vibrational modes at the close-packed surface. From the similarity of the calculated stretch mode of oxygen adsorbed at low coverage at Ag(111) and Pd(111), 50 meV 24 and 60 meV 25 respectively, we furthermore do not expect significant variations in the vibrational free energies of alloy surface structures with different metal composition. What will therefore mostly contribute are the changes of the vibrational modes of oxygen in the gas phase and adsorbed at the surface. To estimate this, we use a simple Einstein model equivalent to the one employed in Ref. 16 and compute the vibrational free energy contribution when the characteristic frequency is changed from ∼ 200 meV (O 2 stretch frequency in the gas phase) to ∼ 55 meV in the adsorbed state. Allowing for O coverages up to 1 ML and varying the characteristic frequency in the adsorbed state by ±50 %, already the total contribution to the surface free energy stayed always below 5 meV/Å 2 for temperatures ranging up to 600 K. What matters for the discussion of O-induced segregation is, however, the difference of this vibrational surface free energy contribution for different surface configurations. For surface configurations with the same O coverage (differing only in the Ag and Pd concentration in the surface layers) this difference will be completely negligible in view of the above cited similarity of the O-metal stretch modes. For surface configurations differing in O coverage, the difference in the surface vibrational free energy will be a bit larger, yet still negligible compared to the total energy terms discussed below which are of the order of tens of meV/Å 2 . Correspondingly, we will neglect the zero point and vibrational free energy contribution to the surface free energy in this study, but emphasize that this is not a general result, and the effect of vibrations needs to be carefully assessed in each application.
This leaves as remaining term the configurational entropy. A full evaluation of this contribution is computationally involved, since it requires a proper sampling of the huge configuration space spanned by all possible bulk and surface structures. Although modern statistical mechanics methods like Monte Carlo simulations 26,27 are particularly designed to efficiently fulfill this purpose, they still require a prohibitively large number of total energy evaluations to be directly linked with electronic structure theories 28 . A way to circumvent this problem is to map the real system onto a simpler, typically discretized model system, the Hamiltonian of which is sufficiently fast to evaluate. Cluster expansions (also termed lattice-gas hamiltonian approach) 6,8,20-23 are a prominent example for such a mapping approach, where the considered system would in practice e.g. be described by a lattice of possible on-surface adsorption sites for the gas phase species 10,41 . Here, we will instead concentrate on screening a number of known (or possibly relevant or representative) ordered surface structures by directly comparing which of them turns out to be most stable under which (T, p O2 )-conditions, i.e. which of them exhibits the lowest surface free energy. For sufficiently low temperatures, the remaining configurational entropy per surface area in G slab is then only due to a limited number of defects in these ordered surface structures and can be estimated to be below 2 meV/Å 2 for any T < 600 K 29 . We will restrict the discussion of the application to Ag 3 Pd(111) below to the temperature range well below this limit, where then also the consideration of an ordered Ag 3 Pd bulk structure with a small amount of anti-sites is justified 30 . Focusing on this temperature range we will neglect the configurational entropy contribution to G slab , as well as the small configurational entropy of the ordered stoichiometric alloy and metal bulk reservoirs (entering µ bulk(stoich) , µ Ag,fcc and µ Pd,fcc ). The only remaining temperature dependence in our model comes then from the entropy of the oxygen gas phase (entering ∆µ O (T, p)), as well as somehow implicitly through the considered bulk alloy chemical potential range, cf. eq. (6). This will not affect the conclusions drawn below, but we note that the additional limitation that comes with the then resulting direct screening approach is that its reliability is restricted to the number of considered configurations, or in other words that only the stability of those structures plugged in can be compared. Its predictive power extends therefore only to those structures that are directly considered, i.e. the existence of unanticipated surface geometries or stoichiometries cannot be predicted. With this in mind, the here employed direct screening approach to first-principles atomistic thermodynamics can still be a particularly valuable tool, since it allows, for example, to rapidly compare the stability of different structural models without the need to have them mapped onto a common lattice. In this sense, we view the present approach as complementary to a cluster expansion study, since it provides first insight into the possible effects of the surrounding gas phase by screening a wide range of possible structures. If more detailed insight is required (e.g. the exact surface concentration profile), this can then selectively be refined by a cluster expansion study, in particular when phase coexistence or order-disorder phenomena at elevated temperatures are of interest.
In the present work, we therefore approximate the difference of Gibbs free energies entering eq. (3) by the corresponding difference of the leading total energy terms, and our final working equation reads What is then required in practice to evaluate the surface energy of a given surface geometry of Ag 3 Pd (111)  of an isolated O 2 molecule. The considered range for the bulk alloy chemical potential difference ∆µ total Ag−Pd is given through eq. (6), for which (with the equivalent approximations) total energy calculations of the parent Ag and Pd bulk metals (E total Ag,fcc , E total Pd,fcc ) are additionally required. If we define the formation energy for the bulk alloy as then this range spanned by ∆µ Ag−Pd is given by

D. Total energy calculations
We compute all total energies within pseudopotential plane wave density-functional theory (DFT) 31 , using the generalized gradient approximation (GGA) 32 for the exchange-correlation functional. The employed Pd ultrasoft pseudopotential (USPP) 33-35 has a 4d 9 5s 1 electronic configuration and is based on a Koelling-Harmon relativistic all-electron calculation 36 . The cutoff radii were set at 2.5, 2.5, 2.0 bohr for the s, p and d channels. The Ag USPP has a 4d 10 5s 1 electronic configuration, with cutoff radii set at 2.5, 2.5, 2.1 bohr for the s, p and d channels. For both elements the f angular momentum channel is used as the local potential, with two s, p and d projectors, and a non-linear core correction 37 is used with r core = 1.2 bohr. The O USPP has a 2s 2 2p 4 electronic configuration, with cutoff radii set at 1.3 bohr for the s and p channels. The d angular momentum channel is used as the local potential, with two s and p projectors, and a non-linear core correction with r core = 0.5 bohr. We found it advantageous to use a less accurate oxygen USPP with a larger partial core radius (r core = 0.7 bohr) to pre-relax our geometries at a lower computational cost. Substituting the more accurate O USPP in a final calculation, we found the maximum forces to be very close to or less than the force criteria set in the initial optimization (< 0.05 eV/Å). In several cases where the maximum forces exceeded the optimization criteria, we reoptimized the geometry, but found that the surface energies changed by less than 1 meV/Å 2 , which is accepted as a negligible error in the remainder of the calculations. A plane wave cutoff of 350 eV was employed in all calculations using the less accurate O USPP, while we employed the double grid method 34 in the calculations using the more accurate O USPP. The higher energy cutoff (700 eV) with which the density is treated in the lat-FIG. 1: Top views of the atomic geometry of five (2×2) surface unit-cell arrangements considered in this study. The shown configurations have a bulk like second layer composition, but a first layer composition that varies from pure Ag (I) to pure Pd (V). White colored spheres correspond to Ag atoms, whereas grey colored spheres correspond to Pd atoms. The larger spheres correspond to the first layer, and the smaller ones correspond to the second layer.
ter method increases the computational time moderately, but resulted in smooth convergence behavior 38 .
In order to ensure maximum compatibility between the bulk metal and bulk alloy calculations, we employed the same bulk unit-cell to describe both the L1 2 Ag 3 Pd and the fcc Pd and fcc Ag structures. (20 × 20 × 20) Monkhorst-Pack k-point grids were then used to optimize the lattice constants to a 0 (Ag 3 Pd) = 4.09Å, a 0 (Pd) = 3.93Å, and a 0 (Ag) = 4.09Å. Using these optimized lattice constants, we obtain a formation energy ∆E f (Ag 3 Pd) as defined in eq. (11) of −0.21 eV per formula unit. The (111) alloy surfaces are modeled in a periodic supercell geometry, employing 7-layer inversionsymmetric slabs with two identical surfaces. The outer two layers on each side of the slab were fully relaxed, keeping the inner three layers fixed. A vacuum region of 14Å ensures the decoupling of consecutive slabs. The Brillouin zone integrations were performed with (8×8×1) Monkhorst-Pack k-point grids for the (2×2) surface unitcells. For the oxygen adsorption calculations, oxygen was placed on each side of the slabs, maintaining inversion symmetry in all cases. For the computation of E total O2(gas) we employed a (12 × 13 × 14) bohr cell and obtain an O 2 binding energy of −5.91 eV. From test calculations with up to twice as many k-points and an increased plane wave cutoff of 500 eV we estimate the numerical accuracy of the reported surface energies to be within 2 meV/Å 2 , which does not affect any of the physical conclusions drawn.

A. O-induced Pd segregation at Ag3Pd(111)
We illustrate our approach by addressing the Oinduced segregation profile at a Ag 3 Pd(111) surface. At the nominal stoichiometry each (111) layer in the assumed L1 2 bulk crystal structure 30 contains 3 Ag atoms and 1 Pd atom in a (2 × 2) arrangement. We account for a possible segregation at the surface by considering (2 × 2) surface unit-cell atomic configurations in which the composition in the topmost two layers is systematically varied from pure Ag to pure Pd. This leads to 25 possible substrate configurations, five of which are shown in Fig. 1. In addition to these clean surface structures, oxygen adsorption was then examined by placing from one up to four O atoms per surface unit-cell into the onsurface fcc sites. This leads to a pool of 473 inequivalent surface structures (including some larger unit cells), in which the oxygen coverage varies from 0 to 1 monolayer (ML) in 1/4 ML steps in the (2 × 2) surface unit-cells.
Evaluating eq. (10) for all configurations we arrive at the surface free energy plots shown in Figs. 2 and 3. Figure 2 corresponds to the limiting case of a Pd-rich Ag 3 Pd bulk reservoir and Fig. 3 to the opposite limit of a Ag-rich Ag 3 Pd bulk reservoir as discussed in Sec- tion IIB, cf. eq. (6). Within our thermodynamic theory, at each oxygen chemical potential the surface configuration yielding the lowest surface free energy (i.e. the lowest line in the figures) results as "most stable". Focusing first on the Pd-rich situation, we observe different stable configurations in different oxygen chemical potential ranges. At the lowest chemical potentials considered, the lowest surface free energies are naturally obtained for the clean surface structures with varying chemical composition in the topmost two layers. Among those, the "most stable" surface configuration under these O-poor environmental conditions (i.e. the lowest horizontal line in Fig. 2) corresponds to a Ag-terminated alloy, with 100 % Ag in the first layer and 100 % Pd in the second layer. As expected 39 , we therefore obtain segregation of the more noble metal to the surface under gas phase conditions that at room temperature are representative of ultra-high vacuum experiments, cf. the pressure scales in Fig. 2.
Due to the last term in eq. (10), increasingly Orich surface structures get more favorable with increasing ∆µ O . The higher the O concentration (N O /A surf ) at the surface, the steeper the slope of the lines in Figs. 2 and 3. In the case of the Pd-rich bulk limit, a first oxygen containing surface structure gets more stable than the clean Ag-segregated surface above ∆µ O = −0.89 eV. This configuration corresponds to an oxygen coverage of a 1/4 ML and contains now already 75% Pd atoms in the first layer. At even higher O chemical potentials further surface structures with higher oxygen coverages get stabilized, each time then containing a first layer composition of only Pd atoms. Induced by the oxygen adsorption is thus a segregation of the more reactive metal to the surface, which effectively reverses the segregation profile corresponding to room temperature ultra-high vacuum gas phase conditions. While this intuitive trend comes out clearly from our theory accounting for the effect of the environment, the exact surface concentrations as reflected by the obtained "most stable" structures have to be considered with care. Our approach considers a finite set of ordered surface structures and neglects the configurational entropy contribution to the solid phase Gibbs free energies as discussed in Section IIC. The closeness of the lowest energy lines in Fig. 2 reflects a number of structures that are just a little bit less stable at a given ∆µ O than the highlighted "most stable" structure. This indicates already that a more exhaustive search of the huge configuration space e.g. by means of a cluster expansion technique would likely reveal new ground state structures with compositions that are intermediate to those considered here and that exhibit a certain degree of disorder at finite temperatures.
While a quantitative determination of the segregation profile is therefore outside the scope of the present work, the direct comparison of Figs. 2 and 3 nicely illustrates the entangled dependencies of the surface composition on both the bulk alloy and the gas phase reservoir. The results for the Ag-rich bulk alloy reservoir summarized in Fig. 3 also indicate an O-induced Pd segregation at increasing ∆µ O . One major difference to the afore discussed Pd-rich bulk alloy reservoir is that the stabilization of these O-containing Pd-enriched surface structures occurs at higher oxygen chemical potentials. This obviously reflects the higher cost at which the additional Pd surface atoms have to be taken out of the Ag-rich alloy bulk reservoir, and correspondingly a higher driving force from the gas phase is required to establish the same Pd content in the surface layers as compared to the Pd-rich bulk alloy reservoir. A closer look at Figs. 2 and 3, and equivalently inspection of eq. (10), reveals, however, that the effect of the different bulk alloy reservoir is not just a mere shift of all surface free energy curves to higher ∆µ O by an amount that corresponds to the range of ∆µ Ag−Pd between our two considered bulk alloy reservoir limits, cf. eq. (12). Within our computational setup this range amounts to 0.28 eV, and this is indeed roughly the order of magnitude by which the lines in Fig. 3 are shifted compared to those in Fig. 2. Still, due to the prefactor (x Ag − x Ag(bulk) ) of the ∆µ Ag−Pd -term in eq. (10) surface energies of different surface configurations are affected differently by the change in the bulk reservoir. Accordingly, we observe also a different sequence of "most stable" structures in Fig. 3 compared to Fig. 2, which underscores the need to explicitly account for both the effect of the environment and the bulk alloy structure to describe the structure and composition of alloy surfaces in realistic applications.
Another intricacy that is nicely illustrated by our data is the effect of compositional changes in the deeper substrate layers. As apparent from Figs. 2 and 3 the firstlayer Pd-segregation induced in more oxygen-rich environments is accompanied by substantial changes of the second layer composition. These changes have quite some impact on the stability of the various surface structures as is revealed when restricting our structure data base to those structures with nominal Ag 3 Pd second layer stoichiometry. When evaluating eq. (10) for this subset of surface structures we obtain as "most stable" clean surface for the case of a Pd-rich bulk reservoir a structure that contains only 75% Ag in the first layer (in contrast to the 100% Ag termination when also considering nonstoichiometries in the second substrate layer). While the remaining sequence of "most stable" structures with increasing Pd content in the first layer is similar to the one shown in Figs. 2 and 3 respectively, the oxygen chemical potentials at which the transitions between the different structures occur are shifted by sometimes more than 0.1 eV compared to the values exhibited in Figs. 2 and 3. As indicated by the drawn pressure scales in the figures this can correspond to several orders of magnitude in pressure in the temperature range discussed here, from which we conclude that a quantitative modeling of the Ag 3 Pd(111) surface requires at least the consideration of non-stoichiometries in the topmost two substrate layers.

B. Contributing factors to adsorbate-induced segregation
The results presented in the previous section clearly show that the adsorption-induced segregation profile at the surface results from the compromise between the tendency to lower the surface free energy by segregating the more noble metal species to the surface, and the tendency to lower the surface free energy by achieving stronger adsorbate bonding at surfaces enriched with the more reactive metal species. In order to further qualify these two competing trends we define as a fundamental quantity describing the adsorbate bonding the average oxygen binding energy E b (N O , x Ag , x Pd ) at a given surface con-figuration containing N O adsorbed O atoms per surface unit-cell and with a detailed surface composition in the surface layers characterized by mole fractions of x Ag Ag atoms and x Pd Pd atoms, . (13) As before this definition holds directly for symmetric slabs with adsorption on both sides, and with the chosen sign convention and zero reference in form of the gas phase O 2 molecule an average binding energy E b < 0 indicates that dissociative adsorption is exothermic.
Similarly we define as a quantity describing the surface segregation the segregation energy E seg (x Ag , x Pd ) as the cost to change the reference clean alloy surface with nominal bulk stoichiometric composition into an alternative clean surface configuration with concentrations of x Ag Ag atoms and x Pd Pd atoms in the topmost layers (either by rearrangement of the surface atoms including atoms in sub-surface layers, or by segregation of atoms to or from the bulk reservoir), The factor 1/2 accounts as before for the fact that we are considering symmetric slabs where the segregation occurs equally at both surfaces, and within our sign convention a segregation energy E seg < 0 indicates that segregation is exothermic. With these definitions, we can now recast our working equation for the surface free energy, eq. (10), into an intuitive structure, where γ bulk−term. is the surface free energy of the reference clean alloy surface with nominal bulk stoichiometric composition. This form of the equation shows clearly that the surface free energy can either be lowered by an energetically favorable segregation of metal atoms (E seg < 0), or by exothermic adsorption of N O oxygen atoms (E b < 0), where increasing oxygen chemical potentials in the gas phase further favor configurations with an increasing number of oxygen atoms at the surface.
Having disentangled these various contributing factors, we now proceed to discuss our detailed results for each one separately. This enables us to carve out in detail the complexity behind each of them, thereby bringing us into the position to scrutinize various approximate treatments that have been suggested in the literature. Figure 4 shows the segregation energies as defined in eq. (14) for the 25 clean Ag 3 Pd(111) surface configurations computed. As expected, the segregation energy decreases overall largely with increasing Ag content in the topmost layer, which is the driving force behind the stabilization of the Ag-terminated clean alloy surface structures at the lowest oxygen chemical potentials in Figs. 2 and 3. In detail, this variation of the segregation energy is, however, not simple, and depends much on the bulk alloy reservoir and the second layer composition. For the 100% Ag-terminated surfaces, the lowest segregation energy is for example obtained for a 0% Pd concentration in the second layer in the case of the Ag-rich bulk reser-voir, but for a 100% Pd concentration in the second layer in the case of the Pd-rich bulk reservoir.

Segregation energy
This complexity needs to be captured in the modeling to obtain quantitatively accurate results for the segregation profile, and thus dictates an explicit consideration of the structure and composition of both alloy bulk and surface. A prevalent alternative approach is instead to obtain estimates for the segregation energy by computing the energy cost to segregate an impurity atom in a metal host to the surface 1,2 . In the present case, this would correspond to computing a diluted Ag impurity in the bulk or at the (111) surface of fcc Pd, or a diluted Pd impurity in the bulk or at the (111) surface of fcc Ag. Within our computational setup, we obtain such impurity segregation energy estimates by calculating the total energy difference when once placing a Ag atom in the central layer and once in the topmost layer of a seven layer (2 × 2) Pd(111) slab (at the optimized Pd lattice constant), as well as respectively doing the same for a Pd atom and a Ag(111) slab (at the optimized Ag lattice constant). The values obtained are −0.31 eV for the case of a Ag atom in Pd, and +0.24 eV for the case of a Pd atom in Ag, which compare well to the equivalent numbers computed by Ruban and coworkers 2 using the local-density approximation as exchange-correlation functional (-0.26 eV and +0.28 eV, respectively).
Using these numbers the total segregation energy of a given surface configuration is then simply obtained by appropriately summing the impurity segregation costs. For a 100% Pd-terminated surface this corresponds e.g. to moving three Ag atoms per (2 × 2) surface unit-cell to the bulk and moving three Pd atoms per (2 × 2) surface unit-cell from the bulk (−3 × (−0.31 eV) + 3 × 0.24 eV = 1.65 eV). Normalized to the area of the (2 × 2) Ag 3 Pd(111) surface unit-cell, this yields E seg /A surf = 57 meV/Å 2 in this simple model. The blue line in Fig.  4 shows the corresponding data for such variations in the topmost layer Ag concentration, which reproduces the overall trend of the calculated DFT data surprisingly well (in particular when comparing to the filled symbols representing the DFT data where also the second layer composition has been kept at the bulk stoichiometry). By construction, the model can however not account for the dependence on the bulk alloy reservoir and leads thus to quantitative errors of the order of several 10 meV/Å 2 in the E seg /A surf contribution to the surface free energy. This is quite sizable considering the very magnitude of E seg /A surf and γ reflected in Figs. 2, 3 and 4, and underscores the necessity for a more refined treatment of the segregation energy contribution which does account for the detailed alloy bulk and surface structure.

Oxygen binding energy
A similar complexity as for the segregation energy is observed, when we turn to the third term in eq. (15) and therewith to the term that is directly responsible for the oxygen-induced segregation. Starting the discussion with the lowest O coverage that can be treated within our (2 × 2) surface unit-cells, i.e. 0.25 ML with 1 O atom per surface unit-cell, Fig. 5 shows the variation of the average binding energy with the number of first-layer Ag atoms to which the adsorbed O atom is directly coordinated. As expected, less exothermic binding energies are obtained with an increasing number of Ag atoms immediately involved in the oxygen bonding. This is in line with the known stronger bond strength of oxygen at the more reactive Pd(111) surface compared to the more noble Ag(111) surface. However, the situation at the alloy surface is not a mere linear composition weighted average between these two limiting cases as also shown in Fig. 5. In fact, we find average binding energies at Pd-rich Ag 3 Pd(111) surface configurations that are even stronger than on Pd(111), and average binding energies at Ag-rich Ag 3 Pd(111) surface configurations that are even weaker than on Ag(111). Judging from the different, known binding properties at the parent metal surfaces would therefore significantly underestimate the effect of oxygen-induced Pd surface segregation at this alloy sur-face, i.e. also in this case a more simplified treatment that does not explicitly account for the alloy surface structure would not be appropriate. The range of average binding energies apparent in Fig 15) of up to ∼ 80 meV/Å 2 , which is then of comparable order, if not larger than the corresponding variation of the E seg /A surf term with nonstoichiometries in the second layer. The non-trivial interplay between these two equally-sized contributions to the surface free energy is thus what dictates at least the consideration of non-stoichiometries in the topmost two layers when aiming at a quantitative description of the adsorbate-induced segregation at this Ag 3 Pd(111) surface.
From our data, we also arrive at a first estimate of the impact that the third layer atoms have on the surface adsorption properties. For surface configurations with pure first and second layer compositions, i.e. either 100 % Ag or 100 % Pd in either of the two layers, the four fcc adsorption sites in our (2 × 2) surface unit-cells are completely equivalent, except for variations in the type of atoms in the third bulk-stoichiometric layer (i.e. with three Ag atoms and one Pd atom per (2 × 2) cell). In these cases, we obtained variations of up to 0.1 eV per O atom in E b , depending on whether the oxygen atom was located directly over a Pd atom in the third layer or not. Compared to the average binding energy variations due to differing compositions in the topmost two layers, these are still noticeable, but already much smaller variations. On this basis, we would therefore conclude that a consideration of non-stoichiometries in the topmost two layers should capture the essential physics of oxygen-induced segregation at Ag 3 Pd(111).
Apart from the dependence on the surface composition the average O binding energy also displays a strong variation as a function of coverage. Figure 6 illustrates this for the subset of computed configurations with either a 100 % Ag or a 100 % Pd concentration in the topmost layer. In both cases, we observe a strong reduction in bond strength with increasing number of oxygen adsorbates at the surface. For the Ag-terminated surfaces, this has the consequence that dissociative adsorption be- comes endothermic at all configurations with coverages θ > 0.25 ML, which is very similar to the findings at the pure Ag(111) surface 40 . Within the scatter in the data shown in Fig. 6, which is due to the varying second layer composition of the surface configurations displayed, the decrease of the average binding energy with coverage is roughly linear for both Ag and Pd terminated surfaces. Again, this is similar to the bonding at the parent metal surfaces Ag(111) and Pd(111) 40 . The amount, by which E b decreases with coverage is, however, not the same for the two differently terminated surfaces. Whereas for the Pd-terminated surfaces the average binding energy decreases by about 1 eV in the computed coverage range 0.25 ≤ θ ≤ 1 ML, it decreases for the Ag-terminated surfaces by about 1.5 eV. In other words, the anyway stronger bonding at Pd-enriched surfaces at low O coverages becomes even more favored at increasing O coverage in comparison to Ag-enriched surfaces with the same oxygen load. This will further enhance the tendency towards Pd surface segregation at increasing O chemical potentials, and shows that a restricted information about low coverage bonding properties will also not be sufficient to correctly model the segregation profile in a reactive envi-ronment like the presently studied oxygen atmospheres.

Gas phase chemical potential
Even though the average binding energies per O atom decrease with increasing oxygen coverage, the total contribution of N O times the smaller E b in eq. (15) may still lead to lowest surface free energies for surface configurations with higher oxygen load. At increasing oxygen chemical potentials such structures are further favored by the last term in eq. (15), i.e. the increasing abundance of oxygen in the gas phase helps to stabilize O-rich surface structures. Within our calculated set of surface configurations, this leads to the sequence of "most stable" structures in Figs. 2 and 3, with increasing concentrations of adsorbed O atoms at higher values of ∆µ O . At sufficiently high oxygen chemical potentials also more O-rich surface configurations beyond the simple on-surface adlayers considered in the present study may get stabilized. Apart from ultra-thin so-called surface oxide films, which have been characterized for both parent metal surfaces 42,43 , this will ultimately be thicker bulk-like oxide films, and due to their higher stability most likely palladium oxides.
Assuming that what forms at the surface is simply the most stable palladium bulk oxide, PdO, we can estimate the oxygen chemical potential at which such a bulk oxidation occurs from the stability condition 44 Within the same approximations as applied in Section IIC, we then arrive in the limit of the Pd-rich bulk reservoir (µ Pd = µ Pd(metal) ) at where the formation energy of bulk PdO is defined as Similarly, one obtains in the limit of the Ag-rich bulk reservoir (µ Ag = µ Ag(metal) ) ∆µ bulkoxid.,Ag−rich with the alloy formation energy as defined in eq. (11).
Within our computational setup, we determine a value of E f (PdO) = −0.86 eV per formula unit, so that the chemical potentials corresponding to the onset of bulk oxide formation are ∆µ bulkoxid.,Pd−rich O = −0.86 eV and ∆µ bulkoxid.,Ag−rich O = −0.65 eV. In both cases, the onset of bulk oxide formation occurs therefore just at slightly higher oxygen chemical potentials than the first stabilization of surface configurations with adsorbed O atoms within the set of structures compared in our study (at ∆µ O = −0.89 eV and −0.68 eV respectively, cf. Figs. 2 and 3). This is largely different compared to the situation at Pd(111), where a p(2 × 2) structure with 0.25 ML on-surface O adatoms in fcc hollow sites gets stabilized at oxygen chemical potentials that are significantly lower than the onset of bulk oxide formation 44 . The reason for this difference is the high cost of segregating the scarce Pd in Ag 3 Pd to the surface that is not necessary at Pd(111). The average O binding energy we compute at 0.25 ML coverage at the Ag-rich surface termination (100% Ag in the first layer, 100% Pd in the second layer) that is most stable at the lowest ∆µ O is only −0.13 eV/O atom, which should be compared to our computed average binding energy of −1.29 eV/O atom at the pure Pd(111) surface at this coverage. Only by segregating Pd atoms to the first layer (and thus involving them into the direct coordination to the O adsorbate), does the E b of corresponding Ag 3 Pd(111) configurations come into the same range as the Pd(111) binding energy, e.g. −1.29 eV/O atom for the configuration with 100% Pd in the first layer and 100% Ag in the second layer, cf. Fig.  5. However, stabilizing such configurations then involves the cost of segregating Pd to the topmost layer and thus delays their stability range up to oxygen chemical potentials, where bulk oxide formation is already about to set in.

IV. SUMMARY
In conclusion, we have presented a first-principles atomistic thermodynamics approach to describe the structure and composition of an alloy surface in contact with a (reactive) environment. Accounting for both the dependency on the bulk alloy and the gas phase reservoir, this approach provides the appropriate framework for a detailed discussion of the driving factors behind adsorbate-induced surface segregation. We illustrate this for the Ag 3 Pd(111) surface exposed to an oxygen atmosphere, where we obtain an inversion of the segregation profile in increasingly oxygen-rich environments. Whereas a minimal segregation energy stabilizes Ag-terminated surface structures in ultra-high vacuum type gas phase conditions, the much stronger oxygen bonding favors increasingly Pd-rich terminations in atmospheres with higher oxygen content. Our analysis shows that the details of this transition are intricately coupled to the bulk reservoir structure, and already small non-stoichiometries in the nominally ordered Ag 3 Pd bulk alloy can not only change the oxygen pressures required for the transition by several orders of magnitude, but also affect the explicit segregation profile itself. This highlights the restricted relevance of detailed segregation data from ultra-high vacuum experiments for applications in realistic environments, as well as the complexity that needs to be captured in the modeling to obtain quantitatively accurate results. In this respect our study shows that at least for the Ag 3 Pd(111) surface simplified treatments based on impurity segregation and adsorption data of the parent metals are not able to correctly describe the adsorbate-induced segregation, nor that restricted information about the low-coverage adsorbate bonding at the alloy surface or the exclusive consideration of non-stoichiometries in the topmost substrate layer would be sufficient. This holds already for the case where the reactive environment leads only to the formation of on-surface adsorbate layers, and is even more pronounced when considering the possibility of a more complex restructuring like here the oxide formation in an O-containing environment. Also in the latter case the ruling factors are the availability of the more reactive metal species in the alloy bulk and the increased oxygen bond strength that this species can provide. In the studied case of the Ag 3 Pd(111) surface the high cost of segregating the scarce Pd to the surface delays the stabilization of Pd-enriched surface structures in fact almost up to gas phase conditions where the formation of thicker palladium oxide films would already become thermodynamically favorable. We speculate that this will be different in more Pd-rich Ag-Pd alloys, or in alloys composed of species with a larger difference in their O binding capabilities.