Enabling Enhanced Lithium Ion Storage Performance of Graphdiyne by Doping with Group-15 Elements: A First-Principles Study.

As a typical two-dimensional material possessing sp and sp2 hybrid orbitals, graphdiyne (GDY) and its derivatives have been proposed as an attractive candidate for high-performance lithium ion batteries (LIBs). In this work, an advanced GDY LIB electrode is designed by doping with group-15 elements. With the aid of first-principles simulations, the geometric properties, electronic structures, theoretical storage capacities, open-circuit voltages, and diffusion path of Li atoms on doped GDY are comprehensively investigated. Specifically, 14 different adsorption sites are proposed, most of which are situated out of plane of the carbon network, resulting from the out of plane Pz orbitals of conduction band minimum and valence band maximum. Among the five doped GDY, phosphorus-doped graphdiyne (P-GDY) exhibits prominent lithium ion storage behavior, i.e., the maximum theoretical capacity is 1949 mA·h·g-1, which is ∼2.6 times higher than that of GDY. Moreover, calculation results in terms of the in-plane migration of lithium ion on P-GDY indicate that Li atoms prefer to diffuse across the carbon network (with a moderate barrier of 0.46 eV) rather than directly through the middle of the hexagonal aperture (with a higher barrier of 1.78 eV). Thus, this approach provides novel insights into the Li ion storage properties of group-15 element-doped GDY from the prospect of theoretical calculations, which would be useful to guide the future design of high-capacity GDY anodes for LIBs.

Introduction

During the pn class="Chemical">ast decades, with the growinclass="Chemical">pan>g market demands for portable energy storage devices and new technologies, ass="Chemical">lithiumn> ion batteries (<ass="Chemical">spn>an class="Chemical">LIBs) have played an increasingly important part in the broad energy storage market.[1] Among the various kinds of anode materials for LIBs, graphite has been the first choice for commercial applications due to its advantages of low lithiation potential and high stability.[2] As a layered material composed of sp2 hybrid orbitals, graphite shows a typical pattern for the storage of lithium atoms. Lithium atoms are embedded between two graphene layers and reside in the center of the phenyl ring, being stored with the occupation pattern of LiC6, which gives rise to a theoretical capacity of 372 mA·h·g–1.[3] Nevertheless, with the increasing desire for higher-capacity batteries, traditional graphite materials with lower capacity and limited promotion space has been gradually phased out. Novel anode materials have been intensively investigated, such as graphite derivatives like carbon nanotubes, graphene, and non-carbon materials like transition metal oxides, sulfides, and silicon monoxide.[4−8]

<n class="Chemical">span class="Chemical">Graphyne (GY), as one of the <ass="Chemical">span class="Chemical">graphite derivatives, is an sp–sp2 hybrid planar carbon sheet with two benzene rings connected via acetylene linkages exhibiting a two-dimensional layered structure.[9] Owing to the special sp–sp2 hybrid structure, GY exhibits excellent characteristics, such as high π-conjunction, tunable electronic properties, two-dimensional plane framework, and uniformly distributed pores,[10,11] making it become a promising anode material with high lithium storage capacity and transport capacity. GY was first synthesized by Li et al. in 2010 and was named graphdiyne (GDY) since its units contain two acetylenic groups between two adjacent benzene rings.[12] These two acetylenic chains make GDY own an 18-C triangular aperture, in contrast with the 12-C triangular aperture of GY. Lithium ions in GY and GDY with large apertures can migrate through these pores vertically to the two-dimensional plane easily,[13] whereas lithium ions can migrate only within graphene interlayers in graphite.

Based on the large aperture of GDY, <span class="Chemical">lithium atoms can be stored in the form of LiC3, thus resulting in a storage capacity of 744 mA·h·g–1, which is twice higher than the theoretical capacity of <ass="Chemical">span class="Chemical">graphite.[14] In the further study of GDY, it is found that due to the unique preparation process of GDY, the bottom-up synthesis method can be adopted for doping various heteroatoms into GDY qualitatively and quantitatively, which is conducive to improve the electrochemical behavior of GDY. To date, a variety of elements and groups have been successfully doped into GDY including H, Cl, F, N, CH3, NH, etc.,[15−20] and greatly changed the electrochemical performance of GDY. Both the heteroatoms substituting the carbon atoms in the benzene ring and the impurity groups replacing the acetylene chain make the morphological characteristics of GDY change from a combination of benzene rings and triangular apertures into a nanosheet consisting of benzene rings and hexagonal apertures. The large-area hexagonal apertures allow GDY to have much more lithium adsorption sites and higher interionic mobility.

