Raman Spectra of Bulk and Few-Layer GeSe From First-Principles Calculations

Raman spectra play a significant role in the study of polar materials. Herein, we report the influence of strain and interlayer shift on vibration responses in bulk and few-layer ferrovalley material GeSe in different polarization states (ferroelectric/FE and antiferroelectric/AFE) based on density functional theory and density functional perturbation theory calculations. We find A g 1 mode shifts by about 10 cm−1 from monolayer to bilayer and trilayer due to the interlayer coupling. The A g 3 mode on behalf of FE mode is observed that is consistent with the experiments in bulk and few-layer GeSe. Meanwhile, in our calculations, with the transition between AFE and FE state in the bilayer and trilayer, the Raman frequency of A g 2 and A g 3 mode decrease obviously whereas that of A g 1 mode increases. Interestingly, the Raman peaks shifted a lot due to the strain effect. We expect these variations in the Raman spectroscopy can be employed to identify the status of GeSe films, e.g., the AFE or FE state, and the number of layers in experiments.


INTRODUCTION
Ferroelectric (FE) materials with a stable spontaneous polarization that could be switched under external electric field have been widely studied and exploited in multifunctional devices such as ferroelectric synapse, field-effect transistors, and ferroelectric tunnel junction (Scott and Paz de Araujo, 1989;Mathews et al., 1997;Velev et al., 2007;Garcia et al., 2009;Huang et al., 2018b;Shen et al., 2019;Tian et al., 2019;Guan et al., 2020) due to their abundant physics. Until now, perovskite ferroelectric oxides are mostly utilized in these devices. However, with the decrease of thickness, ferroelectricity can hardly be maintained because of the existence of depolarization field in ferroelectric films (Junquera and Ghosez, 2003), which severely restrain the application of traditional perovskite ferroelectrics in miniaturized and high-density devices. In addition, the interfacial defect caused by lattice mismatch can also destroy the film ferroelectricity (Duan et al., 2006;Wang et al., 2010). The emergence of twodimensional (2D) ferroelectrics provides an opportunity to resolve these difficulties (Wu and Jena, 2018). Recently, plenty of 2D ferroelectrics have been successfully exfoliated from a bulk structure in experiments, including in-plane ferroelectricity (Chang et al., 2016;Chang et al., 2019;Higashitarumizu et al., 2020), and out-of-plane ferroelectricity Zhou et al., 2017;You et al., 2019;Yuan et al., 2019). Beyond experimental works, theoretical studies have predicted that ferroelectricity can survive in 2D materials (Ding et al., 2017;Huang et al., 2018a), in which some even possess noncollinear ferroelectric ordering Song et al., 2021).
Among these works, group IV chalcogenides (MX, M Sn, Ge; X S, Se) with few-layer have been discovered with intrinsic ferroelectricity and antiferroelectricity in the experiment (Fei et al., 2016), with a fantastic optical selective property as a polarizer (Shen et al., 2018), valley physics (Rodin et al., 2016), a high absorption coefficient as photovoltaic cells (Franzman et al., 2010;Shi and Kioupakis, 2015), robust ferroelectricity as nonvolatile storage (Wang and Qian, 2017) and so on (Yagmurcukardes et al., 2016). On the other hand, FE and antiferroelectric (AFE) phase transition is also predicted by interlayer sliding (Xu et al., 2021), suggesting an AFE tunnel junction can be realized in these materials .
As a most general and powerful tool to study crystal structures and vibration properties, Raman spectroscopy has been diffusely employed due to the advantage of nondestructive to the sample and easy sample preparation. Especially in 2D materials, the Raman spectrum can be used to precisely identify the number of layers because of its high sensitivity of thickness. For instance, in graphene (Ferrari et al., 2006), the G peak down-shifts with the decrease of layers, which is a sign to determine the number of layers. A similar case occurs in the theoretical work of PbI 2 (Yagmurcukardes et al., 2018), the larger number of layers, the more blueshift with the Raman peak of A g 1 mode. On the contrary, a redshift of A g and B 1g modes with the increased layers in few-layer phosphorene (Feng et al., 2015) could be observed. Moreover, the information of interlayer stacking and an effect of external field can be also detected by Raman spectroscopy (Zhang et al., 2015). Therefore, investigation using Raman spectroscopy could help to understand charming properties in 2D materials Saboori et al., 2019).
However, it should be mentioned that few-layer GeSe is rarely explored in the experiment due to the great challenge for sample preparation, leading to difficulty to study their layer numbers and phase transition related Raman spectroscopy. Consequently, systematic study of Raman spectroscopy of layered GeSe is necessary by theoretical calculation, which can provide guidelines for evaluating phonon-related characterization of GeSe based on Raman spectroscopy.
In this paper, we fully investigate the Raman spectroscopy, concerning the bulk and few-layer GeSe from monolayer to trilayer, by first-principles calculations. The van der Waals (vdW) correction is adopted in our calculations by comparing the Raman spectrum and crystal structure with experiment results in bulk GeSe. We find four vibration modes of A g 1 , B 3g , A g 2 and A g 3 in bulk GeSe, in which A g 2 is absent in the experiment but exists in SnSe (Yang et al., 2018). For monolayer GeSe, these vibration peaks are located at 81.84 cm −1 , 97.54 cm −1 , 141.94 cm −1 and 183.14 cm −1 , respectively. With the increase of layer number, A g 1 mode has a redshift but A g 2 and A g 3 exhibit a large blue shift as large as 40 cm −1 and 30 cm −1 , which verify the importance of interlayer interaction and can be used to identify the thickness of GeSe. In addition, we find strain has a remarkable influence on the Raman spectrum in different structures, suggestive of the possibility to probe the strain effect using optical method. More importantly, AFE/FE phase transition could be triggered by some methods such as interlayer sliding or strain. Herein, we use interlayer shift to realize the transition in few-layer as well as study the Raman spectroscopy on transition structures. We find an A g 1 mode increase of about 10 cm −1 whereas A g 2 and A g 3 decrease with the transition from AFE to FE phase, helping us to determine the AFE state or FE state. Our calculation could give information about the number of layers, whether the system under strain state and the ferroelectric phase of layered GeSe.

