p-Type conductivity mechanism and defect structure of nitrogen-doped LiNbO₃ from first-principles calculations

Abstract

Most metal-doped lithium niobates (LiNbO₃, LN) exhibit n-type conductivity. The absence of p-type conductive LiNbO₃ limits its application. Based on the finding that p-type conductive LiNbO₃ can be realized by doping with a non-metallic element N, we investigate the most stable defect configurations and formation energies of LiNbO₃ doped with non-metal nitrogen (LN:N) by first-principles calculations. Nitrogen substitution, interstitial and quasi-substitution point defects in different sites and their effects were explored. The results show that N prefers to occupy the oxygen site with only little lattice distortion. Ab initio molecular dynamics (AIMD) simulations confirm the structural stability of an N ion occupying the O site. The charge-state transition level ε(0/−1) slightly above the valence band maximum (VBM) indicates that N point defects would contribute to p-type conductivity of LiNbO₃. The analysis of the band structure reveals that the partially filled impurity levels can accommodate electrons that jump from valence bands and result in holes to become the main charge carriers. The calculation not only explains the occurrence of p-type conductivity in LN:N but also provides a simple and efficient way to discover p-type conductive candidates in numerous doped LiNbO₃ crystals.

---

Introduction

Lithium niobate (LiNbO₃, LN) is considered as one of the promising platforms for integrated photonics.^1,2 Many photonic devices, such as waveguides, couplers, wavelength converters, and electro-optical modulators, have been fabricated with LN in the past few decades.^3–10 However, the aforementioned devices are passive components and the absence of p-type conductive LN limits its application in active components. Various methods, such as irradiation by ultraviolet light and strong doping with various optical damage resistant ions have been tried, but the photoinduced hole carrier concentration is too low to be detected by the Hall effect. Thermoelectric oxidation to achieve p-type LN is also hard to realize due to heavy out-diffusion of lithium.^11–15 Anyhow, it remains a great challenge to fabricate stable p-type LN of high quality. Recently, Li et al. have successfully fabricated a p-type LN film on a Si substrate by doping with nitrogen atoms,¹⁶ thus introducing non-metallic doping in LN, and promoting the fabrication of integrated devices based on LN thin films. However, the origin of p-type conductivity and the defect structures in LiNbO₃ doped with non-metal nitrogen (LN:N) are not clear yet.

As we know, many efforts have been devoted to the theoretical calculation of LN to explain some experimental phenomena and to provide a predictive guidance for the experiments in the past few years.^17–20 Most of these computational studies concentrated on the replacement of Li or Nb cations in doped LN. Xu et al. reported the structures, energies and site preferences of Er defects in LiNbO₃.¹⁹ And according to the report of Li et al., non-photorefractive ions generally prefer the Li site.¹⁸ Despite the cationic dopants being different from each other, all of them resulted in n-type conductivity.^18–20 Little attention was paid to the relationship between the conductivity type and dopants. For non-metallic dopants, there is no relevant report to show their site occupation in LN, and the exploration of its impact on the properties of the non-metallic doped LN.

Non-metallic dopants have also played an important role in improving the properties of compounds. For example, graphene doped with nitrogen can open the band gap from zero, change the electronic structure of graphene and improve the free carrier density of graphene.^21,22 The N element also plays a role in improving the band structure of titanium dioxide metal compounds.^23,24 Due to the distinctive p-type conductivity performance of LN:N it is of importance to explore its p-type conductivity mechanism and its defect structure.

In this paper, employing first-principles calculations, we investigate the formation energies of nitrogen point defects: substitution of anions, N as interstitials and quasi-substitutions. The local atomic relaxations and their effects on the electronic structure of each defect are explored. The charge states of N dopants corresponding to the analysis of the band structure are also examined. The interaction between the dopants and the inherent atoms has also been studied. Finally, we discuss the p-type conductivity mechanism of LN:N through an analysis of the electronic properties.

Method