Among the ocean of heteroatoms and groups doping in GDY, <span class="Chemical">nitrogen-doped graphdiyne (<ass="Chemical">span class="Chemical">N-GDY) manifests distinctive performance. Nitrogen and carbon are adjacent in atomic number with a similar atomic size; consequently, the substitution of a carbon atom with a nitrogen atom would lead to a slight deformation in the configuration of GDY. Besides, nitrogen doping would provide an extra electron to GDY, making N-GDY an n-type semiconductor and significantly enhance the electronic transport properties of GDY. Since the successful preparation of GDY, a variety of nitrogen-doped GDY with remarkable electrochemical performance have been synthesized like bamboo shoots after a spring rain.[21−24] However, the exact doping location and content of nitrogen doping in these samples are fuzzy and indistinct, and qualitative and quantitative doping in GDY is technically challenging; accordingly, it is difficult to distinguish the influence of N-doping on the macroscopic properties of GDY from a microscopic perspective. In 2019, Yang et al. precisely controlled the location and quantity of nitrogen doping in N-GDY through a bottom-up method.[18] The morphology of the synthesized N-GDY has large-area apertures (∼6 times larger than the apertures in GDY), contributing to the high specific surface area, high lithium storage capacity, and a convenient ion migration path for N-doped GDY.

In this article, we attempt to examine the effect of doping with group-15 elements on Li adsorption and diffusion characteristics of GDY via first-principles simulations. Calculations of the geometric properties, electronic structures, adsorption configurations, maximum storage capacity, and diffusion path are implemented. By directly comparing the geometric and electronic structures of primitive and group-15 element-doped GDY (represented by <span class="Chemical">M-GDY), we obtain a comprehensive understanding of the effects of doping with group-15 elements on the atomic scale. Our simulation results suggest that <ass="Chemical">span class="Chemical">P-GDY reveals outstanding properties in Li storage and mobility, with a maximum theoretical storage capacity of 1949 mA·h·g–1, which is ∼2.6 times the capacity of GDY.

Theoretical Methods

First-principles calculations based on density functional theory (DFT)[25] are carried out to investigate the geometric and electronic properties of group-15 element-doped GDY, as implemented in the Vienna Ab initio Simulation Package (VA<span class="Chemical">SP) using a plane wave basis, the Perdew–Burke–Ernzerhof (PBE) density functional, and the projector augmented wave (PAW) potentials.[26,27] Van der Waals calculations of Grimme’s D3 with the Becke–Johnson damping method (DFT-D3) is employed to account for the diass="Chemical">spersion correction.[28] An energy cutoff of 550 eV is utilized, and a 3 × 3 × 1 Monkhorst–Pack grid centered at the γ point is chosen for sampling the Brillouin zone. The convergence threshold for the energy of electronic self-consistent iteration is set to 10–5 eV, and the residual forces are taken as 0.01 eV/Å until all atomic positions are fully optimized. Along the z-direction, the vacuum region is set to 12 Å to avoid the interaction between the adjacent cells. To verify the thermal and dynamic stability of the <ass="Chemical">span class="Chemical">M-GDY structure, photon spectra are studied with the phonopy code and ab initio molecular dynamics (AIMD) simulations are performed. For ground-state AIMD simulations of the primitive and Li-adsorbed systems, the model systems first go from 0 to 300 K by velocity regulation within 500 fs, equilibrate at 300 K for 1 ps thereafter, and are evolved with a 6 ps MD run in the NVE thermodynamic ensemble.

