Water formation at low temperatures by surface O 2 hydrogenation III: Monte Carlo simulation

Water is the most abundant molecule found in interstellar icy mantles. In space it is thought to be eﬃciently formed on the surfaces of dust grains through successive hydrogenation of O, O 2 and O 3 . The underlying physico-chemical mechanisms have been studied experimentally in the past decade and in this paper we extend this work theoretically, using Continuous-Time Random-Walk Monte Carlo simulations to disentangle the diﬀerent processes at play during hydrogenation of molecular oxygen. CTRW-MC oﬀers a kinetic approach to compare simulated surface abundances of diﬀerent species to the experimental values. For this purpose, the results of four key experiments—sequential hydrogenation as well as co-deposition experiments at 15 and 25 K—are selected that serve as a reference throughout the modeling stage. The aim is to reproduce all four experiments with a single set of parameters. Input for the simulations consists of binding energies as well as reaction barriers (activation energies). In order to understand the influence of the parameters separately, we vary a single process rate at a time. Our main findings are: (i) The key reactions for the hydrogenation route starting from O 2 are H + O 2 , H + HO 2 , OH + OH, H + H 2 O 2 , H + OH. (ii) The relatively high experimental abundance of H 2 O 2 is due to its slow destruction. (iii) The large consumption of O 2 at a temperature of 25 K is due to a high hydrogen diﬀusion rate. (iv) The diﬀusion of radicals plays an important role in the full reaction network. The resulting set of ‘best fit’ parameters is presented and discussed for use in future astrochemical modeling.


Introduction
Water is an important species in molecular astrophysics and a prerequisite for life on Earth; it controls much of the gas-grain chemical interplay in space and is vital to the formation of more complex molecules as star-formation progresses.For this reason, understanding its fundamental properties is vital to disseminating our knowledge of chemical evolution through star and planet formation, and ultimately of the origin of life itself.Yet it has been a perennial question how water forms under the harsh conditions that govern chemistry in the interstellar medium and which physical and chemical processes are most important.
In dense cold regions of the interstellar medium gas-phase formation routes for H 2 O and subsequent freeze-out mechanisms cannot explain the large ice abundances observed.Therefore, it is expected that water is formed on surfaces of cold icy dust grains that act as catalytic sites for molecule formation.Indeed water has been identified as the main component of interstellar ices. 1 The first reaction scheme for a grain surface formation route of H 2 O was theoretically proposed in 1982 by Tielens and Hagen. 2 This scheme focussed on the hydrogenation of atomic oxygen, molecular oxygen and ozone.At that time, only gas-phase data were available to estimate surface reaction rates.In 2007 Cuppen and  Herbst showed in a theoretical study that the solid-state reaction H + OH -H 2 O likely dominates water formation under diffuse cloud conditions, whereas H 2 + OH -H 2 O + H is prominent in dark molecular clouds. 3][6][7][8][9][10][11] These experiments provide a more detailed understanding of the reaction network involved in water formation, extending the originally proposed network 2 to that depicted in Fig. 1.In particular, water formation at low temperatures by surface O 2 hydrogenation has been studied in depth.Ioppolo  Different experiments have been performed and constraints on several reaction rates could be determined.
The WISH key programme of the space telescope Herschel has recently made observations of gas phase water and species related to its chemistry in prestellar cores and young stellar objects at different evolutionary stages.Hot gas containing O(I), 12 O 2 , 13 cold H 2 O, 14 HO 2 15 and H 2 O 2 16 as well as OH, OH + , H 2 O + and H 3 O + have been observed in this context. 17The identification of these molecules is fully consistent with the solid state network (Fig. 1) has been derived experimentally.These experimental results are extended here to Continuous-Time Random-Walk Monte Carlo (CTRW-MC) simulations to disentangle specific reaction mechanisms and to derive more accurate reaction barriers by comparing laboratory surface abundances with those obtained by simulations.These barriers can then be used as an input for astrochemical models in order to meet observational constraints.We show here that through systematical variation of the many input parameters it is possible to extract information concerning key reactions within the full network.
The manuscript is organized in the following way.Section 2 summarizes the key experiments that have been described previously and that are used here as reference for the simulations.Section 3 provides the details of the applied CTRW-MC method and in Section 4 the results are presented and extensively discussed.We conclude with recommendations for future studies and astrochemical considerations in Section 5 as well as with a short summary in Section 6.

Experimental observations
This section briefly summarizes the previous experimental work that is used here as a starting point.For a detailed description of the used setup, procedures and results we refer to ref. 8, 9 and 18.Two types of experiments have been performed: sequential hydrogenation of O 2 ice and co-deposition experiments of O 2 molecules and H atoms.
During sequential hydrogenation experiments, an O 2 ice, several monolayers (ML) thick, is first prepared at 15 K and subsequently exposed to hydrogen atoms at various constant ice temperatures.This allows the study of final and stable products of O 2 hydrogenation, i.e., H 2 O 2 and H 2 O (Fig. 1).Quantitative information concerning the surface abundances is obtained by dividing the integrated absorption of a selected infrared band (cm À1 ) over the so-called band strength (ML cm À1 ) as described in ref. 8.In Part I we showed that during sequential hydrogenation experiments the initial formation rates of H 2 O 2 and H 2 O are temperature and thickness independent.The final yield, however, does depend on these parameters.Furthermore, due to the competition between reaction of atomic hydrogen with solid oxygen and hydrogen diffusion into the ice, the penetration depth of H atoms was found to span up to 16 ML at 25 K. 8 We will elaborate on this topic in Section 4. 3.
During co-deposition experiments, O 2 molecules and H atoms are released onto the cold substrate at the same time and therefore are adsorbed simultaneously.Different stages of the hydrogenation route, i.e., various reactive intermediates, become experimentally accessible by changing the stoichiometric ratios of O 2 and H 9 (Fig. 1).Due to lacking band strengths for the unstable matrix-isolated reactive intermediates, a quantitative study is not trivial.However, quantification of a single species using the ratio between abundances at different temperatures and/or H/O 2 ratio can be easily performed. 9our of the aforementioned experiments are selected for further comparison to simulated results: sequential hydrogenation and co-deposition experiments both at 15 and 25 K, as listed in Table 1.In the following, 15 and 25 K are defined as 'low' and 'high' temperature.These four experiments are considered to be representative of the different experimentally observed features.In Fig. 2 the measured surface abundances and integrated absorbances are depicted for the selected sequential hydrogenation and co-deposition experiments.In the following, we consistently separate the discussion between these two types.Important features that need to be reproduced by the MC simulations are summarized below.
Sequential hydrogenation: evolution of the H 2 O 2 abundance in time with a sharp transition between the initial (T-independent) linear and final steady state (T-dependent) regime, increase of the final H 2 O 2 production by a factor of B7 between 15 and 25 K, and similar behaviour of the H 2 O abundance for both temperatures.It should be noted that the uncertainty of the band strengths is estimated to be within 50% and this affects the accuracy with which the experimental surface abundances can be determined.