The Vienna ab initio Simulation Package (VASP) is employed to perform the density functional theory (DFT) calculation.^25,26 The results are obtained from the projector-augmented-wave (PAW) formalism.^27–29 Based on the standard PAW method, the Li 2s¹, Nb 4p⁶4d⁴5s¹, O 2s²2p⁴, and N 2s²2p³ states are treated as valence electrons. The wave functions are expanded in a plane wave basis set with a cutoff energy of 400 eV. The force convergence criterion to optimize the structure is set as 0.01 eV Å⁻¹. The hexagonal supercells containing 120 atoms with a 4 × 4 × 4 k-points mesh are employed in the calculations of the charge state energies of the N substitution and interstitial point defects within the density functional theory-generalized gradient approximation (DFT-GGA). The Perdew–Burke–Ernzerhof (PBE) functional is used to model the electron exchange–correlation energy within the GGA.²⁹ The atomic concentration of nitrogen is 0.83% to simulate the experimental results where the concentration is 0.87% as obtained from X-ray Photoelectron Spectroscopy (XPS).¹⁶ A larger supercell containing 240 atoms is used to simulate substitution, interstitial and quasi-substitution of N point defects and defect complexes with intrinsic defects. The hexagonal supercell employed here is a = b = 10.308 Å, c = 27.754 Å, in which the length of the c-axis is more than twice the length of the a and b axis; thereby, a 2 × 2 × 1 k-point mesh over the Brillouin zone with half the number of k points in the z direction (c axis) as in x and y directions generated by the Monkhorst–Pack scheme is employed.³⁰ The usually used local functional underestimates the band gap of the semiconductors. Thus the Heyd–Scuseria–Ernzerhof (HSE06) screened hybrid functional was implemented in the electronic structure calculations.^31,32 Although the band gaps calculated with DFT-PBE and HSE06 are 3.28 eV and 4.74 eV, respectively, the positions of the defect levels related to the VBM are nearly the same. Therefore, the partial density of states of N substituting O atoms with intrinsic defect clusters are calculated using DFT-PBE to reflect the main characters of the electronic structures. The test on the pseudopotential selection for the atom Nb and the dispersion effect can be found in the ESI.†

Defect formation energy (DFE) is related to the difficulty of forming defect clusters, and a lower formation energy corresponds to a more stable defect complex.³³ The DFEs of different point defects and defect complexes are calculated to explore the most stable defect complex formed in LN:N. The DFE of defect complex X with charge q is calculated using^34,35

where E^total(X^q) is the total energy of the supercell with defect complex X, while E^total(perfect) is the total energy of the pristine supercell. E_v is the valence band maximum, and E_F is the Fermi level with regard to the VBM. ΔV relates to the reference potential in the defect supercell and the pristine crystal. The values of ΔV are between 0.005 eV and 0.06 eV. Here i is the species of atoms that have been added to or removed from the pristine crystal, n_i is the number of atoms of species i, and μ_i indicates the chemical potential of the corresponding atoms. The chemical potentials of Nb, Li and O atoms are calculated using the DFT-PBE functional (see Fig. 1). The chemical potentials μ_i depend on the preparation conditions, and they vary with the change of different reference phases within the constraint range. To prevent the LiNbO₃ ternary compound decomposition into binary oxides Li₂O and Nb₂O₅ or single components, thermodynamic restriction is required to limit the range of chemical potentials μ_i in LiNbO₃. μ_Nb, μ_Li and μ_O should also satisfy with the formation of enthalpy of their oxides Li₂O, Nb₂O₅, and stable LiNbO₃.^36,37

Fig. 1 Defect stability range under different reference states calculated with the DFT-PBE functional. The purple region CDFG is the thermodynamically allowed range for the chemical potentials, and line CG represents the Li₂O reference state while line DF corresponds to the Nb₂O₅ reference state.

Bringing all the stability criteria together, we plotted the thermodynamically stable region of the LN in Fig. 1 according to the method in ref. 36 and 37. As shown in Fig. 1, the purple region CDFG is an area that meets the limitation of the equilibrium equation of lithium oxide, niobium oxide, and LN, and the line DF shows the Nb-rich condition while the line CG represents the Li-rich condition. Since the as-grown crystals and films are Li-deficient compositions, the chemical potentials of Nb and Li are chosen according to the line DF.^36,37 The chemical potential of N atoms is assumed to be half of that of N₂ (gas). The chemical potentials of Nb, Li, O and N elements are listed in Table 1. Similar results have been reported by other researchers.^33,37

