Effect of binding energies on the encounter desorption

The abundance of interstellar ice constituents is usually expressed with respect to the water ice because, in denser regions, a significant portion of the interstellar grain surface would be covered by water ice. The binding energy (BE), or adsorption energy of the interstellar species regulates the chemical complexity of the interstellar grain mantle. Due to the high abundance of water ice, the BE of surface species with the water is usually provided and widely used in astrochemical modeling. However, the hydrogen molecules would cover some part of the grain mantle in the denser and colder part of the interstellar medium. Even at around ~ 10K, few atoms and simple molecules with lower adsorption energies can migrate through the surface. The BE of the surface species with H2 substrate would be very different from that of a water substrate. However, adequate information regarding these differences is lacking. Here, we employ the quantum chemical calculation to provide the BE of 95 interstellar species with H2 substrate. These are representative of the BEs of species to a H2 overlayer on a grain surface. On average, we notice that the BE with the H2 monomer substrate is almost ten times lower than the BE of these species reported earlier with the H2 O c-tetramer configuration. The encounter desorption of H and H2 was introduced (with ED (H, H2 ) =45 K and ED (H2 , H2 ) =23 K) to have a realistic estimation of the abundances of the surface species in the colder and denser region. Our quantum chemical calculations yield higher adsorption energy of H2 than that of H (ED (H, H2 ) = 23 - 25 K and ED (H2, H2 ) =67 - 79 K). We further implement an astrochemical model to study the effect of encounter desorption with the resent realistic estimation. The encounter desorption of the N atom (calculations yield ED (N, H2 ) =83 K) is introduced to study the differences with its inclusion.


INTRODUCTION
Interstellar grains mainly consist of amorphous silicate and some form of carbonaceous materials (Li, 2004). It is now well established that these grains can significantly constrain the chemical composition of molecular clouds or star-forming regions. In the cloud's denser regions, where the temperature is reasonably low (∼ 10 K), the grain surface would be covered by icy layers. A sizable portion of these icy layers may contain water molecules. Thus, providing the binding energies (BEs) with the water as a substrate is proper. In reality, the surface species would face a bare grain in the diffuse region. Some ice layers would grow on the top of this grain surface in the denser medium and host the incoming species. The composition of this ice layer depends on the initial elemental abundance of the species in that region. It would not necessarily always be H 2 O dominated. There are ample examples of the presence of a notable portion of CO, CO 2 , CH 4 , NH 3 , CH 3 OH, etc., on the ice (Gibb et al., 2004;Das et al., 2010;Das and Chakrabarti, 2011;Das et al., 2016;Gorai et al., 2020a).
Hydrogen molecules are ubiquitous in denser regions of the interstellar medium (ISM). Thus, its accretion rate on the grain is much higher in comparison to the others. However, because of the low adsorption energy, it can be easily desorbed from the grain. Despite this, a significant portion of the grain mantle would be covered by molecular hydrogen, especially around the cold and dense interstellar condition. This inhomogeneous surface coverage could influence the mobility of the other surface species. Initially, the encounter desorption mechanism was introduced by Hincelin et al. (2015) to eliminate the overestimation of the abundance of molecular hydrogen on the grain. This desorption process occurs during surface diffusion and is induced by the presence of repulsive inter-H 2 forces, effectively reducing the BE of H 2 . They considered gH 2 + gH 2 → H 2 + gH 2 , where, "g" designates the grain surface species. They obtained an excellent match within the microscopic Monte Carlo method and the rate equation approach when they implemented this unique approach. The Monte Carlo approach is best suited for monitoring the chemical composition of the grain mantle. However, it is timeconsuming (Chakrabarti et al., 2006a,b;Das et al., 2008aDas et al., , 2010Das and Chakrabarti, 2011;Das et al., 2016;Cuppen and Herbst, 2007). Recently, Chang et al. (2021) considered a similar process and included H's desorption by a similar mechanism. They considered gH + gH 2 → H + gH 2 , which means, whenever the surface H meets one surface H 2 , surface H desorb with a certain probability. They reported a significant difference between the formation of some key surface species with the inclusion of this treatment.
A substantial amount of BE values are available from the temperature-programmed desorption (TPD) studies on various model substrates like graphite, diamond-like carbon, amorphous or crystalline silica, silicates, water, and other ice surfaces (Collings et al., 2004;Ward et al., 2012;Noble et al., 2012;Dulieu et al., 2013). But, the BE of the species with H 2 substrate is yet to be known. Cuppen and Herbst (2007) had estimated the BE of H atom on H 2 substrate ∼ 45 K by following Vidali et al. (1991) A vital impediment in examining the encounter desorption with other species is the shortage of information about the adsorption energy of these species with H 2 molecule. Here, we employ quantum chemical calculations to determine the adsorption energy of these species with H 2 molecule. Obtained BE assessments are executed in our Chemical Model of Molecular Cloud (hereafter CMMC) (Das et al., 2015a,b;Gorai et al., 2017a,b;Sil et al., 2018;Gorai et al., 2020). The encounter desorption effect is vital during the prestellar core phase. Study the formation of stars is one of the essential intricacies of astrophysics. A complete understanding of the star formation process is yet to be established. However, in brief, stars are formed by a long condensation process (Pagani et al., 2013). In the beginning, warm diffuse material (∼ 8000 K) converts into a cold neutral atomic gas (∼ 100 K and ∼ 10 − 100 cm −3 ). After further evolution, it transforms into a more dense region (10 2 − 10 4 cm −3 and ∼ 10 − 20 K). If no other heating source is present, then a dense core (> 10 4 cm −3 ) appears in some places of these turbulent materials. Some of these cores further evolve into prestellar cores (> 10 5 cm −3 ) (Bergin and Tafalla, 2007;Keto and Caselli, 2008). Prestellar cores further continue their evolution for the formation of the protostar. Due to the accretion of atoms and molecules, gas-phase abundance is depleting, whereas the molecular ice mantles form. The chemical composition of the grain mantle is mainly governed by the addition of atomic hydrogen with the atoms or simple molecules. The chemical composition of the bulk ices further varies with the star formation process associated with it. Depending on this, it is expected that the ice composition would be very different in various places. However, from the infrared observations, it was observed that the significant repositories of interstellar hydrogen, oxygen, carbon, and nitrogen are H 2 O, CH 3 OH, H 2 CO, CO, CO 2 , CH 4 , and NH 3 (Gibb et al., 2000;Whittet et al., 2007;Öberg et al., 2008;Boogert et al., 2015). This paper is compiled as follows. In Section 2, we confer computational methodology. Discussion and results are presented in Section 3, and finally, in Section 4, we conclude.

