Niobium oxide dihalides NbOX₂: a new family of two-dimensional van der Waals layered materials with intrinsic ferroelectricity and antiferroelectricity

Abstract

Two-dimensional (2D) ferroelectric (FE) materials displaying spontaneous polarizations are promising candidates for miniaturized electronic and memory devices. However, stable FE orderings are only found in a small number of 2D materials by experiment so far. In the current work, based on high-throughput screening of a 2D van der Waals layered materials database and first-principles calculations, we demonstrate niobium oxide dihalides NbOX₂ (X = Cl, Br and I), a group of experimentally synthesized yet underexplored van der Waals layered compounds, as a new family of 2D materials that simultaneously exhibit intrinsic in-plane ferroelectricity and antiferroelectricity. Similar to FE perovskite oxides, polar displacement of Nb cations relative to the center of the anion octahedral cage can lead to experimentally measurable FE polarizations up to 27 μC cm⁻² in layered NbOX₂. The presence of low-lying antiferroelectric (AFE) phases can effectively reduce the energy barrier associated with polarization switching, suggesting switchable ferroelectricity is experimentally achievable. In addition, the mechanism driving FE phase transitions in NbOX₂ monolayers around Curie temperature T_C is clearly revealed by our finite-temperature simulations. NbOCl₂ monolayer is predicted to be a stable ferroelectric with T_C above room temperature. Moreover, application of NbOBr₂ and NbOI₂ monolayers as 2D dielectric capacitors is further developed, where electrostatic energy storage of nearly 100% efficiency can be achieved in the 2D single-layer regime.

---

New concepts

Two-dimensional (2D) ferroelectric (FE) materials with stable spontaneous polarizations above room temperature hold great promise in miniaturized electronic and memory devices. Here, we demonstrate a new strategy to uncover a wide spectrum of experimentally synthesizable yet unreported 2D FE materials. Combining high-throughput screening of a 2D layered materials database and first-principles calculations, layered NbOX₂ (X = Cl, Br and I), have been identified as a new family of 2D FE materials. The 2D NbOX₂ layers can be readily obtained through exfoliation of bulk NbOX₂ which were synthesized decades ago. In particular, NbOCl₂ monolayer with atomic thickness can exhibit room temperature stable ferroelectricity with sizable spontaneous polarization. Moreover, coexistence of FE and antiferroelectric (AFE) phases and the associated AFE/FE phase transition can enable practical applications of layered NbOX₂ as 2D dielectric capacitors, with the potential to achieve nearly 100% efficiency for electrostatic energy storage.

Introduction

ABO₃ perovskite oxides with broken spatial inversion symmetry are commonly known ferroelectric (FE) materials, which exhibit spontaneous polarizations below Curie temperature T_C. The continuous demands for miniaturized FE devices have motivated extensive explorations of FE perovskite oxide thin films with reduced thickness.¹ In the past decade, the emergence of various two-dimensional (2D) materials with diverse functionalities has offered a new platform to investigate ferroelectricity and cooperative FE properties at the nanoscale. In particular, newly discovered 2D FE materials can have stable out-of-plane (vertical) FE orderings even in atomic-thin layers,^2–4 with interlayer spacing below the critical thickness where ferroelectricity in most FE perovskite oxides would disappear.^5,6 In addition, strong coupling of 2D ferroelectricity with ultra high elastic strain^7,8 or superior optoelectronics reponses⁹ make 2D FE materials display rich coupled physical properties, enabling promising applications such as miniaturized electronic and memory devices.¹⁰

Guided by first-principles calculations, a number of 2D FE materials have been successfully synthesized and demonstrated with stable FE orderings by experiment. Typically, in-plane ferroelectricity and FE phase transitions were predicted in group-IV monochalcogenide (GeS, SnS and SnSe, etc.) monolayers,^11–14 and then experimentally demonstrated in atomic-thick SnTe.^15,16 Computational discovery of 2D FE α-In₂Se₃ layers¹⁷ initiated the extensive experimental studies on intercorrelated in-plane and out-of-plane ferroelectricity within this system.^18–22 Spontaneous symmetry breaking in the 1T phase of MoS₂ monolayers can lead to robust in-plane ferroelectricity,²³ while its structural polarity was detected by second harmonic generation measurement very recently.²⁴ The material discoveries mentioned above largely relied on researchers' experience or chemical intuition, which was usually restricted to identifying 2D FE materials with structure or composition analogous to the known ones. For example, group-IV monochalcogenide monolayers were investigated as structural analogues to phosphorene.¹¹ Therefore, providing rational strategies for a wide spectrum of new 2D FE material discovery will largely accelerate research progress within this field, which is also beneficial for a complete and insightful understanding regarding the unique properties of 2D ferroelectricity at the nanoscale.

In this work, instead of using a traditional intuition guided material design approach, we demonstrate a new strategy for discovering experimentally synthesizable 2D layered FE materials. Recent data-mining studies have provided the complete database covering almost all experimentally reported 2D van der Waals layered materials (including 2D layered FE materials).^25,26 Further screening of these 2D materials databases for acentric layered semiconductors or insulators with switchable polarizations can lead to the discovery of new classes of 2D layered FE materials which are not simply composition or structural analogues to any experimentally reported ones. After searching the 2D weakly bonded layered materials database provided by Cheon et al.²⁵ using basic structural and energetic criteria (Fig. 1), we identified niobium oxide dihalides NbOX₂ (X = Cl, Br and I) as long-sought-after 2D van der Waals layered FE materials.

Fig. 1 High-throughput screening scheme for 2D FE materials. We start by looking for non-metallic 2D monolayers without inversion symmetry, which restrict the search to about 260 types of 2D insulating or semiconducting monolayers with C₁, C₂, C_2v, C_3v, C_4v and C_6v polar point groups. We then analyze polarization for each candidate, and identify two potential FE states with opposite polarizations (±P states). The potential PE phase is obtained after “averaging” the atomic structures of the ±P states. After structural optimization, those 2D materials with unstable FE or PE structures are excluded. More importantly, the energy difference between the FE and PE phases (ΔE) should be small enough so that the polarization can be switchable. The energy criteria for ΔE is chosen to be 40 meV per atom, corresponding to the PE–FE energy difference in PbTiO₃. Based on our material screening scheme, we not only obtain the known 2D FE materials, such as group-IV monochalcogenides, but also identify a new class of layered FE materials – NbOX₂.

Bulk NbOX₂ have been synthesized since the 1960s,^27,28 and they can be potentially exfoliated into 2D layered forms. Despite their structural polarities at room temperature,^28–31 NbOX₂ have not been considered or investigated as 2D layered ferroelectrics. Except for earlier tight-binding band structure calculations,³² accurate predictions regarding electronic and FE properties for NbOX₂ bulk and 2D layers have never been reported. In our work, based on first-principles calculations, model Hamiltonian and molecular dynamics simulations, we identify the presence of stable FE and metastable antiferroelectric (AFE) phases in both 3D bulk and 2D layered forms of NbOX₂. The coexistence of FE and AFE phases can lead to unique polarization switching and structural phase transition features in NbOX₂, distinct from the currently known 2D FE materials. Our findings demonstrate NbOX₂ as a practical system to study the intrinsic ferroelectricity and antiferroelectricity down to the monolayer limit, and to explore the novel functionalities associated with 2D FE/AFE phase transitions.

