Na3V(PO4)2 cathode material for Na ion batteries : Defects, dopants and Na diffusion

Layered Na3V(PO4)2 has been recently identified as a high rate cathode material for Na ion batteries. We use atomistic simulation based on the classical pair potentials to calculate the most favourable intrinsic defect process, Na migration paths and tetravalent dopant incorporation at V and P sites. The Na-V anti-site defect is the most energetically favourable defect process. The Na Frenkel is the second most favourable intrinsic defect but only higher by 0.19 eV than the anti-site. Two dimensional long range Na ion migration with activation energy of 0.59 eV is observed along the ab plane implying that Na3V(PO4)2 could be a promising cathode material for Na ion batteries. The formation of both Na vacancy and interstitial defects can be simultaneously achieved by substituting Ge on the V site and the P site required for vacancy migration and storage capacity respectively. High exoergic solution energy is calculated for La on the V site suggesting that the formation of Na3(VxLa1-x)(PO4)2 composition should be experimentally possible.


Introduction
Growing energy demand arising from non-renewable fossil fuels requires high capacity energy storage systems. Lithium ion batteries (LIBs) are promising candidate devices that have been successfully developed and commercialized for small-scale applications. Low abundance of lithium and its unequal distribution in the earth's crust led to a significant challenge in constructing LIBs needed for large-scale applications such as hybrid electric vehicles. Sodium ion batteries (SIBs) are being actively considered as an alternative to LIBs owing to unlimited resource of Na found worldwide, low cost and similar energy storage mechanism to LIBs [1][2][3]. High natural abundance of Na is an advantage of designing electrode materials for large-scale NIBs.
A novel layered type Na 3 V(PO 4 ) 2 has been recently introduced as a high voltage cathode material for NIBs [18,19]. As this material consists of polyanionic orthophosphate (PO 4 3 -) groups with strong P-O bonds, good structural stability and reasonably high redox potentials are expected. Furthermore, variable oxidation states of vanadium (from +2 to +5) is an advantage of exhibiting good redox reactions by forming different V-based polyhedral structures. Kim et al. [18] reported the synthesis of Na 3 V(PO 4 ) 2 using X-ray diffraction first time. Their study exhibited outstanding power capability (~79% of the theoretical capacity at 15 C) and excellent cyclability. Very recently Kovrugin et al. [19] synthesised Na 3 V(PO 4 ) 2 and showed that as-prepared material can provide high voltage plateaus around 3.6 and 4.0 V vs Na + /Na. Furthermore, a high theoretical capacity of 173 mAh/g with a possible two electron transfer per V is expected. There are only a very few theoretical studies available on this new cathode material in the literature. Density functional theory calculations performed by Kim et al. [18] show that two dimensionally connected Na diffusion path with activation energy of 0.433 eV is present. The bond valence energy landscape method was used by Kovrugin et al. [19] to calculate the activation energies of Na ion diffusion in Na 3 V(PO 4 ) 2 . High activation energies (from 1.6 eV to 5.5 eV) were reported for local Na hops meaning that Na ion migration in this material is expected to be slow.
In the present study, we have used well-established atomistic simulations based on the interatomic potentials to examine the defect processes, Na self-diffusion and dopant properties in Na 3 V(PO 4 ) 2 in order to optimize this material in rechargeable Na ion batteries. The current computational modelling technique has been successfully applied to various battery materials [20][21][22][23][24][25][26][27][28][29][30][31][32][33] and can be useful to future experimental work providing defect energetics, diffusion mechanisms and promising dopants to be considered in Na 3 V(PO 4 ) 2 .

Computational Methods
All calculations were performed using the GULP code [34] which is associated with the classical Born model description of ionic solids. In this method, long-range Coulombic interactions and short-range repulsive interactions representing electron-electron repulsion and van der Waals interactions were used to model interactions between ions. Short range interactions were modelled using Buckingham potentials (Table S1). Geometry optimisation (both atom positions and lattice constants) was performed using the Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm [35]. A supercell containing 603 atoms was used for all defect calculations. In all calculations, a gradient norm of 0.001 eV/Å was used to convergence both non-defect and defect structures. Mott-Littleton method [36] was used to model lattice relaxations around point defects and the migrating ions. The local maximum energy along the diffusion path was considered as the activation energy of migration. As the present model assumes a full charge ionic model within the dilute limit, the defect energies will be overestimated but the trend will be consistent. Thermodynamically the migration and formation energies (i.e. the defect parameters) may be defined with the comparison of the real defective crystal to an isochoric or isobaric (as in the present study as well as elsewhere, e.g., see ref 37) non-defective crystal and be interconnected via thermodynamic relations [38][39][40][41]. Kovrugin et al. [19] in their experiment. Vanadium and phosphorous form VO 6 octahedron and PO 4 tetrahedron with adjacent oxygen atoms respectively (refer to Figure 1). The VO 6 octahedra and PO 4 tetrahedra share their corners to form a VO 6 -PO 4 layer along the ab plane.