Quantum chemical calculations
Here, we have utilized the Gaussian 09 suite of programs (Frisch et al., 2013) for quantum chemical calculations. In a periodic treatment of surface adsorption phenomena, the BE is related to the interaction energy (∆E), as: For a bounded adsorbate, the BE is a positive quantity and is defined as: where E ss is the optimized energy for the complex system where a species is placed at a suitable distance from the grain surface. E surface and E species are the optimized energies of the grain surface and species, respectively.
To find the optimized energy of all structures, we have used a Second-order Møller-Plesset (MP2) method with an aug-cc-pVDZ basis set (Dunning, 1989). We have considered 95 interstellar species for the computation of their BEs with the H 2 substrate. To make the calculation more straightforward, we have considered a monomer configuration of the H 2 molecule as an adsorbent. The adsorbates noted in Table 1 are placed at a suitable distance from the adsorbent with a weak bond so that a Van der Waals interaction occurs during the optimization. All the optimized geometries are provided in the supplementary information. We must confess that the interstellar species considered in this study are often larger than the H 2 . Since the estimated BEs are different in different locations, this may lead to a fallacious estimation. It is recommended to take the average whenever different binding sites are found. Despite these flaws, it can provide us with a general picture and startup initiative to compare the BE of a species with the water and H 2 substrate. Following the BE calculations carried out by Das et al. (2018) (see Tables 2 and 3), here, we have not considered the ZPE and BSSE corrections for our BE calculations. All the obtained BE values are noted in Table 1. It is interesting to note that except for phosphorous, the calculated adsorption energy of most of the abundant atoms (H, C, N, O, and S) with the H 2 substrate is found to be < 100 K.
In Table 1, we have reported the BE values obtained by considering a free-standing H 2 interacting with a species. But in reality, this H 2 would be pre-adsorbed and can feel the surface. It would yield a different BE than the previous case. To check the effect of condensed H 2 O in the ice phase, we also have calculated the interaction energy by considering the molecule embedded in a continuum solvation field.  Sandford and Allamandola (1993).
For this purpose, we have examined the local effects and the integral equation formalism (IEF) variant of the polarizable continuum model (PCM) (Cancès et al., 1997;Tomasi et al., 2005) with water as a solvent (Gorai et al., 2020a). The obtained values for H, H 2 , N with the IEFPCM model are noted in Table 1 (in  parentheses). The two calculations significantly differ. For example, with the free-standing H 2 , we have the BE of H, H 2 , and N ∼ 23 K, 67 K, and 83 K, respectively, whereas with the IEFPCM model, we have obtained ∼ 25 K, 79 K, and 78 K, respectively. So, the free-standing H 2 underestimated the BE of H and H 2 by 2 K and 12 K, whereas it overestimated N's case by 5 K. We also provided in Table 1 the literature BE values (Cuppen and Herbst, 2007;Sandford and Allamandola, 1993) (if available) for the comparison. Das et al. (2018) provided BEs of the ∼ 100 interstellar species considering the c-tetramer configuration of water molecules. Table 1 shows the BE of the roughly same interstellar species with the H 2 monomer. We have noticed that the obtained BEs with H 2 are much smaller than those of the water tetramer configuration. On average, we have received almost ten times lower BEs with the H 2 surface. Table 1 shows the ground state of the species used to calculate the BE. The values of BE are very much sensitive on the chosen ground state spin multiplicity. To evaluate the ground state spin multiplicity of each species, we have taken the help of Gaussian 09 suite of program. The way to check for the ground state spin multiplicity is to run separate calculations (job type "opt+freq"), each with different spin multiplicities, and then compare the results between them. The lowest energy electronic state solution of the chosen spin multiplicity is the ground state noted for the species in Table 1.