Results and discussion

Ground state structures and ferroelectric properties

Niobium oxide dihalides NbOCl₂, NbOBr₂ and NbOI₂ belong to a group of transition metal oxide halide compounds with chemical formula MOX₂ (M = V, Nb, Mo and Ta, and X = Cl, Br and I).^27,28 Bulk NbOX₂ is a typical van der Waals layered material,²⁵ which is formed by stacking of NbOX₂ monolayers along the out-of-plane direction. Due to the weak van der Waals interlayer interactions (the interlayer binding energy E_b = 12–13 meV Å⁻² in Table 1, while the E_b of graphite is predicted to be 20.3 meV Å^−2 33), it should be eminently possible to obtain 2D NbOX₂ multi-layers or monolayers via mechanical exfoliation techniques. Similar to the perovskite crystal structure, the NbOX₂ monolayer is composed of NbO₂X₄ octahedra with mixed edge- and corner-sharing connectivity. Typically, as shown in Fig. 2(a), NbO₂X₄ octahedra build up a 2D structural network through extensive interconnecting X–X edges along one planar direction (crystallographic a axis) and cornered O atoms along the other (b axis).

Table 1 Calculated interlayer binding energy (E_b), Nb polar displacement (d_Nb) relative to the center of the NbO₂X₄ octahedra, FE potential depth ΔE_FE, spontaneous polarization P and nominal energy band gaps E_g for bulk and monolayer FE NbOX₂. We assume the thickness of the NbOX₂ monolayer is same as the interlayer spacing of bulk NbOX₂, when simulating its P. The optical transition allowed E_g are given in parentheses

	NbOCl2		NbOBr2		NbOI2
	Bulk	Monolayer	Bulk	Monolayer	Bulk	Monolayer
E b (meV Å^−2)	12.19	—	12.77	—	13.40	—
d Nb (Å)	0.15	0.14	0.14	0.14	0.13	0.13
ΔEFE (meV f.u.^−1)	13.31	13.63	10.64	11.00	8.12	9.08
P (μC cm^−2)	26.86	26.50	22.28	22.32	17.78	17.66
E g (eV)	1.84 (3.73)	1.88 (3.99)	1.79 (3.07)	1.87 (3.01)	1.69 (1.83)	1.77 (1.84)

---

Fig. 2 Layered structures of PE NbOX₂. (a) Three-dimensional crystal structure for paraelectric (PE) phases of bulk NbOX₂, with Nb in brown, O in red and halogen atoms X in blue. Both side and top views of an isolated PE NbOX₂ monolayer are shown to highlight the one-dimensional Peierls distortion between neighboring Nb atoms along the crystallographic a axis. (b) Calculated phonon spectrum for bulk PE NbOI₂. There are three unstable soft-phonon modes corresponding to FE, AFE and antipolar Nb polar displacement relative to the center of the NbO₂I₄ octahedra along the crystallographic b axis. The unit cell of NbOX₂ is indicated by a black rectangle.

As determined by experimental structural characterizations,^29–31 bulk NbOX₂ can crystalize in a ferroelectric C2 phase with monoclinic symmetry, where Nb cations exhibit non-zero off-center displacement along the Nb–O–Nb atomic chain direction (b axis). In our work, in order to explore all possible NbOX₂ polar configurations, we choose the Nb polar-distortion free, centrosymmetric C2/m phase (ground state structures for MoOCl₂³⁴ and TaOI₂³⁵) as the paraelectric (PE) reference for bulk NbOX₂. Fig. 2(a) displays the crystal structure for the PE bulk phase (crystallographic parameters shown in Table S1 of the ESI†), in which there exits one-dimensional Peierls distortion³⁶ of Nb atoms in each NbO₂X₄ octahedra, leading to an alternation of two unequal Nb–Nb distances along the Nb–X–Nb direction (a axis). As we will discuss later, Nb–Nb Peierls distortion is crucial for determining the electronic structures and band gaps of NbOX₂. After examining the phonon spectra of PE phases, the soft phonon modes and corresponding polar structural instabilities for bulk NbOX₂ can be identified.

We use bulk NbOI₂ as an example and display the calculated phonon spectrum of its PE phase in Fig. 2(b). Three types of soft optical phonon modes with imaginary frequencies at the Γ point (ω = i 83, 82 and 71 cm⁻¹) are identified, corresponding to FE (all Nb atoms have exactly the same polar displacement), antipolar (two neighboring monolayers have anti-parallel Nb polar displacement) and AFE (Nb atoms from neighboring octahedra within each monolayer displace oppositely) Nb polar displacement patterns relative to the center of the NbO₂I₄ octahedra. It is noted that Nb polar displacement in all three modes is restricted along the Nb–O–Nb direction. We next inject each soft mode into PE NbOI₂, followed by structural optimization to obtain lower symmetry polar phases. FE phases with C2 symmetry (Table S2 in the ESI†) are obtained after relaxation of FE-mode related structure. Cooperative displacement of all Nb atoms (polar displacement amplitude d_Nb ≈ 0.14 Å) in the FE phase can break inversion symmetry and generate non-zero spontaneous polarizations along the b axis, which can be qualified using Berry phase calculations. The typical polarization-energy double well curves are obtained for bulk NbOX₂ FE phases (Fig. S1 in the ESI†), after recording the variation of polarization and total energy of the system with respect to the amplitude of the FE mode connecting the FE C2 and PE C2/m phases. Based on our calculations, bulk NbOCl₂ exhibits the largest spontaneous polarization, highest FE potential depth and largest Nb polar displacement amplitude (Table 1), due to the significant electronegativity of the Cl anion. Besides the FE phase, antipolar phases with P2/c symmetry and AFE P1̄ phases (crystallographic parameters shown in Tables S3 and S4 of the ESI†) with anti-parallel interlayer or intralayer Nb polar displacement are also identified. Each Nb cation in the antipolar or AFE phases has non-zero polar displacement, but their contributions to the overall polarization of the system cancel out.

After investigation of bulk NbOX₂, we now turn to 2D NbOX₂ multi-layers and monolayers, where we will demonstrate these NbOX₂ layers as intrinsic 2D FE (AFE) materials. The soft-mode adopted structural optimization scheme is also used to obtain the polar ground state structures for 2D NbOX₂ layers. FE, antipolar and AFE soft modes are found in the PE phase of NbOX₂ multi-layers, while the monolayers only have FE and AFE modes (Fig. S2 and Table S5 in the ESI†). Without a vertical stacking sequence and out-of-plane periodicity, all 2D NbOX₂ layers are optimized into orthorhombic symmetry. Fig. 3(a) summarizes our calculated energetic results for all polar phases in NbOI₂ (results for NbOCl₂ and NbOBr₂ are provided in Fig. S3 of the ESI†). The FE phase is the most stable structure for all NbOX₂ layered systems, independent of the layer numbers. Except for the monolayer, the antipolar phase is the next most stable phase, and is very close to the FE phase in energy. Even the less stable AFE phase has a small energy difference (<2 meV f.u.⁻¹) relative to the FE phase.

