Density functional theory calculations on alkali and the alkaline Ca atoms adsorbed on graphene monolayers

Applied Surface Science 413 (2017) 197–208

Contents lists available at ScienceDirect

Applied Surface Science journal homepage: www.elsevier.com/locate/apsusc

Full Length Article

Density functional theory calculations on alkali and the alkaline Ca atoms adsorbed on graphene monolayers Nicholas Dimakis a,∗ , Danielle Valdez a , Fernando Antonio Flor a , Andres Salgado a , Kolade Adjibi a , Sarah Vargas b , Justin Saenz b a b

Department of Physics, University of Texas Rio Grande Valley, Edinburg, TX, United States Robert Vela High School, Edinburg, TX, United States

a r t i c l e

i n f o

Article history: Received 8 December 2016 Received in revised form 24 March 2017 Accepted 2 April 2017 Available online 3 April 2017 Keywords: Density functional theory Graphene Alkali/alkaline adsorption Overlap population Adatom coverage

a b s t r a c t The adsorption of the alkali Li, K, and Na and the alkaline Ca on graphene is studied using periodic density functional theory (DFT) under various adatom coverages. The charge transfers between the adatom and the graphene sheet and the almost unchanged densities-of-states spectra in the energy region near and below the Fermi level support an ionic bond pattern between the adatom and the graphene atoms. However, the presence of small orbital overlap between the metal and the nearest graphene atom is indicative of small covalent bonding. Van der Waals interactions are examined through a semiempirical correction in the DFT functional and by comparing adatom-graphene calculations between 3% and 1.4% adatom coverages. Optimized adatom-graphene geometries identify the preferred adatom sites, whereas the adatom-graphene strength is correlated with the adsorption energy and the adatom distance from the graphene plane. Calculated electronic properties and structural parameters are obtained using hybrid functionals and a generalized gradient approximation functional paired with basis sets of various sizes. We found that due to long range electrostatic forces between the alkali/alkaline adatoms and the graphene monolayer, the adatom-graphene structural and electronic properties could be well-described by specific DFT functionals paired with high-quality adatom basis sets. For Li, K, and Na adsorbed on graphene, increased adatom surface coverage weakens the adatom-graphene interaction. However, this statement does not apply for Ca adsorbed on graphene. In this case, the Ca adsorption strength, which is stronger at higher coverages, is opposite to increases in the Ca–4s orbital population. © 2017 Elsevier B.V. All rights reserved.

1. Introduction Graphene and carbon-based nanostructures tubes have been extensively studied both experimentally and computationally due to their wide-ranging applications in microelectronics, hydrogen storage, and chemical sensors [1–3]. Graphene has a high mechanical strength and it’s high thermal and electrical conductivity [4–6] made it a favorable choice to serve as flexible electronics for cell phone panel displays and televisions. Graphene emerged through the accidental discovery by Novoselov et al. [7]. Its existence seemed to contradict Mermin-Wagner theorem, which states that perfect 2D crystals are not thermodynamically stable, and therefore do not exist [8]. However, this is not the case for graphene due to the presence of microscopic buckling [9] and intrinsic ripples

∗ Corresponding author. E-mail addresses: [email protected], [email protected] (N. Dimakis). http://dx.doi.org/10.1016/j.apsusc.2017.04.010 0169-4332/© 2017 Elsevier B.V. All rights reserved.

[10] in the monolayer structure that distort its surface symmetry. Graphene structure is hexagonal with carbon atoms forming inplane ␴-bonds via sp2 hybridization, where the out-of-plane carbon pz orbitals contribute to the formation of ␲-bonding and ␲*- antibonding orbitals, the latter being above the Fermi level and thus, unoccupied. Graphene is a zero band-gap semiconductor as the ␲ and ␲* orbitals coincide at the Fermi level [11]. Alkali and alkaline adatom adsorption on graphene could be employed to improve graphene potential use as hydrogen storage device, field effect transistor, and superconducting material. Graphene properties can be tuned by controlling adsorption. Chan et al. [12] studied alkali, alkaline, and transition metal adsorption, including Li, Na, K, and Ca on graphene monolayers, using spin-polarized plane-wave periodic density functional theory (DFT) [13–15], under the PBE generalized gradient approximation (GGA) [16]. Similar calculations were performed by Nakada and Ishii using local density approximation [11] and reported migration energies (i.e., energy required for the adatom to move to a different adsorption site). These authors found that migration energies, which are

198

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

related to the stability of the adatom-graphene systems at room temperature, are the highest for Li on graphene (0.3 eV), whereas these values are much smaller for adsorption of Na, K, and Ca on graphene (0.07–0.013 eV). Var der Waals interactions for singleand double-sided adsorption of eight metals on graphene including Li, Na, and Ca has been studied by Wong et al., [17] in search for a metal-graphene system used as hydrogen storage device. The authors reported plane-wave PBE and vdW-DF2 calculations [18], the latter includes corrections to the PBE functional to account for the adatom-graphene var der Waals interactions. For Li/graphene, Wong et al. reported significant changes in the Li adsorption energies (Eads ) and adatom height h (i.e., the distance between the adatom and the graphene plane) between the two methods for Li adsorption at 16.67% adatom coverage. However, in this work, the inclusion of the correction terms to account for van der Waals interactions for Li/graphene do not significantly shift Eads and h values (vide infra). Calcium adsorption on graphene has been studied by Cazorla et al. [19], where it was found that Ca Eads values are very different between calculations using molecular and periodic DFT and concluded that Ca-C24 H12 is not equivalent to Ca-decorated graphene nanostructures. The adsorption of Ca on graphene is extensively studied in our current work and we found that the Cagraphene interaction critically depends on the method employed (vide infra). Li adatom nucleation has been studied by Zhou et al. [20] due to its applications in hydrogen storage and Li-ion batteries and was found that Li prefer to cluster in a four atom rhombic formation during adsorption on graphene. Moon et al., [21] studied Na adsorption on graphene and graphene oxides using DFT for use in Na-ion based batteries as an alternative to Li-ion batteries. The adsorption of Li, K, Na, and Ca on graphene is examined using periodic DFT under Gaussian-type basis sets centered at the atoms. This method is an alternative to periodic DFT approaches under plane-wave basis sets, which use pseudopotentials. Ulian et al., in their report on modeling talc (i.e., Mg3 Si4 O10 (OH)2 ), compare these two approaches and found that both provide structural and vibrational properties of similar accuracy [22]. When hybrid DFT functionals are used, periodic DFT with Gaussian basis sets require significantly less CPU time relative to plane-wave approaches. The quality of the basis sets affects both calculations. For plane-wave methods, this quality depends on a single parameter, the electronic kinetic energy cutoff, whereas for Gaussian-type approaches it depends on the number of shells used as well as the inclusion of very diffuse functions [22]. In this work, all atoms are described by all-electron basis sets. More specifically, the graphene carbon atoms are described by the contracted triple-␨ (4s3p1d) basis set [23], which is optimized for crystalline calculations, as in our past calculations for Fe/graphene [24]. The Li, K, Na, and Ca adatoms are described by basis sets of various sizes. For example, Li is described by the (4s1p) and (6s3p1d) basis sets. Here, we study the adatom coverage effect on the graphene surface, which affects the electronic and structural properties of the adatom-graphene systems. The adatom-graphene interaction, examined using two hybrid DFT functionals (B3LYP [25] and PBE0 [26,27]) and the GGA PBE functional basis sets of various sizes, is analyzed using changes in the Eads , the adatom height h, the density-of-states (DOS) spectra, and the charge transfers between the adatom and the graphene monolayer. Contrary to past reports, we thoroughly examine the effect of the basis set and the functional employed on the adatom-graphene interaction, with and without corrections that account for van der Waals interactions. Here, we found that the hybrid functionals do not correctly describe the Cagraphene interaction at low coverages (3% and less). For Li, Na, and K, we correlate the adatom-graphene strengths with the charge transferred from the adatom to graphene, whereas the Ca-graphene

