Ab initio study of sodium cointercalation with diglyme molecule into graphite

TL;DR

This study uses ab initio calculations to understand sodium cointercalation with diglyme molecules into graphite, revealing structural, energetic, and electronic properties relevant for sodium-ion battery anodes.

Contribution

It provides detailed atomistic insights into sodium-diglyme cointercalation mechanisms and properties, advancing the development of improved anode materials.

Findings

Na(digl)$_2$C$_n$ has lowest intercalation energy at n≈21

Solvated Na(digl)$_2$ diffuses rapidly in graphite layers

Electronic conductance of cointercalated compounds is enhanced

Abstract

The cointercalation of sodium with the solvent organic molecule into graphite can resolve difficulty of forming the stage-I Na-graphite intercalation compound, which is a predominant anode of Na-ion battery. To clarify the mechanism of such cointercalation, we investigate the atomistic structure, energetics, electrochemical properties, ion and electron conductance, and charge transferring upon de/intercalation of the solvated Na-diglyme ion into graphite with {\it ab initio} calculations. It is found that the Na(digl)$_2$C$_n$ compound has the negatively lowest intercalation energy at $n\approx$21, the solvated Na(digl)$_2$ ion diffuses fast in the interlayer space, and their electronic conductance can be enhanced compared to graphite. The calculations reveal that the diglyme molecules as well as Na atom donates electrons to the graphene layer, resulting in the formation of ionic bonding between the graphene layer and the moiety of diglyme molecule. This work will contribute to the development of innovative anode materials for alkali-ion battery applications.

Tables

Table 1. Table 1: Lattice constants ( a 𝑎 a and c 𝑐 c ) and exfoliation energy ( E exf subscript 𝐸 exf E_{\text{exf}} ) in graphite, cohesive energy in sodium and lithium metals ( E coh subscript 𝐸 coh E_{\text{coh}} ), binding energy in free solvation complex A(digl) x ( E b subscript 𝐸 b E_{\text{b}} ), intercalation energy in A(digl) x C 24 ( E int subscript 𝐸 int E_{\text{int}} ), exfoliation energy in A(digl) x C 20 , and average electrode voltage ( V el subscript 𝑉 el V_{\text{el}} )

		LDA	LDA	PBE	PBE	PBEsol	PBEsol	LYP	LYP	Exp.
			+vdW		+vdW		+vdW		+vdW
Graphite	$a$ (Å)	2.466	2.465	2.465	2.464	2.464	2.464	2.471	2.470	2.461a
	$c$ (Å)	6.217	6.086	6.639	6.509	6.378	6.247	6.816	6.730	6.705a
	$E_{\text{exf}}$ (meV)	$-110$	$-130$	$-60$	$-76$	$-77$	$-95$	$-47$	$-61$	$-35\sim -{52}^b$
$E_{\text{coh}}$ (eV)	Na	$-1.31$	$-1.33$	$-1.13$	$-1.16$	$-1.22$	$-1.25$	$-0.91$	$-0.93$	$-{1.11}^c$
	Li	$-1.72$	$-1.74$	$-1.61$	$-1.62$	$-1.64$	$-1.65$	$-1.28$	$-1.30$	$-{1.63}^c$
$E_{\text{b}}$ (eV)	Na(digl)	$-0.42$	$-0.46$	$-0.29$	$-0.32$	$-0.29$	$-0.32$	$-0.28$	$-0.31$
	Na(digl)2	$-0.65$	$-0.80$	0.27	0.12	0.05	$-0.11$	0.27	0.14
	Li(digl)2	$-1.89$	$-2.02$	$-2.20$	$-2.32$	0.00	$-0.13$	$-0.91$	$-1.03$
$E_{\text{int}}$ (eV)	Na(digl)C20	$-$2.85	$-$3.06	$-$1.68	$-$1.84	$-$1.83	$-$2.16	$-$1.18	$-$1.40
	Na(digl)2C20	$-$1.04	$-$1.36	1.01	0.61	0.38	$-$0.02	1.87	1.44
	Li(digl)2C20	$-3$.11	$-$3.36	$-$0.70	$-$1.18	$-$1.60	$-$1.87	$-$0.12	$-$0.36
$E_{\text{exf}}$ (meV)	Na(digl)C20	$-276$	$-330$	$-171$	$-227$	$-226$	$-265$	$-168$	$-206$
	Na(digl)2C20	$-287$	$-315$	$-202$	$-227$	$-230$	$-255$	$-163$	$-182$
	Li(digl)2C20	$-283$	$-289$	$-211$	$-232$	$-241$	$-259$	$-180$	$-192$
$V_{\text{el}}$ (V)	Na(digl)C24	4.01	4.28	1.23	4.77	2.00	2.00	0.40	1.17
	Na(digl)2C24	2.70	3.80	$-2.51$	$-1.43$	$-0.93$	0.16	$-3.34$	$-3.70$	0.15$\sim$0.45d
	Li(digl)2C20	4.52	5.31	$-$0.89	0.31	0.99	1.75	$-$2.08	$-$1.01	0.75d

