First-principles investigation of the significant anisotropy and ultrahigh thermoelectric efficiency of a novel two-dimensional Ga2I2S2 at room temperature

Two-dimensional (2D) thermoelectric (TE) materials have been widely developed; however, some 2D materials exhibit isotropic phonon, electron transport properties, and poor TE performance, which limit their application scope. Thus, exploring excellent anisotropic and ultrahigh-performance TE materials are very warranted. Herein, we first investigate the phonon thermal and TE properties of a novel 2D-connectivity ternary compound named Ga2I2S2. This paper comprehensively studies the phonon dispersion, phonon anharmonicity, lattice thermal conductivity, electronic structure, carrier mobility, Seebeck coefficient, electrical conductivity, and the dimensionless figure of merit (ZT) versus carrier concentration for 2D Ga2I2S2. We conclude that the in-plane lattice thermal conductivities of Ga2I2S2 at room temperature (300 K) are found to be 1.55 W mK−1 in the X-axis direction (xx-direction) and 3.82 W mK−1 in the Y-axis direction (yy-direction), which means its anisotropy ratio reaches 1.46. Simultaneously, the TE performance of p-type and n-type doping 2D Ga2I2S2 also shows significant anisotropy, giving rise to the ZT peak values of p-type doping in xx- and yy-directions being 0.81 and 1.99, respectively, and those of n-type doping reach ultrahigh values of 7.12 and 2.89 at 300 K, which are obviously higher than the reported values for p-type and n-type doping ternary compound Sn2BiX (ZT∼ 1.70 and ∼2.45 at 300 K) (2020 Nano Energy 67 104283). This work demonstrates that 2D Ga2I2S2 has high anisotropic TE conversion efficiency and can also be used as a new potential room-temperature TE material.


Introduction
Realizing the power generation and cooling function of thermoelectric (TE) materials based on the transmission and interaction of electrons and phonons is a representative energy conversion and clean energy technology [1,2]. At present, the low TE conversion efficiency is the main bottleneck that restricts the development of TE materials. Meanwhile, some advanced fabrication technologies can boost numerous scientific and technological advancements in novel materials and devices for TE applications. Therefore, the discovery of highperformance TE materials is an urgent goal for application. The efficiency of a TE material can be evaluated by the dimensionless figure of merit (ZT = σS 2 T/ (κ e + κ L )), where σ represents the electrical conductivity, S stands for the Seebeck coefficient, and T is the Kelvin temperature. κ e and κ L are the electronic and lattice thermal conductivities, respectively. σS 2 is the power factor (PF) [3]. To obtain high-performance TE materials, previous studies have been devoted to investigating two aspects: enhancing the PF and reducing the total thermal conductivity. As early as 1993, Dresselhaus et al [4,5] predicted that in comparison with bulk materials, the enhancement of the electronic density of states (DOS) caused by the quantum confinement effect in 2D materials would greatly increase the PF of the material, providing an important theoretical guidance. In 2018, Miao et al [6] confirmed the previous prediction of the 2D quantum confinement effect. In their study, the TE transport properties of the 2D layered material γ-InSe were first studied in an experiment combined with first-principles calculations. They revealed that, in 2D materials, the quantum confinement effect will cause a sharper DOS at the conduction band edge, thereby enhancing the PF. At the same time, due to the limited effect size, the layered materials also bring the potential to control the performance of phonon thermal and electrical transport. Therefore, exploring new layered 2D materials with special structures and electronic properties are the key to screen potential candidates for TE applications.
At present, high-throughput calculation methods can efficiently predict new anisotropic TE materials with excellent TE performances, which has application potential in commercialized Peltier cooler and transverse thermal energy collection [7,8]. Sarikurt et al [9] used a density functional theory (DFT) driven solution of the Boltzmann transport equation (BTE) and the electronic fitness function [10] as an index for screening high-performance 2D layered TE materials with high anisotropy. Those related materials are based on the extended orthogonal (space group: Pmmn) layered metal oxyhalide (FeOCl) type structure [11], including Ir 2 Cl 2 O 2 , Al 2 I 2 S 2 [12], and Ga 2 I 2 S 2 , whose layers are combined with van der Waals forces. FeOCl-type compounds usually have relatively weak van der Waals forces, which make them easy to peel off of bulk materials into 2D materials and have the potential for property regulation. Wang et al [13] studied TE properties of FeOCl-type ternary TiNX (X = F Cl and Br). They found that the phonon thermal and electrical transport properties showed significant anisotropy, that is, the transport properties in the yy-direction are usually better than those in the xx-direction. They also found that the maximum ZT values of the three materials along yy-direction at 500 K can reach 1.0, 0.89, and 1.17, respectively. More interestingly, a new family of Janus 2D materials depended on the monolayer MX 2 (M = V, Nb, etc.; or X = S, Se, etc) with halogen (F, Cl, and I) substitution using systematic first-principles calculations were proposed [14]. Meanwhile, Janus 2D materials have been expanded to implement shape transformation under diverse stimuli [15,16]. Therefore, the low-cost and flexible method from a previous work [17] may capture such Ga 2 I 2 S 2 type TE devices.
In summary, the new 2D layered ternary FeOCl-type material may be a potential candidate material in high-performance TE devices in the future. However, explorations of significant anisotropic phonon thermal and TE properties of such materials are still insufficient. To further reveal the intrinsic mechanism of excellent TE materials and accelerate their application in commercial TE fields, we performed first-principles calculations combined with the BTE and relaxation time approximation (RTA) to systematically calculate the phonon thermal and TE properties of the layered FeOCl-type named 2D-connectivity [18] Ga 2 I 2 S 2 . Furthermore, the potential of Ga 2 I 2 S 2 as a high-performance anisotropic TE material is further evaluated for low-to-moderate temperature TE technology module.

