S IMULATIONS OF HYDROGEN OUTGASSING AND STICKING COEFFICIENTS AT A COPPER ELECTRODE SURFACE : D EPENDENCIES ON TEMPERATURE , INCIDENT ANGLE AND ENERGY

Outgassing remains a pertinent issue as it typically is the ﬁrst stage of possible plasma formation, and can lead to effects such as breakdown, surface ﬂashover, and pulse shortening in high power systems. Here two pertinent aspects are probed: (i) a model-based assessment of outgassing and associated temperature-dependent rates from a copper electrode based on molecular dynamics simulations, and (ii) calculations for the sticking coefﬁcients of hydrogen gas atoms as a function of incident energy and angle. Our results of temperature dependent diffusion coefﬁcients for hydrogen in copper, agree well with experimental reports over a wide range from 300 to 1350 K, and show reduction in the presence of vacancies. Results also show low reﬂection coefﬁcients at both high and low energies, with a maxima at around 6.5 eV. A curve ﬁt to the data is predicted to roughly hold for a range of incident angles. Adsorption is predicted to occur for incident energies below 10 eV, with absorption dominating above 10 eV.


I. INTRODUCTION
High power devices operating under large voltages are used in numerous applications including vacuum electronics [1][2][3], particle acceleration [4,5], inertial confinement fusion [6], and microwave generation [2].The strong electric fields that are typically encountered can have two potential outcomes: (a) the creation of strong Maxwell stresses at the surface [7], which can lead to plastic deformations, defect creation, or probable mass ejection [8].Such Maxwell stresses arising at metal tips under large electric fields are also routinely encountered in field-ion microscopy [9,10].Deformation or crack formation due to Maxwell stressing of a material can allow for the enhanced outflow of any gases that might have been trapped or dissolved in accordance with the solid solubility and Henry's law, within the electrode material.(b) Strong field emission of electrons, which contributes to localized heating.Such resistive Joule heating has been shown to occur within ultrashort (∼ns) time scales [8,11].It might even be conjectured that localized heating can lead to an instability at the metal tip, similar to the Rayleigh-plateau [12,13] process for fluids or solid columns [14].Furthermore, the geometric confinement inherent in nanoprotrusions would shorten phonon mean free paths, reduce heat outflow, and hasten the temperature rise to produce strong thermal effects [15].
In the context of diffusion of gases trapped inside a metal, Neuber et al. [16] had established the role of outgassing through time-resolved spectroscopy that revealed the emergence of hydrogen.Detection of hydrogen is the easiest and fastest, given the high mobility associated with the small mass.Hydrogen outgassing has been reported from several metals [17], and known to adversely affect the material properties, degrade mechanical strength, and contribute to embrittlement [18][19][20][21].Outgassing can either be the result of desorption of previously adsorbed molecules near the surface, or arise from bulk diffusion with vaporization following permeation.Transport is initially governed by the rate at which hydrogen atoms diffuse to the surface, where the concentration is depleted compared to the bulk.Dynamics in the regime are described by Fick's diffusion equation with an Arrhenius dependence on energy and temperature.With reductions in the sample hydrogen concentration over time, less hydrogen arrives at the surface per unit time.Eventually, the rate at which hydrogen atoms at the surface combine and desorb as H 2 molecules (referred to as recombination), becomes the dominant factor in the outgassing [22][23][24].Here, we focus more on the short time scales where the diffusion-limited process dominates.
From an electrical standpoint, outgassing can affect operating conditions by leading to plasma formation and degrading the efficiencies of high power devices [25].Information of such outgassing could, for example, be supplied to particle-incell (PIC) simulators for a more comprehensive and realistic analysis of high power microwave (HPM) device operation.Outgassing though, which has received some attention, represents only part of the overall picture.For a complete analysis, it is also important to assess the possible attachment or deposition which can either lead to gaseous entry back into the metal from the near-vacuum region of an HPM device, or physiosorption at the surface.Both are important processes that need to be considered, as they dictate the overall dynamical equilibrium and determine the net concentration of gas within the system.A simple schematic of these processes is shown in Fig. 1.A gas atom incident onto the metal surface at an angle "θ " and energy "E ", will have a probability S(E , θ ) of "sticking" to the surface, or else it would be reflected back.The sticking probability can be expected to depend on the incident energy and angle.Since this aspect has not been studied to the best of our knowledge, it merits quantitative analysis.If the atom is not reflected, then depending on its energy and "angle-of-attack", it could either be physiosorbed or make its way deeper into the metal.This would constitute a reverse mechanism to replenish the gas within the electrodes.Inclusion of such gas trapping or inflow through the surface provides a more complete and comprehensive picture.The associated rates can be used for dynamic analysis that includes both injection for gas build-up in the near-vacuum system and loss terms through "sticking in", via appropriate rate equation analysis.
Calculations for the sticking coefficient have often been based on density functional theory (DFT) [26][27][28][29] that entail calculations of potential energy surfaces for the interactions between molecular hydrogen and copper.Use of the DFT does present a problem, in that many-body effects are not adequately taken into account.Also, since DFT which primarily relies on the use of periodic boundary conditions, simulations of defects in materials require special attention.For example, with defects included, a larger unit cell would be needed to avoid any inadvertent or artificial enhancements in the defect density [30].Also, accurate development of the many-atom potential energy surfaces remains a challenging and complex ongoing problem [31].Two of the commonly used techniques, the Born-Oppenheimer static surface model [32,33] and ab initio molecular dynamics simulations [34] have not yielded very good fits to experimental data [35].Multiscale schemes that combine DFT for the binding energies, grand canonical Monte Carlo (GCMC) for the surface coverage, and molecular dynamics for desorption and sticking coefficients have recently been reported by Lane et al. [36].In that study, desorption of water vapor was carried out based on the LAMMPS software tool.Here, this same software tool is used in the present analysis based on molecular dynamics (MD) simulations.The MD approach has been used previously for evaluations of sticking and reflection of low energy deuterium incident on single crystal tungsten [37].The advantage of the MD scheme is the inclusion of the many-body physics, along with internal defects, or stepped-surface features.
In this contribution, we present simulation results for hydrogen transport in the copper electrode system, based on the MD technique.Gas out-diffusion is probed as a function of lattice temperature.The goal of the present calculations is twofold.First, at the simplest level, the intent is to understand the physics of outgassing, the influence that material defects such as vacancies may have on gas diffusion, the dependency on temperature, and the potential effects that surface imperfections such as local pits or steps might have on transport.Additionally, carrying out calculation at different temperatures, provides the means to predict outgassing during the operation of high power HPM devices since localized heating is a natural outcome.For example, the bombardment of charged particles at the electrodes (e.g., electrons at the anode) would lead to energy loss from inelastic collisions, resulting in large temperature increases.The requisite power densities could be computed, for example, through Monte Carlo simulations that incorporate the energy loss function (ELF) and mean free paths as demonstrated by our group [38].The only piece necessary for such outgassing analysis during device operation with ongoing internal heating, would be the space-and time-dependent thermal calculations from the ELF to obtain the internal temperature profiles.
The technique used here is general and copper has been chosen as a representative metal often used in high voltage applications.In addition to simulations of outgassing from the surface, separate simulations have also been carried out to gauge the probability of re-trapping at the surface or deeper entry into the electrode.The related sticking coefficient is evaluated as a function of both the incident hydrogen energy and the strike angle relative to the metal surface.A curve fit from the results is also derived for practical convenience.It must be emphasized that this represents a modest contribution that could be enhanced to include different crystal orientation, defects, surface oxide layers, and possible Maxwell stress driven explosions of subsurface voids.Other aspects of relevance could be dynamic temperature changes due to charge bombardment of electrode surfaces during device operation, or structural changes due to energy deposition from incident charged particles and local time-dependent heat generation.