The adsorption energy (Eads) of <span class="Chemical">lithium atoms onto <ass="Chemical">span class="Chemical">M-GDY is described by the equation belowwhere n is the total number of Li atoms adsorbed on M-GDY, ELi, Eprimitive, and Etotal represent the free energy of a single lithium atom, M-GDY, and Li-adsorbed M-GDY, respectively. In this definition, when Eads is larger than the cohesive energy of the lithium cluster (1.63 eV), the adsorption process is spontaneous and stable. The larger the value is, the stronger the adsorption is.

The theoretical capacity of <n class="Chemical">span class="Chemical">M-GDY is calculated by the following equationwhere n denotes the total number of Li atoms adsorbed on <ass="Chemical">span class="Chemical">M-GDY, F is the Faraday constant, and M is the molar mass of the selected cell.

The search for transition states is mapped out with the climbing image nudged eln class="Chemical">astic band (<n class="Chemical">span class="Chemical">CI-NEB) method using the VTST tools in VA<ass="Chemical">span class="Chemical">SP.[29−32] When evaluating the electronic band structure, the k-points are generated with the Seekpath module.[33,34] The data analysis and image processing are performed using VASPKIT and VESTA.[35,36]

Results and Discussion

Structural and Electronic Properties of GDY and M-GDY

Structural Properties of GDY and M-GDY

To look into the structural changes in GDY before and after doping with group-15 elements, we have directly compared the optimized geometries of primitive GDY, H-doped GDY, and <span class="Chemical">M-GDY, as shown in Figure . The configurations of doped GDY are constructed based on the experimental results, in which the doping locations of group-15 elements are selected to be consistent with <ass="Chemical">span class="Chemical">PY-GDY in ref (18). In GDY, the benzene ring is linked to six adjacent benzene rings via two acetylene chains. However, when hydrogen atoms and the group-15 elements substitute the carbon atoms (in the meta-position) of the benzene ring, the number of the nearest neighbor of the benzene ring decreases to three, inducing dramatic change in the structure of doped GDY, thus forming a six-times larger aperture in contrast with the primitive GDY.

Figure 1

Top view of optimized geometries of (a) primitive GDY; (b) hydrogen-doped GDY; and (c) the group-15 element-doped GDY. Color scheme: C – dark gray, H – green, group-15 elements (N, P, As, Sb, and Bi) – red.

Top view of optimized geometries of (a) primitive GDY; (b) <n class="Chemical">span class="Chemical">hydrogen-doped GDY; and (c) the group-15 element-doped GDY. Color scheme: C – dark gray, H – green, group-15 elements (N, P, As, Sb, and Bi) – red.

We have further investigated the influence of doping with different elements (varying from H, N, P, As, Sb to Bi) on the geometric properties of GDY. We notice that GDY and H-GDY share similar lattice vectors; however, when group-15 elements were doped into GDY, the variations in the cell size of GDY become noteworthy. Compared to GDY and H-GDY, the length of ass="Chemical">N-GDY in the <ass="Chemical">span class="Species">a-axis direction is somewhat reduced, and the axial angle γ also becomes smaller. When the doped elements change gradually from smaller nitrogen to larger bismuth, the doping-induced in-plane lattice distortion becomes more and more remarkable. As shown in Table , with the increase in the atomic number of doped elements, the lattice vector a and axial angle γ of doped GDY increase (changing from 120 to 125.51°), while the lattice vector b decreases. In addition, the total area of the simulation cell also increases with the growing atomic size of the doped element. These changes make doped GDY become long and narrow in the top view, which may exert a non-negligible influence on lithium storage.

Table 1

Lattice Constants of Primitive GDY and Doped GDY

	GDY	H-GDY	N-GDY	P-GDY	As-GDY	Sb-GDY	Bi-GDY
a (Å)	16.33	16.33	16.10	17.37	17.72	18.37	18.55
b (Å)	16.33	16.38	16.33	16.17	16.14	15.96	15.96
γ (deg)	120.00	120.05	119.49	122.52	123.33	125.10	125.51
s (Å2)	231.03	231.93	229.01	237.02	239.02	240.75	241.81