First-principles calculations
All first-principles calculations in this paper are employed using the Vienna Ab-initio Simulation Package code on the basis of DFT [19]. Electron-ion interactions are described by the plane wave pseudopotential method called Projected Augmented Wave [20,21], and the Perdew-Burke-Ernzerhof (PBE) is selected as the exchange-correlation functional [22]. To ensure that the electrons and ions converge in the structural optimization process, the plane wave cutoff energy of the 2D Ga 2 I 2 S 2 is set to 500 eV, and a Monkhorst-Pack kpoint grid of 15 × 15 × 1 is used to ensure accurate calculation. When optimizing the crystal structure, the residual stress and the maximum force on each atom are less than 0.01 kbars and 10 −6 eV Å −1 , respectively. Meanwhile, the convergence criterion of the electronic self-consistent calculation is that the energy difference between two adjacent iterations is less than 10 −8 eV. The calculation of the second-order interatomic force constants (IFCs) of the harmonic properties of phonons uses the finite displacement method [23], where a 5 × 8 × 1 supercell is used. The third-order IFCs are calculated to depict anharmonic phonon scattering using the same method and supercell size as the second-order IFCs. Among them, the cut-off of the nearest neighbor with the largest thirdorder IFCs is selected as the 14th nearest neighbor, whose corresponding truncation radius is 0.71 nm. To avoid the influence of the interaction between the layers, a vacuum layer of 20 Å is constructed in the out-of-plane direction for all the calculations.

Phonon thermal transport
The BTE is a key method to study the transport process of phonons in solids [24], which is the bridge that connects the DFT with phonon thermal transport problems. The κ L can be calculated by solving the phonon BTE [25]. Based on the obtained second-order and third-order IFCs, the ShengBTE [26] package is used to solve the phonon BTE, and κ L at different temperatures, the three-phonon scattering phase space, and phonon lifetime are acquired. Through the iterative solution, the lattice thermal conductivity κ αβ L can be defined as: where α and β represent Cartesian components in x, y, and z directions. k B stands for the Boltzmann constant, and T represents the absolute temperature. Ω represents the volume of the unit cell, and N is the number of Gamma-centered q-point grids. f 0 depicts the Bose-Einstein distribution, and ℏ stands for the reduced Planck constant. ω λ and υ λ represent the angular frequency and phonon group velocity of the phonon mode λ, respectively. F β λ represents the component of the phonon mean free path (MFP) F λ along the β direction.