Table 2. Table 2: Graphene unit cell, specific capacity ( C 𝐶 C ), average bond length between Na and O ( d̄ Na-O Na-O {}_{\text{Na-O}} ), interlayer distance ( d int subscript 𝑑 int d_{\text{int}} ), relative volume expansion ratio ( r vol subscript 𝑟 vol r_{\text{vol}} ), exfoliation energy ( E exf subscript 𝐸 exf E_{\text{exf}} ), and intercalation energy ( E int subscript 𝐸 int E_{\text{int}} ) in t 𝑡 t -GICs Na(digl) 2 C n with n = 𝑛 absent n= 18, 20, 24, 28, and 32.

	C18	C20	C24	C28	C32
Unit cell	$\left(3\times 3\right)$	$\left(\sqrt{12}\times \sqrt{7}\right)$	$\left(4\times 3\right)$	$\left(\sqrt{12}\times \sqrt{13}\right)$	$\left(4\times 4\right)$
$C$ (mAh g-1)	112.0	101.8	86.1	74.6	65.8
d̄${}_{\text{Na-O}}$ (Å)	2.253	2.315	2.275	2.297	2.277
$d_{\text{int}}$ (Å)	10.30	10.22	9.96	10.19	10.15
$r_{\text{vol}}$ (%)	240	236	223	231	229
$E_{\text{exf}}$ (meV)	$-$245	$-$260	$-$259	$-$288	$-$279
$E_{\text{int}}$ (eV)	1.31	$-$0.02	$-$0.14	0.00	0.24

Equations

$$E_{\text{vdW}}=-s_6\sum _{j<i}\frac{C_6^{ij}}{r_{ij}^6}\frac{1}{1+e^{-d(r_{ij}/R_{ij}^0-1)}}$$

$$x\text{A}+y(\text{digl})\to {\text{A}}_x(\text{digl}{)}_y$$

$$E_{\text{b}}=E_{{\text{A}}_x(\text{digl}{)}_y}-(x{E}_{\text{A}}+y{E}_{\text{digl}})$$

$${\text{A}}_x(\text{digl}{)}_y+{\text{C}}_n\to {\text{A}}_x(\text{digl}{)}_y{\text{C}}_n,$$

$$E_{\text{int}}=E_{{\text{A}}_x(\text{digl}{)}_y{\text{C}}_n}-(E_{{\text{A}}_x(\text{digl}{)}_y}+E_{{\text{C}}_n}).$$

$$E_{\text{exf}}=\frac{1}{n_{\text{C}}}[E(d_{\text{i}}=d_{\text{e}})-E(d_{\text{i}}=\infty )],$$

$$V_{\text{el}}=-\frac{E_{{\text{Na}}_{x_j}(\text{digl}{)}_{2x_j}{\text{C}}_{32}}-E_{{\text{Na}}_{x_i}(\text{digl}{)}_{2x_i}{\text{C}}_{32}}-(x_j-x_i)(2E_{\text{digl}}+E_{{\text{Na}}_{\text{bcc}}})}{(x_j-x_i)e}$$

$$\Delta n(\text{r})=n_{{\text{Na(digl)}}_x{\text{C}}_n}(\text{r})-[n_{\text{Na}}(\text{r})+n_{{\text{(digl)}}_x}(\text{r})+n_{{\text{C}}_n}(\text{r})]$$

Full text

††journal: Nano Energy

Energy and environmental issues become increasingly global, actual and vital to all the nations on the earth. This prompts people to cease from mass consumption of fossil fuels, which are being exhausted and moreover cause greenhouse gas emission, and to create new technologies for utilizing renewable and clean energy sources such as solar and wind power. However, these sources are changeable in time and widespread in space. For solar and wind power to become competitive with fossil fuels in the electricity market, therefore, it is necessary to develop an efficient energy storage system as well as an effective energy harvesting device. Although Li-ion batteries (LIBs) have been widely used for portable electronics during the past two decades due to their high power density, growing concerns for low abundance and high price of lithium resource [1, 2, 3] triggered extensive research on Na-ion batterries (NIBs) recently considering high abundance and low cost of sodium precursor [4, 5, 6, 7]. When compared to Li+ ion, however, Na+ ion has lower ionization potential and larger ionic radius, which make it difficult to find suitable anode materials [8, 9].

Requirements for a good anode material can be (i) high specific capacity and low redox potential for high energy density, (ii) high abundance, low cost and facile fabrication for economics, (iii) low chemical reactivity with electrolyte and small volume change during operation for safety and long cycling life, and (iv) non-toxicity for environmental friendliness. When surveying NIB anodes based on these requirements, metal sodium can be used surely but has a critical problem such as high reactivity with organic electrolyte solvents, while metal oxides and alloy based materials have in general high specific capacities but suffer from large volume change [10, 11, 12, 13]. In spite of relatively low capacity, carbon based materials might be the most widely used anode for NIBs because of good electronic conductivity and excellent electrochemical performance [14, 15]. The research tended to mainly non-graphitic carbons like disordered hard carbon, which has high specific capacity over $\sim$300 mAh g*-1* but still shortcomings such as low coulombic efficiency in the first cycle due to large surface area [16, 17]. During the last few years, the significant progress in using graphitic carbons has been made [18, 19, 20]. From the studies of hard carbons, it was awoken that desirable carbon materials should have long-range order, high interlayer distance, low porosity and small surface area [19]. In this regard, Wen et al. [20] suggested the expanded graphite with a mean interlayer distance of $\sim$0.43 nm (0.334 nm in graphite) and a relatively small surface area of $\sim$30 m2 g*-1*, being able to deliver high specific capacity of $\sim$300 mAh g*-1* at low current density of 20 mA g*-1*. At high current density, however, the modest performance was observed. Using the expanded synthetic graphite materials, Ramos et al. [19] increased the current density up to 100 mA g*-1*, at which the capacity was of $\sim$110 mAh g*-1*.