To confirm the thermal and dynamical stability of <span class="Chemical">M-GDY, the phonon ass="Chemical">spectra are calculated, and ground-state AIMD simulations are performed. Figure S1 shows the calculated phonon ass="Chemical">spectra of <ass="Chemical">span class="Chemical">M-GDY; all vibration modes of M-GDY are real in the entire Brillouin zone, confirming the dynamical stability of M-GDY. The snapshots of top views and side views of M-GDY at t = 0, 3, and 6 ps in AIMD simulations (Figure S2) indicate that the M-GDY is thermodynamically stable during the entire 6 ps MD run. The time evolution of the M–C bond lengths of M-GDY in AIMD simulations is also monitored. As shown in Figure , the M–C bonds (M represents group-15 elements) are very stable, vibrating around their equilibrium values (1.35 Å for N–C bonds, 1.77 Å for P–C bonds, 1.91 Å for As–C bond, 2.11 Å for Sb–C bond, and 2.21 Å for Bi–C bonds, respectively). Therefore, M-GDY is thermally and dynamically stable.

Figure 2

Time evolutions of the M–C bond lengths of M-GDY in AIMD simulations.

Time evolutions of the M–C bond lengths of <span class="Chemical">M-GDY in AIMD simulations.

Electronic Properties of GDY and M-GDY

To further explore the electronic properties of GDY and <span cln class="Chemical">ass="Chemical">M-GDY, the static electronic band structures and projected density of states (PDOS) are calculated and investigated, as shown in Figure . The energy levels near the Fermi level of both primitive GDY and <n class="Chemical">ass="Chemical">spn>an class="Chemical">M-GDY are mainly contributed by the 2Pz orbitals of carbon atoms. From Figure a, the band structure of primitive GDY implies that the primitive GDY is a direct bandgap semiconductor with a bandgap of 0.50 eV at the γ point, which agrees well with the previous DFT calculation result of ∼0.46 eV.[37] When hydrogen atoms are doped on GDY, the conduction band minimum (CBM) and valence band maximum (VBM) of H-GDY separate from each other with a wide energy gap of 2.35 eV (Figure b). We notice that doping with group-15 elements also enlarges the bandgaps of GDY, 1.93 eV for N-GDY, 1.68 eV for P-GDY, 1.37 eV for As-GDY, 1.02 eV for Sb-GDY, and 0.84 eV for Bi-GDY, as shown in Figure c–g. This increase in the energy gap of doped GDY is similar to the situation of GDY doped with 3 nitrogen atoms (N3-GDY),[38] which is supposed to come from the dramatic structural changes induced by element doping. As mentioned before, the benzene ring of the doped GDY is connected with three neighboring benzene rings through two acetylene chains, giving rise to a six-times larger aperture with a weaker π-conjugation, hence leading to a wider bandgap. Compared with the band structure of H-GDY, M-GDY has a smaller bandgap. This is because group-15 elements provide extra electrons to M-GDY, resulting in the down-shift of CBM, thus resulting in a reduced bandgap and improved conductivity. Interestingly, we find that the bandgaps of M-GDY gradually decrease as the atomic number increases when the doping elements vary from nitrogen to bismuth. With the increase of atomic number, the electronegativities of group-15 elements decrease, which means that the overall tendency of group-15 elements to attract electrons from carbon atoms becomes weaker. Consequently, the bonding strength between the doped element and carbon atoms becomes weaker, suggesting that less energy is needed for a valence electron to break free from the binding bond.

Figure 3

The electronic band structure and PDOS of (a) primitive GDY, (b) H-GDY, (c) N-GDY, (d) P-GDY, (e) As-GDY, (f) Sb-GDY, and (g) Bi-GDY.

The electronic band structure and PDOS of (a) primitive GDY, (b) H-GDY, (c) <span class="Chemical">N-GDY, (d) <ass="Chemical">span class="Chemical">P-GDY, (e) As-GDY, (f) Sb-GDY, and (g) Bi-GDY.