Thermoelectric transport
The electrical transport properties of semiconductors are usually governed by their band structures and carrier scattering processes. This paper uses the Boltz Trap software to compute the electrical properties [27]. According to the electron BTE and constant RTA [28], the σ and S can be expressed as: where f 0 depicts the equilibrium Fermi distribution function, µ represents the chemical potential, and υ k and τ k stand for the group velocity and relaxation time in the k state. At this time, the carrier mobility (µ) can be defined as: where ε nk and f nk respectively represent the electron energy of the n band wave vector and the Fermi-Dirac distribution of the electron in the equilibrium state, and υ nk and τ nk represent the electron group velocity and carrier relaxation time, respectively. The effect of phonon scattering is mainly considered in 2D materials, and the carrier mobility (µ 2D ) can be expressed as [29,30]: where α represents charge transfer directions (xx or yy), and m * α represents the effective mass along the α direction, m * d = √ m x m y . E 1α stands for the deformation potential (DP) along the α direction and is usually calculated through the highest occupied state in the valence band (VBM) and the lowest occupied state in the conduction band (CBM). E 1α = ∆E VBM/CBM / (∆α/α 0 ), where α 0 represents the lattice constant in the xx or yy direction in the fully relaxed structure, and ∆α represents the deformation caused by stress. C 2D (computed by(E − E 0 ) /S 0 = C 2D (∆α/α 0 ) 2 /2) stands for the elastic modulus in the transport direction, which can depict the second-order partial derivative of the system energy (E) concerning the deformation, and S 0 represents the area of the studied crystal structure. In addition, the carrier relaxation time can be calculated by τ = µm * /e. Through the Wiedemann-Franz law [31], the total thermal conductivity contributed by electrons is proportional to the electrical conductivity, that is, κ e = LσT, where L represents the Lorentz constant, which is taken as L = 2.45 × 10 −8 WΩK −2 [32,33].

The crystal structure and phonon properties
The 2D FeOCl-type Ga 2 I 2 S 2 has six atoms in the primitive cell, including 2 Ga atoms, 2 I atoms, and 2 S atoms. Each Ga atom has two I atoms and four S atoms as the nearest neighbors. Under PBE functional optimization with no pressure, the lattice constants of Ga 2 I 2 S 2 along the xx-and yydirections are 5.027 and 3.752 Å, respectively, which agree with the stable structure of high-throughput screening of 5.03 and 3.75 Å [9]. The bond length between a Ga atom and an adjacent I atom (L Ga−I ) and the bond length between a Ga and an adjacent S atom (L Ga−S ) are 2.8 and 2.5 Å, respectively. Figures 1(a) and (b) show the top and side views of the 2D Ga 2 I 2 S 2 . The path diagram of the high symmetry point in the BZ is plotted in figure 1(c). To verify the 2D Ga 2 I 2 S 2 structure's stability, we first compute the phonon dispersion (figure 1(d)) along the Γ-X-S-Y-Γ path. From the phonon dispersion, it can be seen that there are 18 phonon branches, including 3 acoustic phonon branches (transverse acoustic (TA), longitudinal acoustic (LA), out-of-plane acoustic flexural modes (ZA)) and 15 optical phonon branches (Optical). For phonon dispersion, there is a phonon softening mode and no band gap of phonon frequency. Besides, a crossover phenomenon between the high-frequency acoustic branches and the low-frequency optical branches can be observed. This phenomenon usually indicates that the acoustic branch and optical branch phonons have strong internal scattering [36,37], resulting in lower phonon lifetimes. For the phonon branch around the Γ point, the dispersion relationship between TA and LA shows a linear trend in the lowfrequency range, while the dispersion relationship of the ZA branch is a quadratic parabola. This is due to the rapid weakening of the lateral restoring force [38,39]. These observed changes are common phonon characteristics of 2D materials. Furthermore, the phonon dispersion relationship shows significant anisotropy in the xx-(Γ-X) and yy-directions (Γ-Y), giving rise to high anisotropy of phonon group velocities in two different directions. The slope of high-frequency optical phonons is smaller along the xx-direction than in the yy-direction. The relatively localized phonon modes show that the coupling for 2D Ga 2 I 2 S 2 is weaker along the xx-direction. The strong phonon anharmonicity caused by chemical bonds in materials often leads to an ideal crystal with a lower κ L . The Grüneisen constant [40] quantitatively describes the strength of phonon anharmonicity, and the square of the Grüneisen constant is approximately proportional to the three-phonon scattering rate [41]. To further describe the phonon anharmonicity of 2D Ga 2 I 2 S 2 in the xx and yy directions, we calculate the dispersion relationship of the Grüneisen constant. As shown in figure 1(e), the Grüneisen constant in the xx-direction is larger than that in the yy-direction, and the corresponding average Grüneisen constants are 2.32 and 1.07, respectively, which can be used to explain why the κ L of 2D Ga 2 I 2 S 2 in the xxdirection is smaller than that in the yy-direction. Since the phonon branches crossover phenomenon can easily occur in figure 1(d), it is not easy to distinguish their respective contributions to κ L in the entire BZ. Figures 1(f) and (g) show the 3D phonon dispersion relationship in the entire BZ at room temperature and the contribution (color map) of different phonons along the xx-and yy-directions to the total κ L . The significant anisotropic of κ L exists in figures 1(f) and (g), and the contributions of the LA, TA, and ZA modes to κ L are 26.2%, 14.8%, and 19.4% along the xx-direction, and 37%, 20%, and 9.5% along the yy-direction.
The participation of the ZA mode to κ L in both directions are far lower for 2D Ga 2 I 2 S 2 than that of graphene (75%) in the same acoustic mode at room temperature [42], and its values are equivalent to anisotropic phosphorene (16% and 8%) [43]. Due to the reflection symmetry perpendicular to the graphene, scattering of the ZA mode can hardly occur. This phenomenon is called the crystal symmetry-based phonon scattering selectivity rule [44][45][46]. Compared with graphene, the phonon scattering selectivity rule of phosphorene and 2D Ga 2 I 2 S 2 is broken by their own structures, leading to a larger scattering rate of the ZA mode, which will reduce its contribution to the κ L .