Then, why not graphite itself for NIB anode? Although the alkali metals are known to reversibly intercalate into graphite, forming binary graphite intercalation compounds ($b$-GICs) in first stage, only sodium can form $b$-GICs in high stages [21, 22]. This leads to a very low specific capacity for cells, e.g., $\sim$35 mAh g*-1* for \ceNaC64 versus $\sim$372 mAh g*-1* for \ceLiC6, hindering practical use of graphite as NIB anode. The reason might be the small interlayer spacing of graphite, i.e., big mismatch between graphite interlayer distance and Na ion size, and most likely, weak binding of Na ion (electron donor) to graphene sheet (acceptor) through $\pi$ interaction [23, 24, 25]. A fortunate side for electrode performance is that the larger radius ions diffuse in graphite more smoothly [23].

To form stage-I Na-GICs, cointercalation process, i.e., intercalation of sodium with solvent organic molecule, has been suggested [26, 27, 28, 29, 30, 31]. This produces donor-type ternary GICs ($t$-GICs), denoting Na${}^{+}_{x}$(solv)yC${}^{-}_{n}$ with Cn the graphene sheet containing $n$ carbon atoms and \ce(solv)_y the $y$ solvent organic molecules [26]. The solvents are typically a series of alkylamines [31, 32, 30], and ether-based electrolytes [26, 27, 28, 29]. Recently, Jache and Adelhelm [26] made use of $t$-GICs containing sodium atom and diglyme (diethylene glycol dimethyl ether: \ceC6H14O3) molecule with an estimated stoichiometry of \ceNa(digl)2C20 as an effective NIB anode, which delivers moderately high reversible capacities close to 100 mAh g*-1* for 1000 cycles and coulomb efficiencies over 99.87%. Hasa et al. [33] demonstrated the NIB full cell consisted of a graphite anode, glyme-based electrolyte (tetraethylene glycol dimethyl ether (TEGDME)-\ceNaClO4), and a layered oxide cathode, i.e., Graphite/TEGDME-\ceNaClO4/\ceNa_0.7CoO2, which offers unique characteristics of long cycle life, high efficiency and high power density, though the lower energy density than LIBs.

In spite of such successful applications, there are still concerns about using $t$-GICs as NIB anode [15]. This is in part due to a lack of sufficient understanding of sodium cointercalation phenomena. Ab initio calculations within the framework of density functional theory (DFT) have proven to give valuable insights into electrochemical performance of NIB anodes [34, 35, 36, 37]. In this work, we perform systematic ab initio DFT calculations of $t$-GICs cointercalated with sodium atom and diglyme molecule, aiming to clarify the mechanism of formation and operation as NIB anodes.

1 Computational methods

Firstly we make models for the crystalline $t$-GICs cointercalated with sodium atom and diglyme molecule, Nax(digl)yCn with $x=1$, $y=2$, and different $n$ numbers. Regarding the intercalant complexes, it was assumed in the previous works [26, 38, 39, 40] that two diglyme molecules are bound to one alkali ion with six coordination number by ether oxygen atoms from the solvent molecules, and thus we make a model of solvated ion as one sodium atom surrounded by two crossing diglyme molecules, i.e., \ceNa(digl)2 (for reference, one Na and one diglyme \ceNa(digl)1, and one Li and two diglymes \ceLi(digl)2 are also considered). Note that different solvated ion models with two alkali atoms and two diglyme molecules \ceNa2(digl)2 were also suggested and studied in other works [27, 28, 29]. Regarding the graphite hosts with different carbon numbers $n$, we use ($3\times 3$), ($4\times 3$), ($4\times 4$), ($\sqrt{12}\times\sqrt{7}$), and ($\sqrt{12}\times\sqrt{13}$) cells, which correspond to $n=18$, 24, 32, 20, and 28, respectively. The solvated ions are assumed to be intercalated into $AA$-stacked graphene layers, as shown in Fig. 1.

The ab initio DFT calculations were carried out using the pseudopotential-pseudo atomic orbital (PAO) method as implemented in the SIESTA package [41]. We have tested different exchange-correlation (XC) functionals, such as the Perdew-Zunger formalism within local density approximation (LDA) [42], and the Perdew-Burke-Ernzerhof (PBE) [43], the PBE for solid (PBEsol) [44] and the Lee-Yang-Parr (LYP) [45] functionals within generalized gradient approximation (GGA). It can be thought that the dispersive van der Waals (vdW) interactions between graphene sheets may play an important role in Na(digl)-cointercalated graphite as well as pristine graphite [46, 47, 48]. To estimate such vdW interactions, we used the semi-empirical Grimme’s approach [49], where the vdW energy is expressed as follows,

[TABLE]