Na3V(PO4)2 crystal structure
Sodium atoms occupy sites in two different interlayers along the ab plane enabling Na ions to migrate via two dimensional channels. Here we employed classical potentials based atomistic simulation to reproduce the experimental crystal structure of Na 3 V(PO 4 ) 2 with well-established pair potentials available in the literature (see Table S1 in the supplementary potential information). There is an excellent agreement between the calculated and experimental structural parameters within a margin of ~1% error (see Table 1).

Intrinsic defects
We calculated Frenkel, Schottky and anti-site defect energies. These defect energies were calculated by combining isolated point defects (vacancy and interstitial) and are useful to understand the electrochemical behaviour of Na 3 V(PO 4 ) 2 . The intrinsic defect reactions (Eqs.1-8) for Frenkel, Schottky and anti-site defects are given below using Kröger-Vink notation [42].
Na Frenkel: Na Na → Na + Na i (1) Intrinsic defect reaction energies were calculated (refer to Figure 2 Schottky-like reaction (relation 5) leading to the formation of further Na and is calculated to be 2.27 eV per defect, implying that this process may also take place at elevated temperatures. The reaction energies for the P Frenkel, the V Frenkel, the O Frenkel and Schottky defects are even higher suggesting that they are unlikely to occur.

Na ion diffusion
Na ion diffusion with low activation energy is a key requirement for a promising cathode material. The current methodology allowed the identification of various possible Na ion diffusion paths together with activation energies in Na 3 V(PO 4 ) 2 . Using experimental techniques it is often difficult to observe Na ion diffusion pathways and calculate activation energies. The current modelling technique has been successfully used to determine various ionic transport paths in different ionic materials including Na and Li ion battery materials.
For example, the lithium ion migration path calculated in LiFePO 4 using classical pair potentials [47] was exactly observed later in the high-temperature powder neutron diffraction and the maximum entropy method [48]. A curved one-dimensional chain for lithium motion with Li-Li separation of 3.01 Å was clearly visualized along [010] direction in the experiment. In our previous study [25], we have shown that three dimensional long range Na ion pathway can be observed in Na 2 MnSiO 4 .
Three different Na diffusion local Na hops (refer to Figure 3) were identified. Table 2 reports the activation energies together with the Na-Na separation. Energy profile diagrams for each Na hops are shown in Figure 4. Possible long-range paths were also constructed by connecting local Na hops. We identified three long-range three dimensional paths (see Figure   3). In the first long range path (along the ab plane), Na ion migrates via a zig-zag pattern (B

Dopant substitution
We considered isovalent dopants on the V site and tetravalent dopants on both V and P sites.
Tetravalent dopants substitution required charge-compensation in the form of vacancies and interstitials. Appropriate lattice energies (calculated using the same potentials) were also calculated and used.
First, tetravalent dopants (Si, Ge, Ti, Sn and Ce) were considered on the V site. The following reaction was used to calculate solution energies: The lowest solution energy is calculated for Ge (see Figure 5). The solution energy calculated for Sn is only higher by 0.05 eV compared that of Ge. Both Ti and Ce exhibit high solution energies, whereas the incorporation of Si 4+ is highly unfavourable.
Next, we considered the incorporation of M 4+ ions on the P site with Na interstitials as charge-compensating defects (equation 10).
X ′ ⦁ 2 MO 2 + 2 P P + Na 2 O → 2M P + 2 Na i + P 2 O 5 Doping of M 4+ ions on the P site introduces Na interstitials enhancing the capacity in Na 3 V(PO 4 ) 2 required for large scale applications. The most favourable dopant for this process is Si 4+ (see Figure 5). The solution energy gradually increases with ionic radius. Doping of Ce 4+ is unlikely to occur due to the high solution energy (9.62 eV/defect). Finally, isovalent dopants (Al, Ga, Sc, In, Y, Gd and La) were considered on the V site. The following reaction was used to calculate the solution energy: Exoergic solution energies are observed for most of the dopants except for Al and Ga (refer to Figure 6). Interestingly, La exhibits high negative solution energy indicating that the experimental investigation of this dopant should be carried out.

Conclusions
The present study aimed to understand the defect chemistry, Na ion migration and the dopant properties in Na 3 V(PO 4 ) 2 using well-established atomistic simulation technique. The Na-V anti-site defect is found to be the most dominant defect. This indicates that a small percentage of Na ions on the V site and V ions on the Na site would be observed at elevated temperatures.
The second most stable defect energy process is the Na Frenkel energy ensuring the number of Na vacancies that are necessary for vacancy assisted Na diffusion. Two dimensional Na ion migration along the ab plane with the activation energy of 0.59 eV indicates that this material would exhibit high ionic conductivity. The tetravalent dopant Ge 4+ on the V site would increase the concentration of Na vacancy needed for Na ion migration while the same dopant on the P site would increase the Na content in the form of Na interstitial needed for capacity.
Isovalent dopant La seems to be promising candidate for the formation of Na 3 (V x La 1-x )(PO 4 ) 2 .  Table 2. Calculated Na-Na separations and activation energies using classical pair-potential method for the Na-ion migration between two adjacent Na sites (refer to Figure 3).