Table 1 Chemical potential of different atomic species in LN according to the defect stability range

Condition	Atom species	Chemical potential (eV)
Li-Deficient	Li	−3.63
	Nb	−18.96
	O	−5.59
	N	−8.33

---

Due to the total internal energies obtained from DFT calculations corresponding to the Helmholtz free energy at zero temperature, there is a difference between the VASP work environment and real conditions. Therefore, the free energy (F = E − TS) should be taken into account. Strain effects can be considered ignorable in the large supercell, and electronic entropy is negligible due to the large band gap of LiNbO₃, too. The free energy is mainly related to the entropy which is related to the configuration contributions of point defects and defect clusters. Boltzmann's entropy can be calculated using^17,37,38

S = k_B ln W (5)

where k_B is the Boltzmann constant, and W is the number of point defect configurations in the LiNbO₃ supercell.^37,38 At a room temperature of 300 K, the entropies of point defects are about 0.11 eV, and for defect clusters the value is about 0.21–0.48 eV.

The ab initio molecular dynamic (AIMD) simulation with a canonical ensemble is performed by using the Nóse algorithm. The temperature of the system is maintained at 300 K.³⁹ The 240-atom supercell with the Γ point is implemented here in the whole AIMD calculation. The time step is 1.0 fs and the total simulation time is 10.0 ps.

Results and discussion

Point defects of LN:N

As N is adjacent to O in the periodic table, they have many similarities in physical and chemical properties, such as the size of atoms, the electronic affinity, and the chemical valence state. This suggests that N may substitute O sites and plays a similar role to O by forming N substitutional or the quasi-substitutional point defects. On the other hand, the atomic radius of N is much smaller than those of Li and Nb, and thus N may exist as interstitial atoms. In LN crystals, the O atoms form three different kinds of octahedron structures: a Li-centered octahedron, an Nb-centered octahedron, and a hollow octahedron. In Fig. 2, three isolated sites of point defects around these three oxygen octahedra are shown: the substitutional N (Fig. 2(a)), the interstitial N (Fig. 2(b)), and the quasi-substitutional N (Fig. 2(c)).^40,41

Fig. 2 Typical ball and stick models for LN with substitutional N (a), interstitial N (b), and quasi-substitutional N (c). The N dopants in the different sites are represented as blue balls and labeled with the numbers 1–7 in (a) and (b). Red balls with numbers 1–7 in (c) represent O interstitials whose original sites have been substituted by N. The blue balls labeled with (a–c) represent the different original O sites that are substituted by N dopants. The pink, green and yellow balls indicate the Li, Nb and O atoms at the normal sites in LN, respectively.

As shown in Fig. 2(a), numbers 1–3 present three different O sites that can be substituted by N (marked as N_O) corresponding to the three different oxygen octahedron. The DFEs of the three N_O point defect models are 5.64, 5.64 and 5.63 eV, respectively, showing a weak difference. This indicates that there is no obvious difference in the position of N_O point defects, and the structure of N_O point defects is fairly stable. The models for N as an interstitial (named as N_i) are shown in Fig. 2(b). Due to the limitation of the N atomic size, the possible interstitial sites distribute between two different O atoms in the same layer (labeled as 2, 4, 7 in Fig. 2(b)) or different O layers (labeled as 1, 3, 6 in Fig. 2(b)), and N interstitials can also exist at the center of oxygen octahedron, for example, number 5 in Fig. 2(b). The formation energies change significantly: the lowest formation energy is 4.68 eV. The most stable structure corresponds to the model where N_i lies between O atoms and it is close to the O hollow octahedron. The results show that there may be a positional priority for N_i point defects. With N as quasi-substitutions (N_q), N occupies the normal site of the O, while the original O atom becomes an interstitial at the same time. In Fig. 2(c), the different combinations of the N substitution and O interstitial are displayed. The distribution of O interstitials is the same as that of N interstitials. For comparison, the example where an O interstitial is far away from the N substitution is given in Table 2. Fig. 2(c) and Table 2 show that the formation energies of N_q point defects change significantly. The case of N as quasi-substitution usually comes along with strong lattice distortions that will seriously destroy the lattice structure of LiNbO₃. Therefore, it is difficult to form N quasi-substitutions in the crystals. As seen from the formation energies, N lying in the specific location as an interstitial point defect has the lowest formation energy, while the distribution of N_O is the most stable.