where $r_{ij}=|\textbf{r}_{i}-\textbf{r}_{j}|$ the distance between the $i$-th and $j$-th atoms, $C_{6}^{ij}=\sqrt{C_{6}^{i}C_{6}^{j}}$ being $C_{6}$ the dispersion coefficient of atom, $R_{ij}^{0}=R_{i}^{0}+R_{j}^{0}$ being $R^{0}$ the vdW radius of atom, and $s_{6}$ and $d$ are the scaling and damping parameters, respectively. We took the values of $C_{6}$ and $R^{0}$ for all atoms from Ref. [49] and used the default values of 1.66 for $s_{6}$ and 20.0 for $d$ given in the SIESTA code. Note that some papers reported negligible vdW effect for $t$-GICs [37, 46], and the localized basis sets like PAO might not properly describe the vdW force, compared to expanded basis sets like plane wave [50].

For all the atoms, standard double-$\zeta$ polarized (DZP) basis sets were generated with 300 meV energy shift for orbital-confining cutoff radii and 0.25 split norm for the split-valence of basis. We constructed the Troullier-Martins [51] norm-conserving pseudopotentials of atoms using the ATOM code provided in the package, where the valence electron configurations are Na-3s13p03d04fs0, Li-2s12p03d04fs0, C-2s22p23d04fs0, O-2s22p43d04fs0, H-1s12p03d04fs0. The non-linear core correction was applied to Na with a core radius of 1.50 Bohr and Li with 0.7 Bohr. The transferability testings were performed for all the generated pseudopotentials. Spin polarization was considered. We have conducted the convergence test of total energy with respect to the major computational parameters$-$plane-wave cutoff energy and $k$-grid cutoff length. After that, they were set to be 300 Ry and 10 Å, respectively, which guarantee an accuracy of total energy of 5 meV per cell (Fig.S1†). The crystalline lattice parameters and internal atomic coordinates were fully relaxed until the forces on each atom converged to within 0.02 eV Å*-1* and the stress on lattice to within 0.01 GPa.

To estimate the possibility of compound formation, we calculate the total energy difference between the compound and its reactants, as the entropic and enthalpic contributions might be negligible. First, upon the formation reaction of solvated sodium or lithium Ax(digl)y ($\text{A}=\text{Na}$ or Li, $x,y=1$ or 2) complexes,

[TABLE]

the binding energy can be calculated using the total energies of isolated Ax(digl)y complex, diglyme molecule and alkali atom,

[TABLE]

The cubic supercells with a lattice constant of 40 Å were used to model the isolated subjects. Then, upon the intercalation of the solvation complex into graphite to form Ax(digl)yCn $t$-GICs,

[TABLE]

the intercalation energy can be calculated as follows [52],

[TABLE]

The negative binding or intercalation energy indicates that the formation of compound is , i.e., spontaneous. The strength of binding between graphene sheets in graphite and GICs can be estimated by the exfoliation energy per carbon atom,

[TABLE]

where $n_{\text{C}}$ is the number of carbon atoms in graphene sheet, $d_{\text{e}}$ the equilibrium interlayer distance and $E(d_{\text{i}}=\infty)$ was replaced by $E(30$ Å$)$ in this work due to little change in the energy beyond $d_{\text{i}}=30$ Å. In these energetic estimations, basis set superposition error (BSSE) was not corrected, because it was reported that the GICs with a strong binding between the intercalant and the graphene sheet are not proper to apply the BSSE correction method [52].

2 Results and discussion

2.1 Assessment of XC functional and vdW inclusion

It was in general accepted that the lattice-related and energetic properties of graphite and GICs depend sensitively on the choice of XC functional and whether the inclusion of vdW correction or not, as evidenced in our previous works [48, 37]. Therefore, it is necessary to evaluate the reliability of various XC functionals and vdW inclusion. We computed the energetic quantities such as binding, intercalation, and exfoliation energies, which are indicative factors whether the compound can be formed spontaneously or not, employing one LDA functional (PZ) and three GGA functionals (PBE, PBEsol and LYP) with and without vdW correction.

We begin with an analysis of binding energy calculation of free Ax(digl)y ($\text{A}=\text{Na}$ or Li, $x,y=1$ or 2) complexes. As it was confirmed from experiments that the alkali ions in the crystalline phases of A-glyme compounds are coordinated by six oxygen atoms [26, 38, 39], two diglyme molecules per alkali atom were suggested to reach a stable solvation shell in dilute solution [26]. Moreover, it was shown from theoretical studies on diglyme-solvated alkali ions that the solvation shell still consists of two diglyme molecules in the vicinity of the electrode surface [40]. Therefore, we considered mainly Na(digl)2 complex in this work and additionally Na(digl) and Li(digl)2 complexes for reference, although Na2(digl)2 complexes are not ruled out [28]. Fig. 2 shows the molecular structures of diglyme molecule and these solvation complexes optimized by PBEsol+vdW method. Note that the initial configurations for geometrical optimizations were selected from conformation searching with the systematic grid scan method. While Na and Li binding energies in Na(digl) and Li(digl)2 are negative from all XC functionals with and without vdW correction as listed in Table 1, it was observed that the binding energy of Na to two diglyme molecules in Na(digl)2 are negative only from LDA, LDA+vdW and PBEsol+vdW functionals. When including the vdW correction, the binding energies become slightly lower than those without vdW energy addition. Regarding the cohesive energy ($E_{\text{coh}}$) of Na and Li metals, PBE and PBEsol could give the closest values to the experiment, with a slightly lowering by including the vdW interaction.