CALCULATION DETAILS
The optimized structures are calculated by employing the PWSCF package of the QUANTUM-ESPRESSO (Giannozzi et al., 2009) within the density functional theory (DFT) (Hohenberg and Kohn, 1964;Kohn and Sham, 1965). We adopt the local density approximation (LDA) with Perdew-Zunger parametrization (Perdew and Zunger, 1981) and the generalized gradient approximation (GGA) with Perdew-Burke-Ernzehof parametrization (Perdew et al., 1996) to evaluate the exchange-correlation energy and consider a modified norm-conserving pseudopotential to describe the valence electron-ion interactions (Gonze et al., 1991). In our calculations, the vdW correction was fully taken into account by using the DFT + D2 method (Grimme, 2006;Barone et al., 2009). To avoid the spurious interactions between periodic images, vacuum spacing of 20 A was set along the c-direction. The energy cutoff was set to 50 Ry and a Brillouin zone (BZ) integration is adopted with a k-grid density of 7 × 7 × 3 for bulk structures and 7 × 7 × 1 for few-layer structures via using the k-points scheme. All structures are optimized until the Hellman-Feynman force is below 10 −6 Ry/Bohr and the convergence of electric energy is of about 10 −4 Ry/atom. The related phonon vibration frequencies are calculated by diagonalizing the force constant matrix within the density functional perturbation theory (DFPT) (Baroni et al., 2001). The BZ integration is adopted with a k-grid density of 14 × 14 × 6 for bulk structures and 14 × 14 × 1 for few-layer structures at gamma point. The force tolerance is set to 10 −10 .
The Raman intensity of Raman activate mode (Supplementary Figures S1,S2) can be written as (Umari et al., 2001;Ceriotti et al., 2006): where e s and e i represent the electric polarization vectors of incident and scattered light, respectively.