Fig. 3 Energy profiles and FE polarization of 2D NbOI₂ layers. (a) Variation of the total energy for different polar configurations in the NbOI₂ system as a function of number of layers, where the energy of the corresponding PE phases are chosen as energy zero. The trend in the energy stability of different NbOI₂ polar configurations is: FE ≈ antipolar > AFE. (b) Nb polar displacement-energy double well plots connecting PE and two FE phases in NbOI₂ monolayer, bilayer and bulk. Symbols are the calculated results, lines are fitted to the data based on the Landau model. Both the magnitude of Nb polar displacement and FE potential depth remain almost unchanged with respect to the number of layers. (c) The calculated energy contour plot (in meV f.u.⁻¹) of an isolated NbOI₂ monolayer as a function of the polar displacement between two Nb atoms in one unit cell. One PE state and two energetically equivalent FE (AFE) states are identified. The energy of the PE state is set as the energy zero.

FE phases of NbOX₂ also show weak layer-dependent in-plane ferroelectricity due to the van der Waals interlayer interactions.³⁷ Shown in Fig. 3(b) are the calculated double-well potential curves for FE NbOI₂, where NbOI₂ monolayer, bilayer and bulk exhibit almost identical Nb polar displacement amplitude and FE potential depth in their FE phases (similar results are also found in NbOCl₂ and NbOBr₂). Based on layer-independent energetic stability and ferroelectricity, it is expected that in-plane FE polarization of almost the same magnitude (Table 1) can be obtained in all structural forms of NbOX₂, including 3D bulk and 2D layers. This unique feature makes 2D NbOX₂ distinct from group-IV monochalcogenides (e.g. GeS, SnS and SnTe), a well studied 2D FE material class with in-plane ferroelectricity. As it is demonstrated by experiment, 2D orthorhombic SnTe layers (γ phase) exhibit strong antipolar interlayer coupling, leading to a vanishing of ferroelectricity in the bulk phase and those multi-layers containing even numbers of layers.¹⁶

NbOX₂ monolayer can be stabilized into two polar phases, either a stable FE or metastable AFE phase, which contain two NbO₂X₄ octahedra in their unit cell. To evaluate the possibility of FE-to-AFE transition, we calculate 2D energy contour plot of NbOI₂ monolayer as a function of polar displacement of two Nb atoms within a single unit cell. As shown in Fig. 3(c), two degenerate FE (AFE) phases locate along the diagonal directions, separated by the PE phase. Transition between two FE phases with opposite polarizations does not need to go through the high-energy PE phase. Instead, it will cross a low-lying AFE phase, where the polar displacement of one Nb atom is reversed, while keeping the other almost unchanged.

Electronic structure and optical absorption properties

Layered NbOX₂, including 3D bulk and 2D layers, have FE ground states with spontaneous polarizations. The coupling between FE ordering and the electronic and optical properties of 2D materials can lead to interesting phenomena, such as the bulk photovoltaic^9,38 and ferrovalley effects.^12,39 We will investigate the electronic structures and optical absorption properties for the FE phases in bulk and monolayer NbOX₂.

Fig. 4 displays our calculated energy band structures and projected density of states (PDOS) for FE NbOX₂, obtained from HSE hybrid density functional calculations. All NbOX₂ systems, including 3D bulk and 2D monolayer, are predicted to be indirect semiconductors with nominal E_g (the energy difference between the filled and empty band edges, shown in Table 1) around 1.8 eV. Due to the hybridization of Nb-4d with anion p orbitals, largely dispersed energy bands around the conduction band minimum (CBM) in NbOX₂ are mainly contributed to by empty Nb-4d states. Below the Fermi level, the highest occupied valence band is almost dispersionless through the entire Brillouin zones of the monoclinic bulk and orthorhombic monolayer. For NbOCl₂ and NbOBr₂, such a flat band is isolated from other valence bands of hybridized anion p orbitals. Further orbital component analysis indicates that the dispersionless valence bands (marked as green circles in the energy bands of the monolayers) originate from localized Nb-d_z² states. Meanwhile, similar weakly dispersed energy bands also appear in conduction states (above CBM), corresponding to empty Nb-d_z² orbitals.

Fig. 4 Electronic structures for FE NbOX₂. HSE calculated energy band structures and projected density of states (PDOS) for FE bulk NbOX₂ with C2 monoclinic symmetry (upper panel) and FE NbOX₂ monolayer with Pmm2 orthorhombic symmetry (lower panel). Fermi energy level is set as energy zero. The choice of k-point path for the monoclinic Brillouin zone follows the convention from ref. 40. For NbOX₂ monolayer, the contribution of Nb-d_z² orbital to energy bands are highlighted by green circles. The variation of the energy separation between the occupied and empty Nb-d_z² orbitals at Γ point (ΔE_Nb-d) as a function of Nb–Nb paring distance (d_Nb–Nb) in FE NbOX₂ monolayers is also presented. NbOCl₂ monolayer with the shortest d_Nb–Nb has the strongest Peierls distortion intensity and therefore the largest ΔE_Nb-d.

Formation of the energy gap between the filled and empty Nb-d_z² states (ΔE_Nb-d) is a direct consequence of 1D Peierls distortion.³⁶ Each Nb⁴⁺ cation in NbOX₂ has one unpaired 4d electron. Paring of two 4d electrons simultaneously occurs in the NbOX₂ single unit cell, after periodic dimerization of neighboring Nb cations created by the 1D Peierls distortion. The energy gap is then formed as two paired electrons are favored to fully occupy the low-energy Nb-d_z² orbital, leaving the high-energy d_z² orbital empty. We double checked the electronic structures for NbOX₂ using spin polarized calculations, which also predict nonmagnetic ground states with paired electrons on the Nb-d_z² orbitals. In fact, the intensity of the 1D Peierls distortion in NbOX₂ can be well qualified by energy separation between the occupied and empty Nb-d_z² orbitals (ΔE_Nb-d). As shown in the last plot of Fig. 4, the strongest Peierls distortion with the shortest Nb–Nb distance in NbOCl₂ monolayer leads to the largest ΔE_Nb-d. The overall extraordinarily large ΔE_Nb-d in all three systems indicates that 1D Peierls distortion in NbOX₂ is ultra-stable, able to survive over thermal fluctuations at high temperature.

In layered NbOX₂, Nb-d_z² orbitals can hardly hybridize with any planar anion p orbitals. As a result, Nb-d_z² states correspond to the dispersionless energy bands that are quite localized in energy. More importantly, optical transitions from the localized Nb-d_z² valence band to other empty Nb-d orbitals from conduction bands are completely forbidden. The effective optical absorption in NbOX₂ comes from the transition between hybridized anion-p orbitals of the valence band and hybridized Nb-4d orbitals (except d_z²) of the conduction bands. Therefore, as far as optical absorption properties are considered, NbOCl₂ and NbOBr₂ are insulators with optical transition allowed E_g > 3.0 eV (Table 1). Due to the weak electronegativity of the I anion, bulk and monolayer NbOI₂ can exhibit semiconducting electronic and optical properties. Typically, as a ferroelectric semiconductor with an indirect E_g, bulk NbOI₂ can effectively absorb visible light with a photon energy above 2.0 eV (Fig. S4 in the ESI†), making it suitable for FE photovoltaic applications.^41,42