Fig. 1. The a) (2 × 2), b) (4 × 4), and c) (6 × 6) graphene sheets. B and H stand for adsorption at the bridge and hcp sites, respectively.

strength, which appears stronger at high coverages, is described by changes in the Ca–4s orbital population (vide infra). 2. Computational methods Monolayer graphene is modeled as two-dimensional 2 × 2, 4 × 4, and 6 × 6 hexagonal lattices with 8, 32, and 72 carbon atoms, respectively. These configurations correspond to about 12.5% 3%, 1.4% adatom coverages, respectively. Fig. 1 shows the 2 × 2, 4 × 4, and 6 × 6 graphene sheets used here, as well as the observed adsorption sites. The alkali atoms of Li, K, and Na as well as the alkaline atom Ca are placed as adatoms to graphene and are free to move during the geometry optimizations of the adatom-graphene systems. Similar to our recent work for metal adsorption on graphene [28], the optimized geometries and the electronic properties of the clean graphene and adatom-graphene systems are obtained using

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

199

Table 1 DFT calculated structural properties (adatom height h and adatom-nearest graphene C distance) and Eads for Li adsorbed on 2 × 2, 4 × 4, and 6 × 6 graphene sheets. All Li adatoms are adsorbed at the hcp site, whereas Li is adsorbed at the bridge site on the 2 × 2 graphene sheet under the B3LYP/(4s1p) calculations. Number in parenthesis refer to past DFT calculations from the literature. Uncorrected Eads values are in square brackets. Adatom

Supercell

Functional/Adatom Basis Set

h (Å)

Metal-C (Å)

Eads (eV)

Li

2×2

B3LYP/(4s1p)a /(6s3p1d) PBE0/(4s1p) /(6s3p1d) PBE/(4s1p) /(6s3p1d) B3LYP/(4s1p) /(6s3p1d) PBE0/(4s1p) /(6s3p1d) /(6s3p1d)b PBE/(4s1p) /(6s3p1d) /(6s3p1d)b PBE0/(6s3p1d)

2.00 1.78 1.68 1.77 1.69 1.75 1.67 1.63 1.67 1.61 1.68 1.66 1.61 1.69 1.64 (5.99, 1.91b [17]) (1.62 [11]) (1.71 [12]) (1.72c , 1.78d [47])

2.15 2.28 2.21 2.28 2.22 2.26 2.20 2.18 2.19 2.17 2.22 2.20 2.17 2.24 2.18

0 [−0.80] −0.26 [−0.78] −0.40 [−1.28] −0.32 [−0.79] −0.50 [−1.34] −0.46 [−0.94] −1.22 [−2.06] −0.98 [−1.65] −1.59 [−2.38] −1.41 [−1.96] −1.65 [−2.23] −1.41 [−2.16] −1.20 [−1.35] −1.54 [−2.15] −1.50 [−2.08] (−0.43, 0.17b [17]) (−1.4 [11]) (−1.096 [12]) (−0.27c , −0.54d [47])

4×4

6×6 1×1 3×3 4×4 Mol. Clusters a b c d

(2.23 [12])

Bridge site adsorption. Grimme correction was used [32]. LiC24 H12 . LiC92 H24 were used.

periodic restricted DFT under the CRYSTAL09 [29] program. CRYSTAL09 computes electronic, structural, and vibrational properties of periodic, polymer, and molecular systems using Gaussian-type basis sets centered at the atoms. All of our calculations have been repeated using unrestricted DFT at low-spin configurations for the adatom-graphene systems. However, no changes were observed in the calculated geometries between these two methods. In this work, we use the B3LYP hybrid semiempirical functional, the PBE0 non-empirical parameter-free functional, and the GGA PBE functional. Here, the B3LYP functional has the same exact and exchange functionals as the original B3LYP functional, whereas the VWN correlation functional is replaced by VWN5 [30]. We must state that the parameter-free PBE0 functional is superior to the semiempirical B3LYP functional for the periodic DFT calculations of solids [31]. Paier et al., reported B3LYP calculated lattice parameters and atomization energies to be overestimated by 1% and underestimated by about 17%, respectively relative to experimental measurements. Therefore, results using the PBE0 and PBE functionals are preferred over corresponding calculations using the B3LYP functional. Some calculations were repeated using the semiempirical correction by Grimme, which improves the description of the DFT functionals for long-range electron correlations that are responsible for van der Waals interactions between the adatom and the graphene monolayer [32]. The graphene lattice parameters obtained by using B3LYP and PBE0 functionals are 2.45 Å and 2.44 Å, respectively. These values are very close to the experimentally obtained graphite lattice parameter of 2.46 Å [33]. Brillouin zone integrations (MonkhorstPack grid) [34], the Fermi energy, and the density matrix calculations (Gilat grid) [35,36] were performed on a 24 × 24 grid. For the systems of this work, a smaller grid of 12 × 12 did not affect the obtained electronic and structural properties. For convergence purposes, the Fermi surface was smeared with a Gaussian of 0.005 Hartrees. Moreover, the SCF energy convergence was achieved by using Anderson quadratic mixing [37], coupled with additional mixing of the occupied with the virtual orbitals. The SCF energy threshold value for our calculations was set at 10−12 Hartrees and 10−9 Hartrees for clean graphene and the adatom-graphene systems, respectively (default value is 10−7 Hartrees). A large integration grid was used (XLGRID keyword): This is a pruned grid with