We then tested graphite unit cell involving 4 carbon atoms, of which lattice constants and exfoliation energy per carbon atom were computed with the different XC functionals. Here, the equilibrium interplanar lattice constant $c$ and the exfoliation energy $E_{\text{exf}}$ were determined by calculating the interlayer potential energies as a function of the length of lattice constant $c$ and then extracting the position where the interlayer potential energy is in minimal. The in-plane lattice constant $a$, governed by the C$-$C covalent ($\sigma$ and $\pi$) bonding, is well reproduced with all XC functionals as listed in Table 1. However, the interlayer distance $c$ along the direction of vdW interactions depends evidently on the choice of XC functional as listed in Table 1 and shown in Fig. 3.

All the XC functionals except LYP, with localized DZP$-$PAO basis sets, yielded distinct underestimation of lattice constant $c$ unlike plane-wave basis sets, with which PBE gave the significant overestimation, LDA the underestimation, and PBE+vdW could give proper value [36, 48, 37]. Inclusion of vdW correction lowered the lattice constant by $\sim$0.12 Å, and the exfoliation energy by 15$\sim$20 meV. It was found that LYP+vdW yields the closest interlayer lattice constant $c$ of 6.730 Å to the experimental value of 6.705 Å [53] and LYP the best exfoliation energy of $-47$ meV to the experimental value of $-35\sim-52$ meV [54, 55]. Tasaki [52] reported the calculated value of $-69$ meV with PBE+vdW using double numerical polarization atomic orbital basis sets, while Dappe et al. [57] obtained $-61\sim-74$ meV through a combination of the DFT approach and the perturbation theory. Their results are close to those from PBE and PBEsol functionals in this work.

For further assessment of the XC functionals (and vdW correction) and getting an insight of $t$-GIC formation (and stability), the intercalation and exfoliation energies were calculated. Among different configurations, the lowest energetic configuration was determined by comparing their total energies (see Fig.S2†). As shown in Figs. 4(a)-(c), all the XC functionals produce the negative intercalation energies for almost Na(digl)Cn and Li(digl)2Cn compounds, while only LDA, LDA+vdW and PBEsol+vdW functionals could give the negative values for Na(digl)2Cn compounds with some $n$ values. It should be emphasized that with PBEsol+vdW functional the lowest intercalation energy in Na(digl)2Cn was found at $n\approx 21$ (i.e., the most stable compound is Na(digl)2C21), corresponding to specific capacity of $\sim$97 mAh g*-1* close to the experimentally identified capacity of $\sim$100 mAh g*-1* [26]. In addition, when including the vdW correction, the intercalation energies become slightly lower in these compounds. We see from Figs. 3 and 4(d)-(f) that the interlayer distance $d_{\text{int}}$ is calculated to become larger and the magnitude of exfoliation energy smaller going from LDA to PBEsol, PBE and LYP. These indicate that LDA induces more binding than GGA and the inclusion of vdW correction leads to stronger binding in Na/Li $t$-GICs.

The electrode voltage represents one of crucial electrochemical properties to determine the suitability of material as an electrode for NIBs. Assuming that during the charge process Na-GICs transform from the compound of higher carbon atom concentration (i.e., lower Na concentration compound) to the compound of lower carbon atom concentration (higher Na concentration), we convert the $t$-GICs from Na(digl)xCn to Na${}_{x_{i}}$(digl)${}_{2x_{i}}$C32 with $x_{i}=32/n$, treating C32 $t$-GIC as the reference compound (starting compound). Then, the electrode potential $V_{\text{el}}$ with respect to Na/Na+ can be calculated as follows [58],

[TABLE]

In Table 1, the electrode voltages calculated with $x_{i}=32/32$, and $x_{j}=32/24$ for Na $t$-GICs and $x_{j}=32/20$ for Li $t$-GICs are presented. It was found that PBEsol+vdW give the closest electrode voltage of 0.16 V for Na(digl)2Cn to the experimental values of 0.15 (50th cycle) $\sim$0.45 V (first cycle) at the capacity of 86 mAh g*-1* corresponding to C24 compound [26].

The most reliable XC functional was determined as a compromise to simultaneously describe the binding properties of graphite and sodium metal, the intercalation energies of Na(digl)xCn compounds and more importantly, the electrode voltage of Na(digl)2Cn that greatly dominates the electrochemical performance of NIB. As a result, PBEsol+vdW was chosen as the most reliable XC functional since it results in excellent agreement with the electrode voltage of Na(digl)2Cn, negative intercalation energy and reasonable binding properties of graphite and Na (Li) metal. This is agreed with that PBEsol can reliably predict material properties of Na-amine $t$-GICs with the plane wave method [37].

2.2 Variation in structural and energetic properties