2D NbOX₂ layers have orthorhombic symmetry with anisotropic planar crystal structures, where Nb atoms bond with X–X edges along one planar direction and O atoms along the other. Similarly to many other 2D FE materials,^12,43 strong anisotropy in crystal structures of 2D NbOX₂ will be reflected in their energy bands and optical absorption properties as well. Using semiconducting monolayer NbOI₂ as an example, we will illustrate its highly anisotropic and strongly coupled electronic and optical absorption properties. Fig. 4(a) shows the band structure for FE NbOI₂ monolayer, where all the energy bands around the Fermi energy level are marked by colored circles whose radii are proportional to the contribution of the corresponding atomic orbitals. I anion has a weaker electronegativity than O, therefore the valence bands around the Fermi level mainly come from I-p orbitals. Except for the localized d_z² bands, most valence and conduction bands show anisotropic dispersion along two planar directions: valence I-p_y bands and low-lying conduction Nb-d_yz bands are highly dispersed along Γ–Y (reciprocal lattice equivalent of the crystallographic b axis), but weakly dispersed along the Γ–X direction.

Besides electronic structures, we have also investigated the anisotropic optical properties of NbOI₂ monolayer by simulating its optical absorption coefficient α, excited by incident light polarized along the x and y directions (corresponding to the a and b axes, respectively). As shown in Fig. 5(b), both α_x and α_y have optical absorption edges around 2.0 eV. Above the absorption edge, a strong optical absorption peak appears at 3.15 eV, when incident light is polarized along the x direction. The strength of this absorption peak in α_x is larger than that of α_y by almost 10 times, indicating a strong direction-dependent photo-absorption in NbOI₂ monolayer. As energy bands around the Fermi level are mainly contributed to by hybridized I-p and Nb-d orbitals, the low-energy optical absorption in NbOI₂ can be assigned to optical transitions from I-p to Nb-d orbitals. When NbOI₂ monolayer is excited by incident light polarized along the x direction, optical transitions between all weakly dispersed Nb-d_yz and I-p_y energy bands along the whole Γ–X direction (indicated by dashed red arrows in Fig. 5(a)) can contribute to optical absorption at 3.15 eV in α_x. While along the y direction, the weaker absorption at 3.15 eV comes from the optical transition (blue dashed arrow) of Nb-d_yz and I-p_y bands around the Γ point only. As a result, when monochromatic light with photon energy of 3.15 eV is incident on FE NbOI₂ monolayer, optical absorption along the nonpolar a axis is significantly larger than that along the polar axis, leading to a nearly linearly polarized optical absorption selectivity in NbOI₂ monolayer.³⁹

Fig. 5 Anisotropic electronic structure and optical absorptions of NbOI₂ monolayer. (a) Calculated band structure for FE NbOI₂ monolayer. The radii of the colored circles are proportional to the contribution of the corresponding atomic orbitals. The Fermi energy level is set as energy zero. (b) The optical absorption spectrum for FE NbOI₂ monolayer under incident light polarized along crystallographic a and b axes (x and y directions), respectively. Blue and red arrows in (a) indicate the optical transitions that contribute to the anisotropic optical absorption peaks around 3.15 eV, which are marked by the dashed line in (b).

Ferroelectric polarization switching and domain-wall properties

It is clearly shown that layered NbOX₂ systems have in-plane ferroelectricity at the monolayer limit. As a fundamental nature of ferroelectricity, FE polarizations in 2D NbOX₂ should be switchable under an applied electric field. Unlike other 2D FE materials with planar ferroelectricity (e.g. SnS and SnSe), FE ordering in 2D NbOX₂ monolayers indeed has a “1D” nature, as the ferroelectricity or antiferroelectricity is restricted along the polar b axis, but forbidden along the other nonpolar axis. Based on its 1D collinear nature, FE polarization switching in 2D NbOX₂ monolayers can be achieved through the direct reversal of the polarization direction by 180°, rather than rotation of polarization across the nonpolar axis.

During the polarization reversal (switching) process, FE domain-walls between oppositely orientated FE mono-domains will be formed in NbOX₂ monolayers. Meanwhile, the polarization reversal process can be carried out through the motion of FE domain-walls.^44,45Fig. 6(a) displays the supercell configuration containing two FE mono-domains with opposite polarization directions, separated by a 180° FE domain-wall located along the I–I lattice plane in NbOI₂ monolayer. Each mono-domain is composed of four NbOI₂ unit cells, stacking along the a axis. To preserve the periodic boundary condition required by first-principles plane wave calculations, we must include two coherent FE domain-walls (twin walls) in one supercell, so that the whole supercell can be repeated periodically along both axes.

Fig. 6 FE domain wall and polarization reversal process. (a) Atomic structure (including side and top views) for a 180° FE domain-wall (marked by the grey area) between two oppositely oriented FE mono-domains in NbOI₂ monolayer. The change of the associated Nb polar displacement as a function of the distance from the domain-wall is also given. Symbols are the calculated results, the solid line is fitted to the data based on domain wall theory and the dashed lines indicate the mono-domain values. (b) Minimum energy pathway (MEP) for polarization reversal and the associated domain wall motion process obtained from solid-state NEB calculations. Reversal of NbOI₂ polarization is achieved by motion of the FE domain wall from one I–I lattice plane to the next. Atomic structures for initial, saddle and final states along the simulated MEP are also presented. Black arrows indicate the polar displacement of Nb cations. All calculations were performed using the supercell configuration containing eight NbOI₂ unit cells. For clarity, only half of the domain-wall supercell is shown.

In order to evaluate how polarization is reversed across the FE domain-wall, we analyze the evolution of Nb polar displacement amplitude (d_Nb) as a function of Nb distance (r) away from the domain-wall. As shown in Fig. 6(a), polarization is reversed as d_Nb changes its sign abruptly across the domain-wall. d_Nb of each Nb cation almost fully recovers to its mono-domain value (indicated by dashed lines). As a result, the 180° FE domain-wall supercell configuration is geometrically similar to a NbOI₂ superlattice, where AFE phases are sandwiched between two FE phases with opposite polarizations. As the energy difference between AFE and FE phases of NbOI₂ monolayer is as low as 1.24 meV f.u.⁻¹, the 180° FE domain-wall has a very small domain wall energy. The d_Nb–r curve shown in Fig. 6(a) can be quantitatively described using a function corresponding to the soliton solution of one-dimensional fourth-order Landau–Ginzburg domain wall theory⁴⁴ as:

where d⁰_Nb specifies the FE mono-domain value and 2ξ_DW refers to the domain-wall width. After least-square fitting to eqn (1), an ultra narrow domain-wall width of only ∼1 Å is predicted for domain-wall configuration in NbOI₂ monolayer. Unlike the 180° FE domain-wall of perovskite PbTiO₃, where there is the intermediate PE phase separating two FE mono-domains,^44,45 the absence of PE phase in the NbOI₂ domain-wall configuration is responsible for its ultra narrow domain-wall width.