75 radial and 974 angular points. Densities-of states (DOS) spectra and crystal orbital overlap populations (COOP) for selected pair of atoms [38] were calculated directly by CRYSTAL09. The XcrySDen graphical package was used for the charge density plots of this work [39]. Charge transfers and adatom orbital populations were calculated using Mulliken population analysis [40]. For comparison reasons, these properties are also obtained via integration of the DOS spectrum calculated from CRYSTAL09 by considering the corresponding bands within the appropriate energy regions, as in our past report [28]. Davidson and Chakravorty reported that the absolute magnitude of the calculated charges through Mulliken population analysis or similar methods have little physical meaning since these calculations may strongly depend on the choice of the basis sets used [41]. However, Segall et al., using periodic planewave DFT on several bulk crystals, showed the importance of the relative charge values derived from the Mulliken population by projecting planes waves onto linear combination of atomic orbitals, which were used to measure the covalency of a bond [42]. The stability of the final geometry conformations was secured via post-geometry optimizations of the final structure (FINALRUN keyword set to 4). The Eads is defined as Eads = EX/graphene − EX − Egraphene , where X is the adatom is a measure of the adatom interaction with graphene. The adatom forms a bond with graphene if Eads < 0. When Gaussian basis sets are used, Eads calculations suffer from basis set superposition errors (BSSE) [43]. The BSSE errors are minimized using either large basis sets or employing counterpoise corrections using “ghost” massless atoms in the fragment energy calculations of the adatom-graphene structures. Calculated Eads values using the counterpoise correction are less negative in energy than corresponding uncorrected values. Hoverer, Eads overcorrection may occur, which leads to less negative or even positive Eads . Here, we report both corrected and uncorrected Eads values. 3. Results and discussion 3.1. Adatom-graphene geometries and Eads Tables 1–3 show the adatom height h, the distance between the adatom and the nearest carbon atom in the graphene mono-

200

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

Table 2 DFT calculated structural properties (adatom height h and adatom-nearest graphene C distance) and Eads for Na, K, adsorbed on 2 × 2, 4 × 4, and 6 × 6 graphene sheets. All adatoms are adsorbed at the hcp site. Number in parenthesis refer to past DFT calculations from the literature. Uncorrected Eads values are in square brackets. Adatom

Supercell

Functional/Adatom Basis Set

h (Å)

Metal-C (Å)

Eads (eV)

Na

2×2

B3LYP/(5s3p1d) /(9s5p1d) PBE0/(9s5p1d) PBE/(9s5p1d) B3LYP/(5s3p1d) /(9s5p1d) PBE0/(9s5p1d) /(9s5p1d)a PBE/(9s5p1d) PBE0/(9s5p1d)

2.34 2.58 2.53 2.50 2.13 2.25 2.19 2.19 2.23 2.09 (2.70, 2.81a [17]) (2.22 [11]) (2.28 [12], 2.285 [21])

2.74 2.95 2.90 2.88 2.56 2.69 2.63 2.64 2.67 2.58

0.13 [−0.69] 0.04 [−0.39] −0.03 [−0.43] −0.09 [−0.50] −0.52 [−1.60] −0.28 [−0.98] −0.58 [−1.24] −0.94 [−1.60] −0.48 [−1.13] −0.79 [−1.45] (−0.01, 0.04a [17]) (−0.7 [11]) (−0.462 [12], −0.504 [21])

B3LYP/(6s4p1d) /(6s6p3d1f) PBE0/(6s6p3d1f) PBE/(6s6p3d1f) B3LYP/(6s4p1d) /(6s6p3d1f) PBE0/(6s6p3d1f) /(6s6p3d1f)a PBE/(6s6p3d1f) PBE0/(6s6p3d1f)

2.78 2.54 2.51 2.51 2.47 2.36 2.33 2.33 2.32 2.36 (2.58 [11]) (2.60 [12], 2.82 [48])

3.11 2.92 2.89 2.89 2.88 2.81 2.79 2.80 2.79 2.83

4×4

6×6 1×1 3×3 4×4 2×2

K

4×4

6×6 3×3 4×4 a

(2.70 [12])

(2.99 [12])

0.07 [−0.78] −0.17 [−1.08] −0.23 [−1.13] −0.32 [−1.18] −1.06 [−2.33] −1.83 [−2.69] −2.36 [−3.01] −2.70 [−3.29] −2.09 [−2.90] −2.73 [−3.34] (−0.8 [11]) (−0.802 [12], −0.44 [48])

Grimme correction was used [32].

Table 3 DFT calculated structural properties (adatom height h and adatom-nearest graphene C distance) and Eads for Ca adsorbed on 2 × 2, 4 × 4, and 6 × 6 graphene sheets. All adatoms are adsorbed at the hcp site. Number in parenthesis refer to past DFT calculations from the literature. Uncorrected Eads values are in square brackets. Adatom

Supercell

Functional/Adatom Basis Set

h (Å)

Metal-C (Å)

Eads (eV)

Ca

2×2

B3LYP/(6s4p1d) /(9s5p3d) PBE0/(9s5p3d) PBE/(9s5p3d) B3LYP/(6s4p1d) /(9s5p3d) PBE0/(6s4p1d) /(9s5p3d) /(9s5p3d)a PBE/(9s5p3d) /(9s5p3d)a PBE0/(9s5p3d) PBE/(9s5p3d)

2.32 2.86 2.60 2.51 2.14 3.90 2.02 3.74 2.78 2.55 2.58 2.89 2.49 (2.33, 2.49a [17]) (2.14 [11]) (2.13 [12], 2.26 [19])

2.72 3.19 2.96 2.89 2.61 4.16 2.53 4.01 3.13 2.94 2.97 2.23 2.90

1.09 [−1.67] 0.31 [−0.19] 0.20 [−0.37] 0.06 [−0.54] 0.87 [−2.62] 0.15 [−0.13] 0.45 [−3.15] 0.12 [−0.16] 0.22 [−0.34] 0.27 [−0.36] −0.04 [−0.66] 0.41 [−0.13] 0.21 [−0.43] (0.52, 0.50a [17]) (−0.5 [11]) (−0.632 [12], −0.443 [19])

4×4

6×6 1×1 3×3 4×4 a

(2.72 [12])

Grimme correction was used [32].

layer, and the BSSE-corrected and uncorrected Eads for Li, Na and K, and Ca, respectively, at various DFT functionals and adatom basis sets for 1.4%, 3%, and 12.5% adatom coverages. We identify possible van der Waals adatom-adatom interactions by comparing the adatom-graphene geometries and Eads calculated at the 3% adatom coverage with corresponding properties using 1) the 1.4% adatom coverage (6 × 6 graphene sheet) and 2) the Grimme correction in the DFT functional [32]. The latter is also used for identification of van der Waals interactions between the adatom and the graphene monolayer. For all cases, with the exception of Li/graphene under the B3LYP functional paired with the Li (4s1p) basis set (i.e., B3LYP/(4s1p)), the adatom preferred site is the hcp, in agreement with past reports [11,12,17,19,21,47]. The Li/graphene optimized geometries (i.e., Li-graphene heights h and Li-nearest C distances) and Eads are obtained using the B3LYP, PBE0, and PBE functionals paired with the Li (4s1p) and (6s3p1d) basis sets at 3% and 12.5% Li coverages by employing all possible functionals and basis sets combinations. Table 1 shows that the Li (6s3p1d) basis set