The Monte Carlo method
This section describes the general Continuous-Time Random-Walk Monte Carlo method used for the simulations, discussing sequentially the simulation of O 2 deposition, the hydrogenation of O 2 ice and the different parameters used for a co-deposition simulation.Subsequently, distinct details of the program are addressed.For a more detailed overview the reader is referred to ref. 19.The present results are obtained with the program described in ref. 3, which has been extended to account for the specific characteristics needed here.The reproducibility of the simulations is monitored by performing each of the standard simulations (see below) three times using different seeds.Standard deviations have typical values of o5% for abundant species (e.g., O 2 , H 2 O 2 ), o10% for less abundant species (e.g., OH) and o40% for O 3 .

Deposition of an O 2 surface
During the deposition phase the O 2 ice is formed starting from a bare surface; the surface used to mimic the experimental gold substrate has a smooth topology and is 2.5% of that depicted in Fig. 1(c) from ref. 20.The simulation then starts by the addition of an O 2 molecule to the surface.This occurs at time where r is the surface site density (1 Â 10 15 sites cm À2 ), F O 2 is the flux of O 2 molecules, X is a random number between 0 and 1, and t current is the current time which is set to 0 at the beginning of the simulation.During the following Monte Carlo cycles a competition between hopping and desorption of the molecules on the surface and deposition of new molecules determines the sequence of events.This sequence is determined by the time at which each event occurs and is given by where X 0 is again a random number, and R hop and R des are, respectively, the rates for hopping and desorption to occur.Both rates are assumed to be thermally activated according to with E a,Y the activation energy (or barrier) for process Y in Kelvin and n the pre-exponential factor, which is approximated by the standard value for physisorbed species, kT h $ 10 12 s À1 .This factor can be seen as a trial frequency for attempting a new event.
The total binding energy, E tot,bind , for each molecule to a site is calculated by additive contributions of its neighbours.The single binding energy value, E, is then obtained by dividing the desorption energy by a factor of 3.8. 21Note that an order of magnitude difference is applied in binding energy to atoms and small molecules (H, H 2 and O) and H 2 O like molecules (all the rest).For example, the binding energy of a hydrogen atom binding to a flat surface of O 2 molecules corresponds to the number of neighbours times the contribution: 1Á2Á66 + 4Á2Á66/8 = 198 K. Binding to a surface of H 2 molecules yields 19.8 K.
The barrier for hopping from site i to j is given by where x is an adjustable parameter to change the diffusion, E the single binding energy value, d is the distance between the sites and DE bind (i,j) is the difference in total binding energy between the two sites.This equation is derived from eqn (10) in ref. 22.The first term represents the diffusion barrier for equal total binding energies of the sites, while the second ensures microscopic reversibility.Only hopping events between nearest neighbours and next-nearest neighbours are considered.
The barrier for desorption is determined only by E tot,bind .Note that the single binding energy value of OH is much smaller than that of H 2 O.This will be briefly discussed in Section 4.2.3,where the low binding energy is associated with a high diffusion rate.
A total of 10 ML of O 2 is eventually simulated to deposit on top of the bare surface.This results in a smooth ice surface according to the definition from ref. 20.In order to compare the results of the different sequential hydrogenation simulations, the same oxygen surface is used for all simulations, except for those that investigate the influence of surface roughness and thickness.

Sequential hydrogenation
After deposition of O 2 molecules, the surface is exposed to H atoms and H 2 molecules in a similar manner as before, except that now also reactions can occur and new species are formed.The flux of H atoms and H 2 molecules is set to SÁ2.5 Â 10 13 particles cm À2 s À1 for both species, where the sticking coefficient S equals 0.2.The reaction rates are calculated by using reaction activation energies (or reaction barriers) in eqn (3) and are listed in Table 3.The reaction rates represent the competition between different reactions and therefore the ratio between the rates is used for future reference.Temperature independent reactions are those that have a low to zero barrier.The temperature dependence for the remaining reactions arises from a combination of thermally activated  processes and tunneling, i.e., reactions always occur at a higher rate for higher temperatures, not entirely according to eqn (3), but rather scaling the rates by an arbitrary factor to account for tunneling.Note that tunneling is thus not explicitly included.Although many reactions leave the products with a large excess energy, this energy is thought to be efficiently dissipated into the ice surface on picosecond timescales, which we conclude from preliminary Molecular Dynamics simulations.Furthermore, in the laboratory the He cryostat provides a suitable dissipation path, while in the interstellar medium the time scales allow thermalization.Since the amount of energy remaining in the molecule does not seem to correlate with the amount of initial excess energy, we include an arbitrary, small excess energy of only 100 K for each reaction product in our model.This excess energy can be applied to overcome barriers for reaction, desorption and/or diffusion.In this way, a chemical desorption mechanism, like the one proposed in Garrod et al., is implicitly included. 23A more accurate implementation of the excess energy, for instance an energy dependence on the exothermicity of a specific reaction, is a subject of future studies.Diffusion reduces the excess energy by an arbitrary factor of 1.6 for each hop.Furthermore, after 10 À8 s, the local temperature of the excited species is set back to the temperature of the surface.This is based on the assumption that a molecule on the surface will be thermalized after 10 ns.This rather subjective choice of time scale does not affect the outcome of the simulations, since we find that hot species either react immediately to form a new molecule or remain in their initial configuration for times much longer than 10 À8 s.
The rates listed in Table 3 are explicitly for thermalized reactions only.Consider a hydrogen atom at room temperature (300 K) that lands atop the O 2 surface with a temperature of 15 K.We assume that half of the energy ( 12 Â (300 + 15)) is immediately dissipated into the surface, leaving the atom with a local temperature of 157.5 K.This energy is dissipated either slowly through hopping or through an immediate step function after 10 À8 s.All reactions that take place before thermalization have a much higher probability of occurrence, according to eqn (3).This is indeed the case for all reactions with a barrier.The thermalization effect is most pronounced at early times during the simulation, when the surface is not yet filled with molecules, preventing fast reactions.Take for example a surface where H 2 O has already formed; hydrogen atoms landing on a site next to H 2 O are forced to hop in order to meet a reaction partner, thereby losing their energy.The effect is much stronger at 15 K, where the competition between the different processes is small, mainly because the rates are less close to the trial frequency.For instance, for the reaction H + H 2 O 2 -H 2 O + OH, the difference between the hot mechanism and the thermalized reaction is a factor of e ÀE a;reaction e ÀE a;reaction 15 % 10 21 .Reactions that occur easily in the experiment are therefore not necessarily accessible under conditions resembling the darker regions of the interstellar medium.
Finally, in this paper we only consider physisorbed hydrogen atoms, both on the surface as well as trapped in pores.
The solid phase interactions are incorporated indirectly through positioning of the hydrogen atoms on the surface and accounting for the binding energy at a specific lattice site.If a species has more neighbours, the cumulative binding energy increases.Furthermore, the effect of the surroundings on the reaction itself is included in the barrier.In this way, a stabilizing effect through for instance hydrogen bonding of a neighboring species is not included explicitly as is the case when real potentials are used.Rather, a change to the activation energy is applied.In Section 4.2.7 we will briefly talk about the possibility of explicitly incorporating solid state effects involving H atoms.