We then explore the kinetic process for domain-wall motion by simulating the corresponding minimum energy pathway (MEP) trajectory and the associated energy barrier height, using a generalized solid-state NEB method.⁴⁶ Such a method allows both atomic positions and lattice parameters to relax along the pathway. As shown in Fig. 6(b), transition from the initial to the final domain-wall configurations resembles the motion of the domain wall from one I–I lattice plane to the next, accompanied by the reversal of polarization direction of one group of Nb cations by 180°. At the saddle point configuration along the MEP, the domain-wall moves to the lattice plane where Nb cation is located. This very Nb cation is then constrained to have zero d_Nb, while d_Nb of other Nb cations remain unchanged. As a result, the MEP associated with domain-wall motion is similar to the energy pathway for transition from single FE to AFE phase, as we obtained earlier from the energy contour plot in Fig. 3(c). The overall barrier height for polarization reversal of the Nb cation in the NbOI₂ monolayer is predicted to be 14.9 meV, which is much smaller than the polarization rotation barrier (∼40 meV) in FE perovskite PbTiO₃.⁴⁷ Therefore, the presence of the low-lying AFE phase can lead to a small energy barrier, which is beneficial for easy domain-wall motion and polarization switching in NbOI₂ monolayer.

Using the simulated energy barrier, we can estimate how large the electric field needs to reverse FE polarization in NbOI₂ monolayer. The critical electric field E_C is given by: E_C ≃ E_barrier/P·V, where P is the spontaneous polarization and V is the normalized cell volume for NbOI₂ monolayer. The estimated critical electric field E_C ≃ 0.63 MV cm⁻¹, which is readily accessible by experiment.⁴⁸ Therefore, electric field induced switching of FE polarization can be achievable under laboratory conditions.

Structural phase transition at finite temperature

Up to now, we have simulated the structural, electronic and polarization switching properties for FE NbOX₂ layers using zero-temperature DFT calculations. For practical applications of ferroelectrics, 2D NbOX₂ should have relatively high Curie temperature T_C, so that their polarizations can still persist above room temperature. In the following, we will investigate the ferroelectricity and structural phase-transition of NbOX₂ monolayers at finite temperature using both Monte Carlo (MC) and ab initio molecular dynamics (MD) simulations.

The geometry of 2D NbOX₂ monolayers can be represented by a 2D lattice grid containing a number of Nb cations, where each Nb cation at 2D grid point i [i is a collapsed index of its 2D position (m, n)] has a unique polar displacement d_i. Due to the 1D nature of polarization, d_i of Nb cations are collinearly arranged into columns along the polar axis. Using polar displacement d_i as the order parameter, the configuration of NbOX₂ monolayer at any condition can be unambiguously specified. Therefore, we can express the free energy of NbOX₂ monolayer as Landau–Ginzburg expansion of order parameter d_i as:¹³

where the first two terms describe the FE energy for each Nb cation and parameters A and B can be obtained after fitting a double-well potential for FE NbOX₂ monolayers. The last two terms correspond to the interactions of d_i and d_j between two nearest neighboring Nb cations. Here, we include the parameters C_x and C_y, which measure the interactions along the non-polar and polar axes respectively (Fig. S5 in the ESI†). Due to the high energy cost to form a head-to-head/tail-to-tail polar configuration (Nb cations within the same column have antiparallel d_i) in NbOX₂, parameter C_y ≫ C_x. As shown in Table 2, the absolute values for most of the parameters A–C follow the trend: NbOCl₂ > NbOBr₂ > NbOI₂, consistent with the fact that FE NbOCl₂ has the strongest ferroelectricity.

Table 2 Parameters used for MC simulation in eqn (2). Except the fourth-ordered parameter B (in meV Å⁻⁴ per f.u.), parameters A, C_x and C_y are all in units of meV Å⁻² per f.u. Curie temperatures T_C are obtained via fitting to MC data based on eqn (3)

	A	B	C x	C y	T C (K)
NbOCl2	−1.29 × 10^3	3.03 × 10^4	48.23	6.08 × 10^3	396
NbOBr2	−1.11 × 10^3	2.74 × 10^4	38.68	4.28 × 10^3	283
NbOI2	−1.01 × 10^3	2.85 × 10^4	19.14	4.82 × 10^3	242

---

Based on the effective Hamiltonian we developed for NbOX₂ monolayers, we can investigate temperature induced structural transition using MC simulations (computational details in the ESI†). For comparison, we also perform parameter-free ab initio MD simulations, where the time evolution of the instantaneous temperature, total energies and cation polar displacement of NbOX₂ can be obtained (Fig. S6 in the ESI†). Due to the large computational cost, we restrict our MD simulations of NbOX₂ monolayers to selected target temperature. Fig. 7 shows our simulated temperature dependent macroscopic average polar displacement 〈d(T)〉 in NbOX₂ monolayers. The structural transition from FE to PE phase is identified, as the simulated 〈d(T)〉 drops abruptly to zero around T_C. Moreover, at the selected temperature, ab initio MD predicted ensemble average 〈d(T)〉 are in overall good agreement with MC results (the largest derivation occurs in NbOI₂, with ΔT ∼ 60 K), which validates the effective Hamiltonian and parameters we used for MC simulations. Curie temperature T_C can be quantitatively determined after fitting MC simulated 〈d(T)〉 as follows:¹³

where δ is critical exponent for the FE–PE phase transition. As shown in Table 2, the tendency of T_C for FE NbOX₂ monolayers also follows: NbOCl₂ > NbOBr₂ > NbOI₂. NbOCl₂ monolayer is predicted to have the highest T_C of 396 K, and can exhibit considerable polarization at room temperature.

Fig. 7 Structural phase-transitions and Curie temperatures of NbOX₂ monolayers. Upper panel: Temperature dependence of the average Nb polar displacement obtained from MC and ab initio MD simulations for NbOCl₂, NbOBr₂ and NbOI₂ monolayers. Average polar displacement 〈d(T)〉 is normalized with respect to the corresponding 0 K value. Black lines are fitted to MC data based on eqn (3). The error bars of Nb polar displacement highlight the thermal fluctuations during the MD simulations. Lower panel: MD simulated polarization reversal process in NbOI₂ monolayer. The normalized polar displacements (d_i(T)/|d(0 K)|) for each Nb cation in a 8 × 8 lattice grid are presented when the temperature is below, approaching and above Curie temperature T_C. Polarization reversal occurs when T ≈ T_C, and is initiated by switching d_i of a single Nb from one cation column along the polar axis (y direction).