provides similar Li/graphene geometries per 3% and 12.5% Li coverages for all functionals used here. However, the B3LYP/(6s3p1d) BSSE-corrected Eads are less negative relative to corresponding values from the PBE0 and PBE calculations, which is indicative of weaker Li-graphene interaction. The use of the smaller Li (4s1p) basis set shows mixed results in the Li/graphene geometries and Eads , when compared with corresponding calculations with the larger Li (6s3p1d) basis set. For example, at 12.5% coverage, the B3LYP/(4s1p) calculations describe Li-graphene as a non-bonding interaction (Eads ≥ 0 eV) contrary to calculations with other functionals and Li basis sets. However, the same calculations at 3% coverage show only small differences on the Li height h and the BSSE-corrected Eads , when compared with corresponding calculations with the larger Li (6s3p1d) basis set (Table 1). We must state that Li/graphene calculations performed with the smaller Li basis set show larger differences between the BSSE-corrected and uncorrected Eads relative to corresponding calculations with the larger Li

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

basis set under the same DFT functional. This statement applies to the majority of the calculations here. At 3% coverage, Li/graphene optimal geometries using the PBE0 and PBE functionals paired with the Li (6s3p1d) basis set are almost identical. This statement applies to corresponding calculations for K/graphene, whereas for Na/graphene we observe small differences in the geometry between the two calculations (|h| = 0.04 Å, vide infra). However, at this same coverage, the PBE0 calculations shows more negative BSSE-corrected Eads values, when compared to corresponding PBE calculations, whereas the exact opposite is observed at 12.5% coverage. At 3% coverage, we repeated the PBE0 and PBE calculations paired with the larger Li (6s3p1d) basis set using the Grimme correction, which accounts for van der Waals interactions. The Grimme correction was not used with the B3LYP functional since this functional appears inferior to PBE0 and PBE functionals, as stated before [31]. Grimme-corrected PBE0 and PBE calculations systematically lengthen the Li height h and shift Eads towards more negative values (Table 1). At 3% coverage, although the Li-Li distance appears about 9.8 Å, long range Li-Li interactions are still evident. The PBE0/(6s3p1d) calculations with the large 6 × 6 graphene sheet show small but noticeable differences in the Li/graphene optimal geometry and Eads , relative to corresponding calculations at 3% coverage. Our calculated Li heights h and Li-nearest graphene C distances under the PBE0 and PBE functionals paired with the Li (6s3p1d) basis set are within the ranges reported by other researchers using periodic DFT [11,12,17]. Cluster DFT calculations by Kheirabadi and Shafiekhani [47] using the B3LYP functional show that Eads shifts to less negative values as coverage increases, in agreement with our current report. Our BSSE-corrected Eads using the PBE0 and PBE functionals and the Li (6s3p1d) basis set are in agreement with the corresponding nonGrimme corrected values of Wong et al. [17] for calculations at 12.5% Li coverage and with Nakada and Ishii [11] for calculations at lower coverages. Chan et al. [12] reported significantly less negative Eads values at 3% coverage. Some disagreements of our calculated Eads values with past reports are also observed for calculations of K and Ca adsorbed on graphene (vide infra). Zhou et al. [20] reported that the energy landscape of atom adsorption on graphene contains multiple local minima and determining the global minima may depend on the initial guess. Moreover, most of similar past work was performed using periodic DFT methods by employing a single functional typically under plane-waves approaches and using a single-size graphene sheet. Therefore, it is important to examine the behavior of more DFT functionals in the electronic and structural properties of adatom-graphene systems as well as the effect of the adatom coverage in these properties. For DFT approaches that use Gaussian-type basis sets it is advantageous to also examine the effect of the basis set size in these calculations. We must also state that although a larger basis set will improve the description of the adatom-graphene interaction, the use of a large Gaussian basis sets in periodic calculations may not be always feasible or possible. Large basis sets, which contain small exponents, may lead to poor energy convergences and to linear dependencies within the crystalline basis sets. We now examine the Na and K adsorption on graphene in a similar fashion as in the Li/graphene case. Here, we study the effect of the Na and K basis set size only for the B3LYP functional, whereas for the calculations with the PBE0 and PBE functionals we use the larger Na (9s5p1d) and K (6s6p3d1f) basis sets. Table 2 shows that the adatom basis set size significantly affects the Naand K-graphene geometries, as well as the corresponding adatom Eads values. At 12.5% adatom coverage, the Na-graphene interaction is non-bonding to weak bonding (BSSE-corrected Eads ∼0 eV and uncorrected Eads ≤ −0.5 eV, for Na (9s5p1d) basis set, Table 2), in agreement with calculation of Wong et al. [17] at the higher 16.67% coverage. At 3% coverage, the BSSE-corrected Eads values

201

from the PBE0 calculations are more negative relative to the corresponding values from the PBE calculations for both Na- and K-graphene systems, as in the Li/graphene case. Contrary to the Li/graphene case, the Grimme correction in the PBE0 calculations at 3% coverage does not affect the optimal geometries of the Naand K-graphene systems. However, the Grimme correction here shows stronger adatom-graphene interaction by the shifting of the BSSE-corrected Eads to more negative values. Although this correction does not affect the Na/graphene geometry, var der Waals interactions are still evident via comparison between calculations at 3% and 1.4% coverages (Table 2). Therefore, since alkali and alkaline adsorptions on graphene are governed by ionic-type forces and by long range electrostatic interactions [44], it is important to examine adsorption using larger graphene sheets in addition to the use of the Grimme correction. The Na/graphene geometries and non-Grimme corrected Eads are within the ranges reported in the literature [11,12,17,21], whereas our calculated K height h values are shorter than the ones previously reported [11,12,48], leading to significantly more negative Eads values (Table 2). For K/graphene at 3% coverage, the BSSE-corrected Eads using the PBE/(6s6p3d1f) is closer to the value of −2.00 eV reported by Lou et al. [45] using the C54 H18 K system under the local density approximation (Table 2). Lugo-Solis and Vasiliev [44] reported that discrepancies in the calculated Eads for K/graphene could be attributed to possible electrostatic long range interactions that are not correctly described by DFT. Fig. 2 shows the effect of the adatom coverage on the adatom height h and the BSSE-corrected Eads for Li, Na, and K and Ca adsorbed on graphene under the PBE0 and PBE calculations, respectively, using the larger of the adatoms basis sets. We observe that for Li, Na, and K adsorption on graphene, increased adatom coverage weakens the adatom-graphene interaction, as verified by the shift of the BSSE-corrected Eads to less negative values. A similar effect has been observed during carbon monoxide adsorption on Pt and Ru surfaces, where increased CO coverage lead to less negative CO Eads [46]. For Li, Na, and K adsorbed on graphene, the adatom height h is increased along with increased coverage (up to 3%), in agreement with changes in the Eads , whereas there is no general trend for h values between 1.4% and 3% coverages. The Ca/graphene is a special case, as there is no trend drawn between h and Eads with the Ca coverage. Table 3 shows that Ca interaction with graphene is non-bonding to weak bonding for all coverages examined here irrespective of the method employed (i.e., BSSE-corrected Eads > 0 and uncorrected Eads ≤ −0.63 eV for Ca (9s5p3d) calculations). Past reports showed this interaction as weakly bonding [11,12,19] at 3–6.25% coverages and non-bonding at higher coverages [17]. As it was stated in the last section, counterpoise Eads corrections may lead to overestimated Eads . Therefore, it is important to report both BSSE-corrected and uncorrected Eads , when non-bonding and/or weak bonding interactions are studied. At 3% coverage, the B3LYP and the PBE0 calculations paired with the larger Ca (9s5p3d) basis set show a relatively large Ca height h (e.g., PBE0/(9s5p3d) with h = 3.74 Å) comparted to corresponding PBE calculations (PBE/(9s5p3d) with h = 2.55 Å). The corresponding PBE0 calculations using the Grimme correction showed a substantially different Ca/graphene geometry relative to calculations without this correction, thus verifying that the PBE0 functional is inadequate to correctly describe the Ca-graphene interaction. This is due to the strong Ca-Ca and Ca-graphene van der Waals interactions that cannot be correctly described by the PBE0 functional paired under the Ca (9s5p3d) basis set for coverages ≤3%. This statement is also verified by the substantial changes in the Ca/graphene geometry between the 1.4% and 3% coverages using the PBE0/(9s5p3d) method (Table 3). Computational calculations on Ca/graphene strongly depend on the method employed (e.g., cluster or periodic DFT). Cazorla et al. [19] using the PBE func-