Co-deposition
During a co-deposition simulation, H, H 2 and O 2 all settle on the bare ('gold') surface simultaneously, after which their possible events are determined and evaluated again using the Arrhenius behaviour.Reactions can thus take place from the start.The individual deposition rate of O 2 , H and H 2 is set to be 2.5 Â 10 13 , 5 Â 10 12 and 5 Â 10 12 particles cm À2 s À1 , respectively, assuming again a 1 : 1 ratio between H and H 2 .The sticking coefficient S is thus already accounted for in the fluxes:

Size and ice morphology
We use a lattice-gas model, where the lattice surface consists of 50 Â 50 sites.This is large enough to overcome the finite size limit, which is set by the rate of atomic hydrogen desorption, as outlined in ref. 19.Larger lattices would result in computationally too expensive simulations.We performed one simulation using standard parameters (see below) on a surface of 100 Â 100 sites as a check and found no difference with respect to the smaller surface.
Although experimental surfaces are probably amorphous, the adsorption sites will be clearly defined and likely distributed with some kind of order in terms of the number of neighbours.A lattice-gas model can thus be seen as a grid of potential wells, and the process results in a change of the occupancy of the lattice sites.The largest advantage of using a lattice-gas model is the reduction of computational costs as result of working with a predefined event table.Since the molecules are confined to a lattice, only a fixed amount of processes can occur and large time scales can be simulated.The disadvantage is that the level of molecular detail is limited and some mechanisms are not included.These types of models have, however, demonstrated how powerful they can be in mimicking ice chemistry, covering hours of simulated time. 24,25o account for the penetration of H atoms observed experimentally, small species (i.e., H, H 2 , OH and O) are allowed to diffuse to interstitial sites in the oxygen ice.For this purpose each monolayer in the ice is represented by two fields in the matrix that holds the ice; one can hold all species present in the simulations, the other contains mainly small species.The diffusion takes place only when an O 2 or a HO 2 molecule is atop the final position, since penetration in water-like structures has not been observed experimentally.Larger species can be present in the intermediate layer, though, as a result of the positioning of reaction products.

Results and discussion
In this section we present and discuss Monte Carlo simulations of the four experiments listed in Table 1.The reaction input parameters for these simulations are given in Table 3.The second column indicates the settings used for the so-called standard simulations (R).These are based on ab initio barriers, experimental constraints or other models and represent the starting point for optimization.Subsequently, for each reaction with a barrier, this activation energy E a is either increased or decreased, after which the effect of the reaction on the final result is studied.The reaction rates are listed in the third and fourth columns of Table 3.Note that while changing a rate, all other parameters are kept constant.Decreasing or increasing the reaction rates in a systematic way forms the core of this paper and is described in Section 4.2.The final goal is to see whether these variations bring us closer to the experimental results, while simultaneously providing insight into the underlying molecular mechanism.Furthermore, we hope to obtain information on the sensitivity of the system on a specific reaction, i.e., information on the error bar of the barrier.Besides studying the chemistry itself, we also investigate the effect of the diffusion barrier of H and O 2 as well as the effect of the interstitial positioning.
We first discuss the standard simulations and how the different input parameters affect the production of various species.The simulations are then compared with the experimental results by using the time evolution of the surface abundances of the molecular species, i.e., comparing the simulated trends to the experimental findings shown in Fig. 2. For the sequential hydrogenation experiments, the observations are not only based on the absolute abundances, but also take into account the ratio between H 2 O and H 2 O 2 as well as trends in the time evolution.For the co-deposition experiments, accurate band strengths are not available.Therefore, the approximation is made that each hydroxyl group contributes equally strongly to the OH-stretching and bending vibrations, which allows for an estimation of the relative band strengths and therefore relative abundances.To strengthen this analysis, we performed simulations for selected sets of parameters for a H/O 2 ratio of 2. In this way, we can compare the molecular abundances between different experiments using infrared features.Furthermore, the OH/H 2 O 2 ratio is studied, following the conclusion from ref. 9 that their abundances are correlated for    different conditions using the standard settings.The colour coding indicates the different molecular species.Each square represents one molecule or atom.The black squares on the bottom represent the initial substrate on which the O 2 ice (red) is deposited.H and H 2 can diffuse into the ice at interstitial positions (intermediate white rows) and form new species there.The sequential O 2 hydrogenation cross sections (Fig. 3) show that the penetration of the H atoms at 15 K is less than at 25 K where some hydrogenation has occurred even at the deepest layers.Indeed a larger penetration depth is required to reproduce the 25 K experimental results.The final ice structure at high temperature is more irregular with respect to that at low temperature.Pores (or empty sites) are formed upon reaction in the original O 2 structure.This, combined with the temperature effect of eqn (3), allows for easier hydrogen diffusion.Moreover, at 15 K the top of the surface has been hydrogenated to water almost entirely whereas at 25 K some H 2 O 2 in the top layer is still available for hydrogenation.
For the co-deposition simulations (Fig. 4), the resulting ice mantle consists mainly of O 2 with other species embedded in the O 2 matrix.The top 5 ML of the matrix have a different composition than the lower layers because the hydrogen atoms landing on the surface can still penetrate some of the ice and induce further reactions.The main difference between the 15 and 25 K ice mantles is in the amount of HO 2 and OH versus H 2 O 2 .At 25 K most of the small species have reacted further as a result of the higher mobility of H and OH.Finally, the total mantle thickness for a higher temperature is smaller as a result of the higher O 2 desorption probability.The thermal desorption value of pure O 2 ice has been experimentally determined as 31 K. 26,27 4.1.2Time evolution of the surface abundance.The solid curves in each of the panels of Fig. 5 and 6 represent the time evolution of a specific molecular species with respect to the initial conditions before hydrogenation for the standard simulations.In Fig. 5 the sequential hydrogenation and in Fig. 6 the co-deposition simulations are shown.The O 2 molecules are consumed and the time evolutions in Fig. 5 (red) are therefore negative.O 2 cannot be spectroscopically observed unless it is abundantly mixed with other species, since it is infrared inactive.All other species are formed in the process.Note the differences in the vertical scales in Fig. 5 between 15 and 25 K as well as between Fig. 5 and 6.
From these figures we can conclude that the trends in the time evolution of water and hydrogen peroxide surface abundances are in good qualitative agreement with the experimental data.Like for the experimental results (Fig. 2), the initial increase in H 2 O 2 abundance is similar at low and high temperature and there is a clear transition between the linear increase and final steady state behaviour for sequential hydrogenation simulations.The slight decrease before reaching the steady state for the H 2 O 2 signal at 15 K is also reproduced by the standard model.The experimental transition is sharper than the simulated one and the level of the steady state at 25 K is much higher (14 ML vs. 4 ML).The H 2 O production at 15 and 25 K follows roughly the same behaviour and the curve shape agrees with the experimental one as well.The H 2 O abundance is slightly higher for lower temperature, and as a result the H 2 O/H 2 O 2 ratio is overestimated.During experiments, the surface abundance of O 3 stays below the detection limit, and this is in agreement with the low value of the simulated abundances.
The molecular time evolution for the co-deposition simulations is plotted as a solid curve in Fig. 6.The final O 2 thickness is of the order of 80 ML, which is beyond the dynamic range of the panels (see also Fig. 4).At later times some O 2 is observed experimentally only at 15 K for H-atom rich conditions.The simulations confirm this trend, since the final O 2 abundance at 15 K is higher than at 25 K.The abundance of most species increases linearly over time except for HO 2 which slowly levels off.These abundances are compared to the linear experimental regime (t o 60 minutes) where H 2 O 2 and OH have not yet reached steady state.H 2 O and O 3 remain below the IR detection limit in the experiments and simulated abundances are correspondingly low.From the comparison of the bending and stretching vibration modes it seems that the simulations overproduce the intermediate species HO 2 for both temperatures.This can partially be explained by the hydrogen diffusion rate incorporated in the standard simulations.A higher mobility would allow hydrogen to scan both a larger surface area and bulk volume, hence allowing for the reaction H + HO 2 to take place more often.The dependence of the molecular surface abundances on both temperature and H/O 2 ratio is summarized in Table 4. Comparison with the experimental data is also shown and for H 2 O 2 a good agreement is found.The arrows show the influence of increasing temperature or H/O 2 ratio on the amount of species formed.Indeed, the HO 2 and OH production decreases at higher temperature, whereas for H 2 O 2 it increases.The OH/H 2 O 2 ratio stays equal throughout a single simulation, in accordance with the experiments.From the above, it appears that the main inconsistencies are caused by radical species.The dynamics of these species is not fully understood, but one should also keep in mind the experimental error bars involved in this analysis, as stated in Section 2.
Arrows indicate the effect on the molecular surface abundances for an increase in temperature (column two and three) or H/O 2 ratio (last two columns).A double arrow represents a large increase (up) or decrease (down), a single arrow a small increase (decrease) and a hyphen no change.