In order to understand the microscopic mechanism governing FE–PE phase transition, we further examine temperature dependent cation displacement evolution. The lower panel of Fig. 7 plots the real space mapping of d_i from the whole NbOI₂ lattice grid during MD simulations after the system reaches thermal equilibrium at target temperature below, approaching and above T_C, respectively. Below T_C, the magnitude of d_i starts to decrease from its zero temperature value as temperature increases. The polarization reversal event (changing sign of d_i) does not occur until T ≈ T_C, where |d_i| of some Nb cations become small enough so that thermal energy is comparable with the energy cost to form a head-to-head/tail-to-tail polar configuration. The polarization reversal process around T_C is initiated by switching d_i of a single Nb cation from one cation column. Then the whole column simultaneously switches its polarization to the opposite direction, to avoid the unfavorable head-to-head/tail-to-tail configuration. The FE–PE phase transition is achieved through dynamic reversal of column polarization as mentioned above. As AFE phase in NbOX₂ is more stable than the distortion-free PE phase. Even when T > T_C, the system stays in a thermal-equilibrium PE state, where the overall polarization is almost zero, but individual columns can still have non-zero and randomly oriented polarization along the polar axis. As a result, PE state above T_C is more like a disordered AFE phase. Moreover, the energy cost associated with “collective” reversal of whole column polarization is negligible compared to the “isolated” cation dipole switching energy (C_y ≫ C_x, Table 2), thus T_C of NbOX₂ monolayer is mainly determined by “isolated” rather than “collective” cation dipole switching events.⁴⁹ Therefore, T_C is not simply proportional to the energy difference associated with “collective” polarization reversal (FE-to-AFE phase transition). Even though AFE phase is close to FE phase in energy, NbOX₂ monolayers can still have T_C around or even above room temperature.

Experimental outlook and FE device by design

We have provided strong evidence that 2D NbOX₂ possess intrinsic in-plane ferroelectricity and antiferroelectricity. Experimentally, bulk NbOX₂ have already been prepared using solid-state sintering methods, while their structure polarity at room temperature was verified by XRD measurement.^28–31 Direct identification of ferroelectricity or antiferroelectricity in bulk NbOX₂ requires further measurement of their polarization–electric field (P–E) hysteresis loops. 2D NbOX₂ layers, including multilayers and monolayers can be obtained by direct exfoliation of bulk samples or chemical vapor deposition. As long as the temperature is below T_C, 2D NbOX₂ layers can crystalize in a ground state FE phase with polarization comparable with that of bulk NbOX₂, after poling by electric field (above critical field E_C) applied along the polar axis. Based on our simulation, 2D NbOX₂ monolayers exhibit variable FE polarizations and T_C. Chemical alloys made of individual NbOX₂ end members can offer a wide range of material choice to achieve optimal T_C and FE properties. For example, a 2D NbOCl_xI_2−x monolayer is expected to possess T_C above room temperature and semiconducting properties comparable with NbOI₂. More specifically, increasing concentration x can lead to a higher T_C and larger optical transition allowed E_g in the NbOCl_xI_2−x monolayer. Alternatively, one can also use experimentally accessible tensile strain to engineer the FE properties of NbOX₂ monolayers. As shown in Fig. S7 of the ESI,† more than 6 times enhancement of the FE potential depth and twice of the Nb polar displacement can be achieved with 3% tensile strain applied along the polar axis of NbOX₂ monolayers. As a result, tensile strained NbOBr₂ and NbOI₂ monolayers are expected to have elevated T_C above room temperature. Due to linearly polarized optical selectivity of the semiconducting NbOI₂ monolayer, its in-plane polar and non-polar axes can be easily distinguished using incident polarized light. As NbOX₂ have the interconnected octahedral structural framework similar to perovskite oxides, after imaging of Nb off-center polar displacement (0.13–0.14 Å) using annular bright-field scanning transmission electron microscopy (ABF-STEM),^50,51 one can directly observe and quantify FE polarization appearing in 2D NbOX₂ layers. Moreover, measuring the shear piezoelectric effect arising from in-plane polarization by lateral piezoelectric force microscopy (PFM)³⁷ can also provide a direct experimental proof of FE polarization in 2D NbOX₂.

The 1D polarization nature and coexistence of FE/AFE phases make 2D NbOX₂ exhibit the unique polarization reversal and structural phase-transition features. In experiment, AFE perovskite oxides have been widely used as dielectric capacitors for electrostatic energy storage through electric charging and discharging processes.⁵² In particular, high energy densities are achieved in perovskite solid solutions near the phase boundaries where FE and AFE phases coexist.⁵³ We therefore propose the application of NbOX₂ monolayers as 2D dielectric capacitors.

Fig. 8(a) displays a schematic diagram for a 2D single-layered capacitor obtained by lateral growth of a heterostructure between a NbOX₂ monolayer and another 2D metallic material, where an electric field can be applied along the polar axis through a metal electrode. In order to evaluate the performance of a 2D NbOX₂ based capacitor, we choose NbOI₂ monolayer as an example, and investigate the electric-field induced phase transition by simulating its P–E loop around room temperature using MC simulations. When electric field E is applied along the polar axis, the additional energy term −(E·P)V will be incorporated into the effective Hamiltonian in eqn (2). The effective polarization P can be computed from Nb polar displacement d_i and Born effective charge Z_i* by: . Therefore, the electric field related energy term can be simplified as: −E(Σ_id_i·Z_i*). Based on our predictions, NbOI₂ monolayer adopts a disordered AFE (PE) phase at room temperature. The simulated P–E loop in Fig. 8(b) shows that under the charging process (increasing E) P increases abruptly around 0.15 MV cm⁻¹, indicating that a small critical field E_C can trigger the AFE-to-FE transition in NbOI₂ monolayer at room temperature. Above E_C, P increases linearly within FE phase upon further charging. Remarkably, an almost identical P–E curve associated with the FE-to-AFE transition is obtained under the discharging process (decreasing E). As a result, under the applied electric field, NbOI₂ monolayer possesses P–E loops with almost no hysteresis, which is consistent with its small FE/AFE energy difference and low polarization switching energy barrier at room temperature.⁵⁴ Similarly to NbOI₂, NbOBr₂ monolayer also crystalizes in a disordered AFE phase at room temperature and exhibits nearly hysteresis-free P–E loops. While NbOCl₂ monolayer has a non-zero FE polarization at room temperature, it features P–E double hysteresis loops of typical ferroelectrics. Due to the energy loss associated with the hysteresis, NbOCl₂ monolayer is less suitable for electrostatic energy storage than NbOBr₂ or NbOI₂. A 2D capacitor based on NbOBr₂ or NbOI₂ monolayers can exhibit nearly zero energy loss during the charging and discharging processes, yet has the advantage of operating at a low electric field, which can potentially enable nearly 100% efficiency for electrostatic energy storage in the 2D monolayer.

Fig. 8 Electrostatic energy storage in a 2D FE capacitor. (a) Schematic model for a 2D FE capacitor device based on NbOX₂ monolayers. (b) MC simulated P–E loops for NbOX₂ monolayers at room temperature. Electric field E is applied along the polar axis. Results corresponding to the electric charging and discharging processes are marked in red and blue colors, respectively.

Conclusions