Table 2 Formation energies of different typical models corresponding to Fig. 2

Model	Substitutional			Interstitial
	①	②	③	①	②	③	④	⑤	⑥	⑦
Formation energy (eV)	5.75	5.75	5.74	5.16	5.81	7.26	5.71	7.53	4.79	10.82

---

Quasi-substitutional
Model	A + ①	A + ②	A + ③	B + ③	B + ④	B + ⑥	C + ⑤	C + ⑦	C + far away
Formation energy (eV)	5.16	7.52	5.20	5.16	7.68	5.21	5.24	5.20	9.71

---

Besides knowing their formation energies, it is also important to know how N point defects influence the crystal lattice. The knowledge of the lattice distortion around point defects greatly contributes to the understanding of their interaction with inherent atoms. The local structures of point defect N_O before and after structure relaxation are shown in Fig. 3(a) and (b), respectively. As shown in Fig. 3(b), N_O defects cause only a quite slight local distortion of the LN lattice. In detail, the distances from the upper O layer decrease while the distances from the lower O layer increase. The largest distortion is only 1.15% compared to that in the pristine crystals. In addition to the change in the distances of O atoms between different layers, the distances between the same layers of O atoms are increased, which indicates an increase of the O octahedral volume. As a whole, N_O point defects are fairly stable. This fact can be ascribed to the similarity between N and O atoms.

Fig. 3 Local structure of N substitution at the O site (a) before and (b) after structure relaxation. The local structure of an N interstitial (c) before and (d) after structure relaxation. The pink, green, yellow and blue balls represent Li, Nb, O and N atoms respectively. The numbers on the N–O bond represent the distances of N and O atoms. The distance is given in the unit Å.

Fig. 3(c) and (d), clearly show the lattice relaxation of N_i point defects. In Fig. 3(d), the N_i point defect moves against the z-direction due to the ferroelectric distortion and the distance between N_i and the nearest-neighboring O atoms along the z-direction increases. Due to the Coulomb repulsion of the electrons around N and O atoms, the O plane undergoes serious changes, in which some O atoms get out of the original O plane. Compared with lattice distortions of N_O, N_i point defects have a stronger influence. This effect on the lattice structure is mainly reflected in the influence of O atoms, which means the influence on crystal electronic properties. We ignore the situation of N_q because of its strong distortion of the lattice structure.

Based on the above analysis, N_O defects are superior to interstitial and quasi-substitution point defects because of their lower and more stable formation energies, and substitutions have the least influence on LN crystals compared with the other two forms of point defects. The charge states of the N_O and N_i point defects are studied further to investigate their changes.

AIMD simulations were done to investigate the stability of the N_O defect structure. On the basis of first-principles calculations and AIMD simulation, Fig. 4 shows the stability of the N_O defect structure at a temperature of 300 K. The oscillation of energy (red line) tends to be stable, and the similarity of the two structures before and after 10 ps AIMD simulation indicates that it is fairly suitable for N ions to substitute O atoms.

Fig. 4 Total energy fluctuation during AIMD simulations of LN:N at a temperature of 300 K. The structural snapshots of the LN:N crystal at times of 0 and 10 ps during AIMD simulation.

Charge states of N point defects

Differently from metal dopants replacing Nb or Li ions with a positive charge state, N prefers to occupy O atoms or to be an isolated interstitial. The charge states of the further two possible point defects N_O and N_i were explored, too. Fig. 5 shows the formation energies of N_O and N_i point defects in all charge states at the most stable charge states. We can find that the charge states of the point defects vary with the increase of the Fermi energy from the VBM towards the CBM, a range that corresponds to the band gap of LN with an experimental value of 3.8 eV.^42,43 This value is somewhat larger than the calculated result of 3.28 eV, but the difference is in allowable error scope.^44,45 The multiple charge states of the point defects N_O and N_i indicate that they can form acceptors in LN crystals.