Phonon thermal transport properties
In this section, we calculate the temperature dependence of the κ L of 2D Ga 2 I 2 S 2 , as plotted in figure 2(a). Due to the limitations of ShengBTE's solution to the phonon BTE under low temperatures, the results of solving the equation will diverge. Hence, in this work, we choose 300-800 K as the temperature range. Meanwhile, the κ L of 2D Ga 2 I 2 S 2 , with respect to the tensors of the IFC and q-mesh, has been tested and illustrated in figure S1 (available online at stacks.iop.org/IJEM/4/025001/mmedia) (supporting information). To compare the contribution of the N-Process and U-Process to phonon scattering, we use RTA and the iterative solution (ITE) to obtain κ L . The detailed calculations of N-and U-Processes at 300 K are plotted in figure S2 (supporting information). From figure 2(a), it appears that κ L gradually decreases as the temperature increases, mainly due to scattering between the acoustic branch and optical branch phonons. The phonon lifetime decreases as the temperature rises, and κ ITE L > κ RTA L , especially along the yy-direction, indicating that the N-Process is dominant in the phonon scattering of 2D Ga 2 I 2 S 2 . We note that the κ xx L calculated by the RTA in figure 2(a) is nearly the same as those of the ITE method at 300-800 K, while κ yy L largely deviates from the ITE results. This is mainly because the group velocities of all phonon modes along the two directions show large anisotropy, and the velocities along the yy-direction contribute more heavily to κ yy L than those along the xx-direction in 2D Ga 2 I 2 S 2 system. At 300 K, the κ L of Ga 2 I 2 S 2 along the xx-and yydirections are 1.55 and 3.82 W mK −1 , respectively. Herein, dimensionless parameter η (η = (κ yy L − κ xx L ) /κ xx L ) is used to evaluate the anisotropy of κ L . The calculated η values under the two methods are shown in figure 2(a). The η values calculated from the two methods converge to 1.44 and 1.18, respectively. These strongly anisotropic κ L characteristics can be compared with phosphorene (2.46) [43] and borophene (1.30) [39], which further confirms the strong directional dependence of κ L .
The phonon group velocity and lifetime describe the size of the phonon MFP for an ideal crystal material. To deeply explore the influence of phonons with different MFP on the anisotropy of κ L , figure 2(b) shows the κ L of 2D Ga 2 I 2 S 2 along two different directions at 300-500 K. Traditionally, with the variation of the phonons' MFP, the cumulative κ L is the summation of the contributions of phonons with MFPs smaller than the specific value. Overall, the MFP of phonons in 2D Ga 2 I 2 S 2 is relatively small, and phonon MFP in two directions shows significant anisotropy. The maximum MFP of phonons at 300 K is 586 and 335 nm in the xx-and  yy-directions, respectively, which are close to the values of excellent TE materials Bi 2 Te 3 [47] and PbTe [48]. As temperature increases, the maximum phonon MFP in different directions gradually decreases, and the shorter phonon MFP further restricts the size effect to regulate thermal transport. We calculate the cumulative κ L of 2D Ga 2 I 2 S 2 along two directions at 300-500 K versus the phonon frequency and plot them in figure 2(c). The orange shaded area in figure 2(c) represents the contribution of all acoustic and low-frequency optical branches to the κ L . This area contributes 75.1% and 84.5% to the κ L along the xx-and yy-directions at 300 K, which further confirms the dominance of acoustic phonon branches to phonon thermal transport. We also calculated the frequencydependent κ L of 2D Ga 2 I 2 S 2 at 300 K, which is presented in figure S3 (supporting information). Due to the limitation of technology in material synthesis, 2D materials are usually synthesized on a micro-nano scale. Therefore, the influence of phonon boundary scattering on the κ L of a material must be considered. When the system size of 2D Ga 2 I 2 S 2 is in the scope of 0.1-100 µm, the κ L varies with temperature, which is illustrated in figure 2(d). When the system size is close to ∞, it refers to the default value in the ShengBTE [26] software (infinite system). Compared with other sizes at 300-800 K, when the system size is 0.1 µm, the κ L of 2D Ga 2 I 2 S 2 along both directions has significantly reduces to 30% and 28%, respectively. The boundary regulation provides some guidance for the actual size design of 2D Ga 2 I 2 S 2 TE devices.
To reveal the intrinsic physical mechanism, which clarifies the highly anisotropic κ L of 2D Ga 2 I 2 S 2 , the phonon group velocity, scattering mechanism, and anharmonicity of the entire BZ at room temperature are further calculated, as illustrated in figure 3. The group velocity of 2D Ga 2 I 2 S 2 along the xx-and yy-directions shows significant anisotropy, as plotted in figure 3(a). The group velocity along the xx-direction is greater than that along the yy-direction in the low-frequency region, and the squared average phonon group velocity in the two directions are 95.59 (Å/ps) 2 and 129.67 (Å/ps) 2 , respectively, which confirms that the significant anisotropy of the group velocity is the main factor for high anisotropy of the κ L . Meanwhile, we show the phonon group velocity projected phonon dispersion curves of 2D Ga 2 I 2 S 2 , as illustrated in figure S4 (supporting information). The phonon lifetime is another important parameter to determine the κ L of 2D Ga 2 I 2 S 2 . It is largely governed by the phonon scattering mechanism and anharmonicity. The phonon scattering rate includes the anharmonic, isotope, and boundary scattering rate. Among them, boundary scattering plays an important role in analyzing the low temperature dependence and the size effect of the thermal conductivities of the materials. Here, we only study and compare the first two scattering processes, as depicted in figure 3(b). It appears that the three-phonon scattering rate is far higher than the isotope scattering rate, especially for the acoustic branches, whose values are almost four orders of magnitude larger than the isotope scattering rate. This indicates that the three-phonon scattering process may play a key role in determining a phonon's lifetime. As is illustrated in figure 3(c), the lifetimes of the ZA modes are equivalent to those of the TA modes (∼100 ps), and the lifetimes of the optical phonons are obviously shorter than those of the three acoustic phonon branches. To further identify and understand the phonon transport mechanism, we calculate the phase space of the total three-phonon scattering process (P 3 ) [49,50]: where Θ represents the normalization factor, D ± λ (q) represents two types of three-phonon scattering channels: (1) D + λ (q) corresponds to the absorption process P + The threephonon scattering phase space P 3 contains two processes, P + 3 and P − 3 (see figure S5 in the supporting information), which are used to evaluate the number of phonon scattering channels [50]. Our calculation results for P 3 of 2D Ga 2 I 2 S 2 are plotted in figure 3(d). It clear that the P 3 of the ZA branch in the low-frequency scope is far larger than the P 3 of the TA, LA, and Optical branches. Therefore, ZA + ZA/LA/TA ↔ LA/TA/Optical are the most significant scattering channels for the three-phonon scattering process in 2D Ga 2 I 2 S 2 .