From the PDOS of <n class="Chemical">span class="Chemical">N-GDYn>, it is seen that the valence band maximum (VBM) and conduction band minimum (CBM) are mainly contributed by the <ass="Chemical">spn>an class="Chemical">carbon atoms. However, for P-, As-, Sb-, and Bi-doped GDY, the density of states of CBM and VBM consist of contributions from both the carbon atoms and the doped element. Notably, the contributions from the doped element are increasing with increasing atomic number. Figure shows the wave functions of CBM and VBM frontier orbitals of M-GDY. Consistent with the band structure, the CBM and VBM of M-GDY are primarily constituted by the Pz orbitals, which are spreading out of plane of monolayer GDY. These Pz orbitals would offer several delocalized adsorbed locations, which is of great significance to determine the adsorption sites of lithium atoms. Interestingly, as the doping element changes from N to Bi, the distribution of wave functions over the benzene ring becomes more delocalized, especially at the doped site, suggesting that the bond strength between the doped atom and carbon atoms becomes weaker, which is the main reason for the decreased bandgap. However, the bandgap of doped GDY increases rapidly from 0.5 to 2.35 eV due to the dramatic structural change, when impurity elements vary from N to Bi, the bandgap decrease to an acceptable value of 0.83 eV and become suitable for electronic transmission gradually.

Figure 4

Wave functions of CBM and VBM frontier orbitals of M-GDY.

Wave functions of CBM and VBM frontier orbitals of <span class="Chemical">M-GDY.

Lithium Storage Properties of M-GDY

Adsorption Configurations of a Single Lithium Atom on M-GDY

To examine the adsorption properties of a single <span class="Chemical">lithium atom on <n class="Chemical">ass="Chemical">span class="Chemical">M-GDY, 14 individual lithium adsorption sites are considered, which can be classified into three categories: (i) lithium atoms adsorbed above the acetylene chain, with two hydrogen atoms on both sides (Figure a–e); (ii) lithium atoms adsorbed above the acetylene chain surrounded by two group-15 elements and two hydrogen atoms (Figure f–j); (iii) lithium atoms absorbed in and out of the monolayer M-GDY plane by bonding to group-15 elements (Figure k–n). From the side view of the optimized adsorption structures, most of the adsorption sites are situated outside the plane of M-GDY, which are determined by the out of plane Pz orbitals of the CBM and VBM.

Figure 5

Top and side views of fourteen individual lithium adsorption configurations of M-GDY.

Top and side views of fourteen individual <n class="Chemical">span cln class="Chemical">ass="Chemical">lithium adsorption configurations of <n class="Chemical">ass="Chemical">span class="Chemical">M-GDY.

The adsorption energies of the 14 different adsorption configurations are also calculated and collated in Table . Note that because some bonding sites are not stable, we have fixed the ass="Chemical">lithium atom and the surrounding <ass="Chemical">span class="Chemical">carbon atoms in calculating the free energies of 14 different adsorption configurations, leading to overestimated free energies, while all atoms are allowed to move freely in three dimensions for calculating the free energies of primitive M-GDY. Consequently, the adsorption energies, which are defined as Eads = (nELi + Eprimitive – Etotal)/n, are underestimated (by 0.3–0.4 eV). Although the adsorption energies at some sites are lower than 1.63 eV, which is the cohesive energy of bulk Li, a Li atom can still adsorb steadily at these sites. The adsorption energies of the 14 sites on N-GDY, P-GDY, and As-GDY are all larger than 1.1 eV, indicating stable adsorption of Li on these sites in N-GDY, P-GDY, and As-GDY. Intriguingly, we discover that as the doped elements change from N to As, the adsorption energies of locations above the acetylene chain (a–j) increase, while the adsorption energies of the sites to which doped elements (k–n) are bound decrease, resulting from the gradually reduced electronegativities of doped elements. In N-GDY, sites k–n where Li directly bonded to N exhibit larger adsorption energies in contrast to sites a–j where Li is far from N, on account of the stronger electronegativity of N in comparison with C. However, there are four sites (b, d, g, and i) in Sb-GDY and 7 sites (a, b, d, f, g, i, and j) in Bi-GDY with negative adsorption energies, suggesting that these locations are energetically incapable of adsorbing lithium atoms. This can be attributed to the lattice distortion induced by doping with elements having a large atomic size, that is, the lattice vector a and axial angle γ increase while the lattice vector b decreases, leading to an increased distance between Li and doped atoms, and shortened distance between Li and H. Therefore, N-GDY, P-GDY, and As-GDY all reveal overwhelming performance in Li adsorption.

