Li-decorated Pmmn8 phase of borophene for hydrogen storage. A van der Waals corrected density-functional theory study

Li-decorated Pmmn8 phase of borophene for hydrogen storage. A van der Waals corrected density-functional theory study. A. Lebona,∗,1, R.H Aguilera-del-Torob, L.J. Gallegoc, A. Vegab aLaboratoire de Chimie et ´ Electrochimie Mol´eculaire et Analytique, UEB/UBO, UMR CNRS 6521, 29238 Brest Cedex, France bDepartamento de F´ısica Te´orica, At´omica y ´ Optica, Universidad de Valladolid, E–47011 Valladolid, Spain cDepartamento de F´ısica de la Materia Condensada, Facultad de F´ısica, Universidad de Santiago de Compostela, E–15782 Santiago de Compostela, Spain Abstract We performed standard and van der Waals-corrected density functional theory calculations to investigate the hydrogen storage capacity of a phase of borophene with Pmmn symmetry and nonzero thickness. This borophene sheet (Pmmn8) has 8 atoms in its unit cell and is more stable than the planar αsheet and that the corrugated Pmmn2 sheet (2 atoms in the unit cell). Our results show that, in pristine form, the Pmmn8 sheet is not suited for hydrogen storage applications. However, decoration with Li atoms and strain increase the hydrogen storage ability of the sheet. We performed also a detailed quantum chemical topological analysis that shows that the B-Li interaction in the hydrogenated Li-decorated Pmmn8 sheet is ionic. Our results for the adsorption of H2on the Li-decorated Pmmn8 sheet are compared with those obtained for the adsorption of H2on Ti-decorated zigzag graphene nanoribbons. Key words: Hydrogen storage on Li-decorated borophene, van der Waals-corrected DFT, quantum chemical topology analysis ∗A. Lebon Email addresses: [email protected] (A. Lebon) Preprint submitted to International Journal of Hydrogen Energy September 5, 2018 1. Introduction Hydrogen is an abundant and clean energy source that is considered to be a firm candidate to replace fossil fuels in the car industry in the future [1]. One of the main problems for that purpose is to find safe and practical means for storing hydrogen. In principle, liquid H2could be stored in cryogenic tanks at a temperature below 20 K. However, this method requires the use of large and heavy containers with strict insulation conditions to avoid boil-off losses. In recent years, a wide range of solid-state materials have been investigated to find out their hydrogen-storage ability. In particular, special attention has been devoted to carbon nanostructures, such as carbon nanotubes and graphene-based systems because they are light-weight, inexpensive and have large surface area, high thermal stability and exceptional mechanical properties (see, e.g., Ref. [2] and references cited therein). The studies performed have shown that in pristine form these nanostructures do not fulfill the H2-binding energy criterion specified by the U.S. Department of Energy (DOE) (see, e.g., Ref. [3]), a value within the range ∼0.160.26 eV, which will allows storage and discharge to be cycled at ambient conditions [4, 5]. However, it has been predicted that these nanostructures could be used as hydrogen storage media if they are doped or decorated with atoms of appropriated species such as early transition metals and alkali and alkali-earth metals, which resist cluster formation on such substrates. For instance, Yildirim and Ciraci [6] have performed density-functional-theory (DFT) calculations with the generalized gradient approximation (GGA) to show that Ti-decorated single-walled carbon nanotubes (SWCNTs) are potential high-capacity hydrogen-storage media. Similar conclusion has been obtained by Lebon et al. [7] for (hydrogen-passivated) Ti-decorated zigzag graphene nanoribbons (ZGNRs) by means of van der Waals (vdW)-corrected DFT calculations. Experimental studies on SWCNTs-Ti metal composites have confirmed that they have an efficient hydrogen uptake of up to 4.74 wt% (the DOE target is 5.5 wt%) and a hydrogen average binding energy of 0.4 eV [8]. After the first isolation of graphene [9], growing interest has been devoted to other single-layer or few-layer two-dimensional (2D) materials, observing new phenomena that occur when charge or heat transport is confined to a plane. Some recent examples are MoS2[10] and black phosphorene [11, 12], which are particularly promising for electronics and optoelectronics. The hydrogen-storage capacities of Li-decorated phosphorene have been investi2 gated by DFT calculations [13, 14]. A 2D nanostructure that has attracted great attention of late is borophene. The possibility of a boron-based 2D nanostructure was predicted some time ago, but the experimental realization in the forms of sheets on a Ag(111) substrate has been reported recently by Mannix et al. [15]. Early DFT calculations suggested that single-atomiclayer boron sheets composed of triangular and hexagonal motifs were locally stable, the more stable structure of this kind being called the αsheet [16]. However, subsequent calculations predicted two novel 2D boron phases with nonzero thickness that are considerably more stable than the αsheet [17]. One of these phases shows a Dirac cone in its electronic band structure, which renders it especially interesting from electronic perspectives. It has a buckled structure containing 8 atoms in the unit cell, being considerably more complex than other 2D materials such as graphene and phosphorene; it will be henceforth referred to as the Pmmn8 phase, or simply the βsheet. The structure of the boron monolayer obtained by Mannix et al. on the Ag(111) substrate has Pmmn symmetry with lattice constants a and b equal to 2.89 ˚ A and 5.00 ˚ A, respectively [15]. However, free standing relaxation of this structure removes the slight corrugation along the a direction, preserving the buckling along the b direction [15]. The resulting Pmmn structure of the freestanding borophene has 2 B atoms per unit cell and lattice constants a = 2.865 ˚ A and b = 1.67 ˚ A; this structure will henceforth be denoted as the Pmmn2 or the γphase. Recent calculations of the phonon spectra of the βand γphases have shown than both of them consist exclusively of real frequencies, meaning that both structures are mechanically stable [18, 19]. However, the phonon spectra of the γphase presents a low-frequency valley in the ZA branch between Γ and X, suggesting that the energetic cost of inducing a transition to the βphase is small. By contrast, the βphase has not corresponding soft modes, and it is thus a better candidate for mechanical stability in practice. Experimental measurements and theoretical calculations show that both supported and freestanding Pmmn2 borophene are highly anisotropic 2D metals [15]. On the other hand, the βphase of borophene has been predicted to have a rather large conductivity comparable to that of MoS2, and a significant in plane anisotropy comparable in order of magnitude with that of black phosphorene [18]. Additionally, recent DFT calculations have shown that boron nanoribbons constructed from the βphase are more stable and more interesting for potential technological applications than those derived from the γphase [19]. The questions arising in the context of the present paper is whether 3 borophene is a material appropriate for hydrogen storage. Early DFT investigations in this area were performed by Er et al. [20], who predicted that the Li-decorated αsheet is a good candidate for hydrogen storage on the basis of the GGA. In the fully loaded case, this nanostructure can contain up to 10.7 wt% of molecular hydrogen with an average binding energy of 0.15 eV, thus being a better material for hydrogen storage than its graphene-based counterpart. More recent DFT calculations have shown that the Li-decorated γsheet can reach an hydrogen-storage capacity of 13.7 wt%, with H2binding energies of 0.142 and 0.176 eV, as predicted using the vdW-corrected DFT functional of Grimme and the local density approximation (LDA), respectively [21]. Since, as indicated above, both the αand γsheets of borophene are much less stable than the βsheet, it is natural to inquire what is the hydrogen-storage capacity of borophene in its most stable phase, something that, to our knowledge, has not been performed yet. In the present work, we investigate the H2storage capacity of freestanding (or suspended) pristine and Li-decorated βsheets. For this purpose, we use the Perdew-Burke-Ernzerhof (PBE) form of the GGA [22] and the optB88vdW functional proposed by Klimeˇs et al. [23]. We note that in a recent paper on H2adsorption of Al-doped ZGNRs it was shown that the GGA yielded H2binding energies similar to those obtained using the optB88-vdW functional, although the GGA underestimated binding energies on pristine ZGNRs [24]. This means that for systems in which vdW forces are critical, results obtained with standard DFT may in some cases be not reliable. Apart from the investigation of the hydrogen-storage capacity of the βsheet, this paper aims at understanding the nature of the chemical bonding between the different species of the hydrogenated Li-decorated βsheet. For this purpose, we perform an extensive quantum chemical topology analysis by coupling results based on the electronic density and the electron localization function. Our results are compared with those obtained recently for Ti-decorated ZGNRs [7]. Details of the computational method used are given in Sec. 2. In Sec. 3 we present and discuss the structural and energetic results obtained for the H2uptake on pristine and Li-decorated βsheets. Our analysis of the chemical bonding is performed in Sec. 4. Finally, in Section 5, we summarize our main conclusions. 4 2. Computational Details First-principles calculations were performed with the Vienna ab initio simulation package (VASP) [25, 26], using the PBE form of the GGA for exchange and correlation [22]. Calculations employing vdW forces were performed using the optB88-vdW functional proposed by Klimeˇs et al.[23]. VASP is a plane-wave code that solves the spin-polarized Kohn-Sham equations within the projector-augmented wave (PAW) approach [27]. The following electrons were treated as valence electrons: 1s (H), 2s2p (B) and 1s2s (Li). A cut-off energy of 500 eV was used for the plane-wave basis. Brillouin zones of all the structures were integrated with k-point spacing of ∼ 0.02 ˚ A−1. Gaussian smearing was used with a smearing width of 0.01 eV. Periodic images of the sheet along the surface normal were well separated by a distance of 22 ˚ A, so that the interactions between them are negligible. Structural relaxations were carried out by employing the conjugated-gradient algorithm. The relaxations were assumed to be complete when the total force remaining on each atom and the stress tensor components acting on the cells were smaller than 0.01 eV/˚ A and 2.5 kbars, respectively. Total energies then were converged to within an energy of 10−6eV/unit cell. In our calculations on the βsheet, we considered a 2 ×2 supercell made of 32 B atoms; in what follows, we shall refer to this supercell simply as the βsheet or β-borophene. Tensile strains were tested in the x and y directions for the uniaxial strains and in both x and y directions for the biaxial strain. Prior to the investigation of the H2storage capacity of the Li-decorated βsheet, molecular hydrogen was put on the pristine βsheet, and we examined 38 input geometries. Subsequently, we decorated the β sheet with 4 Li atoms on both sides. This system was then loaded with an increasing number of H2molecules per Li atom, 1, 1.5 and 2, and in each case we examined a reasonable number of starting configurations. Specifically, for the Li-decorated βsheet with one adsorbed H2molecule per Li atom, we tested ten input geometries. The two more stable geometries obtained for this case (denoted as D1and D2) were later used to add progressively higher amounts of H2molecules. The gravimetric hydrogen densities for the Lidecorated βsheet with 1, 1.5 and 2 H2adsorbed molecules per Li atom are 3.84, 5.64 and 7.39 wt%, respectively. We also performed DFT calculations on the γsheet that showed that, after full relaxation, this sheet remains stable. However, calculations of the cohesive energies of the βand γsheets, defined as the difference between 5 the energy of an isolated B atom and the average energy per atom in the nanostructures, show that the βsheet is more stable (6.33 eV against 6.19 eV with the GGA, and 6.12 eV against 5.99 eV with the optB88-vdW functional). This result is consistent with the greater mechanical stability of the β sheet, recently demonstrated by comparing the phonon spectra of both sheets [18, 19]. Additionally, we performed PBE-GGA calculations that showed that the βsheet decorated with 4 lithium on both sides is more stable, by 2.43 eV, than the similarly decorated γsheet. The charge on atoms was evaluated according to the analysis provided by Bader’s method [28], which is based upon partitioning the structures into atomic volumes by locating the zero-flux surfaces of the electron density gradient. This method is implemented efficiently in an open source code due to the Henkelman group [29]. The evaluation of the charges is carried out with the total charge density available from the VASP output. To perform an analysis of bonding it is necessary to complete the charge analysis by a characterization with two topological approaches, the electron localization function (ELF) [30] and the Quantum Theory of Atoms In Molecules (QTAIM) [28, 31]. The quantum chemical topology part will focus on population or value assessment of two scalar fields, the ELF and the electron density. The output files of VASP enable to plot the ELF or the charge density and this is meaningful as far as valence electrons are concerned. However, a technical problem impedes the use of these grid files to analyze the population of the ELF basins or derive the position and values at the bond critical point (BCP). In fact, PAW or pseudopotentials usually introduce a smoothing of the electron density at the cut-off radius. As the analysis of these two scalar fields relies on the use of the mathematical gradient, it might add some spurious minima around a sphere whose radius is the cut-off radius [32]. Hence, it is recommended to use scalar fields produced by all-electron codes. The optimized geometries obtained with VASP were then used to perform single-point energy calculations to generate wave-function files with Gaussian09 Rev C01 [33] at the same level of theory with a 3-21G* basis set. The ELF and QTAIM calculations have been done after post-processing the wavefunction files through the TOPMOD [34] and CRITIC2 [32] programs, respectively. The ELF streamlines and specifies the chemical bond, in line with the concepts introduced by Lewis. This interpretative tool allows for a partition of space in terms of basins. The number of electrons involved in each basin (i.e. basin’s population) can be evaluated by integrating the electron density over the basin volume. In the QTAIM framework developed 6 by Bader [28], the identification of BCPs as well as a bond path between two atoms proves that the two atoms are chemically bonded. The electron density (ρ), the Laplacian of ρ(∇2ρ) and other quantities can be evaluated at the BCP and help to elucidate the nature of bonding [31, 35]. For the sake of clarity, a few words may be stated about the ELF function. The ELF was introduced by Becke and Edgecombe within the context of the Hartree-Fock calculations [30]. In the VASP code and after post-processing the output wavefunction of Gaussian09 Rev C01 with TOPMOD [34], it is expressed according to an alternative but equivalent formulation given by Silvi and Savin [36] that is valid in the DFT framework. The ELF that is a dimensionless function is an inverse and normalized local measurement of the amount of Pauli repulsion at a given point of space: the function value is close to 1 the lower the Pauli repulsion, and closer to 0 in the opposite case. This is why ELF values are associated with the spatial localization of electrons. The ELF is generally presented with isosurface plots and in solidstate physics it is used to infer the ionic or metallic character of bonding between fragments of a system. The purpose here is to go beyond this first approach and provide more details when discussing the population analysis in the QCT part. 3. Structures and energetics As indicated above, we performed calculations on the β-borophene, a 2 ×2 supercell containing 32 B atoms. Top and side views of the structure are shown in the upper panels of Fig. 1. We first investigated the adsorption of a single H2molecule on the pristine βsheet by considering 38 distinct configurations, by locating the H2molecule in different sites such as a furrow of the sheet, on the top of a B atom and at a bridge site over a B-B bond. Parallel and perpendicular initial orientations of the H2molecule with respect to the sheet were considered. The binding energy of the H2molecule was computed as Eb=E(β) + E(H2)−E(β+ H2) (1) where E(β), E(H2) and E(β+ H2) are the optimized energy of the (supercell) βsheet, that of the free H2molecule and that of the final optimized structure upon H2adsorption. Our results show that there is a general tendency of the B atoms to dissociate the hydrogen atoms of the H2molecule. The most stable configuration in the case of non-dissociative, molecular hydrogen adsorption on the βsheet is shown in the middle panels of Fig. 1. 7 Figure 1: Front and side views of the 2 ×2 supercell of the β(Pmmn8) sheet of boron in pristine form (upper panels), with a H2molecule adsorbed in the most stable nondissociative configuration (middle panels), and decorated with Li atoms (purple balls) at furrow sites on both sides of the sheet (lower panels). 8 In that position, the H2molecule is adsorbed on the βsheet with a binding energy of 8 meV with the PBE-GGA and 47 meV with the optB88-vdW functional (see Table 1, where we also show the B-H and H-H distances). However, this configuration is less stable than other in which the H atoms of the H2molecule are dissociated. The very small H2binding energy and the dissociative character of the hydrogen adsorption on the pristine βsheet make this system not suitable for hydrogen storage. We note that very weak H2binding energies have been predicted for the adsorption of H2on the α and γsheets of borophene [20, 21], without observing a dissociative tendency. After this preliminary study, we decorate the (supercell) βsheet with 4 Li atoms on both sides and relax this hybrid nanostructure. The preferential adsorption sites of the Li atoms on the βsheet are the furrow sites (see the lower panels of Fig. 1). The average binding energy of the Li atoms on the βsheet was computed using the expression Eb= [E(β) + nE(Li) −E(β+nLi)]/n (2) where E(Li) and E(β+ n Li) are the energies of a single Li atom and that of the optimized structure of the sheet upon adsorption of n = 8 Li atoms. The average binding energies of the Li atoms on the βsheet were 2.36 eV and 2.37 eV with the GGA and the optB88-vdW functionals, respectively. Therefore, the stabilizing effect of vdW forces is in this case negligible. We note that the computed binding energies of the Li atom on the βsheet are somewhat lower than the value obtained by Jiang et al. [37] for the binding energy of a Li atom on the γsheet using the GGA, 2.75 eV. The differences can be attributed to the different geometries of the βand γsheets. In fact, in the Li-decorated βsheet each Li atom has four neighboring B atoms at a distance of 2.42 ˚ A, while in the most stable configuration of the Li-decorated γsheet each Li atom has six B atoms at a distance of 2.39 ˚ A [37]. Using now as supporting material the βsheet decorated with Li atoms on both sides, we investigate the adsorption of H2molecules. Figure 2 shows some of the input geometries used in our calculations. The upper panels show two of the ten input geometries used in our calculations when one H2molecule is adsorbed per each Li atom; the middle panels show two of the input geometries used when 1.5 H2molecules were adsorbed per each Li atom, which were obtained by adding more H2molecules to the most stable configuration (D1) obtained in the case of the adsorption with a single H2molecule per Li atom; and the lower panels correspond to two input 9 Table 2: Bader charges obtained with VASP for Li-decorated β-borophene (Bfor short) and Ti-decorated ZGNR (Zfor short) with and without adsorbed H2molecules. The charge expressed as the number of the electronic charge e is defined as δQ =Q−Qval, where Qval is the valence charge. First rows of data correspond to results obtained with the optB88-vdW functional, and last rows to those obtained with the GGA. Numbers in parenthesis are standard deviations. Hpas refers to passivating H atoms. B B+1H2B+1.5H2B+2H2Z Z+4H2 B/C 0.2(1) 0.2(1) 0.2(1) 0.2(1) 0.1(2) 0.1(1) Li/Ti -0.86(1) -0.86(3) -0.86(3) -0.86(1) -0.88(2) -1.22(1) H 0.01(2) 0.01(3) 0.02(2) 0.06(4) Hpas -0.06(4) -0.07(2) B/C 0.2(1) 0.2(1) 0.2(1) 0.2(1) 0.1(2) 0.1(1) Li/Ti -0.87(1) -0.87(1) -0.87(4) -0.87(1) -0.91(2) -1.24(1) H 0.01(2) 0.01(3) 0.01(2) 0.06(4) Hpas -0.05(4) -0.08(2) Ti-decorated ZGNR. In fact, for this system the passivating H atoms are distinguished from the saturating H2molecules. The H2molecules display a trend to accept 0.05 electron charge, an effect that was explained in Ref. [7] within the framework of the Dewar interaction, where a σ∗state is partially occupied. The adsorption of H2molecules on the Ti-decorated ZGNR deeply modifies the charge balance between the Ti atom and the C atoms of the underlying graphene ribbon: Ti atoms loose more charge after the hybridation process that takes place with the H2molecules. However, in the Li-decorated βsheet, it is clear that no such a behaviour occurs when the system is loaded with H2molecules: the interaction between the Li atoms and the B atoms of the βsheet has ionic character, and the subsequent hydrogenation does not modify the ionic character of the B-Li bonding. We complete below our analysis based on the Bader charges with a topological analysis of pertinent scalar fields that allows deeper insights into the interactions that exist in the H2/Li-decorated β-borophene systems. The electron localization function (ELF) is a scalar field derived from the wave function that provides a three dimensional picture of the formation of pairing electrons. We have carried out calculations for the Li-decorated β sheet with n H2(n = 1, 1.5, 2) adsorbed molecules per Li atom. Figure 4 shows the one-electron density and the ELF for the Li-decorated βsheet 16 with 1 H2molecule per Li atom; the isosurfaces of the one-electron density are set at ρ= 0.05 e and for the ELF at η= 0.8. The upper panel of Fig. 4 corroborates the ionic character of the B-Li interaction according to the repartition of the isosurface of the one-electron density. In fact, the B atoms are gathered within an electronic isosurface that illustrates the metallic character of the borophene layer. The Li atoms seem deposited on top of the layer on both sides with no electronic connection to the underlying borophene surface. H2molecules are then in equilibrium above the Li atoms. As far as the lower panel of Fig. 4 is concerned, no basins of the ELF connects the B atoms with the Li atoms. Valence basins connect the B atoms of the layer together. The connection between B atoms illustrates the metallic character of the borophene layer. On the other hand, the absence of valence basin between boron and lithium atoms is another signature of the ionic character of the B-Li bond, where the charge of the Li atom flows to the borophene layer. The interaction between the Li atoms is also ionic as no basin connects two Li atoms on one side of the βsheet, preventing them from clustering. The relaxed structures obtained with VASP were employed as input structures for Gaussian calculations. The output wave functions enabled the calculation of the ELF structure with the help of the TOPMOD programs [34]. We remind that potentials such as PAW or pseudopotentials introduce a smoothing of the density in the vicinity of a sphere whose radius equals the cut-off radius. The construction feature of the PAW leads to a bad assessment of the position of the basins and may even add spurious basins. The handling of the effects of the pseudopotential transformation stems from the missing core electrons. The valence pseudo-density is mostly coincident with the all-electron density in the bonding and interstitial regions, but smoother close to the nuclei, particularly for elements with few valence electrons or strongly cationic atoms. In these cases, maxima of the electron density may appear at positions different from the nucleus positions [32]. With the ELF derived from the Gaussian09, we are able to count the number of electrons within the basins. There exist three types of basins: bonding and non-bonding valence basins, which correspond to bonds and lone pairs, respectively, and core basins, which are attached to a nucleus. Attractors are special points in basins where the value of the ELF reaches a maximum. Basins are defined around attractors and the gradient of the ELF within a basin converges toward the attractor. The core basins are organized around nuclei (with Z ≥2) and the valence basins occupy the remaining space. The core basins are labeled by C(X), X being the atomic 17 Figure 4: Isosurface plots for Li-decorated β-borophene loaded with 1 H2molecule per Li atom. The isosurface values are ρ=0.05 for the one-electron density (top) and η= 0.8 for the electronic localization function (bottom). 18 Table 3: Population of the ELF basins in number of electronic charge for β-borophene and ZGNR in pristine forms, and for Li-decorated β-borophene (Bfor short) and Ti-decorated ZGNR (Zfor short) without and with adsorbed H2molecules. Numbers in parenthesis are standard deviations. C(B) C(Li) V(B,Li) V(Li,H) V(H,H) V(H,H,H) β2.06(0) B2.06(1) 2.03 0.00 B+1H22.09(1) 2.03 0.00 0.02 2.02 B+1.5H22.09(1) 2.03 0.00 0.03(1) 2.03(5) 0.06(3) B+2H22.09(1) 2.03 0.00 0.03(3) 2.0(1) 0.11(3) C(C) C(Ti) V(C,Ti) V(Ti,H) V(H,H) V(H,H,H) ZGNR 2.08(1) Z2.12(2) 19.9(2) 0.7(1) Z+4H22.12(2) 19.92(3) 0.8(2) 0.5(3) 1.96(3) 0.4(1) label of the atom to which the core belongs. The core basins closely match the inner shell structure of the atoms. A valence basin is characterized by its synaptic order, which is the number of core basins with which it shares a common boundary. Monosynaptic basins (labeled V(A)) usually correspond to lone pair regions, whereas disynaptic and polysynaptic basins (labeled V(A,B,C,...)) characterize the covalent bonds. Note that there is no need to associate synapticity to a core basin. Overall, the spatial distribution of the valence basins closely matches the non-bonding and bonding basins of the valence shell electron pair repulsion (VSEPR) model. To substantiate the previous analysis that leads us to infer the ionic character of the B-Li interaction, the core basins populations of the Li and B atoms were determined for the βsheet, the Li-decorated βsheet, and the Li-decorated βsheet with different amounts of H2molecules per Li atom. Calculations were also performed for the Ti-decorated ZGNR and the H2/Tidecorated ZGNR. The results listed in Table 3 indicate that the population of the core basin on Li is approximately 2.03 electrons, which is consistent with the fact that only 1s2electrons are bond to the Li nucleus. The Li atoms adopt therefore a cationic Li+oxidation state. The Li atoms share a valence basin with hydrogen atom V(Li, H) but it is noted that its population is small as it amounts at most to 0.03 electrons. The B atoms present more or less the same population for the core basins C(B) than the Li atoms, which is 19 2.06 to 2.09 electrons. In contrast to the Li atoms, disynaptic and trisynaptic basins, V(B, B) and V(B, B, B), connect the B atoms together and ensure the metallic character of the borophene sheet. Hence, there exists a degree of delocalization corresponding to a metallic state. For the H2/Li-decorated βsheet, H atoms within the H2molecules are connected by a bond whose population amounts to about 2 electrons, which is typical of a covalent bond. The ELF picture for the H2/Ti-decorated ZGNR is at odd with that of the H2/Li-decorated βsheet. In fact, Table 3 illustrates an exchange of electron in the Ti-C bonding ( V(C,Ti) ∈[0.6,1.0] e), whereas the Li-B bonding stems from an electrostatic interaction (V(B,Li) = 0.0 e). The analysis of the ELF basins population suggests a noticeable degree of Ti-H hybridization as well, since 0.5 electron on average are involved in the bonding between the H2 molecules and the Ti atom, as it is shown for the V(Ti,H) basin in the fourth column of Table 3. The V(Ti,H) population exceeds the V(Li,H) population by more than one order of magnitude. This quantum topological quantity confirms the existence and the efficiency of the σ∗hybridization termed as the Dewar interaction for the H2saturated Ti-decorated ZGNR [7]. This last finding can obviously be related to the H-H separation difference in the Li-decorated βsheet (0.76˚ A) and in the Ti-decorated ZGNR (0.82˚ A). In addition to that, the partial charge transfer implied by the σ∗hybridization partially stretches the H-H bond giving way to a length increase. The almost negligible V(Li, H) population of 0.03 electrons at most can be put in parallel with the tiny elongation of molecular hydrogen in Li-decorated β-borophene for which H-H separation increases by hardly 0.01 ˚ A. Table 3 provides a quantitative analysis of the population of the V(H,H) basins. For both boron and carbon-based systems, the population amounts to nearly 2 electrons. This feature indicates that, when attached to Li or Ti atoms, the H2molecules are not dissociated. At last according to the analysis of the ELF function, we note the presence of V(H,H,H) basins. These trisynaptic basins connect together two hydrogen atoms of a H2molecule with one hydrogen atom of another H2molecule. The population related to this kind of basins is small for the H2loaded Li-decorated βsheet. However, it increases with the number of H2molecules from 0.0 electrons at n = 1 to 0.11 electrons at n = 2 via 0.06 electrons at n = 1.5 (last line of Table 3). Nevertheless, for the H2saturated Ti-decorated ZGNR the average V(H,H,H) population reaches 0.4 electrons. Thus, the Ti atom favours interactions between H2molecules. Although of weak intensity, this feature is observed with the Li-decorated βsheet. The difference in V(H,H,H) populations can be related with the 20 Table 4: Values of the electron density and its Laplacian at the BCPs for the Li-decorated βsheet (Bfor short) loaded with H2molecules and the Ti-decorated ZGNR (Zfor short) saturated with 4 H2molecules. All quantities are in atomic units. Numbers in parenthesis are standard deviations. Field H-H B-Li B-H Li-H H2-H2 B+ 1H2 ρ0.234(1) 0.018(2) 0.0058(2) 0.0066(6) 0.0019(1) ∆ρ-0.64(1) 0.06(5) 0.010(2) 0.04(1) 0.0056(2) B+ 1.5H2 ρ0.233(1) 0.017(2) 0.005(1) 0.006(1) 0.002(1) ∆ρ-0.65(3) 0.11(4) 0.009(3) 0.04(1) 0.007(2) B+ 2H2 ρ0.235(2) 0.018(3) 0.005(1) 0.0071(6) 0.002(1) ∆ρ-0.65(3) 0.10(2) 0.008(3) 0.052(5) 0.007(3) Field H-H C-Ti C-H Ti-H H2-H2 Z+ 4H2 ρ0.243(3) 0.037(2) 0.0064(3) 0.036(1) 0.0017(5) 0.21(1) ∆ρ-1.0(6) 0.10(1) 0.022(4) 0.02(1) 0.004(2) -0.6(1) difference of dext (last line of Table 1) in both systems. The analysis of the one-electron density allows to get the BCPs. The value of the one-electron density at these specific points and the value of the Laplacian at these BCPs give information about the type of bonding between the chemical species. Some of the most common methods of QTAIM [31] have been implemented within the CRITIC2 program [32]. The method was applied here to the Li-decorated βsheet with different amounts of H2of coverage (1, 1.5 and 2 H2per Li atom). As previously, the obtained values were compared to those of the H2saturated Ti-decorated ZGNR. All the relevant data are gathered in Table 4. To better illustrate our results, Figure 5 displays three views of the positions of the BCPs, namely, a top and a side view, together with a zoom on a region around a Li atom marked with an asterisk. The small pink balls in the figure indicate the position of the BCPs. It is evident that within the D1hydrogenated system selected for display - the one with 245(131) meV binding energy per Li atom - there exist different kinds of chemical bonding. The BCPs can be classified into covalent bonds, ionic bonds and non covalent bonds. In the first category enter the H-H bonds of the hydrogen molecule, with an average number of 0.238(4) atomic unit at the BCP and a negative value of the Laplacian of -0.6. The 21 amount of density at the BCPs and the value of the Laplacian are suggestive of a covalent bond [35]. Even though it is not listed in Table 4, the B-B bond in the βsheet can be considered as well as covalent with an average number of 0.14(2) atomic unit at the BCPs and a negative Laplacian of -0.2. Nevertheless, this B-B bond seems soft in comparison to the C-C bond where 0.27 atomic unit are found at the BCPs and a large negative Laplacian at the BCP of -0.7 is noted. This point might be related to the softness of the boron structure in comparison to the carbon structure. This is illustrated in our study by the non-negligible distortion of the Li-decorated βsheet. The softness of the βsheet agrees with the vision exposed by De et al. [38], who proposed that the energy landscape of B clusters is glass-like. Four types of BCPs are isolated through the analysis of the electronic density. They correspond to B-Li, Li-H, B-H and H2-H2bonds following the decreasing order of intensity at the BCP. For the three aforementioned BCPs, the respective values are 0.018(2), 0.007(1) and 0.0055(4), and the positive values of the Laplacian at the BCP discard the covalent nature of these bondings. The B-Li bond has an evident ionic character if one refers to the 0.1 positive value of the Laplacian at the BCP. B and Li atoms sit some 2.3 ˚ A from each other. The value of the Laplacian at the BCP for the Li-H bond is approximately 0.04 and cannot be regarded as a hydrogen bond since the BCP connects the Li atom to the center of the H2molecule instead of the directional bonding expected for a hydrogen bond. The Li-H distance amounts to about 2.2 ˚ A. It is very likely that the electron cloud of the hydrogen molecule feels the interaction of the neighbouring Li+ion. The B-H bond is by far the weakest, as the value of the Laplacian at the BCP (0.01) suggests. B-H bonds are long-range interactions that cannot be neglected as they seem to contribute to the stabilization of the H2molecules. The inclination of the molecules that are located close to the Li atom might be accounted for this long-range interaction. As already pointed out, the boron sheet has a metallic character. The electron charge mobility within the boron sheet might provoke this long-range interaction with the electrons that are shared by the two hydrogen atoms within the molecule. The H2 molecule has a trend to orient one of its H atoms towards a B atom. This feature suggests the existence of a weak hydrogen bond. This is illustrated on the lower panel of Fig. 5, where one of the H atoms of the H2molecule points towards a B atom. Finally, we note the existence of BCPs between two H2 molecules, for which the BCPs is around 0.002 atomic unit and the Laplacian at these BCPs is the smallest one reported in Table 4. This interaction 22 Figure 5: Top and side views of the positions of the BCPs for Li-decorated β-borophene with one H2molecule per Li atom. Zoom around one Li atom (bottom). 23 can be viewed as a vdW interaction since the H2molecules are connected within a network where the molecules are separated on average by 3 ˚ A. This last point is an argument for incorporating long-range forces in simulation studies on the stability of H2molecules on possible storage media. Table I confirms the stabilization effect of the vdW forces since a few tens of meV are gained for the binding energies of the H2molecules as compared to the PBE calculations. If we now turn to the H2saturated Ti-decorated ZGNR, all the interactions are more pronounced. For instance, the value at the BCPs for the Ti-H interaction are 3 to 5 times larger than the Li-H interaction in the Li-decorated βsheet. The possibility of the occupancy of a σ∗orbital termed as Dewar interaction finds some additional support with the present analysis. A possibility that is ruled out by the cationic state of the Li atom sitting at the furrow of the βsheet. 5. Conclusions We have performed an extensive DFT study of the H2storage capacity of a recently proposed phase of borophene with Pmmn symmetry and nonzero thickness, the βsheet [17]. This sheet, generated by two inter-penetrating lattices, has 8 atoms in its cell unit, and is most stable than the αsheet [16], composed of triangular and hexagonal motifs, and that the Pmmn2 sheet, with only 2 atoms in the unit cell, the γsheet [15]. Our calculations were performed using both standard DFT, with the PBE form of the GGA [22], and vdW-corrected DFT, using the optB88-vdW functional [23]. Our results show that, in pristine form, the βsheet is not suitable for H2storage: the binding energy of molecular hydrogen on that substrate is too weak, and the B atoms tend to dissociate the hydrogen atoms of the H2molecule. Then, we decorated the βsheet with Li atoms on both sides. The Li atoms do not form clusters on the borophene surface and prevent the dissociation of adsorbed H2molecules. We performed calculations for the Li-decorated βsheet with 1, 1.5 and 2 adsorbed H2molecules per Li atom. In the first case, we found H2binding energies that do satisfy the criterion specified by the DOE, but the hydrogen gravimetric density is 3.84 wt%, somewhat lower than the DOE target 5.5 wt%. However, in the other cases we found H2binding energies which can be considered as acceptable for hydrogen storage since they are similar to those obtained for Li-decorated αand γsheets, which have been recently proposed as potential media for H2storage [20, 21], and the hydrogen gravimetric densities are greater than the threshold value of 5.5 wt%. We 24 also showed that strain can be a way to increase the H2storage capacity of the Li-decorated βsheet. On the other hand, our results show that the inclusion of vdW forces is essential to appropriately describe the adsorption of H2molecules on the Li-decorated βsheet. In order to have information of the interactions in the H2/Li-decorated β sheet, we have performed a detailed quantum chemical topological analysis using two scalar fields, the ELF and the one-electron density. The results of this analysis show that the B-Li interaction in the H2/Li-decorated βsheet is ionic, in strong contrast with the covalent character of the C-Ti interactions in H2saturated Ti-decorated ZGNRs, which have been recently proposed as possible templates for H2storage [7], and that most of the chemical interactions in that nanostructure are of weak intensity. 6. Acknowledgments This work was supported by the Spanish Ministry of Science and Innovation (Projects No. FIS2014-59279-P), the Junta de Castilla y Le´on (Project VA124G18) and the Xunta de Galicia (AGRUP2015/11 and GRC ED431C). These research projects were partially supported by FEDER. R.A. acknowledges the financial support and the kind hospitality from the Universit de Brest (UBO), France. Facilities provided by the Pole de Calcul Intensif pour la Mer (DATARMOR, Brest) are also acknowledged. We also thank Amador Garc´ıa-Fuente and Jes´us Carrete Monta˜na for useful information on the Pmmn8 structure of borophene. 25