202

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

Fig. 2. Adatom height h and BSSE-corrected Eads for a) Li, b) Na, c) K and d) Ca adsorbed on graphene relative to the adatom coverage.

tional found that Ca/graphene is described very differently between molecular and periodic DFT. More specifically, in the former case, graphene is mimicked by the coronene molecule (i.e., C12 H24 ) and Ca-graphene interaction is weak (i.e., h = 3.84 Å and Eads = −0.06 eV), whereas corresponding periodic calculations showed a substantially stronger interaction (i.e., h = 2.26 Å and Eads = −0.443 eV). Here, the Ca interaction with graphene will be examined further using the PBE functional. 3.2. DOS spectra and charge transfers Tables 4 and 5 show the Li, K, Na, and Ca orbital populations of the adatom-graphene systems and the charge transfers between these adatoms and graphene using Mulliken population analysis and DOS integration, respectively. Table 4 also shows the adatomnearest graphene carbon COOP, which accounts for orbital overlap, whereas the differences in the values obtained by these two calculations are shown in Table 5. For Li, K, and Na calculations, the PBE0 functional was used, whereas for Ca calculations the PBE functional is used for the reasons explained in the last subsection. For all of these calculations, the Grimme correction does not affect the electronic properties of the adatom and the graphene sheet (Tables 4 and 5). Although COOP cannot describe the strength of a primarily ionic bond, our calculated COOP values show that some orbital overlap occurs between the adatom and the nearest graphene C atom. This effect indicates that the ionic alkali/alkalinegraphene interaction also contains small covalent parts (e.g., 3% Li coverage, COOP = 0.040, Table 4). Past reports at 6.25% [11] and 3% adatom coverages [12] show that for the systems of this work charge is transferred from the adatom towards the graphene sheet. This statement is in agreement with our Mulliken population calculations with the exception of Na/graphene at the high 12.5% coverage, where minor charge flow

is observed in the opposite direction (i.e., from graphene to Na). The DOS integrations calculations show the same charge flow between the adatoms and graphene as the Milliken population, with the exception of Na/graphene. For this case, the DOS integration calculations show that the charge flow from Na to graphene is decreased as Na coverage increases, and reverses at low coverage. This trend is opposite to corresponding calculations from Mulliken population and will be explained below using the information from the graphs of the DOS spectra. Charge values reported in the literature may vary due to the different methods employed for these calculations [11]. However, as we mentioned before, the absolute charge values are of minimal importance, and our focus is on possible trends that may be drawn from relative changes of these values due to method employed, adatom coverage, etc. We now compare adatom orbital populations and charge transfers obtained via Milliken population and DOS integration. For Li and K adsorbed on graphene, there is an excellent agreement between these two methods (Tables 4 and 5). The only exceptions to this statement are the Li–2s and K-2p orbital populations at 3% adatom coverages that are overestimated by DOS integration relative to Mulliken population, leading to underestimation by the former method in the charge transferred to graphene. However, for Na/graphene and Ca/graphene, some differences between the charges obtained by these two methods are noteworthy. As mentioned above, these mismatches can be explained by examining the DOS spectra for these two cases. Figs. 3 and 4 show the total DOS spectra and the DOS per adatom orbital population and per atom of the spd graphene for Li and Na and K and Ca, respectively, adsorbed on graphene, at coverages examined in this work. For Na/graphene and Ca/graphene, some energy bands that correspond to DOS peaks in the vicinity of the Fermi level are of partial occupancy and thus, their occupancies are difficult to be obtained correctly by CRYSTAL09. For example, for Ca/graphene, the Ca–4s orbital appears as

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

203

Table 4 Adatom orbital populations, charge transferred between the adatom and graphene using Mulliken population, and the Metal-C COOP for Li, Na, K, and Ca adsorbed on 2 × 2, 4 × 4, and 6 × 6 graphene sheets. The values in parenthesis correspond to calculations using the Grimme correction. The values in square brackets refer to past reports from the literature. Adatom

Supercell

Functional/Adatom Basis Set

Metal charge per orbital

Charge transferred to graphene

Metal-C COOP

Li

2×2

PBE0/(9s5p3d)

Li–2s = 0.43 e Li-2p = 0.46 e Li–2s = 0.21 e (0.21 e) Li-2p = 0.43 e (0.42 e) Li–2s = 0.21 e Li-2p = 0.41 e

0.10 e

0.031

0.34 e (0.36 e) [0.86a , 0.3912 ] 0.37 e

0.040 (0.039)

Na–3s = 0.47 e Na-3p = 0.47 e Na–3s = 0.29 e (0.30 e) Na-3p = 0.40 e (0.40 e) Na–3s = 0.18 e Na-3p = 0.34 e

−0.02 e

0.031

0.20 e (0.19 e) [0.62a , 0.32b , 0.657c ] 0.34

0.042 (0.041)

K–4s = 0.05 e K-4p = 0.60 e K-3d = 0.23 e K–4s = 0.04 e (0.04 e) K-4p = 0.13 e (0.14 e) K-3d = 0.27 e (0.27 e) K–4s = 0.03 e K-4p = 0.12 e K-3d = 0.24 e

0.17 e

−0.007

0.59 e (0.59 e) [0.63a , 0.41b ]

0.018 (0.018)

0.65 e

0.020

Ca–4s = 1.24 e Ca-4p = 0.13 e Ca-3d = 0.48 e Ca–4s = 1.39 e (1.41 e) Ca-4p = 0.11 e (0.11 e) Ca-3d = 0.27 e (0.27 e) Ca–4s = 1.29 e Ca-4p = 0.12 e Ca-3d = 0.29 e

0.14 e

−0.004

0.22 e (0.21 e) [0.87a , 0.18b ]

−0.012 (−0.011)

0.30 e

−0.010

4×4 6×6 Na

2×2

PBE0/(9s5p1d)

4×4 6×6 2×2