Electronic band structure and carrier mobility
The electrical transport properties mainly depend on results from the corresponding band structure. Figure 4(a) demonstrates the electronic structure of 2D Ga 2 I 2 S 2 under the PBE and HSE functionals. The calculation results show that the band gap of 2D Ga 2 I 2 S 2 produced by the PBE functional is 0.95 eV, while the band gap produced by the HSE functional is 1.67 eV. Therefore, the HSE functional is considered to revise the band structure, especially for the width of the band gap correction [51][52][53]. Thus, the following TE-related parameters depend on the HSE calculation results. The 2D Ga 2 I 2 S 2 is a semiconductor with a direct band gap. Figures 4(b)-(d) show the projected band structure in different atomic orbits. The size of the solid sphere depicts the weight of the atom projected to orbit. It is clear that the I atom and Ga atom represent the main contribution to CBM and VBM, respectively, indicating that I atoms have a greater contribution to the crystal electrical conductivity compared to the other two element atoms. Simultaneously, the doping concentration of the S atom site can adjust the carrier concentration without tuning the band structure, which will help regulation of the TE performance of 2D Ga 2 I 2 S 2 . Since the Seebeck coefficient is usually governed by Table 1. Calculated electronic relaxation time τ , effective mass m * , DP constant E l , elastic constant C 2D , and carrier mobility µ 2D for electrons and holes in the xx-and yy-directions for 2D Ga 2 I 2 S 2 at room temperature.