Fig. 5 The defect formation energies of the most stable charge state of N_O and N_i point defects and intrinsic defects with the growth of Fermi energy from the VBM to the CBM. For No point defect, the transition level ε(0/−1) lies in 0.34 eV shown as the pink line.

In general, the formation energy of the point defect N_O is lower than that of point defect N_i throughout the whole Fermi energy range. Therefore, the existence of N_O is more plausible. Due to the limitation of the doping concentration, the concentration of the N_i point defect can be considered negligible. The transition levels ε(q₁/q₂) are defined as the points q₁ and q₂ where the formation energy and Fermi energy become equal.

where E_f(X^q₁) and E_f(X^q₂) are the formation energies of the defect X in the charge state q₁ and q₂. A defect is a typical shallow acceptor when the transition from negative to the neutral charge state is near or below the VBM.⁴⁶ For the N_O point defect, the transition level ε(0/−1) lies at 0.34 eV, which is close to the VBM. The result illustrates that N_O is a shallow acceptor and would contribute to p-type conductivity of LN. In the all possible Fermi energy range, the stable charge state is negative. N_O with −1 charge state is the dominated charge state for the Fermi level of LN:N lies in the half lower part of the band gap. The stable negative charge state of the point defect N_O supports the p-type behavior of the conductivity. Combined with the former part, as an isolated point defect, N dopants prefer to substitute non-metallic O atoms forming the point defect N_O with a negative charge state of −1. Thus the point defect N_O contributes to p-type conductivity of LN as well.

Defect clusters

In the growth process of LN, the preparation conditions result in a large number of intrinsic point defects and defect complexes. Under Nb₂O₅ rich conditions, intrinsic defect niobium antisite (Nb_Li⁴⁺) is supposed to be considered, and it is usually compensated by four lithium vacancies (V_Li⁻) to maintain a neutral environment in LN crystals.^47,48 In LN:N, Nb_Li⁴⁺ can also be complemented by point defects N_O⁻. Table 3 lists the different clustering patterns, in which Nb_Li⁴⁺ co-exists with V_Li⁻ and N_O⁻. Similar to point defects, the defect formation energy is the main criterion considered here to indicate structural stability.

Table 3 The per formation energies of N_O point defect co-existence with intrinsic defects in LN

Clusters	The formation energy (eV per defect)
NbLi^4+ + 4VLi^−	2.60
NbLi^4+ + 3VLi^− + 1NO^−	2.95
NbLi^4+ + 2VLi^− + 2NO^−	3.46
NbLi^4+ + VLi^− + 3NO^−	3.94
NbLi^4+ + 4NO^−	4.35

---

The formation energy of the intrinsic defect cluster Nb_Li⁴⁺ + V_Li⁻ is 2.60 eV, remaining to be the most stable defect cluster compared with the defect clusters including the N point defects. With the increase of N_O⁻ and the decrease of V_Li⁻, the formation energy becomes higher and higher, indicating a more unstable structure. The results show that it is difficult for N dopants to form a defect cluster with intrinsic defects, and N dopants cannot improve the situation of lithium deficiency.

Corresponding to the defect clusters above, the partial density of states (PDOSs) of the main single atomic states including O 2p, Nb 4d, and N 2p that can reflect the contributions of different atoms in electron density are shown in Fig. 6. Based on the PDOS, we can discriminate whether LN contains only intrinsic defects or also N point defects. The magnified area shows that the introduction of N atoms not only leads to the formation of the new energy level in the band gap, but also causes a strong interaction between N 2p and O 2p electrons. Most N 2p electrons and some O 2p electrons contribute to the defect energy level. With the increase in the number of N point defects, some 4d electrons of Nb ions take part in the formation of defect levels that get closer and closer to the bottom of the conduction band. Although the area of defect levels expands, there is not much overlap between these defect levels as O atoms. This may indicate a preference for isolated N point defects.