Astrochemical model
We have included the encounter desorption phenomenon in our CMMC code (Das et al., 2015a(Das et al., ,b, 2016Gorai et al., 2017a,b;Sil et al., 2018;Gorai et al., 2020) to study its effect. The surface chemistry network of our model is mostly adopted from Ruaud et al. (2015); Das et al. (2015b); Gorai et al. (2020). The gas-phase network of the CMMC model is mainly adopted from the UMIST database (McElroy et al., 2013). Additionally, we have also included the deuterated gas-phase chemical network from the UMIST. A cosmic ray rate of 1.3×10 −17 s −1 is considered in all our models. Cosmic ray-induced desorption and nonthermal desorption rate with a fiducial parameter of 0.01 is considered. For all the grain surface species, we have adopted a photodesorption rate of 1 × 10 −4 per incident UV photon (Ruaud et al., 2015). A sticking coefficient of 1.0 is considered for the neutral species except for the H and H 2 . The sticking coefficients of H and H 2 are considered by following the relation proposed by Chaabouni et al. (2012). Following Garrod and Pauly (2011), here, we have implemented the competition between diffusion, desorption, and reaction. For the diffusion energy (E b ), we have considered R× adsorption energy (E D ). Here, R is a scaling factor that can vary between 0.35 and 0.8 (Garrod et al., 2007). The BE of the species is mostly considered from Wakelam et al. (2017), and a few from Das et al. (2018). Table 2 refers to the adopted initial abundances concerning the total hydrogen nuclei in all forms. Except for HD's value in Table 2, elemental abundances are taken from Semenov et al. (2010). We considered the initial abundances of HD from Roberts and Millar (2000).
The encounter desorption effect was first introduced by Hincelin et al. (2015). The rate of encounter desorption of H 2 on the surface of H 2 is defined as: where gH 2 is the surface concentration of H 2 molecules in cm −3 , P(H 2 , H 2 ) defines the probability of desorption over the diffusion, and k H 2 ,H 2 is the diffusion rate coefficient over the H 2 O substrate. k H 2 ,H 2 Table 2. Initial elemental abundances considered in this study.

Frontiers
Species Abundances H 2 5.00 × 10 −1 He 9.00 × 10 −2 N 7.60 × 10 −5 O 2.56 × 10 −4 C + 1.20 × 10 −4 S + 8.00 × 10 −8 Si + 8.00 × 10 −9 Fe + 3.00 × 10 −9 Na + 2.00 × 10 −9 Mg + 7.00 × 10 −9 Cl + 2.00 × 10 −10 HD 1.60 × 10 −5 is defined as follows (Hasegawa et al., 1992): In the above equation, n d is the dust-grain number density, κ is the probability for the reaction to happen (unity for the exothermic reaction without activation energy), and R dif f is the diffusion of the species. P(H 2 , H 2 ) in Equation 3 is defined as:  Chang et al. (2021), the encounter desorption of species X is defined as: where h X,H 2 is the hopping rate over H 2 O surface (ν exp(−E b (X, H 2 )/T)), P(X,H 2 ) is the desorption probability of gX while encountering with gH 2 , and P X denotes the probability of gX to migrate at the location of gH 2 over the H 2 O substrate. P(X,H 2 ) and P X are defined as, P(X, H 2 ) = Desorption rate of X on H 2 substrate Desorption rate of X on H 2 substrate + Hopping rate of X on H 2 substrate ,

Encounter desorption of H 2
First of all, we have benchmarked our model with Hincelin et al. (2015). In Figure 1, we have compared our results with those obtained in Hincelin et al. (2015). For this comparison, following Hincelin et al. (2015), we have used T= 10 K, E D (H 2 , H 2 O) = 440 K, E D (H, H 2 O) = 450 K, E D (H 2 , H 2 ) = 23 K, and R=0.5. Solid curves in Figure 1 represent the cases obtained here, and the rest are extracted from Hincelin et al. (2015) by using the online tool of Rohatgi (2020). Our results with and without encounter desorption show an excellent match with Hincelin et al. (2015). Presently in the KIDA database (kida.astrophy.u-bordeaux.fr), more updated BE values were listed. It suggests that E D (H, H 2 O) = 650 K. The results obtained from our quantum chemical calculations shown in Table 1 represent the estimated BE values with the H 2 substrate. In the following section, we have used these updated energy values, and the effects of their changes are discussed. Figure 2 shows the time evolution of gH 2 by considering n H = 10 7 cm −3 , T = 10 K and R = 0.35 − 0.80. Interestingly the abundance of gH 2 seems to be invariant with R's changes, whereas it strongly depends on R in encounter desorption. R's lower value means a quicker hopping rate, whereas a higher value represents a delayed hopping rate. With the increase in R, gH 2 abundance raises for the encounter desorption case. It means that as we rise R's value, the encounter desorption effect depreciates. The left panel of Figure 5 exposes that with the increase in R's value, a steady decrease in the ratio between the gH 2 abundance with no encounter desorption case (NE) and with encounter desorption case (EN) is obtained. The probability of the encounter desorption is inversely proportional to the rate of diffusion Table 3. The obtained abundance of gH 2 , gH, gH 2 O and gCH 3 OH for the effect of encounter desorption of H 2 under various situation with R = 0.35, n H = 10 7 cm −3 , and T = 10 K.

Case
Case specification Abundance at 10 6 years with n H = 10 7 cm −3 No. gH No encounter desorption 2.0011749 × 10 −4 2.311976 × 10 −24 8.8375039 × 10 −5 8.8416179 × 10 −6 (0.00) (0.00) 2 E D (H 2 , H 2 ) = 23 K 1.183836 × 10 −11 3.2799079 × 10 −24 8.927643 × 10 −5 6.3629209 × 10 −6 (Hincelin et al., 2015) (  (Equation 7). Since the increase in the value of R induces faster diffusion and hopping, it is lowering the encounter desorption probability of H 2 as expected. Figure 3 shows the time evolution of gH 2 with NE and EN when we have used R = 0.35, T = 10 K and n H = 10 4 − 10 7 cm −3 . In both cases, abundances of gH 2 increase with the density. The middle panel of Figure 5 shows that the gH 2 abundance ratio between NE and EN with density. It depicts that the effect of encounter desorption is more pronounced for the higher density. Figure 4 shows that the gH 2 abundances when we have used n H = 10 7 , R = 0.35, and T= 5 − 20 K. In the right panel of Figure 5, we have shown the gH 2 abundance ratio obtained between NE and EN with the temperature changes. From the figures, it is seen that the effect of encounter desorption is maximum toward the lower temperature (∼ 10 K), and it ceases around 20 K. The curve is similar to the H 2 formation efficiency discussed in Chakrabarti et al. (2006a,b) for olivine grain. With the decrease in temperature, H atoms' mobility decreases. Thus, the formation rate decreases. With the increase in temperature, the hopping rate increases, which can increase the formation efficiency, but at the same time, the residence time of H atoms decreases which affects the H 2 formation efficiency. As a result, the H 2 formation efficiency is maximum at around ∼ 10 K, and the encounter desorption effect is pronounced at the peak hydrogen formation efficiency.
For a better illustration, the obtained abundances with R = 0.35, T=10 K, and n H = 10 7 cm −3 are noted in Table 3 at the end of the total simulation time (∼ 10 6 years). Chang et al. (2021) considered the competition between hopping rate and desorption rate of H 2 (Equation 7), whereas Hincelin et al. (2015) considered the battle between the diffusion and desorption rate of H 2 (Equation 5). This difference in consideration resulting ∼ two times higher abundance of gH 2 with the consideration of Chang et al. (2021) compared to Hincelin et al. (2015) (see case 2 and 3 of Table 3 and Figure 2). Our quantum chemical calculation yields E D (H 2 , H 2 ) = 67 K, which is higher than it was used in the earlier literature value of 23 K (Cuppen and Herbst, 2007;Hincelin et al., 2015;Chang et al., 2021). The computed adsorption energy is further increased to 79 K when we have considered the IEFPCM model. Table 3 shows that increase in the BE (E D (H 2 , H 2 ) = 67 K, and 79 K, case 4 and 5 of Table 3) results in sequentially higher . We have seen significant differences when we have used different energy barriers and different methods (Hincelin et al., 2015;Chang et al., 2021). Obtained values of gH 2 are further noted in Table 3 for better understanding.
surface coverage of gH 2 than it was with E D (H 2 , H 2 ) = 23 K (case 3 of Table 3). In case 5 of Table  3, we have noted the abundance of gH 2 when no encounter desorption effect is considered, but a higher adsorption energy of H atom is used (E D (H, H 2 O) = 650 K). Case 6 of Table 3 also considered this adsorption energy of H atom along with E D (H 2 , H 2 ) = 67 K, and the method of Chang et al. (2021) is used. A comparison between the abundance of gH 2 of case 4 and case 6 (the difference between these two cases are in consideration of the adsorption energy of gH) yields a marginal decrease in the abundance of gH 2 when higher adsorption energy of gH is used. It depicts that the effect of encounter desorption increases with the increase in density.

gH
The obtained abundance of gH is noted in Table 3. The gH abundance is marginally decreased in Chang et al. (2021) compared to Hincelin et al. (2015). The use of higher E D (H 2 , H 2 ) (∼ 67 K and 79 K) lowers the value of gH compared to case 2. However, the use of the H atom's higher adsorption energy (650 K) can increase the gH abundance by a couple of orders of magnitude (see case 7 of Table 3).

gH 2 O, gCH 3 OH
The effect of the encounter desorption on the other major surface species (gH 2 O and gCH 3 OH) is also shown in Table 3. In the bracketed term, we have noted the percentage increase in their abundances from the case where no encounter desorption was considered (for E D (H, H 2 O) = 450 K and 650 K, respectively). Table 3 depicts that the consideration of encounter desorption of H 2 can significantly change (decrease by ∼ 27 − 30%) the methanol abundance (case 3 and case 7) from that was obtained with the no encounter desorption (case 1 and case 6). However, the changes in the surface abundance of water are minimal (∼ ±1%) for the addition of the encounter desorption of H 2 . These changes (increase or decrease) are highly dependent on the adsorption energy of H, temperature, density, and the value of R (∼ 0.35 noted in Table 3). The changes in E D (H 2 , H 2 ) from 23 K to 67 K can influence the surface abundance of methanol and water. For example, in between case 3 and case 4 of Table 3, we can see that there is a significant increase (∼ 15%) in the abundance of gH 2 O when higher adsorption energy (E D (H 2 , H 2 ) = 67 K) is used. However, this higher adsorption energy can marginally underproduce the methanol on the grain. In brief, from Table 3, it is clear that the encounter desorption can significantly change the abundances of surface species. Still, these changes are highly dependent on the adopted adsorption energy with the water and H 2 ice and adopted physical parameters (n H , R, T ).

Encounter desorption of other species
The idea of encounter desorption (Hincelin et al., 2015) primarily arose to eliminate the enhanced surface coverage of H 2 in the relatively denser and colder medium. Since H 2 has lower adsorption energy with the water surface (∼ 440 K), it could move on the surface very fast and occupy a position on the top of another H 2 molecule. Comparatively, in the denser and colder region, the chances of this occurrence enhance. Since the H 2 molecule on H 2 has negligible BE (23 K used in Cuppen and Herbst, 2007;Hincelin et al., 2015), it could easily desorb back to the gas phase. Other surface species can, of course, meet with H 2 , but the idea of this encounter desorption arises when the species can occupy a position on the top of the H 2 molecule. For example, a carbon atom is having a BE of 10000 K (Wakelam et al., 2017). H 2 could quickly meet one C atom on the grain surface, but due to the lower mobility of atomic   carbon at a low temperature, every time H 2 will be on the top of the carbon atom. Since the whole C-H 2 system is attached to the water substrate; this will not satisfy the encounter desorption probability. Among the various key elements considered in this study, gH, gN, and gF have the BE of 650 K (Wakelam et al., 2017), 720 K (Wakelam et al., 2017), and 800 K (listed in the original OSU gas-grain code from Eric Herbst group in 2006), respectively, with the water ice. It yields a reasonable timescale for hopping even at a low grain temperature (∼ 10K). Since the initial elemental abundance of F is negligible, we can neglect its contribution. The hopping time scale is heavily dependent on the assumed value of R. For example, by considering R = 0.35, at 10 K, the hopping timescale for gH and gN is 1.12×10 4 years (with E D (H, H 2 O) = 650 K) and 4.61×10 −3 years (with E D (N, H 2 O) = 720 K), respectively. It changes to 1.9 years and 226 years for H and N atoms, respectively, for R = 0.5. Since the typical lifetime of a dark cloud is ∼ 10 6 years, the criterion related to the encounter desorption is often satisfied. Among the di-atomic species, H 2 is only having a faster-swapping rate (having BE 440 K, which corresponds to a hopping time scale of ∼ 1.24 × 10 −7 years and 9 × 10 −5 years, respectively with R=0.35 and R=0.5). Looking at the faster hopping rate and their abundances on the grain surface, we have extended the consideration of the encounter desorption of these species. We have considered gX + gH 2 → X + gH 2 , where X refers to H 2 , H, and N.
In Figure 6, we have shown the time evolution of the abundances of gH, gH 2 , gN, gD, and gHD with n H = 10 7 cm −3 , T = 10 K, and R = 0.35. The encounter desorption of H 2 and without the encounter desorption effect are shown to show the differences. Figure 6 depicts that the abundances of gN, gH, and gH 2 , have a reasonably high surface coverage. Since these species have a reasonable hopping rate at the low temperature, encounter desorption of these species need to be considered in the chemical model. Here, we have included the encounter desorption of these species sequentially to check their effect on the final abundances of some of the key surface species (gH 2 O, gCH 3 OH, and gNH 3 ). To check the effect of encounter desorption of the other species, we have sequentially included the encounter desorption of H 2 , H, and N. Figure 7 shows the time evolution of the encounter desorption of gH 2 O, gCH 3 OH, and gNH 3 . We have already discussed the encounter desorption of gH 2 in section 3.1. Figure 7 shows that when we have included the encounter desorption of the H atom and N atom, the time evolution of the abundances shows significant changes in abundance. It depicts that considering the effect of encounter desorption of N atom can substantially increase the abundances of gH 2 O, gCH 3 OH, and gNH 3 for the physical condition considered here (n H = 10 7 cm −3 , T = 10 K, and R = 0.35). We further have included the encounter desorption of D and HD by considering the same BE as it was obtained for H and H 2 with the H 2 substrate. The cumulative effect (by considering the encounter desorption of H, H 2 , N, D, and HD