We compared the binding properties of different solvated ions Na(digl), Na(digl)2 and Li(digl)2, which are formed in ether-based electrolyte solution, intercalate into graphite and diffuse in the interlayer space between graphene sheets. From the calculated binding energies as listed in Table 1, it can be said that the Na atom is less binding to two diglyme molecules in Na(digl)2 ($E_{\text{b}}=-$0.11 eV) than Na binding to one diglyme molecule in Na(digl) ($-$0.32 eV), and also slightly less than Li binding to two diglyme molecules in Li(digl)2 ($-$0.13 eV). Considering that Na atom is coordinated by three or six oxygen atoms, the average Na$-$O bond length in Na(digl)2 (2.617 Å) is estimated to be much shorter than that in Na(digl) (2.931 Å), and also slightly shorter than that in Li(digl)2 (2.697 Å), as depicted in Fig. 2. This indicates that stronger binding of Na to diglyme molecule does not induce smaller Na$-$O bond length, implying neither ionic nor covalent bonding. The bonding is actually coordinate bonding.

When the solvated A(digl)x ion intercalates into graphite, the binding between graphene sheets becomes stronger due to the additional interaction between the intercalants and the graphene sheet, as evidenced by that the exfoliation energies in these $t$-GICs are calculated to be negatively larger than that in graphite ($-$95 meV per carbon atom). Of the three different $t$-GICs at $n=20$, as listed in Table 1, Na(digl)2C20 ($-$255 meV) has the smallest exfoliation energy in magnitude, while Na(digl)C20 ($-$265 meV) the largest, and Li(digl)2C20 ($-$259 meV) in-between, indicating that the binding between graphene sheets in Na(digl)2Cn is weaker than in Na(digl)Cn and also Li(digl)2Cn. However, we should note that the differences between binding energies are not very significant, and moreover, due to an uncertainty of selecting $d_{\text{i}}=\infty$ ($\sim$30 Å in this work) in Eq. 6 and some degrees of errors in the calculations, the differences are almost negligible. On the other hand, as can be seen in Fig. 5, when changing the number of carbon atoms of graphene sheet in these $t$-GICs, the exfoliation energy changes wavily as the interlayer distance dose so. In the cases of intercalants with two diglymes, the interlayer distance is minimal at $n=24$, while the exfoliation energy is also minimal in magnitude, indicating an unrelation between binding strength and interlayer distance.

For further insight into the structural variation upon the intercalation of the solvated Na(digl)2 ion into the graphite, we considered the average Na$-$O bond length and the relative volume expansion ratio, $r_{\text{vol}}=(V-V_{0})/V_{0}\times 100$ (%), together with the interlayer distance. As listed in Table 2 and shown in Fig. 1, the average Na$-$O bond lengths (2.253$\sim$2.315 Å) in Na(digl)2 intercalants that are constrained between graphene sheets are remarkably contracted compared to that (2.617 Å) in free solvated complexes, indicating an enhancement of Na$-$diglyme binding upon its intercalation into graphite (see Fig.S3 for Na(digl)Cn†). This might be due to somewhat complicated interactions between Na atom, diglyme molecules and graphene sheets. It is worthy noting that such enhancement of Na binding to the solvent can prevent decomposition of solvated ion, resulting in formation of a negligible solid-electrolyte interphase (SEI) film on the graphite surface, enabling Na+-solvent migration into the graphite lattice [27]. On the other hand, the interlayer distances $d_{\text{int}}$ are almost over 10 Å, much larger than that of graphite (3.353 Å). Compared with experiments, Jache and Adelhelm reported a 10% expansion for Li intercalation and a 15% expansion for sodium [26], while Kim et al. [28] reported 11.62 Å for Na2(digl)2 cointercalated compound. Therefore, our calculated interlayer distances are reasonable compared with the latter experiment, though the number of intercalated sodium atoms is different. At this point, we note that the interlayer distances are more or less constant with the maximum relative change of 3% as increasing the number of carbon atoms of graphene sheet. This implies that, once $t$-GIC is formed in the first cycle, there would be no serious change of interlayer distance and volume as shown in Table 2 (Table S1† for Na(digl)Cn and Li(digl)2Cn), during further cycles.

As listed in Table 2 and discussed above, the intercalation energies for Na(digl)2Cn compounds are negative only for $n=20\sim 28$, as far as we rely on the PBEsol+vdW functional calculations. This is agreed well with the experiment that Na(digl)2C20 is preferably formed [26]. When one diglyme molecule instead of two molecules or Li atom instead of Na atom intercalates into graphite to form Na(digl)Cn or Li(digl)2Cn, the intercalation energies become negative with larger magnitude for almost $n$ values by all the XC functionals. This indicates that it is easier to form Na-diglyme intercalated graphite with one diglyme molecule than with two diglyme molecules, and to form Li-GICs than Na-GICs. Nevertheless, the change tendency of intercalation energies in Li(digl)2Cn as increasing the number of carbon atoms is similar to Na(digl)2Cn. To check the suitability of these compounds for electrode, we computed the electrode potential as a function of capacity by applying Eq. 7, as shown in Fig. 6. It is found that Na(digl)2Cn has the lowest electrode potential below 0.5 V for the range of capacities considered in this work, indicating it’s the most suitable anode material among three kinds of $t$-GICs.

2.3 Ionic diffusion and electronic conduction