Fig. 6 The PDOS of (a) ; (b) ; (c) ; (d) ; (e) . The magnified area shows the influence of N atoms.

Based on the defect formation energies and PDOS analysis, we found that N point defects do not form so easily defect clusters with intrinsic defects or other N point defects. As we know, the formation of defect clusters has a close relationship with the concentration of N atoms, and the results provide an answer to the question why it is difficult to introduce N atoms into LN and why N dopants were experimentally observed to maintain at a relatively low concentration.

Mechanism of p-type conductivity in LN:N

N ions have significantly improved the properties of LN, and its most distinctive property is p-type conductivity, which will now be explored. In the part above analyzing the charge state of N point defects, the charge state transition of the N point defect from 0 to −1 already indicates that the conductive type of LN:N is p-type. In order to find out the mechanism of p-type conductivity, we study the band structure and the PDOS of LN:N in Fig. 7(a) and (b). In the band structure of LN:N, the bandgap is about 4.74 eV. The new energy level in the band gap is 3.75 eV away from the bottom of the conduction band.

Fig. 7 (a) Band structure of LN:N calculated with HSE06. The energies from the CBM to the VBM and the new energy levels are 4.74 eV and 3.75 eV respectively. (b) PDOS of the main single atomic states including O 2p, Nb 4d, and N 2p with an amplified region of impurity states in the band gap of LN:N.

The new energy level displays an important role in improving the proprieties of LN. Seen from the PDOS, we can find that the contribution of N 2p electrons and O 2p electrons to the defect level of LN:N is almost the same. The interaction between N 2p and O 2p electrons is obvious. N ions influence the electronic distribution of the valence band, especially the electrons near the VBM, and the electrons of N enhance the activity of valence electrons. The existence of new energy levels is associated with the electronic transition in LN:N. The defect levels are approximately 0.39 eV above the VBM, which provides an opportunity for the O 2p electrons to transfer from the valence band with little energy which can easily be obtained from thermal activation. The result of the valence band electron transition leaves a large number of holes as the main carriers for LN:N. In addition, the occupation of defect energy levels has been explored. The occupation states of the valence band and the conduction band are 2 and 0, thus representing the full occupation state and the empty occupation state, respectively. In the band structure, we can clearly see that there are two energy lines. The occupation of the upper one changes from 0.27 to 1.95, and the other does not also have a fully occupied energy level either. Obviously, the holes in the energy level have a good chance to be thermally activated from the valence band.

Based on the analysis above, we infer that N dopants have a strong influence on the distribution of electrons in LN:N. Therefore, the electronic transformation of N_O point defects has been studied. The charge density distributions of an O plane that include N atoms and the electronic charge difference⁴⁹ between the LN:N and pristine system are given in Fig. 8.

Fig. 8 (a) The charge density distributions within an O plane that includes a N atom in LN:N. The color bar stands for the continuous variation of the charge density from 0.00 to 0.60 e Bohr⁻³. (b) Electronic charge difference between LN:N and LN. The purple and blue ellipsoids represent electron depletion and accumulation, respectively. The pink, green, yellow and blue balls denote the Li, Nb, O and N atoms respectively.

For the N_O point defect, N takes the responsibility of the O atom as an electron-withdrawing center. In a sense, N_O guarantees the basic charge-transfer in crystals. Fig. 8(a) shows that there is an electron cloud gathering around N_O, but it is smaller than the charge density of O atoms. Hence, compared with normal O sites, there is a clear depletion in the O site occupied by N. This is consistent with the conclusion we obtained above that the N point defect shows a −1 charge state. As shown in Fig. 8(b), the yellow ellipsoid around N ions presents the electron depletion compared to the original O site. While the yellow area is surrounded by a blue ellipsoid, it indicates aggregation of electrons between N and O atoms, which enables easy formation of covalent bonds.

In general, the charge density and the electronic charge difference prove the electronic interaction between N and O atoms and confirm that N plays a similar role to O with a negative charge state. Therefore, some of the O 2p electrons have the chance to transit from the valence band to the defect levels, leaving holes in the valence band as the main carriers.