Key reactions
The investigation of the reaction rates is done by varying one reaction barrier at a time with respect to the standard values as indicated in Table 3. Fig. 5 and 6 show the resulting time evolution of the surface abundances.The procedure adopted is as follows: if the difference between a given simulation and the standard simulation is larger than the derived standard deviation, it is considered to be a strong effect.The reactions that exhibit the largest effect upon changing their corresponding barriers are: The implications of these reactions, along with those of reactions with H 2 , are discussed in the sections below.For reasons of brevity, we will comment only on features that distinctly differ from those in the standard simulations.

H + O 2 .
Since the reaction H + O 2 is the start of the hydrogenation pathway, changing its rate affects the full reaction scheme.In the gas phase the barrier of the reaction depends on the incoming angle of the molecule 28 and a change in the rate can thus be interpreted as a change in the fraction of successful approaches.
In general, an increase in the rate results in a different competition between the first reaction in the scheme (Fig. 1), H + O 2 , and the follow-up reactions.For a higher rate, more HO 2 and subsequently more H 2 O 2 is formed.Therefore, less hydrogen atoms are available for other reactions and the formation yields of those products indeed decrease.The lower rate results in more reactions that occur 'deeper' in the reaction network, enhancing for instance H + OH -H 2 O.
Since the O 2 molecules in the pre-deposited oxygen layers probably have some preferred orientation (local crystallinity), the H atoms approach them under the same angle when penetrating the ice.This is reflected by a certain barrier.For co-deposition experiments, the angle dependence is of less physical importance since the oxygen beam provides O 2 molecules with a range of different orientations at the surface.This means that there are more molecules available for a reaction pathway with a lower barrier.Since in the Monte Carlo simulations the ice geometry is not explicitly taken into account, the barrier can be lowered for co-deposition simulations.
The For the analogous gas-phase reaction, H 2 O 2 is typically not found as a final product since there is no third body to remove the excess energy. 29,30Since experimentally OH and H 2 O 2 were found to behave similarly, a common formation route was suggested. 9In the standard case, (R2a): The latter branching ratio leads to a more open ice structure for the sequential simulations; the regularity of the O 2 matrix is destroyed.This is because HO 2 will split into two species which will take up interstitial space and distort the ice structure.Two OH radicals can then react again with H 2 O 2 , which forms additional pores.A more open ice structure allows H atoms to penetrate deeper in the ice layers and thus more reactions can take place.
In the co-deposition experiments the changed branching ratio is reflected in the OH/H 2 O 2 ratio at 15 K.At 25 K most of the formed OH immediately reacts with H 2 O 2 and little to no difference is observed in the final abundances.In the experiments approximately 3 times more OH is formed than H 2 O 2 at low temperature.In the simulations, this is a factor 2 for the standard simulations and 9 for the changed branching ratio.The real branching ratio is therefore probably in between the two adopted values.This value can only be quantitatively constrained in combination with the reaction rate of OH + OH -H 2 O 2 , which is discussed in the following subsection.have reaction enthalpies of À210 and À70 kJ mol À1 , respectively, and have small to zero gas-phase reaction barriers. 31he results from the simulations are trivial in the sense that an increase in the fraction of reaction (R3b) results in an increase of the H 2 O and O 3 and a decrease in the H 2 O 2 abundance for all parameter settings.
Using a non-zero contribution of reaction (R3b), we are able to correctly reproduce the trend in H 2 O abundance at 25 K for sequential hydrogenation, in particular the initial rise.Based on the small abundance of O 3 found in the experiments, we expect that reaction (R3a) is strongly favoured, which is reflected by a branching ratio of (R3a) : (R3b) = 90 : 10.