In summary, we have identified niobium oxide dihalides NbOX₂ as a new family of 2D van der Waals layered materials, exhibiting intrinsic in-plane ferroelectricity and antiferroelectricity. Comprehensive investigations regarding structural, ferroelectric, electronic, optical absorption and domain-wall properties of 2D NbOX₂ have been carried out using first-principles calculations. NbOX₂ adopts a 2D structural framework composed of interconnected NbO₂X₄ octahedra, where considerable cation off-center polar displacement along one planar axis can lead to 1D ferroelectricity or antiferroelectricity. Along the other planar axis, the ultra-stable Peierls distortion between Nb⁴⁺ cations gives rise to the semiconducting or insulating electronic properties, which are preferable for stable polarizations. In-plane ferroelectricity or antiferroelectricity in 2D NbOX₂ can be sustained down to the monolayer limit, and the energetic stabilities and polar magnitudes are almost independent of the number of layers. Coexistence of FE and AFE phases makes 2D NbOX₂ exhibit easily switchable polarizations. Moreover, the strong coupling between FE order and optical absorption leads to a nearly linearly polarized optical selectivity in the semiconducting NbOI₂ monolayer. Besides zero temperature DFT calculations, we also develop an effective Hamiltonian to simulate the structural phase transitions of NbOX₂ monolayers under finite temperature and an applied electric field. The NbOCl₂ monolayer is predicted to be a robust FE material with room temperature stable ferroelectricity, while NbOBr₂ and NbOI₂ monolayers adopt a PE state composed of disordered AFE phases at room temperature, but will undergo a structural transition into a stable FE phase under a small applied electric field. Based on our simulations, P–E loops characterizing AFE–FE transitions have nearly zero hysteresis, therefore a capacitor composed of NbOBr₂ or NbOI₂ has a great advantage for highly efficient electrostatic energy storage. We hope our work will motivate future experimental exploration of 2D ferroelectricity and antiferroelectricity in layered NbOX₂, and realization of superior energy storage performance that our simulations predicted.

Computational methods

Our first-principles calculations are performed based on density functional theory (DFT) as implemented in the Vienna ab initio Simulation Package (VASP-5.4.1).^55,56 A plane-wave basis set within the projector augmented-wave method⁵⁷ is employed, using a 600 eV plane-wave energy cutoff. We simulate NbOX₂ bulk and 2D layered structures using the strongly constrained and appropriately normed (SCAN) meta-GGA functional.⁵⁸ After inclusion of the revised Vydrov–van Voorhis nonlocal correlation (rVV10),^33,59 SCAN + rVV10 functional can accurately predict the interlayer spacings as well as intralayer lattice parameters for bulk NbOX₂ (Table S2 in the ESI†). 2D NbOX₂ layers are modeled as slabs with a vacuum region of more than 20 Å along the out-of-plane direction. A Γ-centered Monkhorst–Pack k-point grid of about 40 k-points per Å⁻¹ spacing is used for k-point sampling. NbOX₂ bulk and layers are fully optimized until the residual Hellmann–Feynman forces are smaller than 1 meV Å⁻¹ and the stresses less than 0.1 kbar. Phonon frequencies and eigenvectors are calculated based on the finite difference method implemented in the Phonopy package.⁶⁰ The electronic structures and optical absorption spectra of NbOX₂ systems are further calculated using the Heyd–Scuseria–Ernzerhof (HSE06) hybrid functional.^61,62 As band gap changes due to spin–orbital coupling (SOC) effects are quite small (ΔE_g < 0.02 eV for bulk NbOI₂), the SOC effect is excluded during HSE calculations. The electronic contribution to the polarization is calculated following the Berry phase formalism.⁶³ The generalized solid-state nudged elastic band (ss-NEB) method^46,64,65 is used to simulate the minimum energy path (MEP) and associated energy barriers for various kinetic processes.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

Work in XJTU is supported by funding from the National Science Foundation of China under Contract No. 11574244, 51320105014 and 51621063, the State Key Laboratory for Mechanical Behavior of Materials, as well as computational support from the National supercomputer center (NSCC) in Tianjin. X. C. Z. is supported by the National Science Foundation (NSF) through the Nebraska Materials Research Science and Engineering Center (MRSEC) (Grant No. DMR-1420645), and UNL Holland Computing Center. Y. J. thanks Dr Liang Ma, Dr Chongqin Zhu, Dr Ruixiang Fei and Prof. Yunhao Lu for helpful discussions.

References