Raman Active Modes of Bulk and Few-Layer GeSe
Bulk GeSe belongs to a layered structure in AB stacking with vdW interactions, as plotted in Figure 1A, with the space group of Pnma. We noted the longer axis in the x-y plane is the armchair direction (y) and the other axis is the zigzag (x) direction. Figures  1B,C shows the monolayer GeSe that retains the symmetry of the bulk structure. The relative displacement between the Ge atom and Se atom indicates the system is in the ferroelectric phase.
We first optimize the lattice parameters of the bulk structure to compare to experimental lattices, shown in Table1. As we can see, there is a large lattice difference of the b axis (b axis is the armchair direction) under various pseudopotentials and methods for bulk structure. DFT-D2 (vdW-D2) and DFT-D3 (vdW-D3) method of Grimme is applied as vdW correction (Grimme, 2006). The optimized structure with the vdW-D2 method is in good agreement with experimental results compared to others. As to LDA and GGA methods, the lattice parameters is underestimated by 4 and 2.5% of the b axis, respectively. Based on these calculated structures, we investigate the Raman frequencies, and we find vdW-D2 gives a credible result in comparison to other methods, consistent with previous work (Park et al., 2019). Hence the interlayer interactions will be considered in the following calculations of bulk and multilayer structures. We should note that the A g 2 mode obtained in theoretical calculations was not observed in experiments on GeSe  but appears in the similar system SnSe (Yang et al., 2018), which is due to the Raman tensor of A g 2 mode is too small for its Raman peak to be observed in GeSe.
According to the above comparison and analysis, we use the vdW-D2 method to relax our structure and calculate the optical phonon frequencies. To understand the Raman frequencies of bulk GeSe, we analyze the crystal structure by combining the irreducible representations of Γ points. The primitive unit cell of bulk GeSe include eight atoms, resulting in 24 vibrational modes as following: ( 1.2) and which are three acoustical (B 1u , B 2u , B 3u ), seven infrared active modes (3B 1u , B 2u , 3B 3u ), two silent modes (2A u ), as well as twelve Raman modes (4A g , 2B 1g , 4B 2g , 2B 3g ). Four active Raman modes for bulk GeSe as 3A g at 68.21 cm −1 , 179.05 cm −1 , and 189.6 cm −1 and B 3g at 153.02 cm −1 are listed in Table 2, and we define 3A g as A g 1 , A g 2 and A g 3 , respectively. To study the four active phonon vibration modes, we projected the eigenvectors of the dynamical matrix on Ge and Se atoms. Figure 2 gives the outline of the four vibration active modes of bulk GeSe. We find that A g 1 , A g 2 and A g 3 couple the in-plane and out-of-plane vibrations with various contribution. A g 1 and A g 2 modes are mainly dominated by out-of-plane motion, in which Ge and Se atoms (Ge1 and Se2) move to the opposite direction of A g 1 mode whereas in A g 2 mode they move to the same direction. A g 3 mode is most contributed by armchair direction, and only the Ge atom has a small out-of-plane contribution. For A g 3 mode, it  forms by an in-plane vibration along the armchair direction, indicating an intrinsic ferroelectric vibration mode in one layer of bulk GeSe, which agrees with previous experiment work . Monolayer GeSe has been predicted a member of ferrovalley materials with plenty of fantastic physical properties (Morales-Ferreiro et al., 2017;Shen et al., 2018;Liu et al., 2021). Thus, it is necessary to investigate the Raman spectrum dependence of fewlayer. For monolayer GeSe, the four active Raman peaks of A g 1 , B 3g , A g 2, and A g 3 modes are located at 81.84 cm −1 , 97.54 cm −1 , 141.94 cm −1, and 183.14 cm −1 illustrated in Table 3. Even though, both of paraelectric and ferroelectric phases exist in monolayer GeSe, we investigate the FE phase in the current work due to its lower energy state and stable structure. With the increase of thickness from monolayer to trilayer, the Raman shift changes a lot. Indeed, A g 1 mode generates a redshift from 81.84 cm −1 in the monolayer to 69.26 cm −1 in the bilayer. As to the remanent modes, they all increase. However, few shifts are induced between bilayer and trilayer. Our calculation suggests that GeSe has a large Raman shift with increasing layer thickness compared to other 2D materials (Lee et al., 2010;Tan et al., 2012;Yagmurcukardes et al., 2018;Kong et al., 2021). It provides a valuable approach for distinguishing the structure between monolayer and few-layer GeSe. The reason for such a large discrepancy of Raman peaks between monolayer and multilayer is attributed to the interlayer vdW interactions in multilayer GeSe, giving rise to structures difference. As a result, the lattice parameters that is listed in Table 3 along the armchair direction of monolayer GeSe are largely smaller than that of bilayer and trilayer.   The Influence of Strain on Raman Peaks Considering materials always are in strain in heterostructure or different temperatures, the influence of strain on Raman active modes is further explored in bulk and few-layer GeSe by changing the lattice parameter on both armchair and zigzag direction, as depicted in Figure 3. In monolayer GeSe, with the varying strain from −2 to 2%, the vibration mode of A g 1 linearly decreases in frequency. On the contrary, A g 2 and A g 3 modes exhibit a clear increase with strain. And for B 3g mode, an inconspicuous change has been observed. Raman peaks between bilayer and trilayer GeSe possess the same trend in tensile and compressive strains. It should be mentioned that Raman frequency of A g 1 mode is most steady and rarely changed compared with the other three modes. Whereas, B 3g and A g 3 modes increase by compressive strain and decrease along with tensile strain in Figures 3B,C. Moreover, as the strain increases,  the gap of Raman shift among A g 2 and A g 3 modes is shrunk. Our results suggest that the stain has a great effect on the Raman shift, which can be detected in the experiment, demonstrating the stress dependence of Raman spectra can also be used for the determination of crystal direction.