Table 2

Adsorption Energies (eV) of 14 Individual Adsorption Sites on M-GDY

sites	N	P	As	Sb	Bi
a	1.20	1.25	1.36	1.50	–0.84
b	1.25	1.28	1.39	–0.20	–0.37
c	1.19	1.23	1.35	1.47	1.51
d	1.28	1.30	1.39	–0.21	–0.38
e	1.27	1.30	1.39	1.51	1.54
f	1.37	1.48	1.60	1.72	–0.73
g	1.31	1.43	1.54	–0.07	–0.24
h	1.26	1.31	1.44	1.56	1.58
i	1.33	1.34	1.42	–0.14	–0.32
j	1.31	1.34	1.42	1.55	–0.77
k	1.45	1.64	1.79	1.91	1.96
l	2.20	1.48	1.48	1.47	1.35
m	2.00	1.25	1.01	0.86	0.56
n	1.73	1.60	1.68	1.79	1.71

Theoretical Storage Capacity and Open-Circuit Voltage of M-GDY

Bn class="Chemical">ased on the calculation results concerning adsorption energies of different sites, and the overall structural stability, the preferred adsorption configurations of multiple ass="Chemical">lithium ions on <ass="Chemical">span class="Chemical">M-GDY and the corresponding maximum lithium storage capacity can be easily determined.[39] All atoms are allowed to move freely during the structural optimizations of Li-C22M2H4, and more than 10 distinct initial geometries have been tested for heavy Li loadings to determine the most favorable adsorption sites. After the structural optimization of M-GDY with the maximum number of lithium atoms, the structural stability is judged on the basis of the structural integrity and the average adsorption energy. The average adsorption energies of Li on M-GDY with different loadings are shown in Figure S3, which are basically larger than 1.63 eV, suggesting stable Li bonding.

Figure shows the optimized geometries of Li-C22M2H4 (x = 18, 20, 22, and 24) with M representing group-15 elements. From Figure a,b, the maximum number of <span class="Chemical">lithium atoms stored in <ass="Chemical">span class="Chemical">N-GDY is 20 with average Eads = 1.64 eV/Li, corresponding to a theoretical capacity of 1811 mA·h·g–1. When more Li atoms (x = 22) are added to N-GDY, the average adsorption energy is down to 1.62 eV, which is lower than the cohesive energy of the lithium cluster, meaning that Li atoms will be clustered on the surface of doped GDY.

Figure 6

Top and side views of optimized geometries of Li-C22M2H4 complexes.

In the meantime, for <span class="Chemical">P-GDY, the simulation cell can hold up to 24 Li atoms and maintain its structural stability. We have also arranged a <ass="Chemical">span class="Chemical">Li26-C22P2H4 adsorption configuration; however, the optimized structure of Li26-C22P2H4 shown in Figure d lost its integrity and is significantly distorted, indicating that the maximum number of lithium reserved in P-GDY is 24, related to its theoretical capacity of 1949 mA·h·g–1. The averaged adsorption energy for 24 Li atoms in P-GDY is 1.65 eV/Li. The lithium storage properties of As-GDY is similar to P-GDY, which can hold up to 24 Li atoms with a theoretical capacity of 1538 mA·h·g–1. When 26 Li atoms are added to As-GDY, although the average adsorption energy of lithium in this complex is 1.71 eV, one of the As–C bonds is broken, and the structure is no more stable (see Figure f). Although As-GDY can hold the same number of lithium atoms as P-GDY, the greater molar mass made the theoretical capacity of As-GDY lower than that of P-GDY.