Conclusions

In summary, we have performed the detailed first-principles calculations exploring the possible site preferences of isolated defects and their corresponding charge state in LN:N. We discussed the mechanism leading to p-type conductivity in LN:N based on the analysis of its band structure and of the density of states. The formation energy calculation indicates that the most plausible structure to be formed is N occupying the site of O with only a tiny lattice distortion. In addition, N as the isolated point defect is more reasonable in LN:N for lower formation energies. We found that the transition level ε(0/−1) is near the VBM, which proves that LN:N is a candidate p-type material, and the analysis of the band structure confirms that the defect energy level is close to the VBM. This leads to p-type conductivity in LN:N. Owing to the electron overlap of N and O atoms, the defect level includes a contribution of O. Therefore not all occupation states can be thermally activated from the valence band electrons, leaving a large number of positive holes in the valence band. These are then the dominant charge carriers. N doping opens a reliable pathway for producing LN crystals with p-type conductivity, which may find important applications in integrated optics. Our study on N-doped LN presents a powerful example on how to discover good candidates for p-type conductivity in LN among yet unknown dopants.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

We gratefully acknowledge financial support from the National Natural Science Foundation of China with grant number [11674179] and [61705116], and the Program for Changjiang Scholars and the Innovative Research Team in University with grant number [IRT_13R29], and the Slovenian research agency (ARRS, grant P1-0192).

References