The fundamentals of sodium ion transport such as ionic diffusion pathways and their activation barriers are an important aspect for charge-discharge kinetics, and in some cases this is an critical issue for the cycle stability. In fact, the activation energy for sodium ion diffusion during cycling is a determinant factor how fast charge process takes place, and diffusion pathways should be short to enable fast migration of Na ions during battery operation. In this sense, we inquire into the diffusion process of the solvation ions with a calculation of activation energies by applying the nudged elastic band (NEB) method [59] as implemented in Python script Pastafarian111This code was originally developed by J. M. Knaup, and we modified the code to debug some minor errors and allow parallel running with a permission. in conjunction with SIESTA code. We typically choose A(digl)xC20 compounds, and utilized doubled unit cells in the $b$-direction (doubled unit cell) to allow migrations of only alkali ions and in both $a$ and $b$-directions (quadrupled unit cell) for A(digl)x ions from one position to another identical position. The number of NEB images to discretize the path are 15$\sim$17, and atomic relaxations were allowed during the NEB search.

In Fig. 7(a) and (c), we show the calculated activation energies for migrations of the Na(digl)2 solvated ions along the shortest pathways, with a depiction of geometries of some transition states. The activation energy of this migration in Na(digl)2C20 compound was calculated to be about 0.40 eV, which is slightly higher than that in Na(digl)C20 (0.28 eV) and that in Li(digl)2C20 (0.25 eV) (see Fig.S4 for Na(digl)C20 and Li(digl)2C20†). The latter comparison indicates that Li solvated ions diffuse more smoothly than Na solvated ion, which is contrast to the previous ab initio calculations for different $t$-GICs [37] as well as for alkali atom $b$-GICs [23]. It can be said in safe that the $t$-GICs investigated in this work can offer reasonably fast charge-discharge rate since the activation energies are comparable with Li diffusion barrier in graphite (0.40 eV) [23] and Na diffusion barriers in layered anode materials like MoS2 (0.28 eV) [35] and cathode materials like NaCoO2 (0.3 eV) [60]. We also looked over migration of only sodium ion, which might be helpful for drawing an inference whether the formation of SEI film is easy or not. As shown in Fig. 7(b) and (d), the activation energy for this migration was calculated to be much high over 3 eV, indicating that the SEI film formation is not easy.

We then turn to a question how electrons are conducted inside materials during charge-discharge process. This is an important issue for understanding the battery operation since the electrons are transmitted to/from the collector the moment the solvated sodium ions de/intercalates out from/into the graphite. To qualitatively estimate the electron conductance of these materials, we analyse the density of states (DOS) of the relevant compounds. In Fig. 8, we show the total DOS (TDOS) and atomic resolved partial DOS (PDOS) of major spin component (those of minor spin are almost identical to those of major spin) of graphite, isolated diglyme molecule, and Na(digl)2Cn compounds. It was well known that graphite has excellent in-plane electric conductance (no inter-plane conductance), whose charge carrier is just $p_{z}$ valence electron of carbon [21, 37]. When the solvated sodium ion intercalates into graphite, the in-plane electronic conductance is retained well, even better than graphite, as the density of electronic states of carbon around Fermi level in $t$-GICs is relatively higher than in graphite. In Fig. 8, we see the unoccupied state of sodium around 3 eV over Fermi level, indicating that sodium atom donates its $s$ electron to become a cation (electron donor) while carbon atoms in graphene layer accept it to occupy anti-bonding $\pi^{*}$ ($p_{z}$) orbital (electron acceptor). We also observe significant orbital overlaps between graphene layer and diglyme below and above Fermi level, indicating hybridization between $\sigma$ bonding orbitals of carbon atoms in graphene layer and molecular orbitals of diglyme molecule, resulting strong interaction between the electrons. Similar characteristics is observed in Na(digl)Cn and Li(digl)2Cn compounds (Fig.S5(a) and (b)†).

2.4 Formation mechanism of compounds

It is generally accepted that charge transferring is occurred in the formation of GIC from its components. Therefore, consideration of charge transferring is indispensable to make it clear the nature of chemical bonding and thus the formation mechanism. As discussed above, this could be done by analysing the electronic structure such as DOS. To get more intuitive insight into the charge transfer, we investigate a distribution of the charge density difference in Na-solvent cointercalated graphite compounds as calculated by using the following equation,

[TABLE]

Fig. 9 shows the electronic density accumulation (brown colour) and depletion (green colour) in Na(digl)C20 and Na(digl)2Cn with $n=$18, 20, 24, 28, 32 (isosurface of $\pm$0.0016 $|e|$ Å*-3*). In overall, the Na atom has lost its valence electronic density, and the solvent molecules have parts of electron loss and gain together, while the graphene layer has gained the electrons. Observing the distribution of charge density difference does not tell whether the solvent molecule in the aggregate has lost electrons or not. In the case of Na(digl)xCn, some part of electrons has transferred from Na atom to the bonding carbon rings, and this induces the formation of ionic bond between Na ions and carbons, with electrons donated from C$-$C bonds in the graphene layers (Fig.S6(a) and (b)†). It is found that as increasing the number of carbon atoms in graphene layer the amount of accepted electrons becomes smaller. It should be noted that in the case of two diglyme molecules intercalation the Na atom could not interact directly with graphene layer due to a screening by atoms of solvent molecule.