Interlayer Shift-Induced Raman Shift
Besides FE phase, AFE phase GeSe is also considered to be the lowest energy phase in an even number of layers. However, in some ways, the FE phase is more desired by researchers to utilize in multifunctional devices. To date, only the theoretical studies of AFE/FE transition for GeSe multilayers have been reported rather than any experimental works (Xu et al., 2021). Therefore, it is remarkable things to investigate the interlayer shift-induced AFE/FE transition and relevant Raman shift in few-layer GeSe. To realize the phase transition, we consider the AFE phase as a start structure. The structure evolution of the AFE/FE transition is illustrated in Figure 4. In bilayer GeSe, the bottom layer GeSe is fixed and the top layer artificially shifted along-y-direction in Figure 4A. When the top layer was shifted about 0.4b (b is the lattice parameter of y-direction), the system translates to an FE state. The same phenomenon appears in trilayer GeSe, in which we move interlayer between the top and bottom layer that is in fixation, showing in Figure 4B. Each shifted structure has been fully relaxed so that the Raman spectrum is obtained with a stable structure. Figures 5A,B represents the Raman shift in the intermediate phase. A g 1 mode of FE phase increases by of about 10 cm −1 compared with AFE phase in both bilayer and trilayer GeSe. Meanwhile. B 3g , A g 2 , and A g 3 mode show a redshift, in which the same variation rule in A g 2 mode is observed in both of bilayer and trilayer GeSe. In bilayer GeSe, A g 3 mode decreases with a nearly linear slope, yet it firstly increases and then diminishes with the interlayer shift of the trilayer. These large discrepancies in the Raman shift could be taken as an explicit evidence to identify the induced AFE/FE transition in the experiment.

CONCLUSION
To summarize, the Raman frequency of few-layer and bulk GeSe is systematically investigated using first-principles calculations. A large difference of Raman peaks between the monolayer and the multilayer, demonstrating the significant effect of vdW interlayer interaction. We then find the strain influence on these systems, in which A g 2 linearly increases and B 3g mode linearly diminishes from −2 to 2%. More important, the interlayer shift could also induce a phase transition from AFE to FE. As a result, the Raman frequency discrepancy of four vibration modes is natural. We truly expect our calculations can pave a way to verify the number of layers, strain influence, and polarity in few-layer GeSe.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
C-GD and NZ conceived the idea and supervised the work. Y-FZ carried out theory calculations. All authors analyzed the results and contributed to writing the manuscript.