1. M. Dawber, K. M. Rabe and J. F. Scott, Rev. Mod. Phys., 2005, 77, 1083–1130.
2. A. Belianinov, Q. He, A. Dziaugys, P. Maksymovych, E. Eliseev, A. Borisevich, A. Morozovska, J. Banys, Y. Vysochanskii and S. V. Kalinin, Nano Lett., 2015, 15, 3808–3814.
3. F. Liu, L. You, K. L. Seyler, X. Li, P. Yu, J. Lin, X. Wang, J. Zhou, H. Wang, H. He, S. T. Pantelides, W. Zhou, P. Sharma, X. Xu, P. M. Ajayan, J. Wang and Z. Liu, Nat. Commun., 2016, 7, 12357.
4. L. Niu, F. Liu, Q. Zeng, X. Zhu, Y. Wang, P. Yu, J. Shi, J. Lin, J. Zhou, Q. Fu, W. Zhou, T. Yu, X. Liu and Z. Liu, Nano Energy, 2019, 58, 596–603.
5. J. Junquera and P. Ghosez, Nature, 2003, 422, 506–509.
6. D. D. Fong, G. B. Stephenson, S. K. Streiffer, J. A. Eastman, O. Auciello, P. H. Fuoss and C. Thompson, Science, 2004, 304, 1650–1653.
7. J. Feng, X. Qian, C.-W. Huang and J. Li, Nat. Photonics, 2012, 6, 866–872.
8. A. Castellanos-Gomez, R. Roldán, E. Cappelluti, M. Buscema, F. Guinea, H. S. J. van der Zant and G. A. Steele, Nano Lett., 2013, 13, 5361–5366.
9. T. Rangel, B. M. Fregoso, B. S. Mendoza, T. Morimoto, J. E. Moore and J. B. Neaton, Phys. Rev. Lett., 2017, 119, 067402.
10. P. Luo, F. Zhuge, Q. Zhang, Y. Chen, L. Lv, Y. Huang, H. Li and T. Zhai, Nanoscale Horiz., 2019, 4, 26–51.
11. M. Wu and X. C. Zeng, Nano Lett., 2016, 16, 3236–3241.
12. H. Wang and X. Qian, 2D Mater., 2017, 4, 015042.
13. R. Fei, W. Kang and L. Yang, Phys. Rev. Lett., 2016, 117, 097601.
14. M. Mehboudi, B. M. Fregoso, Y. Yang, W. Zhu, A. van der Zande, J. Ferrer, L. Bellaiche, P. Kumar and S. Barraza-Lopez, Phys. Rev. Lett., 2016, 117, 246802.
15. K. Chang, J. Liu, H. Lin, N. Wang, K. Zhao, A. Zhang, F. Jin, Y. Zhong, X. Hu, W. Duan, Q. Zhang, L. Fu, Q.-K. Xue, X. Chen and S.-H. Ji, Science, 2016, 353, 274–278.
16. K. Chang, T. P. Kaloni, H. Lin, A. Bedoya-Pinto, A. K. Pandeya, I. Kostanovskiy, K. Zhao, Y. Zhong, X. Hu, Q.-K. Xue, X. Chen, S.-H. Ji, S. Barraza-Lopez and S. S. P. Parkin, Adv. Mater., 2018, 31, 1804428.
17. W. Ding, J. Zhu, Z. Wang, Y. Gao, D. Xiao, Y. Gu, Z. Zhang and W. Zhu, Nat. Commun., 2017, 8, 14956.
18. Y. Zhou, D. Wu, Y. Zhu, Y. Cho, Q. He, X. Yang, K. Herrera, Z. Chu, Y. Han, M. C. Downer, H. Peng and K. Lai, Nano Lett., 2017, 17, 5508–5513.
19. C. Cui, W.-J. Hu, X. Yan, C. Addiego, W. Gao, Y. Wang, Z. Wang, L. Li, Y. Cheng, P. Li, X. Zhang, H. N. Alshareef, T. Wu, W. Zhu, X. Pan and L.-J. Li, Nano Lett., 2018, 18, 1253–1258.
20. J. Xiao, H. Zhu, Y. Wang, W. Feng, Y. Hu, A. Dasgupta, Y. Han, Y. Wang, D. A. Muller, L. W. Martin, P. Hu and X. Zhang, Phys. Rev. Lett., 2018, 120, 227601.
21. S. Wan, Y. Li, W. Li, X. Mao, W. Zhu and H. Zeng, Nanoscale, 2018, 10, 14885–14892.
22. F. Xue, W. Hu, K.-C. Lee, L.-S. Lu, J. Zhang, H.-L. Tang, A. Han, W.-T. Hsu, S. Tu, W.-H. Chang, C.-H. Lien, J.-H. He, Z. Zhang, L.-J. Li and X. Zhang, Adv. Funct. Mater., 2018, 28, 1803738.
23. S. N. Shirodkar and U. V. Waghmare, Phys. Rev. Lett., 2014, 112, 157601.
24. Y. Fang, X. Hu, W. Zhao, J. Pan, D. Wang, K. Bu, Y. Mao, S. Chu, P. Liu, T. Zhai and F. Huang, J. Am. Chem. Soc., 2019, 141, 790–793.
25. G. Cheon, K.-A. N. Duerloo, A. D. Sendek, C. Porter, Y. Chen and E. J. Reed, Nano Lett., 2017, 17, 1915–1923.
26. N. Mounet, M. Gibertini, P. Schwaller, D. Campi, A. Merkys, A. Marrazzo, T. Sohier, I. E. Castelli, A. Cepellotti, G. Pizzi and N. Marzari, Nat. Nanotechnol., 2018, 13, 246.
27. H. Schnering and H. Wöhrle, Angew. Chem., 1963, 75, 684.
28. H. Hillebrecht, P. Schmidt, H. Rotter, G. Thiele, P. Zönnchen, H. Bengel, H.-J. Cantow, S. Magonov and M.-H. Whangbo, J. Alloys Compd., 1997, 246, 70–79.
29. D. Drobot and E. Pisarev, Russ. J. Inorg. Chem., 1984, 29, 1561–1563.
30. J. Beck and C. Kusterer, Z. Anorg. Allg. Chem., 2006, 632, 2193–2194.
31. J. Rijnsdorp and F. Jellinek, J. Less-Common Met., 1978, 61, 79–82.
32. M.-H. Whangbo, Inorg. Chem., 1982, 21, 1721–1723.
33. H. Peng, Z.-H. Yang, J. P. Perdew and J. Sun, Phys. Rev. X, 2016, 6, 041005.
34. P. Zonnchen, G. Thiele, C. Hess, C. Schlenker, H. Bengel, H. Cantow, S. Magonov, D. Seo and M. Whangbo, New J. Chem., 1996, 20, 295–300.
35. M. Ruck, Acta Crystallogr., Sect. C: Cryst. Struct. Commun., 1995, 51, 1960–1962.
36. S. van Smaalen, Acta Crystallogr., Sect. A: Found. Crystallogr., 2004, 61, 51–61.
37. L. You, F. Liu, H. Li, Y. Hu, S. Zhou, L. Chang, Y. Zhou, Q. Fu, G. Yuan, S. Dong, H. Fan, A. Gruverman, Z. Liu and J. Wang, Adv. Mater., 2018, 30, 1803249.
38. S. M. Young and A. M. Rappe, Phys. Rev. Lett., 2012, 109, 116601.
39. X.-W. Shen, W.-Y. Tong, S.-J. Gong and C.-G. Duan, 2D Mater., 2018, 5, 011001.
40. W. Setyawan and S. Curtarolo, Comput. Mater. Sci., 2010, 49, 299–312.
41. I. Grinberg, D. V. West, M. Torres, G. Gou, D. M. Stein, L. Wu, G. Chen, E. M. Gallo, A. R. Akbashev, P. K. Davies, J. E. Spanier and A. M. Rappe, Nature, 2013, 503, 509–512.
42. H. Wang, G. Gou and J. Li, Nano Energy, 2016, 22, 507–513.
43. S. Guan, C. Liu, Y. Lu, Y. Yao and S. A. Yang, Phys. Rev. B, 2018, 97, 144104.
44. B. Meyer and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 2002, 65, 104111.
45. S. Liu, I. Grinberg and A. M. Rappe, Nature, 2016, 534, 360.
46. D. Sheppard, P. Xiao, W. Chemelewski, D. D. Johnson and G. Henkelman, J. Chem. Phys., 2012, 136, 074103.
47. J. Hong and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 84, 115107.
48. J. Zhang, X. Ke, G. Gou, J. Seidel, B. Xiang, P. Yu, W.-I. Liang, A. M. Minor, Y.-h. Chu, G. Van Tendeloo, X. R. Ren and R. Ramesh, Nat. Commun., 2013, 4, 2768.
49. Q. Yang, M. Wu and J. Li, J. Phys. Chem. Lett., 2018, 9, 7160–7164.
50. S. Van Aert, S. Turner, R. Delville, D. Schryvers, G. Van Tendeloo and E. K. H. Salje, Adv. Mater., 2012, 24, 523–527.
51. C. Ma, Y. Lin, H. Yang, H. Tian, L. Shi, J. Zeng and J. Li, Adv. Mater., 2015, 27, 6328–6332.
52. H. Wang, Y. Liu, T. Yang and S. Zhang, Adv. Funct. Mater., 2019, 29, 1807321.
53. B. Peng, Q. Zhang, X. Li, T. Sun, H. Fan, S. Ke, M. Ye, Y. Wang, W. Lu, H. Niu, J. F. S. Scott, X. Zeng and H. Huang, Adv. Electron. Mater., 2015, 1, 1500052.
54. B. Xu, J. Íñiguez and L. Bellaiche, Nat. Commun., 2017, 8, 15682.
55. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169–11186.
56. G. Kresse and D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758–1775.
57. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953–17979.
58. J. Sun, A. Ruzsinszky and J. P. Perdew, Phys. Rev. Lett., 2015, 115, 036402.
59. R. Sabatini, T. Gorni and S. de Gironcoli, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 87, 041108.
60. L. Chaput, A. Togo, I. Tanaka and G. Hug, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 84, 094302.
61. J. Heyd, G. E. Scuseria and M. Ernzerhof, J. Chem. Phys., 2003, 118, 8207–8215.
62. J. Heyd, G. E. Scuseria and M. Ernzerhof, J. Chem. Phys., 2006, 124, 219906.
63. R. D. King-Smith and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 1651–1654.
64. G. Henkelman, B. P. Uberuaga and H. Jónsson, J. Chem. Phys., 2000, 113, 9901–9904.
65. G. Henkelman and H. Jónsson, J. Chem. Phys., 2000, 113, 9978–9985.

---

Footnote

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

---