1. M. Zhang, C. Wang, R. Cheng, A. Shams-Ansari and M. Lončar, Optica, 2017, 4, 1536.
2. A. Guarino, G. Poberaj, D. Rezzonico, R. Degl'Innocenti and P. Guenter, Nat. Photonics, 2007, 1, 407.
3. Y. Park, B. Min, K. J. Vahala and H. A. Atwater, Adv. Mater., 2006, 18, 1533.
4. M. F. Volk, S. Suntsov, C. E. Rüter and D. Kip, Opt. Express, 2016, 24, 1386.
5. R. Wu, J. Zhang, N. Yao, W. Fang, L. Qiao, Z. Chai, J. Lin and Y. Cheng, Opt. Lett., 2018, 43, 4116.
6. Z. Chen, Y. Wang, Y. Jiang, R. Kong and H. Hu, Opt. Mater., 2017, 72, 136.
7. C. Wang, M. Zhang, X. Chen, M. Bertrand, A. Shams-Ansari, S. Chandrasekhar, P. Winzer and M. Lončar, Nature, 2018, 562, 101.
8. R. Wolf, Y. Jia, S. Bonaus, C. S. Werner, S. J. Herr, I. Breunig, K. Buse and H. Zappe, Optica, 2018, 5, 872.
9. A. Ashkin, G. D. Boyd, J. M. Dziedzic, R. G. Smith, A. A. Ballman, J. J. Levinstein and K. Nassau, Appl. Phys. Lett., 1966, 9, 72.
10. C. Wang, M. Zhang, B. Stern, M. Lipson and M. Lončar, Opt. Express, 2018, 26, 1547.
11. D. L. Staebler and J. J. Amodei, J. Appl. Phys., 1972, 43, 1042.
12. Y. Kong, S. Wu, S. Liu, S. Chen and J. Xu, Appl. Phys. Lett., 2008, 92, 251107.
13. F. Liu, Y. Kong, W. Li, H. Liu, S. Liu, S. Chen, X. Zhang, R. Rupp and J. Xu, Opt. Lett., 2010, 35, 10.
14. J. Deng, J. Wen, Z. Wu and H. Wang, Appl. Phys. Lett., 1994, 64, 2622.
15. Z. Pei, Q. Hu, Y. Kong, S. Liu, S. Chen and J. Xu, AIP Adv., 2011, 1, 032171.
16. W. Li, J. Cui, W. Wang, D. Zheng, L. Jia, S. Saeed, H. Liu, R. Rupp, Y. Kong and J. Xu, Materials, 2019, 12, 819.
17. A. F. Kohan, G. Ceder, D. Morgan and C. G. Van de Walle, Phys. Rev. B: Condens. Matter Mater. Phys., 2000, 61, 15019.
18. L. Li, Y. Li and X. Zhao, Phys. Chem. Chem. Phys., 2018, 20, 17477.
19. H. Xu, D. Lee, S. B. Sinnott, V. Gopalan, V. Dierolf and S. R. Phillpot, Phys. Rev. B: Condens. Matter Mater. Phys., 2009, 80, 144104.
20. F. Bridges, C. Mackeen and L. Kovács, Phys. Rev. B, 2016, 94, 014101.
21. D. Wei, Y. Liu, Y. Wang, H. Zhang, L. Huang and G. Yu, Nano Lett., 2016, 9, 1752.
22. C. Zhou, J. Kong, E. Yenilmez and H. Dai, Science, 2000, 290, 1552.
23. R. Asahi, T. Mikawa, T. Ohwaki, K. Aoki and Y. Taga, Science, 2001, 293, 269.
24. M. Batzill, E. H. Morales and U. Diebold, Phys. Rev. Lett., 2006, 96, 026103.
25. G. Kresse and J. Furthmüller, Comput. Mater. Sci., 1996, 6, 15.
26. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169.
27. P. Pulay, Chem. Phys. Lett., 1980, 73, 393.
28. Q. Li, B. Wang, C. H. Woo, H. Wang and R. Wang, J. Phys. Chem. Solids, 2007, 68, 1336.
29. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865.
30. H. J. Monkhorst and J. D. Pack, Phys. Rev. B: Solid State, 1976, 13, 5188.
31. J. Heyd, G. E. Scuseria and M. Ernzerhof, J. Chem. Phys., 2003, 118, 8207.
32. A. V. Krukau, O. A. Vydrov, A. F. Izmaylov and G. E. Scuseria, J. Chem. Phys., 2006, 125, 224106.
33. H. Xu, A. Chernatynskiy, D. Lee, S. B. Sinnott, V. Gopalan, V. Dierolf and S. R. Phillpot, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 82, 184109.
34. C. G. Van de Walle and J. Neugebauer, J. Appl. Phys., 2004, 95, 3851.
35. S. Lany and A. Zunger, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 78, 235104.
36. H. Xu, D. Lee, J. He, S. B. Sinnott, V. Gopalan, V. Dierolf and S. R. Phillpot, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 78, 174103.
37. Y. Li, W. G. Schmidt and S. Sanna, Phys. Rev. B: Condens. Matter Mater. Phys., 2014, 89, 094111.
38. J. He, R. K. Behera, M. W. Finnis, X. Li, E. C. Dickey, S. R. Phillpot and S. B. Sinnott, Acta Mater., 2007, 55, 4325.
39. D. Bucher, L. C. T. Pierce, J. A. McCammon and P. R. L. Markwick, J. Chem. Theory Comput., 2011, 7, 890.
40. Y. Bar-Yam and J. D. Joannopoulos, Phys. Rev. Lett., 1984, 52, 1129.
41. R. Car, P. Kelly and A. Oshiyama, Phys. Rev. Lett., 1984, 52, 1814.
42. D. Redfield and W. J. Burke, J. Appl. Phys., 1974, 45, 4566.
43. A. Dhar and A. Mansingh, J. Appl. Phys., 1990, 68, 5804.
44. M. S. Hybertsen and S. G. Louie, Phys. Rev. Lett., 1985, 55, 1418.
45. R. W. Godby, M. Schlüter and L. J. Sham, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 10159.
46. P. Rinke, A. Janotti, M. Scheffler and C. G. Van de Walle, Phys. Rev. Lett., 2009, 102, 026402.
47. D. M. Smyth, Prog. Solid State Chem., 1984, 15, 145.
48. H. Donnerberg, S. M. Tomlinson, C. R. A. Catlow and O. F. Schirmer, Phys. Rev. B: Condens. Matter Mater. Phys., 1989, 40, 11909.
49. X. Li and S. Li, J. Mater. Chem. C, 2019, 7, 1630.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c9cp05019a

---