Direction
Carrier type  [54], we find that the band of VBM is flatter than that of CBM, indicating that the DOS of VBM is steeper than that of CBM. The S of the p-type doping 2D Ga 2 I 2 S 2 is expected to be greater than that of the n-type doping system. At the same time, the definition of the effective mass of electrons indicates that a wider band structure will lead to a smaller effective mass. To further reveal the effect of different atomic orbitals near the Fermi level on electric transport, we calculate the DOS of the 2D Ga 2 I 2 S 2 and the crystal orbital Hamilton population (COHP). As is shown in figures 4(e) and (f), COHP can quantitatively describe the strength of one or several chemical bonds in a crystal structure, which can be calculated by the LOBSTER program [55]. Meanwhile, COHP is related to energy, where '-COHP' represents the contribution of the chemical bond to the bonding part. Figure 4(c) shows that VBM is mainly derived from the hybridization of the p orbitals of I and S atoms, and CBM is mainly derived from the hybridization of the s orbitals of Ga atoms and the p orbitals of S and I atoms. COHP partitions the band structure into orbital-pair interactions. Generally speaking, it is a bond-weighted DOS between a pair of adjacent atoms. The Ga-S bonding energy in the crystal structure is greater than that of Ga-I. The bond lengths of Ga-S and Ga-I in the illustration are 2.5 and 2.8 Å, respectively, indicating that the interaction force between Ga and S atoms is stronger than that between Ga and I atoms. The carrier mobility of 2D Ga 2 I 2 S 2 is calculated according to the commonly used DP theory [56]. Using the DP theory requires some parameters, including elastic modulus C 2D , DP constant E 1α , and effective mass m * α ; Detailed effective mass calculations are presented in figure S6 in supporting information. These parameters can all be calculated using first-principles calculations. Table 1 shows some calculated parameters required by the DP theory at room temperature. Compared with common TE materials, the electron relaxation time of 2D Ga 2 I 2 S 2 is longer, which improves the material's TE conversion efficiency. A smaller E l (1.85 eV) usually appears in the xx-direction, and the electron-phonon coupling strength is weak.