Due to the smaller electronegativities of n class="Chemical">Sb and n class="Chemical">Bi, <n class="Chemical">span class="Chemical">Sb-GDY and <ass="Chemical">span class="Chemical">Bi-GDY have fewer adsorption sites compared with N-GDY, P-GDY, and As-GDY. Besides, their shortened lattice vector in the b-axis also conduce to limited Li adsorption. As a result, the maximum number of lithium in Sb-GDY is 18 with a theoretical capacity of 942 mA·h·g–1 (Figure g), and the maximum number of lithium in Bi-GDY is 16 with a capacity of 625 mA·h·g–1 (Figure i). The corresponding adsorption energies are 1.91 eV/Li for Sb-GDY and 1.70 eV/Li for Bi-GDY, respectively. Generally speaking, compared with the theoretical capacity of 744 mA·h·g–1 of primitive GDY, group-15 element-doped GDY have higher lithium storage capacity owing to their larger apertures. Among the five group-15 element-doped GDY, P-GDY has the highest capacity of 1949 mA·h·g–1, ∼2.6 times the capacity of intact GDY, making P-GDY a promising candidate for anode materials in LIBs.

AIMD simulations of Li-adsorbed LiC22M2H4 complexes (at 300 K) are implemented to examine the thermodynamic stability of Li storage in ass="Chemical">M-GDY. As shown in Figure S4, the <ass="Chemical">span class="Chemical">lithium atoms adsorb steadily onto M-GDY and are not clustered on the M-GDY surface during the entire 6 ps simulation time. For N-GDY and P-GDY, the overall structures preserve their integrities and remain almost unchanged, demonstrating excellent lithium storage performance at room temperature. Nevertheless, in As-GDY, Sb-GDY, and Bi-GDY, heavy Li loadings destabilize the M–C bonds. The snapshots of LiC22M2H4 (M = As, Sb, and Bi) complexes at t = 3 and 6 ps in AIMD simulations clearly show the breaking of one M–C bond in the benzene ring and more distorted carbon network. This can also be demonstrated by the large fluctuations (>1 Å) of M–C bonds in LiC22M2H4 (M = As, Sb, and Bi) complexes (Figure S5), which is in great contrast to the case of M–C bonds in primitive M-GDY systems (Figure ). Therefore, N-GDY and P-GDY display excellent thermodynamic stability with heavy Li loadings at room temperature, while the As-GDY, Sb-GDY, and Bi-GDY undergo severe structural distortion and manifest poor stability with a large amount of adsorbed-Li, resulting from the weaker bonding strength between the doped element and carbon.

To explore, in-depth, the electrochemical performance of <span cln class="Chemical">ass="Chemical">LIBs utilizing <n class="Chemical">ass="Chemical">spn>an class="Chemical">M-GDY as the anode material, the open-circuit voltage (OCV) of lithium intercalation in M-GDY are calculated, which is defined aswhere ELi, ELi, and ELi represent the free energy of LiC22M2H4, LiC22M2H4, and bulk lithium, respectively. M represents group-15 elements, and x1 and x2 are the number of adsorbed lithium atoms. Figure shows the variation of OCV of M-GDY with different Li loadings. The OCV of all M-GDY materials decreases as the number of adsorbed Li atoms increased. When x1 = 0 and x2 = 2, the calculated average open-circuit voltages are high in these five configurations. However, with lithium atoms being constantly filled into doped GDY x1 = 2 and x2 = 4, and the average open-circuit voltages decrease sharply and maintains a relatively flat low voltage in the following volume. The overall average voltages are 0.298, 0.31, 0.404, 0.381, and 0.379 V for N-GDY, P-GDY, As-GDY, Sb-GDY, and Bi-GDY, respectively. N-GDY employs the highest OCV, while P-GDY has the lowest OCV. The calculated voltages of M-GDY are in the required voltage window, which is 0–2 V with bulk Li used as the electrode.[40] Moreover, as anode materials, the low open-circuit voltage allows group-15 element-doped GDY to provide a higher operation voltage for lithium ion batteries.

Figure 7

Open-circuit voltage as a function of x Li in LiC22M2H4.

Open-circuit voltage n class="Chemical">as a function of x Li inpan> LiC22M2H4.

Lithium Diffusion on P-GDY