We performed Mulliken population analysis of these compounds to quantitatively estimate the amount of transferred charge, although it is said that such analysis could not yield accurate values for atomic or molecular charges. It turns out that the Na atom in Na(digl)xCn compounds donates a relatively large amount of electrons (0.142$\sim$0.191 $|e|$), but contrastingly the Li atom in Li(digl)2Cn receives electrons (0.056$\sim$0.074 $|e|$) due to probably weaker electropositivity of Li than Na atom (Table S2†). In the entirety, diglyme molecules just release electrons, of which the amount increases in general as increasing the number of carbon atoms in graphene layer, indicating that diglyme molecule is an electron doner. The amount of released electrons from diglyme in Na(digl)2Cn (0.138$\sim$0.233) is larger than in Na(digl)Cn (0.058$\sim$0.122 $|e|$), but smaller than in Li(digl)2Cn (0.342$\sim$0.421 $|e|$). Inside the diglyme molecules, carbon atom releases electrons of $\sim$0.60 $|e|$ in average, hydrogen atom accepts majority of electrons of $\sim$0.25 $|e|$, and oxygen atom accepts minority of electrons of $\sim$0.04 $|e|$, indicating the formation of covalent bond between C and H, and ionic bond between C and O atoms inside the molecule. On the other hand, the amount of accepted electrons per carbon atom of graphene layer decreases in general as increasing the number of carbon atoms, and that in Na(digl)2Cn (0.017$\sim$0.012 $|e|$) is larger than in Na(digl)Cn (0.009$\sim$0.010 $|e|$) but similar to Li(digl)2Cn (0.017$\sim$0.011 $|e|$).

Finally we try to describe the formation mechanism of these $t$-GICs. When forming graphite from carbon atoms, the four valence electrons of carbon make three $sp^{2}$ hybrid orbitals and one $p_{z}$ orbital. The three $sp^{2}$ electrons form strong $\sigma$ bonding in the carbon hexagon on the plane of graphene sheet, while the $p_{z}$ electron forms weak delocalized $\pi$ bonding with those from other carbon atoms in vertical direction to the hexagon plane. Meanwhile, the graphene layers are weakly binding through the vdW interactions, and thus atoms, molecules or ions can be easily intercalated between the graphene sheets. Depending on the electronegativity of the intercalant, the graphene layer can release the $\pi$ electrons to the intercalant or accept electrons from the intercalant, forming ionic bonding through $\pi$-conjugate interaction between graphene layer and intercalant (covalent bonding can also be formed). When only Na atom intercalates into graphite, the electrons move from the Na atom (strong electropositive element) to the graphene layer but the amount of transferred electrons is relatively so small that the $\pi$ interaction is not enough strong to form a stage-I $b$-GIC. When the diglyme molecules intercalate with Na atom into graphite, the diglyme molecule as well as Na atom donates electron to the graphene layer as discussed above, and therefore the $\pi$-conjugate interaction becomes enough strong to form a stable stage-I GIC. In the case of Li atom, the amount of released electrons from diglyme is so much that the Li atom as well as graphene layer should receive electrons in contrast to the Na case. On the other hand, the graphene layer affects the binding of Na atom to diglyme molecule, such that the Na atom could not easily escape but diffuse in the aggregate fast into the graphite inside. In this context, the structure of cointercalants surrounded by two diglyme molecules could be more favourable than one diglyme molecule. Consequently, the interactions between diglyme molecule and graphene layer play a major role in the formation and operation of these $t$-GICs, while the interaction between sodium atom and graphene layer is rather indirect.

3 Conclusions

We investigated the cointercalation of sodium and diglyme molecule into graphite aiming to clarify the atomistic structure, energetics and electrochemical properties for Na-ion battery application with ab initio DFT calculations. The ternary graphite intercalation compounds, Na(digl)2Cn with $n=$18, 20, 24, 28, and 32, were considered with the graphene unit cells of $(3\times 3)$, $(\sqrt{12}\times\sqrt{7})$, $(4\times 3)$, $(\sqrt{12}\times\sqrt{13})$, and $(4\times 4)$. In addition, Na(digl)Cn and Li(digl)2Cn compounds were also considered for reference. We first determined the reliable exchange-correlation functional in connection with vdW inclusion that can reproduce the experimental results for structural and binding properties of graphite and Na/Li metal and for electrode potential of NIB, and that can give the negative intercalation energy of $t$-GICs, as PBEsol+vdW functional. With this functional, the variations in structural properties such as interlayer distance, average Na$-$O bond length and relative volume change, and energetic properties such as exfoliation, intercalation energies and electrode potentials, have been studied systematically, finding out that Na(digl)2Cn has the lowest intercalation energy at $n\approx 21$ corresponding to the capacity of $\sim 100$ mAh g*-1*. We found that the solvated ion can diffuse relatively fast with the migration barrier of $\sim$0.40 eV and the electron conductance can be enhanced after cointercalation when compared with pristine graphite. The calculation results clarify the mechanism of the formation of these $t$-GICS and reveal new prospects for tailoring innovative anode materials of NIBs based on graphitic materials.

Acknowledgments

This work was supported partially from the State Committee of Science and Technology, Democratic People’s Republic of Korea, under the fundamental research project “Design of Innovative Functional Materials for Energy and Environmental Application” (No. 2016-20). The calculations in this work were carried out on the HP Blade System C7000 (HP BL460c) that is owned and managed by Faculty of Materials Science,

Kim Il Sung University.

Appendix A. Supplementary data

Supplementary data related to this article can be found at URL.

Notes

The authors declare no competing financial interest.