Electric transport and TE properties
In the section, we first discuss the relationship between the S, σ, and κ e of 2D Ga 2 I 2 S 2 along different directions with carrier concentration. Then, we use the κ L and relaxation time obtained above to calculate the ZT values corresponding to different carrier concentrations of p-type and n-type doping at 300-500 K. Under the rigid energy band approximation, the chemical potential represents the doping type and carrier concentration of the studied system. Here we only consider the chemical potential near the Fermi surface. In this case, the doping concentration is generally low, and the carrier concentration is in the scope of 10 8 -10 13 cm −2 . Thus, the rigid band model is applicable. Figure 5 shows the changes in the S, σ, σS 2 and κ e of p-and n-type doped 2D Ga 2 I 2 S 2 along different directions with various carrier concentrations. From figures 5(a) and (b), both doping types have large Seebeck coefficients, and the maximum value can reach 1.2 mV K −1 at 300 K, which is comparable to traditional, excellent TE materials. At the same time, although the effective mass calculated in table 1 shows anisotropy along different directions, the S of the doped p-type and n-type systems at different temperatures are very similar in the two directions. However, the S for p-type and n-type doping systems is approximately isotropic, which is a common characteristic of semiconducting materials. Compared with the results shown in figures 5(c) and (d), when the S peaks, the corresponding σ is very small. Obviously, there must be a trade-off between S and σ, which drives the PF of a system with a specific carrier concentration to its maximum value. Furthermore, the σ of 2D Ga 2 I 2 S 2 under the two doping types shows strong anisotropy. Using the results above, we calculate the PF value of 2D Ga 2 I 2 S 2 for the two doping types. The change in PF along different directions with the carrier concentration is shown in figures 5(e) and (f).
The PF of 2D Ga 2 I 2 S 2 under p-type and n-type doping exhibits strong anisotropy, where the orange region represents the ratio of PF anisotropy of the two doping types along different directions. The inset in figure 5(e) shows that the PF under p-type doping suddenly changes at a carrier concentration of 10 12 cm −2 . A similar phenomenon also exists for n-type doping. In addition, it can be easily observed that the PF under the two doping types both reach the peak value near the carrier concentration of 10 12 cm −2 . The PF under the doped n-type system is larger than that of the doped p-type, which is mainly owing to the higher the electrical conductivity of n-type doping comparing p-type doping. Taking into account the contribution of the κ e , we calculate the changes of the κ e along different directions with the carrier concentration under the two doping types. Figures 5(g) and (h) show the variation of κ e , which is in agreement with the variation of σ. Furthermore, as the temperature increases, the total thermal conductivity increases.
Based on the theoretical data above, we calculate the ZT value along different orientations with the carrier concentration under the two doping types, as illustrated in directions under the n-type doping are 7.12 and 2.89, respectively. We find that the ZT value of doped n-type 2D Ga 2 I 2 S 2 along different directions is greater than that of the p-type Calculated ZT of (a) p-type and (b) n-type doped 2D Ga 2 I 2 S 2 versus carrier concentration in the xx-and yy-directions at 300-500 K. Contour map of ZT versus both temperature and carrier concentration for 2D Ga 2 I 2 S 2 : (c) p-type in the xx-direction, (d) p-type in the yy-direction, (e) n-type in the xx-direction, and (f) n-type in the yy-direction. The four-color bars are independent, and the white dotted curves correspond to ZT values.
doping. To better describe the relationship between temperature and carrier concentration parameters and ZT, we calculate the contour map of the ZT value of p (n)-type doped 2D Ga 2 I 2 S 2 versus temperature and carrier concentration, as shown in figures 6(c)-(f), where the white dotted curves are the corresponding ZT value contour curves. It can be seen from figures 6(e) and (f) that the maximum ZT values along different directions at 500 K under n-type doping are 11.5 and 5.42, respectively, and the corresponding carrier concentrations are 2.37 × 10 10 cm −2 and 1.72 × 10 11 cm −2 . For doped p-type system, the maximum ZT values in the two directions are 1.55 and 3.69, respectively, and the corresponding carrier concentrations are 1.51 × 10 12 cm −2 and 8.94 × 10 11 cm −2 . The carrier concentration corresponding to the maximum ZT value shows anisotropy. In general, our theoretical prediction results suggest that 2D Ga 2 I 2 S 2 can be customized by different doping and carrier concentrations, providing some guidance for future experimental work.

Conclusions
In summary, we utilized advanced quantum computation techniques to predict and study phonon thermal/TE transport properties of FeOCl-type 2D-connectivity Ga 2 I 2 S 2 . We found that the in-plane κ L along the xx-and yy-directions at 300 K were 1.55 and 3.82 W mK −1 , respectively. The anisotropy ratio of the κ L was 1.46, which shows strong anisotropy. Meanwhile, the VBM values of 2D Ga 2 I 2 S 2 along the two directions were close to each other, which enhanced the S and ultrahigh PF values. Under the action of the above characteristic, the peak ZT values of p-type doped 2D Ga 2 I 2 S 2 along the xx-and yy-directions at room temperature were 0.81 and 1.99, respectively. Meanwhile, the ultrahigh ZT values along the two directions for n-type doping were 7.12 and 2.89, respectively. Obviously, both directions of the doping types exhibited excellent TE conversion efficiency. Our findings provide some guidance and motivation for future experimental manufacturing efforts on ultrahigh-performance TE materials subject to lowand moderate-temperature ranges, especially for the recently developed reliable manufacturing method of producing such Ga 2 I 2 S 2 type materials from a single crystal material [17]. This work also paves the way toward manufacturing novel materials and devices and governing thermal transport properties by phonons, which expand in TE energy conversion and thermal management fields.

Data availability statement
The data that support the findings of this study are available from the corresponding author upon request.