K

PBE0/(6s6p3d1f)

4×4

6×6

Ca

2×2

PBE/(9s5p3d)

4×4

6×6

a b c

0.040

0.046

Ref. [11]. Ref. [12]. Ref. [18].

Table 5 Adatom orbital populations and charge transferred to graphene by integrating DOS in energies up to EFermi for Li, Na, K, and Ca adsorbed on 2 × 2, 4 × 4, and 6 × 6 graphene sheets. The values in square brackets show the difference between these values and the corresponding values from Mulliken population analysis. The values in parenthesis correspond to calculations using the Grimme correction. Adatom

Supercell

Functional/Adatom Basis Set

Metal charge per orbital

Charge transferred to graphene

Li

2×2

PBE0/(9s5p3d)

Li–2s = 0.43 e [<0.01e] Li-2p = 0.48 e [0.02 e] Li–2s = 0.30 e [0.09 e] (0.30 e [0.09 e]) Li-2p = 0.47 e [0.04 e] (0.46 e [0.04 e]) Li–2s = 0.21 e [<0.01 e] Li-2p = 0.41 e [<0.01 e]

0.08 e [−0.02 e]

Na–3s = 0.39 e [-0.08 e] Na-3p = 0.47 e [<0.01 e] Na–3s = 0.41 e [0.12 e] (0.42 e [0.12 e]) Na-3p = 0.47 e [0.07 e] (0.47 e) [0.07 e] Na–3s = 0.46 e [0.28 e] Na-3p = 0.48 e [0.14 e]

0.07 e [0.09 e]

K–4s = 0.04 e [−0.01 e] K-4p = 0.57 e [−0.03 e] K-3d = 0.25 e [0.02 e] K–4s = 0.06 e [0.02 e] (0.06 e) [0.02 e] K-4p = 0.39 e [0.28 e] (0.39 e) [0.28 e] K-3d = 0.27 e [<0.01 e] (0.26 e [−0.01 e]) K–4s = 0.03 e [<0.01 e] K-4p = 0.12 e [<0.01 e] K-3d = 0.24 e [<0.01 e]

0.18 e [0.01 e]

Ca–4s = 1.31 e [0.07 e] Ca-4p = 0.13 e [<0.01 e] Ca-3d = 0.25 e [-0.23 e] Ca–4s = 1.52 e [0.13 e] (1.54 e [0.13 e]) Ca-4p = 0.12 e [0.01 e] (0.11 e [<0.01 e]) Ca-3d = 0.20 e [0.07 e] (0.19 e [0.08 e]) Ca–4s = 1.06 e [−0.23 e] Ca-4p = 0.12 e [<0.01 e] Ca-3d = 0.68 e [0.39 e]

0.07 e [−0.07 e]

4×4 6×6 Na

2×2

PBE0/(9s5p1d)

4×4 6×6 K

2×2

PBE0/(6s6p3d1f)

4×4

6×6

Ca

2×2

4×4

6×6

PBE/(9s5p3d)

0.22 e [−0.12 e] (0.23 e [−0.13 e]) 0.37 e [<0.01 e]

0.00 e [−0.20 e] (0.0 e [−0.19 e]) −0.06 e [−0.40 e]

0.34 e [−0.25 e] (0.35 e [−0.24 e]) 0.62 e [−0.03 e]

0.10 e [−0.12 e] (0.09 e[−0.12 e])

0.15 e [−0.15 e]

204

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

Fig. 3. Partial DOS spectra per adatom orbital and spd graphene and the total DOS for (a) Li and b) Na adsorbed on graphene using the PBE0 functional. The vertical line is the Fermi level.

a sharp peak on the Fermi level at both 1.5 and 3% adatom coverages and its calculated population differs between the two methods employed here (Fig. 4b). Therefore, orbital populations and charge transfers obtained via DOS integration should be used with caution, when sharp peaks appear at the Fermi level region. Hence, in this work, charges obtained via Mulliken population are preferred over corresponding values obtained via DOS integration. Fig. 5 shows the charge transferred from the adatom to graphene as calculated by Mulliken population relative to Eads for Li, Na,

K, and Ca adsorption on graphene. For Li, Na, and K adsorbed on graphene, a clear trend is drawn between the charge transferred from the adatom to graphene and the BSSE-corrected Eads , when we consider only calculations with large adatom basis sets. Here, the increased charge transfer to graphene is along with the more negative Eads (i.e., stronger adatom-graphene interaction), in agreement with the ionic nature of the alkali-graphene interaction. However, no corresponding trend is drawn for Ca/graphene. The Ca-graphene interaction appears stronger at the higher coverages of 12.5%, in

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

205

Fig. 4. Partial DOS spectra per adatom orbital and spd graphene and total DOS for (a) K and b) Ca adsorbed on graphene. The PBE0 and PBE functionals are used for K/graphene and Ca/graphene calculations, respectively. The vertical line is the Fermi level.

contrast to the other cases examined here (Table 3). We will show that changes in the Ca 4s population govern the Ca-graphene interaction (vide infra). We now discuss the location of the Dirac point (i.e., the energy point in the DOS spectra where DOS is zero) relative to the Fermi level. Here the Dirac point, is clearly visible and appears at values negative to the Fermi level (Figs. 3 and 4). Chan et al. [12], stated that in the case of alkali/alkaline adsorption on graphene surfaces, integrating the DOS spectra region from the Dirac point to the Fermi level could provide information about the charge transferred

from the adatom to graphene. However, these authors reported that charge transfers from the adatoms to graphene obtained via integration of the DOS spectra in the region between the Dirac point and the Fermi level are significantly larger than corresponding values obtained from electron density calculations. In this work, since the Dirac point is to the left of the Fermi level, all adatoms of this work must donate change to graphene. This statement is in agreement with corresponding Mulliken population analysis calculations, with the exception of the Na/graphene at 12.5% coverage using the PBE0/(9s5p1d) method (Table 4 and Fig. 3). This single

206

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

Fig. 5. Charge transferred from the adatom to graphene relative to the BSSE-corrected Eads for a) Li, b) Na, c) K, and d) Ca adsorbed on graphene. Charges have been computed using Mulliken population analysis. Shapes and colors refer to the functional/basis set used. Linear fits are obtained by using only the larger basis sets in the calculations. No trend is drawn for Ca/graphene.

Fig. 6. Charge density 3D plots for Na/graphene using a) the 2 × 2and b) the 4 × 4 graphene supercell relative to corresponding superposition of atomic orbitals. The PBE0 functional and the Na (9s5p1d) basis set were used. Blue areas denote change deficient regions, while red areas denote charge surplus regions. (i.e., charge flows from the blue areas to the red areas). These results agree with charge changes obtained from the Mulliken population analysis (Table 4). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