Diffusion and reaction.
Since OH is formed through the same reaction as H 2 O 2 but can still be detected in the ice, there must exist a mechanism that prohibits recombination to peroxide.
Experimentally, OH is only observed during co-deposition where it eventually freezes out in the O 2 matrix.When two OH radicals are formed from the reaction H + HO 2 the excess energy of the reaction allows the radicals to move away from each other on the surface just before thermalizing and subsequently freezing out in the matrix.Diffusion of thermalized OH at a later stage is not expected.If (R3) takes place, it immediately follows reaction (R2).This makes the two reactions hard to constrain separately.
To model the possibility of two products separating upon reaction is non-trivial, due to the lack of directionality in our simulations, i.e., a high diffusion rate allows fast hopping in all directions, both back and forth.This would allow two species to separate quickly, followed by a fast recombination as a result of hopping back.To prevent recombination, we introduce a relatively high diffusion rate and, contrary to current gas-phase studies, an effective barrier for the reaction OH + OH -H 2 O 2 .Different barrier heights are used for the reactions at 15 and 25 K in order to reproduce the experimental features.
The rate of this reaction has a strong effect on the appearance of the ice in the sequential simulations.For a high rate, H 2 O 2 is formed both at the surface and deep in the ice and a small amount of OH resides in the top layers.For a low rate, a layer of H 2 O covers the surface and OH occupies a large fraction of interstitial positions in the top layers of the matrix.For both the sequential and co-deposition simulations the OH/H 2 O 2 ratio decreases with increasing reaction rate.The decrease in OH also has an effect on the final H 2 O yield which decreases as well.
The intermediate reaction rate employed in the standard simulations reproduces the experiments best, since the high rate reproduces the sequential simulations but overproduces H 2 O 2 in the co-deposition simulations and the low rate overproduces OH in the sequential simulations but is in better agreement with the co-deposition experiment.As mentioned above, the exact rate can only be determined considering multiple input parameters, of which the OH diffusion and the H + HO 2 branching ratio are the most important.
Note, however, that OH diffusion inside the pre-deposited O 2 matrix is not likely.Therefore, in reality the OH fragments formed during sequential hydrogenation will react and the H 2 O 2 abundance correspondingly increases whereas the H  32 The barriers used in our simulations are considerably lower since we include a tunneling contribution as suggested by Miyauchi et al. 6 In general the H 2 O/H 2 O 2 ratio increases with increasing rate.In the sequential simulations this mainly occurs at later times, at the expense of H 2 O 2 , of which the abundance declines in time.By increasing the H + H 2 O 2 reaction rate, more H atoms are consumed in this reaction on top of the surface; less of the initial O 2 is thus affected upon hydrogenation.A decrease in the rate allows for the build-up of hydrogen peroxide and therefore a more correct ratio with respect to H 2 O at 15 K.However, the leveling off of the peroxide signal occurs at later times, contrary to the sharp transition observed experimentally.
For the co-deposition simulations, H 2 O is predominantly formed through (R2d) and (R3b) and the change in H + H 2 O 2 reaction rate only has an effect upon significant decrease of the barrier.The latter, however, results in an overestimation of H 2 O.
Additional simulations are performed to focus only on the H + H 2 O 2 reaction.A surface consisting only of H 2 O 2 molecules (20 K) is used, on top of which H and H 2 (300 K) are deposited.Using a reaction barrier of 1000 K indeed results in a reasonable correspondence to the hydrogenation experiment performed in ref. 9.However, on closer inspection of the dynamics, it becomes clear that all reactions between hydrogen and peroxide occur through the hot mechanism.We therefore expect this reaction to be rather inefficient in the dense interstellar medium.
4.2.5 H + OH.The reaction H + OH is barrierless, but with a very large excess energy of more than 5 eV.The dissipation of this energy with only one reaction product, H 2 O, might be problematic, although some of the excess energy is probably absorbed by the ice surface.Changing the rate of this barrierless reaction should therefore be interpreted as changing the efficiency by which the excess energy of this reaction can be dissipated.For a more in-depth discussion see ref. 11.
An increase in rate leads to a higher water production.At the same time, the production of H 2 O 2 goes down since less OH is available for its formation.The low and intermediate rate results appear to be in better agreement with the selected experiments since the H 2 O/H 2 O 2 ratio is better reproduced.
4.2.6 H 2 + HO 2 and H 2 + OH.Reactions involving molecular, rather than atomic hydrogen are of great importance in an astrochemical context due to the large abundance of H 2 in the interstellar medium.Molecular hydrogen that has been formed on the surface via is formed with an excess energy of 4.5 eV and we assume that it therefore desorbs upon formation, i.e. the fraction that stays on the surface, m, is 0. Only a very small fraction of gas-phase H 2 is hot, both experimentally and in the interstellar medium, and in our simulations it is implemented without substantial excess energy.Runs with m = 0.5 resulted in abundances within the limits of the standard deviations.
The reaction H 2 + HO 2 has a large barrier 33 and is even endothermic, therefore it is not expected to proceed at these low temperatures.However, co-deposition experiments with the same H fluence but with ten times higher In the simulations this reaction is included with a high barrier, and as a result the reaction does not contribute to the H 2 O 2 production.If the activation energy is substantially lowered, a small increase in the H 2 O 2 is observed for the co-deposition simulation at 25 K.The current simulations cannot constrain the value of the reaction barrier.
The reaction H 2 + OH has a gas-phase barrier of 2700 K. 34 Gas-phase experiments indicate that tunneling starts to play a role for OH + H 2 for temperatures below 250 K. 35 Recently, Oba et al. performed solid state experiments, confirming a dependency of the reaction rate on effective mass. 36The reaction rate at 15-25 K should therefore be substantially raised through tunneling.In the present simulations the reaction only plays a minor role in the water production regardless of the barrier for reaction.Even in the low barrier case, the H 2 O formation rate is only slightly changed by this reaction.However, upon closer inspection of the different contributions to the H 2 O formation rate, one notices that the fraction of water that is produced through H 2 + OH increases with reaction rate.The contribution of H + OH declines at the same time, resulting in a zero net effect.
4.2.7 Potentially important reactions.To conclude this section, we discuss the reactions that have been included in the standard and optimization simulation runs, namely Where possible, we performed several test runs including one of these aforementioned reactions in the standard simulation scheme.
Adding the first reaction, (R8), to the scheme with zero barrier has a negligible effect, since it accounts at most for 3% of the total HO 2 budget as a result of the low O production.
Concerning reaction (R9), Vakhtin et al. (2003) reported an empirical relation with a negative barrier in the gas-phase at relatively low temperatures (96-296 K). 37 This would render the reaction effectively barrierless for our purpose.However, a kinetic isotope effect was previously found, suggesting in fact the presence of an activation energy. 38Liquid-phase studies of this reaction indicate that the reaction rate can be decreased by as much as a factor of 1840 with respect to the gas phase. 39,40he presence of water can thus have a large effect on the probability for this reaction to occur.Moreover, a mechanism involving a reorientation into a five-membered prereactive complex was predicted, 39,40 but efficient reorientation is not likely in the solid state.Therefore, we exclude reaction (R9) from our simulations.Test runs were performed by using various settings and unequivocally poor reproductions of the experimental results were found through overproduction of water.
Reactions (R10) to (R13) should be considered part of the same train of thought, namely that hydrogen atoms are not all physisorbed but some are trapped in a chemical bond with for instance O 2 .These complexes or molecules can then continue to react with other species.Only incorporating reaction (R10) (with E a = 0 K 41 ) has little effect, except for co-deposition at 15 K, where the OH/H 2 O 2 ratio changes to a large extent.This can be easily overcome through changing the branching ratio of (R2).Incorporating also reaction (R11) (with E a = 0 K) results in a heavy overproduction of water.Reactions (R12) and (R13) allow for a type of proton exchange and may be considered as an alternative to hydrogen diffusion.Since these are reactions that go back and forth, they strongly increase the computational costs for barriers below 2500 K. Preliminary results show a decreased penetration depth for both temperatures as well as underproduction of H 2 O 2 .We postpone this concept for the time being and continue using our assumption of physisorbed hydrogen atoms diffusing on and through the ice.