II. MODEL DETAILS
Two sets of calculations were performed and are reported here.First, MD simulations were performed to evaluate the diffusion coefficient of hydrogen in copper (ideal as well as in material containing vacancies) as a function of temperature.In addition, the diffusion in the presence of internal defects and vacancies were also carried out.Next, the sticking coefficient for hydrogen incident on a copper surface was evaluated as a function of kinetic energy and incident angle.The model details pertaining to these two sets of calculations are given next.

A. Transport within copper
Diffusion of hydrogen atoms in a copper lattice was simulated using the classical MD scheme.Such simulations provide the inclusion of the many-body effects with multiple collisions, and can be justified where the temperatures or energies are not very low [39], as in the present situation.Though details have been given elsewhere [40], a brief outline is presented here for clarity.In this technique, forces on any atom at a lattice site due to its many neighbors, are evaluated based on relevant interatomic potentials.It thus allows "numerical experiments" with atomic-level resolution.Limitations include the inability to include electrons, treatment based on potential energies (i.e., force fields) that are developed through experimental data fits, the neglect of quantum mechanics calculations, errors due to the finite number of entities included in actual simulations, and inability to simulate large systems or long times due to the computational complexity.Nonetheless, the MD does provide useful and fairly realistic physics-based model analyses at the atomistic level.
Here in the present MD calculations, the interactions between Cu-Cu, H-H and Cu-H atoms were treated on the basis of the bond-order potential (BOP).Traditionally, atomic interactions in metals were based on potentials derived from the embedded-atom method (EAM) [41].Numerous Cu-H potentials had been developed based on the EAM formulation [42][43][44], which generally yields good results for face-centered-cubic (fcc)/transition metals combinations.In theory, compared to pair potentials, EAM provides a better description on metallic bonds and includes both pair-potential energy term and the embedding energy from the surrounding atoms.However, none of the EAM potentials could capture the H 2 gas phase, or appropriately model Cu-H clusters [45].A more refined potential form based on the so-called modified embedded-atom method (MEAM) was subsequently created by EAM to consider the directionality of bonding [46].Here, however, the BOP was chosen for analysis of the metallic copper system as it is superior to the MEAM.The BOPs are a class of empirical, analytical interatomic potentials that are often used, and have the advantage of describing several different bonding states of an atom, with the same parameters.
Early versions of the BOP potential [47] were improved to include distinct bond orders for σ and π bonds [48].Other developments include the adaptive intermolecular reactive empirical bond-order potential [49] obtained from modifications to the Tersoff potential for covalent bonding.Analytic formulations have been developed by Pettifore et al. [50] for Cu, H and Cu-H configurations.Finally, the Ziegler-Biersack-Littmark repulsive potential was not used, as it applies to screened ion-metal interactions, or to describe high-energy collisions between atoms [51].
The MD technique involves solving the equations of motion for interacting particles numerically, subject to appropriate initial and boundary conditions [52].The open source LAMMPS code [53] with its library of potentials was used in the present calculations.A fcc structure was used for copper with a lattice constant of 3.61 Å.This was thus effectively a rectangular box with dimensions of 50 × 68.5 × 500 Å 3 representing a slab of copper with atoms lying within it at the lattice points.The top surface above the box was taken to be vacuum.Atoms could escape the top surface, while periodic boundary conditions were applied at the four lateral faces.Specular reflection was applied at the bottom surface.For the outgassing simulations, a finite number of hydrogen atoms were initially placed at interstitial sites of the copper structure.
In all simulations for outgassing in ideal copper, the starting system configuration was initially relaxed by running the LAMMPS simulations up to 10 ns to achieve thermodynamic equilibrium and energy minimization.A 1 fs time step used [8].During thermalization, the NPT (isobaric-isothermal) ensemble conditions were used with 104 250 copper atoms.Following relaxation, the operating temperature was maintained at a constant value for each of the simulations, and the NVT (canonical) ensemble conditions used.There was no scalar externally applied pressure in these simulations.The Langevin algorithm [54] was used to control pressure, with a damping parameter (P damp ) taken to be 0.5 ps.Finally, the Nose-Hoover thermostat [55,56] was applied for temperature control.
For simulations of diffusion in copper containing vacancies, the possibility of hydrogen trapping needs to be taken into account.Hydrogen in metals tends to gather at vacancies, and hence has reduced diffusive motion [57][58][59].From an energy standpoint, a minima exists at a vacancy site in a copper lattice.This was confirmed through (DFT calculations as discussed in the Appendix.

B. MD for sticking coefficients
For simulations of the sticking coefficient, simulations were carried out to gauge the reflection or absorption of hydrogen atoms from the surface of a copper slab.The reflection coefficient was simply obtained from the fraction of incident particles that re-emerged back after interacting with the copper after a finite time.The sticking coefficient was the remaining fraction and consisted of both those that were trapped at or remained close by the surface, plus those that penetrated deeper into the copper material.Figure 2 shows a representative figure of the geometry, with an atom located in the space above an ordered copper lattice.In the MD simulations carried out, each solitary atom would be incident on to the copper (seen in the bottom portion of the figure) starting from a chosen angle and initial energy.The dimensions of the copper box were taken to be 50, 85, and 100 Å along the x, y, and z directions, respectively.Thus, the incident hydrogen atoms could penetrate and go inside the copper box, or emerge back out after having gone inside the box and being scattered by the copper atoms, or get reflected without any penetration.At each energy five to seven separate simulation sets were carried out to obtain better statistics.For the simulation cases, ten hydrogen atoms were projected onto the copper.Thus, the starting atom would see a pure copper slab.The atom would either go into the copper, or hover close to the surface, or be reflected away.In our simulations, the position of any incident atom that either ended up inside the copper or remained within two lattice constants of the surface, was recorded and used during the dynamics of subsequent incident atoms.If, however, an atom happened to be further than two lattice constants from the FIG. 2. Representative simulation geometry showing a hydrogen atom located in the space above an ordered copper lattice.In the MD simulations, the solitary hydrogen atom would be incident onto the copper, seen in the bottom portion of the figure starting, from a chosen angle and initial energy.
surface, it was excluded from consideration.In this manner, the possible presence and effects of all preceding atoms were taken into account for the dynamical trajectory of subsequent atoms.In this procedure, the number of reflected atoms were counted, and the mean (as well as standard deviation) over all simulations for a given energy were obtained.
The diffusion coefficient (D) of hydrogen gas trapped in the copper lattice was obtained by tracking the usual mean square displacement |r(t )| 2 of all the diffusing hydrogen atoms over a given period of time at each temperature.The mean square displacement (MSD) is related to time-dependent variations in position for an N-atom system as [60] where r i (t ) is the position vector of the i th atom at time t.The MSD connects the diffusion coefficient through the relation [60]: Over longer times the diffusion coefficient saturates to a constant value at the operating temperature.

III. RESULTS AND DISCUSSION
The results obtained from the implementation of the model described in Sec.II are presented and discussed next.

A. Transport within copper
As a validity check of our MD implementation, results of the temperature-dependent diffusion coefficients for hydrogen in pure copper were obtained, and transport in the presence of a discrete set of vacancies was also probed.Our calculations, FIG. 3. Comparison between the temperature dependent diffusion coefficients for hydrogen obtained from the simulations, and various data points reported in the literature.
given in the Appendix, confirmed that a single H atom would be absorbed or trapped at a vacancy.
Results for temperature-dependent diffusion coefficients for hydrogen in pure copper and in the presence of a discrete set of vacancies, are shown in Fig. 3.The diffusion coefficient was obtained from the slope of the mean square displacement versus time.Though some of the data shown in Fig. 3 for the ideal lattice had been reported previously [40], the temperature scale was extended slightly to higher values.For comparison, data reported previous by various groups [61][62][63] are also shown in Fig. 3 for ideal copper.A fairly good agreement between data and the calculations is evident.An Arrhenius-type behavior naturally emerged from the calculations for pure copper, with an activation energy of 0.43 eV and a pre-exponential factor of 2.1 × 10 -7 m 2 /s, which are comparable to other reports [21].Since the melting point of copper is 1358 K, simulations for outgassing in Fig. 3 were carried out up to 1350 K, which is close to this value.
For simulations in copper containing vacancies, the motion of the hydrogen was frozen if its distance from any vacancy got to be less than the separation between the energy barrier and the minima as shown in Fig. 12.In the case of vacancies, the copper system was chosen to have a 20 Å × 20 Å × 40 Å volume containing 8 vacancies at random sites and 20 hydrogen atoms.The behavior in Fig. 3 with the presence of vacancy sites shows lower diffusion as expected and reported elsewhere [57][58][59].It may be pointed out that the reductions on diffusion would depend on the density of vacancies, and this aspect will be probed elsewhere.Due to the presence of a barrier, the difference in hydrogen diffusion coefficient in copper with and without vacancies is not large at the lowest temperature.This is the result of fewer atoms having sufficient energy to transcend the energy barrier and get trapped.From the standpoint of mitigating out-gassing, the presence of discrete vacancies would then be a desirable objective.
Based on the outgassing dynamics, the rate constants for the outflow were computed for hydrogen gas.For these calculations, 500 hydrogen atoms were randomly placed inside the simulation box containing copper, and the system was thermalized at different temperatures.Figure 4 in this regard, represents the simulation system with hydrogen atoms randomly placed in a copper lattice.To obtain outflow rates, the number of atoms remaining inside the simulation box were tracked as a function of time.Periodic boundary conditions were imposed on the four lateral faces, with reflection at the bottom surface.Atoms reaching the top face were allowed to escape and were removed from the system once their distance from the surface exceeded twice the lattice constant.This allowed for time-dependent calculation of the population N (t ) remaining inside the main simulation region.Results obtained from MD simulations are shown in Fig. 5 at different lattice temperatures.The trends seem to qualitatively follow the experimental report by Mamun et al. [64].However, the data by Mamun et al. was for stainless steel and not copper, with a much lower temperature range.An Arrhenius-type behavior though was seen to hold.
Approximating the time dependence of the number of atoms [=N (t )] inside the copper by a first-order process such as N (t ) = N 0 exp(−t/τ ), with τ a suitable time constant, the value of this parameter could be obtained.The time constant extracted from the simulation data of Fig. 5, is presented in the plot of 1/τ versus temperature (T ) thus obtained is presented in Fig. 6.To ensure statistical robustness, at least five simulations were carried out at each temperature.The associated results at a 95% confidence interval have been shown.A curve fit then yielded the following temperature dependence of the outgassing rate: Next, a set of simulations were performed to probe the outgassing in the presence of a surface depression or void.It attempts to roughly mimic the occurrence of a dent on the surface during device operation due to possible material ejection or surface damage.Figures 7(a) and 7(b) presents the side view and the perspective, respectively, of the geometry simulated.Figure 8 shows our MD results obtained for outgassing for the two geometries, one being ideal copper and the other containing a surface void (Fig. 7), as a function of time.A temperature of 300 K was used.Slightly slower outgassing is predicted for the case with the crater/void.The difference arises from the broken lateral symmetry which changes the net molecular forces near the boundary with the void.Hydrogen atoms in pure copper, whose trajectories puts them directly on a path towards the void, can escape out of the metal via the metal-void interface.Similarly, those hydrogen atoms that are far removed from any void and are on a net upward path, would also be able to escape and emerge out of the metal.However, the hydrogen atoms in the metal that pass close to the void, face asymmetric molecular forces, which helps steer them away from the void.
Such asymmetry in forces is also responsible for the surface tension and wettability of droplets on surfaces.The net outcome from the collective forces is an effective time delay in emergence of hydrogen gas atoms from regions near the presence of such surface craters/voids.This delay on outgassing has practical implications, though other influences such as the sticking dynamics on patterned or stepped surfaces would also need to be considered, but are beyond the present scope.

B. Sticking coefficients
Next, simulations for the reflection and absorption of hydrogen incident from the vacuum onto the copper surface were performed at different energies and angles.At each energy 5-7 separate simulation sets were carried out to obtain better statistics.For each simulation case, ten entities were projected onto the copper at definite time intervals.Also, to avoid any lattice distortion, the atomic system was continually thermalized after each run to eliminate any memory effects.The number of reflections were counted, and the mean (as well as standard deviation) over all simulations for given energy were obtained.Results for normal incidence are shown in Fig. 9. Temperature of the copper was taken to be 800 K.The main feature is that at low incident energy, the ability to bounce back easily is low.This leads to low reflection coefficients.At high energies, the atoms are mainly transmitted through, again leading to a lower reflection.Hence, a local maxima in the reflection curve is predicted in Fig. 9.For practical convenience, a curve fit to the data was obtained, and is FIG.9. Energy-dependent reflection coefficient at normal incidence of hydrogen impinging on the surface of a copper electrode, calculated from molecular dynamics simulations.
given as Results at an incident angle of 30°for incoming hydrogen are given in Fig. 10.The overall characteristics are similar to those for the normal incidence with a local maxima.A slightly higher reflection at 30°is predicted compared to the 0°incidence, and is associated with the shallower impact.The simulations yielded adsorption and absorption data which is presented in Fig. 10 for normal incidence.As evident from the figure, adsorption is predicted to occur for incident energies below 10 eV, while absorption would dominate above 10 eV.The overall plot of Fig. 11 exhibits a minima in energy which is in keeping with trends shown by Kaufmann et al. [35] FIG. 10.MD results for the energy-dependent reflection coefficient incident at a 30-degree angle on the surface of a copper electrode.for copper surfaces, and by Maya for deuterium incident on tungsten [37].Results were also obtained for a 60°incident angle, and MD simulations were also carried out at a higher temperature of 1200 K.However, the results were very similar to those discussed in Figs. 9 through 11, and so have not been shown for brevity.Finally, it is likely that the presence of grain boundaries, dislocations or surface adsorbates could change the energy-dependent sticking coefficients, but it is beyond the present scope.
The above calculations and results could form the basis for predictive analysis of gas emissions for comprehensively modeling the electrical response and overall plasma behavior in pulsed power devices.The net temperature-dependent gas emission rates could be used to account for pressure changes, dynamically alter scattering rates in PIC simulations, and probe contributions to plasma growth and field distortions.Currently for example, PIC simulations assume some constant level of background gas, while ignoring gas release rates, or the possible local inhomogeneities in the concentration due to outgassing.However, evaluations of the system efficiency and its space charge could change if more realistic and physically based input to the PIC formulations were made available.The MD simulations reported here are a step in the direction and could be extended to include the role of adsorbates, oxides, Maxwell stress induced blowout of subsurface voids, or sputtering induced defects.A practical validity test for such outgassing predictions might be through comparisons to measured increases of pressure during high power device operation in test chambers.
An aspect not addressed here, but perhaps important in the context of charged particle incidence upon electrodes in pulsed power systems, could be the time dependent changes in temperature and associated spatial gradients.The latter could add a thermodiffusive contribution to gas transport, and detailed analysis would require a three-step process of (i) Initial calculations of the deposited energy from the incident particles, possibly through Monte Carlo calculations with energy-dependent stopping power and mean free paths [38,65], which have been reported elsewhere.(ii) Solutions of heat transport to obtain the evolutionary spatially variable temperature profiles, and (iii) Use of the spatially dependent temperatures thus obtained as inputs in the MD technique to naturally include thermal gradients.The temperature effects could also arise in the context of controlled experiments employing short pulse duration lasers which would give rise to local heating and puffs of emergent gas from the irradiated portions, in addition to features such as possible surface annealing and reconstruction [66].In this context, it may be mentioned that MD simulations at high temperatures, for example, a ramp from 300 to 1300 K to mimic experimental heating conditions in recent pulsed-power experiments have been reported already [36].
Finally, for completeness, it may be mentioned that dissociation of water vapor and the resulting increases of partial pressures in vacuum systems is another practical consideration that arises in pulsed power and high power microwave devices.Accurate simulations of water desorption from metal electrodes would require multiscale approaches that combine, for example, DFT for the binding energies, GCMC for the surface water coverage, and MD for desorption [36].The desorption of water vapor in the above-mentioned study had been carried out based on the same LAMMPS software tool used in the present analysis.So in theory, it would be possible to probe such effects for outgassing using multiscale approaches.However, such calculations are beyond the present scope.

IV. SUMMARY AND CONCLUSIONS
Outgassing is typically the first stage of possible plasma formation in a vacuum, and can lead to various effects such as breakdown, surface flashover, anode-cathode gap closure [67], and pulse shortening [68] in high power systems.It is therefore important to understand and quantify the outgassing process.Here, preliminary simulation results of hydrogen diffusion and outgassing in a representative metal (such as copper) were obtained based on the MD technique.Such calculations will aid in the design assessment of cathodes for high power microwave devices.The calculated temperature dependent diffusion coefficient values were seen to agree with experimental reports in the literature over a range spanning from 300 to 1350 K.An analytic curve fit for the temperaturedependent outgassing time constant was obtained.It has also been shown that outgassing in the presence of surface voids or craters would likely be slower.In addition, calculations for the sticking coefficients of hydrogen gas atoms as a function of incident energy and angle were carried out for a copper electrode.The results showed low reflection coefficients at both high and low energies, with a maxima at around 6.5 eV.A curve fit to the data was obtained and predicted to roughly hold for a range of incident angles.Finally, adsorption and absorption at normal incidence predicted adsorption to occur for incident energies below 10 eV, while absorption is predicted to dominate above 10 eV.The minima in energy is in keeping with trends shown by Kaufmann et al. [35] for copper surfaces.

APPENDIX: DENSITY FUNCTIONAL THEORY BASED ENERGY CALCULATIONS FOR MONOVACANCY
The Vienna ab initio simulation package, a densityfunctional pseudopotential software tool [69][70][71] that relies on the projector augmented-wave method [72,73] was used.The calculations were performed within the generalized gradient approximation to the exchange-correlation potential as parameterized by Perdew-Burke-Ernzerhof [74].The plane wave kinetic energy cutoff was set to 460 eV and a 3 × 3 × 3 Monkhorst-Pack [75] mesh is used to sample the Brillouin zone.
The crystal was represented by a periodic supercell containing 64 host copper atoms obtained by a 4 × 4 × 4 multiplication of a fcc unit cell.A lattice constant of 3.63 Å was used for all structures and was calculated from the fcc Cu bulk structure.The lattice constant was determined by fitting the Murnaghan equation of state with Cu bulk total energy and volume data [76].The monovacancy was created by removing a Cu atom from the center of a supercell.For all calculations the internal atomic positions were fully relaxed while the external parameters (volume and shape of supercell) were fixed.
In order to verify the stability of defect and have a better insight into hydrogen atom diffusion into the bulk, we calculated the energy barrier between the most stable configuration and the next local minima (0.135 eV higher in energy).For the global and local minima energy configuration, the H atom was positioned at 1.21 Å and 3.08 Å from the center of vacancy along the 111 direction, respectively.We used the nudged elastic band (NEB) method for finding the minimum energy path between the two configurations.The result of NEB calculations is shown in Fig. 12.According to the NEB results, H atoms must be displaced along the 111 direction for the minimum energy pathway.Two energy minimas are seen to be separated by a barrier height of about 0.38 eV, corresponding to H atom positioned at 2.26 Å from the center of vacancy.To our knowledge there is no experimental data for the diffusion energy barrier of the system studied here.However, the binding energy of hydrogen to a vacancy has been measured accurately by Besenbacher et al. [77].The binding energy (E b ) between hydrogen atom and vacancy can be calculated as where E (Cu 64 ), E (Cu 64 H ), E (Cu 63 V ), and E (Cu 63 V H ), are the total energies of defect free fcc super cell with 64 Cu atoms, with an interstitial hydrogen atom, one mono-vacancy, and one interstitial hydrogen atom with one mono-vacancy, respectively.Our calculated binding energy (not shown) between H atom and the vacancy is 0.43 eV, which is in excellent agreement with estimated 0.42 eV [77].The result of our binding energy calculation could provide validation of energies predicted by the NEB method.
Hydrogen has slightly higher electronegativity (Pauling value ∼2.2) compared to copper (Pauling value of ∼1.9).Hence, a hydrogen atom in a copper lattice would have a slightly negative charge, while the copper in the neighborhood would be somewhat positive.In a perfect lattice, a hydrogen atom at any site, would on average experience symmetric forces from the surrounding copper atoms.However, if the hydrogen were to approach a vacancy, there would be no positive charge around the vacancy, resulting in an asymmetric net force that would pull it away or impede its motion towards the vacancy.This is reflected in the "energy barrier" in Fig. 12 as a hydrogen atom approaches a vacancy.Once in the energy well with lower energy, the hydrogen would effectively be "trapped".The process can alternatively be thought to arise from reductions in charge density in the vicinity of a vacancy that provide an isosurface for collective hydrogen binding and segregation [78].Apart from calculations of the DFT, ab initio GCMC simulations [79] have also confirmed the vacancybased retardation of hydrogen diffusion in copper containing vacancies.

FIG. 1 .
FIG. 1. Schematic showing the possibility of gas atoms either re-enter the metal surface from the near-vacuum region of a HPM device, or undergo physiosorption.

FIG. 4 .
FIG. 4. LAMMPS simulation system for the outgassing.The copper atoms are shown in red, while the hydrogen atoms (shown in blue) are randomly placed.Hydrogen size not to scale for clarity.

FIG. 5 .
FIG. 5. Simulated result showing number of hydrogen atoms in the copper electrode versus time.

FIG. 7 .
FIG. 7. Simulation setup for copper with a void at the top surface.(a) Side view, and (b) perspective.

FIG. 11 .
FIG. 11.Molecular dynamics calculations of energy-dependent adsorption and absorption of hydrogen at normal incidence on the surface of a copper electrode.