case, which shows minor charge transfer of 0.02 e from graphene to Na, is also examined by the use of charge density 3D plots at the 12.5% and 3% coverages relative to corresponding superposition of atomic orbitals (Fig. 6). We recall that for this case DOS integration show minor charge transfer of 0.05 e towards the graphene substrate. At 12.5% coverage, the Na atoms are colored red (i.e., charge excess) and at 3% coverage are colored blue (i.e., charge deficiency), in agreement with Mulliken population analysis. Therefore, for small charge transfers between the adatom and graphene the charge density 3D plots and the Mulliken population analysis are preferred over calculations using shifts in the Dirac

point and DOS integration calculations. Moreover, Fig. 6 shows no significant orbital overlap between the Na adatom and graphene, thus verifying that Na-graphene interaction is primarily ionic. Increased adatom adsorption from 1.4–3% to 12.5% broadens the adatom orbital DOS as well as the graphene DOS (Figs. 3 and 4), in agreement with our past work on Fe/graphene [24]. Moreover, Li, Na, K, and Ca adsorption does not significantly alter DOS graphene spectra near and below the Fermi level, in agreement with past reports [11,12]. At 12.5% coverage, for Li/graphene, the Li 2p orbital broadens and hybridizes with the 2s orbital, both appearing as partially populated (Fig. 3a and Tables 4 and 5). As Li coverage decreases, its 2s and 2p orbitals are further pushed above the Fermi level and thus, more depopulated. Similar effects appear for Na/graphene, whereas the situation is slightly different for K/graphene. In this latter case, decreased K coverage is along with the decrease of the K 4s and 4p populations, whereas a “tail” of the K-3d orbital appears below the Fermi level, and thus partially populated. The K-3d orbital population remains about constant at all coverages (Tables 4 and 5). As we mentioned before, the Ca/graphene case is a special case and Ca-interaction strength is not described by changes in the overall charge transferred from Ca to graphene. At 12.5% Ca coverage, the PBE/(9s5p3d) calculations show the Ca–4s orbital is broadened around the Fermi level and partially populated (Fig. 4b). Decreased Ca coverage from 12.5% to 3% pushes the Ca–4s orbital below the Fermi level and thus, more populated. However, further deceases in the Ca surface coverage down to 1.4% narrows the Ca–4s orbital spectra and pushes this orbital in the opposite direction towards higher energies. Therefore, the strength of the Ca interaction with graphene is opposite to increases in the Ca 4s population. Moreover, the higher Ca-3d population at the 12.5% coverage further strengthens the Ca-graphene interaction, which explains why at this coverage Ca adsorbs stronger on graphene than any other coverages examine here. Fig. 7 shows that the Grimme correction strongly affects the DOS spectra of the

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

207

the understanding of the atomistic interactions of adatoms with graphene.

References

Fig. 7. Partial DOS spectra per Ca-s and -d orbitals and spd graphene and total DOS for Ca/graphene under the (a) PBE0 and (b) PBE functional calculations with and without the Grimme correction. The vertical line is the Fermi level.

PBE0 calculations, whereas its effect on PBE calculation is minimal. This observation verifies that the PBE0 functional is insufficient to describe the Ca-Ca and Ca-graphene van der Waals interactions at coverages of 3% and less, thus producing bogus results.

4. Conclusions DFT calculations on graphene monolayers with Li, Na, K, and Ca adatoms have been used to analyze alkali and the alkaline Ca adsorption on graphene at 1.4%, 3%, and 12.5% adatom surface coverages. The adatom Eads and geometry optimizations for the adatom-graphene structures were obtained using the hybrid B3LYP and PBE0 functionals and the GGA PBE functional at various adatom basis sets. We found that PBE0 and PBE functionals paired with a sufficiently large adatom basis sets predict structural and electronic information accurately for the Li, Na, and K adsorbed on graphene. For Ca/graphene, the PBE0 functional cannot be used due to poor description of the Ca-Ca and Ca-graphene van der Waals interactions. Mulliken populations analysis for the adatom-graphene systems of this work revealed that alkali/alkaline-graphene interaction is mostly ionic, as expected. This is verified by the calculated charge transfers between the adatom and graphene and the almost unchanged DOS in the region near and below the Fermi level. However, the presence of small adatom-nearest graphene carbon COOP values is indicative of small covalent interaction between these adatoms and graphene, in addition to the dominant ionic interaction, which is described above. We observed that for Li, Na, and K adsorbed on graphene, increased adatom surface coverage weakens the adatom-graphene interaction, as verified by the less negative Eads and increased adatom height h. In these cases, the adatom-graphene interaction strength is opposite to the amount of charge transferred from the adatom to graphene. For Ca/graphene, no trend is drawn between the Ca coverage and Eads , the latter used as a measurement of the adatom-graphene interaction strength. Moreover, the adsorption of Ca on graphene is the strongest at higher coverages, in contrast to the other adsorptions examined here. This unusual behavior is explained by correlating the Ca–4s orbital population with the Cagraphene adsorption strength in an opposite relationship. In the future, we will examine the co-adsorption of these adatoms and water molecules on graphene, with and without defects in the monolayer using periodic DFT. These calculations, complimented with information from experimental spectroscopy, will elucidate