Penetration depth
In Section 4.2.4 we mentioned that the standard simulations cannot reproduce the large experimental surface abundance of H 2 O 2 , 14 ML, at 25 K for sequential hydrogenation.This large abundance can accumulate in time through different mechanisms, for instance: (1) replenishing of the top layers by O 2 mobility and (2) competition between reaction and diffusion of H atoms. Several control experiments have been performed and it was established in Part I that the most likely scenario is the second mechanism, where hydrogen atoms can either react with a partner (initially O 2 ) or diffuse deeper into the ice.The extent to which the atoms enter the icy mantle is called penetration depth.The reaction H + O 2 is (close to) barrierless and temperature independent.Diffusion, however, is enhanced with increasing temperature and hydrogen atoms are thus able to reach deeper layers in the ice at higher temperature.In Fig. 7 a summary of the study on diffusive effects and penetration depth is presented, depicting the influence of initial ice thickness, availability of interstitial positions, the diffusion parameter x and the diffusion and positioning of H and O 2 in the ice.
The effect of the surface itself is discussed first.Increasing the initial surface roughness allows a higher production of both water and hydrogen peroxide, simply as a result of a larger surface area available for hydrogenation.Starting from a more porous ice, i.e., 10% pores in the standard surface, there is no significant effect on the penetration.The most probable explanation is that during the course of the simulation, pores are created as a result of products that continue to react further.A possibility for the low simulated use-up of O 2 could be the initial ice thickness, since the final experimental H 2 O 2 abundance is 14 ML which is more than our initial surface thickness (10 ML).This is illustrated in Fig. 7(a) and (f), and no dependence on ice thickness can be concluded.The fact that the ice thickness does not change abundances significantly is a sign of a lacking penetration mechanism.After all, if hydrogen atoms penetrate deeply into the matrix, a large O 2 reservoir would allow for more reactions to occur.The lack of this observation urges further investigation of diffusive effects with the purpose of understanding how we can increase the consumption of O 2 .
First of all, in the standard simulations there is already a mechanism included to account for some form of penetration, namely the availability of interstitial positions for small species, in particular hydrogen atoms.If this penetration mechanism is switched off, hydrogenation only occurs in the top layer of the ice and the final abundances drop to unrealistically low values as shown in Fig. 7(b) and (g).Furthermore, the shape of the H 2 O 2 curves changes and the sharp transition from linear growth to steady state is no longer present.Also, cannot be properly produced since it needs a second lattice site upon production.
The parameter x in eqn ( 4) is used to determine the ease with which molecules can diffuse, which regulates the diffusion both on the surface and into the bulk.A value higher than for the standard simulations decreases the diffusion rate, which leads to a low penetration depth since molecules are not formed deep in the ice and the top layers therefore easily block further reactions.A value of x lower than for the standard simulations increases the diffusion rate; this allows for more HO 2 production deeper in the ice.The abundance of stable species increases and that of the intermediates decreases as they are more accessible.A high diffusion rate, however, overestimates the abundance of OH and H 2 O for the sequential simulations and changes the OH evolution in the co-deposition simulations (not shown) so that it no longer follows the same trend as H 2 O 2 .Therefore, we conclude that the intermediate diffusion rate reproduces the experiments best.This brings us back to the two possible scenarios suggested in the experimental paper. 8A higher O 2 mobility is obtained through a lower barrier of O 2 diffusion or by easier access to interstitial positions.From Fig. 7(d) and (i) it is clear that both have no effect.Therefore, the second scenario was tested by varying the H mobility as depicted in Fig. 7(e) and (j).A lower barrier results in a steeper increase in the abundances of H 2 O and H 2 O 2 at 25 K and reproduces the transition between the initial rise and subsequent steady state behaviour at 15 K. Since this leads to a higher penetration depth at 25 K, we conclude that fast hydrogen diffusion could be a key feature in reproducing the competition between the reaction H + O 2 and diffusion into the ice, confirming the suggestions made in ref. 8.

Recommendations for future studies
We have gained insight into the surface processes linked to the formation of water ice by simulating previous experimental results with Continuous-Time Random-Walk Monte Carlo simulations and explicitly taking into account relevant surface effects.A systematic approach to varying reaction barriers is applied here to obtain a best fit model and characterize the sensitivity of the full reaction network to the various reactions.With this model, experimental trends are reproduced and specific reaction barriers have been obtained as given in Table 5.These rates are compared to literature values and recommended for future use in astrochemical models.Below, we describe the best fit found by using our simulations and we highlight some astrochemical considerations.