In consideration of the crucial influence of Li mobility onpan> the charge/discharge rate of <span class="Chemical">LIBs, the diffusion of <ass="Chemical">span class="Chemical">lithium atoms on M-GDY has been investigated using the CI-NEB method. Figure shows the energy barrier of the in-plane diffusion path of P-GDY with the migration path of a single lithium atom moving among some preferred adjacent adsorption sites. Migration paths such as 1–4, 4–5, and 8–1 are calculated, while other paths are structurally symmetric with these paths. A quintessential example should be cited that paths 5–8 can be viewed as the mirror image of path 1–4. Sites 1, 4, 5, and 8 lies above the benzene ring center, corresponding to the adsorption position in Figure k. Sites 2, 3, 6, 7, and 9 are equivalent to the adsorption location in Figure n, which is the most favorable adsorption site. When the migration path and energy barrier of these paths are obtained, a large cyclic lithium migration path is formed, and the energy barriers can be found in the diagram below.

Figure 8

In-plane diffusion of Li on P-GDY. The inset shows the diffusion path for a single lithium atom (path 1–4–5–8–1). Energy barrier is displayed in eV units.

In-plane diffusion of Li on <span class="Chemical">P-GDY. The inset shows the diffusion path for a single <ass="Chemical">span class="Chemical">lithium atom (path 1–4–5–8–1). Energy barrier is displayed in eV units.

The <n class="Chemical">span class="Chemical">lithium atom needs to overcome an energy barrier of 0.203 or 0.205 eV to move from the <ass="Chemical">span class="Chemical">benzene center to the diacetylene chain (path 1–2, 4–3, 5–6, 8–7, 8–9). Migration of Li atoms on the diacetylene chain surrounded by 2 P atoms and 2 H atoms requires an energy of at least 0.300 eV (path 2–3, 6–7). Greater energies are required for Li diffusion on the diacetylene chain surrounded by 4 H atoms, which is 0.461 eV (path 4–5). Hence, the energy barrier for in-plane Li diffusion on P-GDY along the carbon network would be 0.461 eV, which is smaller than that of the in-plane diffusion in GDY (up to 0.70 eV for diffusion between neighboring 18-C hexagons), indicating improved carrier transport properties by P-doping. Besides, due to the adsorption of lithium by impurity atoms, the energy barrier to be overcome for lithium to migrate away from the benzene ring is greatly reduced (from 0.84 to 0.2 eV). Unexpectedly, it is discovered that it should pass a large barrier of 1.78 eV for the Li atom to directly move through the middle of the hexagonal hole (path 1–8). This phenomenon is also different from the diffusion path of Li in GDY as lithium atoms prefer migrating in the large aperture of GDY, in great contrast, the lithium atom prefers migrating along the carbon chains in P-GDY. The peculiar in-plane migrating behavior of Li on P-GDY can be attributed to the Pz orbitals of the carbon network, which provide an easier transport pathway for lithium atoms like a highway.

Conclusions

In this study, the storage and diffusion properties of Li on <n class="Chemical">span cln class="Chemical">ass="Chemical">M-GDY have been systematically investigated based on DFT calculations, including their geometric properties, electronic structures, adsorption modes, maximum storage capacity, thermodynamic stability, theoretical <n class="Chemical">ass="Chemical">span class="Chemical">OCV, and migration behavior. It is found that doping with group-15 elements would result in significant changes in GDY configuration, in which the aperture between the adjacent benzene rings varies from a triangle to a 6 times larger hexagon, providing more adsorption sites for Li. This doping-induced dramatic deformation in the geometry also gives rise to the reduction of π-conjugation of the carbon network, leading to an increased energy gap for M-GDY. In addition, the major contribution of the Pz orbital in the CBM and VBM indicates that most adsorption sites of Li are located out of plane. Considering both electronegativity and pore size, P-GDY and As-GDY can hold the largest amount of lithium atoms during the formation of Li24-C22M2H4, with storage capacities of 1949 and 1538 mA·h·g–1, respectively. Moreover, simulations of the in-plane diffusion of lithium on P-GDY suggest that lithium prefers moving along the carbon chains (with an energy barrier of 0.461 eV) rather than moving directly across the middle of the hexagonal aperture (with an energy barrier of 1.78 eV). The outstanding properties of P-GDY, including large storage capacity, remarkable thermodynamic stability, low OCV, and energy barrier for Li migration, imply that P-GDY is a promising anode material for LIBs.