[1] F. Schedin, A.K. Geim, S.V. Morozov, E.W. Hill, P. Blake, M.I. Katsnelson, K.S. Novoselov, Detection of individual gas molecules adsorbed on graphene, Nat. Mater. 6 (2007) 652–655. [2] J.O. Sofo, A.S. Chaudhari, G.B. Barber, Graphane: a two-dimensional hydrocarbon, Phys. Rev. B 75 (15) (2007) 153401–153404. [3] Z. Luo, T. Yu, K. Kim, Z. Ni, Y. Yu, Y. You, S. Lim, Z. Shen, S. Wang, J. Lin, Thickness-dependent reversible hydrogenation of graphene layers, ACS Nano 3 (2009) 1781–1788. [4] Y. Liu, B. Xie, Z. Zhang, Q. Zheng, Z. Xu, Mechanical properties of graphene papers, J. Mech. Phys. Solids 60 (2012) 591–605. [5] E. Pop, V. Varshney, A. Roy, Thermal properties of graphene: fundamentals and applications, MRS Bull. 37 (2012) 1273–1281. [6] A.H. Castro Neto, F. Guinea, N.M.R. Peres, K.S. Novoselov, A.K. Geim, The electronic properties of graphene, Rev. Mod. Phys. 81 (2009) 110–155. [7] K.S. Novoselov, A.K. Geim, S.V. Morozov, D. Jiang, M.I. Katsnelson, I.V. Grigorieva, A.A. Firsov, Electric field effect in anatomically thin carbon films, Science 306 (2004) 666–669. [8] N.D. Mermin, Crystalline order in two dimensions, Phys. Rev. 176 (1968) 250–254. [9] A. O’Hare, F.V. Kusmartsev, K.I. Kugel, A stable flat form of two-dimensional crystals: could graphene, silicene, germanene Be minigap semiconductors? Nano Lett. 12 (2012) 1045–1052. [10] A. Fasolino, J.H. Los, M. Katsnelson, Intrinsic ripples in graphene, Nat. Mater. 6 (11) (2007) 858–861. [11] K. Nakada, A. Ishii, DFT Calculation for adatom adsorption on Graphene, in: J.R. Gong (Ed.), Graphene Simulation, Tech Inc., Croatia, 2011, pp. 3–20. [12] K.T. Chan, J.B. Neaton, M.L. Cohen, First-principles study of metal adatom adsorption on graphene, Phys. Rev. B. 77 (2008) 235430–235432. [13] P. Hohenberg, W. Kohn, Inhomogeneous electron gas, Phys. Rev. B 136 (1964) 864–871. [14] W. Kohn, L.J. Sham, Self-consistent equations including exchange and correlation effects, Phys. Rev. A 140 (1965) 1133–1138. [15] R.G. Parr, W. Yang, Density Functional Theory of Atoms and Molecules, Oxford University Press, New York, 1989. [16] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett. 77 (1996) 3865–3868. [17] J. Wong, S. Yadav, J. Tam, C. Veer Singh, A van der Waals density functional theory comparison of metal decorated graphene systems for hydrogen adsorption, J. Appl. Phys. 115 (2014) 224301–224311. [18] K. Lee, É.D. Murray, L. Kong, B.I. Lundqvist, D.C. Langreth, A higher-accuracy van der Waals density functional, Phys. Rev. B 82 (2010), 081101(R)-081104(R). [19] C. Cazorla, S.A. Shevlin, Z.X. Guo, First-principles study of the stability of calcium-decorated carbon nanostructures, Phys. Rev. B 82 (2010), 155454–12. [20] J. Zhou, Q. Sun, Q. Wang, P. Jena, Tailoring Li adsorption on graphene, Phys. Rev. B 90 (2014), 205427–7. [21] H.S. Moon, J.H. Lee, S. Kwon, I.T. Kim, S.G. Lee, Mechanisms of Na adsorption on graphene and graphene oxide: density functional theory approach, Carbon Lett. 16 (2) (2015) 116–120. [22] G. Ulian, S. Tosoni, G. Valdrè, Comparison between Gaussian-type orbitals and plane wave ab initio density functional theory modeling of layer silicates: talc [Mg3 Si4 O10 (OH)2 ] as model system, J. Chem. Phys. 139 (2013) 204101–204110. [23] M.F. Peintinger, D.V. Oliveira, T. Bredow, Consistent Gaussian basis sets of triple-zeta valence with polarization quality for solid-state calculations, J. Comput. Chem. 34 (2013) 451–459. [24] N. Dimakis, N.E. Navarro, J. Velazquez, A. Salgado, Electronic and vibrational properties of graphene monolayers with iron adatoms: a density functional theory study, Appl. Surf. Sci. 334 (2015) 2–6. [25] A.D. Becke, Density-Functional thermochemistry. III. The role of exact exchange, J. Phys. Chem. 98 (1993) 5648–5652. [26] M. Ernzerhof, G.E. Scuseria, Assessment of the Perdew-Burke-Ernzerhof exchange-correlation functional, J. Chem. Phys. 110 (1999) 5029–5036. [27] C. Adamo, V. Barone, Toward reliable density functional methods without adjustable parameters: the PBE0 model, J. Chem. Phys. 110 (1999) 6158–6170. [28] N. Dimakis, F.A. Flor, A. Salgado, K. Adjibi, S. Vargas, J. Saenz, Density functional theory calculations on transition metal atoms adsorbed on graphene monolayers, Appl. Surf. Sci. (2016), http://dx.doi.org/10.1016/j. apsusc.2016.10.144. [29] R. Dovesi, V.R. Saunders, C. Roetti, R. Orlando, C.M. Zicovich-Wilson, F. Pascale, B. Civalleri, K. Doll, N.M. Harrison, I.J. Bush, P. D’Arco, M. Llunell, CRYSTAL09. http://www.crystal.unito.it/. [30] S.H. Vosko, L. Wilk, M. Nusair, Accurate spin-dependent electron liquid correlation energies for local spin density calculations: a critical analysis, Can. J. Phys. 58 (1980) 1200–1211. [31] J. Paier, M. Marsman, G. Kresse, Why does the B3LYP hybrid functional fail for metals? J. Chem. Phys. 127 (2007) 024103–024113. [32] S. Grimme, Semiempirical gga-type density functional constructed with a long range dispersion correction, J. Comput. Chem. 27 (2006) 1787–1799.

208

N. Dimakis et al. / Applied Surface Science 413 (2017) 197–208

[33] W.S. Wesselowski, K.W. Wassiliew, Zwillingsbildungen bei graphitkristallen, Zeitschrift für Kristallographie 89 (1934) 494–496. [34] H.J. Monkhorst, J.D. Pack, Special points for brillouin-zone integrations, Phys. Rev. B 13 (1976) 5188–5192. [35] G. Gilat, L.J. Raubenheimer, Accurate numerical method for calculating frequency-distribution functions in solids, Phys. Rev. 144 (1966) 390–395. [36] G. Gilat, Analysis of methods for calculating spectral properties in solids, J. Comput. Phys. 10 (1972) 432–465. [37] D.G. Anderson, Iterative procedures for nonlinear integral equations, J. Assoc. Comput. Mach. 12 (1965) 547–560. [38] R. Hoffmann, A chemical and theoretical way to look at bonding on surfaces, Rev. Mod. Phys. 60 (1988) 601–628. [39] A. Kokalj, Computer graphics and graphical user interfaces as tools and simulations of matter at the atomic scale, Comput. Mater. Sci. 28 (2003) 155–168, Code available from http://www.xcrysden.org/. [40] R.S. Mulliken, Electronic population analysis on LCAO-MO molecular wave functions I, J. Chem. Phys. 23 (1955) 1833–1840. [41] E.R. Davidson, S. Chakravorty, A test of the Hirshfeld definition of atomic charges and moments, Theor. Chim. Acta 83 (1992) 319–330.

[42] M.D. Segall, R. Shah, C.J. Pickard, M.C. Payne, Population analysis of plane-wave electronic structure calculations of bulk materials, Phys. Rev. B 54 (1996) 16317–16320. [43] S.F. Boys, F. Bernardi, The calculation of small molecular interactions by the differences of separate total energies. Some procedures with reduced errors, Mol. Phys. 19 (4) (1970) 553–566. [44] A. Lugo-Solis, I. Vasiliev, Ab initio study of K adsorption on graphene and carbon nanotubes: role of long-range ionic forces, Phys. Rev. B 76 (2007) 235431–235438. [45] L. Lou, L. Österlund, B. Hellsing, Electronic structure and kinetics of K on graphite, J. Chem. Phys. 112 (2000) 4788–4796. [46] N. Dimakis, N. Navaro, T. Mion, E.S. Smotkin, Carbon monoxide adsorption coverage effects on platinum and ruthenium surfaces, J. Phys. Chem. C 118 (22) (2014) 11711–11722. [47] N. Kheirabadi, A. Shafiekhani, Graphene/Li-ion battery, J. Appl. Phys. 112 (2012) 124323–124325. [48] M. Caragiu, S. Finberg, Alkali metal adsorption on graphite: a review, J. Phys. Condens. Matter 17 (2005) R995–R1024.