Best fit
The activation energies given in Table 5 are largely similar to those of the standard simulations, except for: the reaction barrier for the reaction H + O 2 in the case of co-deposition simulations, the branching ratio of the reaction OH + OH, which is set to (R3a) : (R3b) = 90 : 10, the binding energy of hydrogen atoms to the surface, E, lowered to 53 K.
Key reactions are indicated in boldface in Table 5.In Fig. 8 the surface abundances as obtained by the best fit simulation are plotted versus time and can be compared to both the standard simulations and the experimental trends.
In general, a clear agreement between experimental and simulated data can be observed.Prominent is the reproduction of the H 2 O 2 abundance evolution where the shape of the curves for sequential hydrogenation is correct as well as the overall behaviour under influence of a change in temperature.Moreover, with the 'best fit' settings the H 2 O production is roughly similar at 15 and 25 K, as it should be.Note that water does remain slightly overproduced at low temperature, likely due to the dynamics of the OH radical (diffusion versus a barrier for (R3)).Furthermore, the best fit settings are used to test once more the influence of the initial ice thickness.With these improved settings, indeed a small increase in H 2 O 2 abundance is found at 25 K using a 16 ML thick surface, indicating an enhancement of hydrogen penetration.
Considering the co-deposition simulations, the 'best fit' runs are performed with parameters equal to those for sequential hydrogenation, except for the rate of the reaction of H + O 2 .As mentioned in Section 4.2.1, the angle dependence of this reaction is less of a restriction during co-deposition and a lower barrier for this reaction is better at reproducing the HO 2 abundances.The overproduction of HO 2 decreases slightly and the production of H 2 O 2 increases with respect to the standard simulations.Moreover, the production of O 3 is enhanced, which is in agreement with the experimental results, judging by Fig. 5.3 from ref. 9. Furthermore, the branching ratio of the reaction H + HO 2 is found to be within the following range:

Astrochemical considerations
Monte Carlo simulations have the potential to bridge the difference in conditions (time constraints, abundances and fluxes) between experiments and the interstellar medium by investigating surface reaction mechanisms under first laboratory, and then interstellar medium conditions.
For instance, a high hydrogen diffusion rate in O 2 ice is necessary to reproduce the large penetration depth observed experimentally.However, the penetration depth is a bulk effect, which is largely absent under astrochemical conditions.A pure O 2 ice is not likely to be present in space and diffusion will be restricted mainly to the surface layers.In this work the goal is to understand the physico-chemical principles, and therefore the penetration depth is an effect that does need to be investigated.
Penetration depth plays a somehow less important role in the co-deposition simulations and since most features are correctly reproduced, we are confident that this set of parameters can be used in astrochemically relevant simulations to reproduce observations of water in different regions of the interstellar medium.This will also be the topic of a future study in which the influence of directionality in a hot diffusion mechanism will be incorporated.
Another  for the surface reaction (R4). 46They found a very high gasphase H 2 O 2 abundance for a range of physical conditions which can be explained by the lack of gas-phase destructive mechanisms.Moreover, the abundance of surface H 2 O 2 as obtained in their model is affected by the height of the barrier of reaction (R4). 47Increasing the barrier for reaction (R4) to our values would result in higher solid H 2 O 2 abundances, further strengthening their and our suggestion to investigate the surface formation and destruction of H 2 O 2 under interstellar relevant conditions in more depth.
The unknown diffusive behaviour of radical species, such as OH and H, is key to predict water formation under conditions relevant to the interstellar medium.However, experiments investigating diffusion are scarce due to the non-trivial competition with the high reaction probabilities involved.In the present study the two most relevant species in this respect are OH and H.The first because the experimental detection in co-deposition experiments poses a challenge as reaction (R4) was initially thought to occur barrierless.Furthermore, the diffusion of hydrogen on a surface was studied experimentally only by Katz et al. in 1999 and barriers found on carbon and silicate surfaces are not necessarily applicable to icy grains as well. 48,49Although the H diffusion is not fully constrained, it is a key parameter for our simulations as we extensively discussed in Section 4.3.
Finally, for astrochemical applications, the reactions with H 2 , here discussed in Section 4.2.6, are actually relevant to dark cloud conditions. 3Since molecular hydrogen possesses a rather strong intramolecular bond, barriers are typically high and tunneling is required to increase the rate of reaction with other ice species.This tunneling should scale with the so-called effective mass of the total system, since the system needs to be considered as a whole as outlined in ref. 36.This is by no means trivial and depends strongly on the binding of the molecule with the surface as well as on the evolution of the effective mass with reaction coordinate. 50Also here further dedicated experimental research as well as (theoretical) modeling to study reactions with H 2 is needed for a deeper understanding of the full reaction scheme.

Practical use of the best fit parameters
In our Monte Carlo simulations the competition between different processes for a single species is explicitly taken into account.This is not easy to implement in classic rate equation models.However, in the interstellar medium, ice chemistry is limited by diffusion of surface species.Especially in dark molecular clouds, species on the surface are thermalized before engaging in reactions.Therefore, reactions with high rates (either due to a low activation energy or the result of tunneling) will dominate.Here we find a number of key reactions for the O 2 + H route, as summarized in Table 5.When implementing these reaction rates in (astrochemical) models, the rates and activation energies given in Table 5 need to be considered in such a way that (i) the ratio between two reaction rates corresponds to the ratio found here and (ii) the rates are scaled according to the prefactor being used that amounts here to 1 Â 10 12 s À1 .

Conclusions
Solid state hydrogenation reactions of O 2 ice have been simulated with a Continuous-Time Random-Walk Monte Carlo method, explicitly taking into account the recent findings of a number of experimental studies.A strategy of systematically varying the different processes at play in the simulations, and comparing the with experimental values, allows us to derive the key processes.In this way, we can gain insight into the actual processes taking place.Our model consists of a combination of reaction barriers (or activation energies) and binding energies for all species present in the reaction network (Fig. 1).In order to probe the large parameter space, we selected four previously performed different experiments.These are sequential hydrogenation and co-deposition experiments at 15 and 25 K that are all reproduced with a single set of parameters.The best fit model reproduces experimentally observed trends using the binding energies given in Table 2, complemented with the reaction barriers given in Table 5.From an extensive set of simulations of the selected laboratory data, we conclude the following: the key reactions for the reaction route starting from O 2 are (R1) H + O 2 , (R2) H + HO 2 , (R3) OH + OH, (R4) H + H 2 O 2 and (R5) H + OH, a relatively slow destruction of H 2 O 2 explains the high accumulation of this species, a high hydrogen diffusion rate is necessary to reproduce the large penetration depth of H into the O 2 ice observed experimentally, and the diffusive behaviour of radical species, such as OH and H, is a key parameter.
et al. (part I, ref. 8) and Cuppen et al. (part II, ref. 9) explored the dependency on a variety of experimental conditions, such as temperature, thickness, H-atom flux, etc.

Fig. 1 A
Fig. 1 A schematic representation of the extended reaction network as obtained in Cuppen et al. 9 Black, straight lines indicate reactions for which the influence of the barrier is studied here.Grey, dashed lines indicate reactions that are kept constant for all simulations.

Fig. 2
Fig. 2 Evolution of the experimental surface abundance of H 2 O 2 (light blue), H 2 O (dark blue), HO 2 (green) and OH (orange) as a function of time for a sequential hydrogenation simulation (left) and a co-deposition simulation (right) at 15 and 25 K.

Fig. 3 A
Fig. 3 A cross section of the simulated ice mantle for a hydrogenation simulation of predeposited solid O 2 at 15 K (left) and 25 K (right).Standard values are used for the reaction rates.O 2 is represented by red, H 2 O 2 by light blue, H 2 O by dark blue, HO 2 by green, OH by orange and O 3 by violet.On the vertical axis the number of monolayers (ML) is given.

Fig. 4
Fig. 4 Similar to Fig. 3 for co-deposition simulations.For reasons of clarity, the center 50 ML of the matrix are omitted from the surface.
Fig. 5 Evolution of the surface abundance of O 2 (red), H 2 O 2 (light blue), H 2 O (dark blue), HO 2 (green), OH (orange) and O 3 (violet) as a function of time for a hydrogenation simulation of predeposited solid O 2 for the five key reactions at 15 K ((a)-(e)) and 25 K ((f)-(j)).Solid curves represent simulation results with the standard value, dashed curves with lower barrier (higher rate), and dash-dotted curves with higher barrier (lower rate).

4 . 2 . 3
OH + OH 4.2.3.1 Branching ratio.The competing reactions, OH + OH -H 2 O 2 (R3a) OH + OH -H 2 O + O, (R3b) H 2 fluence show an increase in H 2 O 2 and decrease in HO 2 in agreement with the proposed reaction.Also in ref. 11 an observed isotope effect between hydrogenation and deuteration of ozone suggested this reaction to occur under laboratory conditions.It was therefore hypothesized that the excess energy from the reaction H + O 2 can be used to overcome both the reaction barrier and the endothermicity of the H 2 + HO 2 reaction since the total reaction H + H 2 + O 2 -H 2 O 2 + H is exothermic.This journal is c the Owner Societies 2013 Phys.Chem.Chem.Phys., 2013, 15, 8287--8302 8297

Fig. 7
Fig. 7 As for Fig. 5, but investigating different parameters relating to penetration depth.Solid curves represent simulation results with standard values.Dashed curves in panel (a) represent standard values applied to a 16 ML surface, and in panel (b) easy access to interstitial positions was switched off.In panel (c) dashed curves represent a low x value (eqn (4)) and dashed-dotted curves a high value.In panels (d) and (e) dashed curves represent a lower barrier of diffusion and dashed-dotted curves represent easier access to interstitial positions for O 2 and H, respectively.Panels (f)-(j) give the 25 K situation.
(R2a):(R2b):(R2c):(R2d) = (56-65):(35-26):(2):(7).The effect of changing this translates into a changed OH/H 2 O 2 ratio for co-deposition at 15 K while no change is observed for co-deposition at 25 K or sequential hydrogenation simulations.Note that the uncertainty in the dynamics of the OH radical can also affect the time evolution of the abundance of this species.The comparison made in Table 4 between different temperatures and H/O 2 ratio remains unchanged.Finally, we tried to reproduce the sequential hydrogenation study from Miyauchi et al. (2008) using the best fit parameters combined with experimental conditions mentioned in ref. 6.We corrected for the different band strengths used as well as the systematic shift in temperature reading of 2 K (compare Fig. 2.6 in ref. 24 and Fig. 3 in ref. 42).Our simulations confirm their experimentally observed trends to within 50%, in terms of H 2 O and H 2 O 2 abundance as well as evolving trends in time.

Fig. 8
Fig. 8 Evolution of the surface abundance of O 2 (red), H 2 O 2 (light blue), H 2 O (dark blue), HO 2 (green), OH (orange) and O 3 (violet) as a function of time for a sequential hydrogenation simulation (left) and a co-deposition simulation (right) at 15 and 25 K compared to the experimental values.The thin lines represent the standard simulations, whereas the thick lines represent the simulation connected to the best fit parameters of Table5.Note that the experimental surface abundances for the co-deposition simulations are in units of integrated absorbance, due to a lack of available band strengths for matrix isolated species.

Table 1
Experimental conditions for the four selected experiments; two sequential hydrogenation and two co-deposition experiments at 15 and 25 K Co-deposition: linearly increasing behaviour of all abundances with time at t exp o 60 min, decrease of HO 2 and OH abundances between 15 and 25 K, increase of H 2 O 2 abundance between 15 and 25 K, low to zero production of H 2 O, and constant ratio of OH/H 2 O 2 for 15 K and for 25 K.

Table 2
shows the single binding energy values, E, used throughout this paper.Each species in the lattice has 6 neighbours and 12 nearest neighbours.Nearest neighbours add a contribution of E and next-nearest neighbours of E/8.

Table 3
List of surface reactions with the used rates for thermalized reactions Standard values used throughout the paper.b Values used to test the effect of the reaction barrier on the overall performance.c Total rate of 1.0 Â 10 12 s À1 for H + HO 2 .Individual channels are corrected for their branching ratios.d Total rate of 4.2 Â 10 À6 or 1.2 Â 10 2 s À1 for 15 and 25 K, respectively.Individual channels are corrected for their branching ratios.
NOTE: rates indicated in ITALICS should be used with care, see text.a

Table 2
Single binding energy values of a species to one surface site, E, as implemented in the simulations

Table 4
Experimental and simulated relative surface abundances For the reaction H + HO 2 , not the rate, but the branching ratio for the product channels OH + OH and H 2 O 2 is varied.
is also allowed in the simulations with a rate of 2 Â 10 À4 s À1 .4.2.2H + HO 2 branching ratio.
2 O abundance decreases.Indeed, this corresponds with a scenario where the H 2 O/H 2 O 2 fraction gets closer to the experimental value.4.2.4H + H 2 O 2 .H 2 O production can proceed through various reactions, of which H + H 2 O 2 is a special case as a result of the high barrier associated with it.Gas-phase calculations suggest a barrier of approximately 2000 K.
important difference between experiments and astronomical observations is the abundance of formed solid H 2 O 2 with respect to H 2 O. Starting from O 2 , a high abundance of H 2 O 2 is found.Here we show that a relatively slow destruction of H 2 O 2 explains the high accumulation of this species.Fig. 2 in ref. 9 indicates that a high and more interstellar

Table 5
Parameters that reproduce the selected experiments the best.A comparison to (mean) gas-phase literature values is also included.Key reactions are listed in bold face NOTE: the hydrogen diffusion barrier used was 53 K instead of the value listed in Table 3. a This work.b A value of 8.3 Â 10 5 is used for the co-deposition simulations.c Reaction practically does not take place in our simulations.d See Section 4.2.3.relevant H/O 2 ratio in co-deposition experiments is in favor of water formation.A more detailed study of the role of the H 2 O 2 surface destruction channel (R4) within the whole reaction network under interstellar conditions is highly needed.Du et al. used a gas-grain model to reproduce the observed gasphase abundance of H 2 O 2 in r Oph A using a low barrier of 92 K