Introduction

Renewable energy is a very important field due to the insufficiency of natural energy source and global warming1. One of the best renewable energy sources is waste heat, which can be converted into electricity via the Seebeck effect2, 3. The performance of thermoelectric materials is measured by a dimensionless quantity ZT called the figure of merit4, 5:

$$ZT=\frac{\sigma {S}^{2}}{\kappa }T$$
(1)

where σ, S, κ, and T are electrical conductivity, Seebeck coefficient, thermal conductivity, and temperature, respectively. The lattice thermal conductivity κ l and the electronic thermal conductivity κ e are included in the thermal conductivity κ in eq. 1. A large Seebeck coefficient, large electrical conductivity, and low thermal conductivity are needed for high thermoelectric performance, but a low amount of charge carrier is required to improve the Seebeck coefficient, which reduces the electrical conductivity6.

Improving ZT has been a big challenge, and different approaches have been used. Semiconductors composed of heavy elements such as Zn4Sb3, PbTe and BiSb have been used to reduce the thermal conductivity7,8,9. Point defects (R1−y Fe4−x Co x Sb12 and Ce y Fe x Co4−x Sb12) have been produced to decrease the lattice thermal conductivity and the optimized electrical conductivity10, 11. Some bulk complex materials also show very good thermoelectric performance such as filled skutterudites (La0.9 Fe3 CoSb12), half-Heusler alloys (ZrCoSn x Sb1−x ), and clathrates (Sr8 Ga16 Ge30) because of their low thermal conductivity and high periodicity in the crystal structure12,13,14,15,16. Zhao et al. recently reported that bulk SnSe is a very good thermoelectric material with a ZT of 2.6 at 973 K17. It was theoretically predicted that bulk SnS, GeSe, and GeS would also show very good thermoelectric performance18.

One of the efficient methods to increase ZT is reducing the dimensionality of the material, which increases the Seebeck coefficient due to the increased density of states near the Fermi level19,20,21. It is reported that the reduction in dimensionality enhances the energy storage and conversion22, 23, the ZT of bulk Bi2 Te3 is improved 13 times by converting into the quantum well. Fei et al. and Cheng et al. reported that a bismuth monolayer and phosphorene showed very promising thermoelectric properties24, 25.

We studied two-dimensional SnSe, SnS, GeSe, and GeS materials for thermoelectric applications. Monolayers of these materials have already been experimentally synthesized, and they have band gaps less than 2 eV26,27,28,29. They have been recently reported to have low lattice thermal conductivity as well30, which is a requirement for efficient thermoelectric materials. Group IV–VI compounds in bulk form have very good thermoelectric efficiency and a simple orthorhombic SnSe crystal was reported to have outstanding thermoelectricity17, 31, 32. It was recently discovered that even a monolayer of SnSe shows optimal thermoelectric properties33, which motivateed us to study the thermoelectric properties of monolayer IV–VI compounds SnSe, SnS, GeSe, and GeS.

Results and Discussions

Bulk SnSe, SnS, GeSe, and GeS have an orthorhombic crystal structure with the space group Pnma(62), while their monolayer analogs have the space group Pmn21(31) (see Fig. 1). The structures are optimized with a large vacuum space of 15 Åin the z-direction until the forces on each atom become zero. The optimized lattice parameters are given in Table 1, and they are in good agreement with previous reports34, 35.

Figure 1
figure 1

Crystal structure of two-dimensional monochalcogenides SnSe, SnS, GeSe, and GeS. (a) Side view perpendicular to the zigzag direction, (b) top view, (c) side view perpendicular to the armchair direction.

Full size image
Table 1 Calculated lattice parmeters and band gaps of SnSe, SnS, GeSe, and GeS. The values in the parentheses are from the refs 18 and 34.
Full size table

Electronic structures are very important for understanding the thermoelectric behavior of materials. The band gaps of SnSe, SnS, GeSe, and GeS are calculated using the exchange-correlation functional within a generalized gradient approximation (GGA), as shown in Table 1. The GGA functional calculations show indirect band gaps for SnSe, SnS, and GeS and a direct band gap for GeSe, as shown in Fig. 2. Density of states (DOS) of the SnS and GeS monolayers has sharp peaks near conduction band minima and valence band maxima as shown in Fig. 2(b,d), which may enhance the Seebeck coefficient. All these monolayers have band gaps less than 2 eV, which suggests that they can be used as thermoelectric materials. It is very difficult to get the optimal value of the ZT for wide band gap materials because heavy doping is required.

Figure 2
figure 2

Band structures along the high-symmetry k-points Γ, X, S, and Y and density of states of (a) SnSe, (b) SnS, (c) GeSe, and (d) GeS.

Full size image

The carrier mobility (μ) of the group IV–VI monolayers is calculated in order to get relaxation time (τ). Our method to calculate the mobility is based on deformation potential theory used extensively to calculate carrier mobility andrelaxation time of two-dimensional materials36,37,38,39, the expression to calculate the mobility is given by36, 39, 40:

$$\mu =\frac{e{\hslash }^{3}{C}^{2D}}{{k}_{B}T{m}^{\ast }{m}_{d}{E}_{1}^{2}}$$
(2)

where C 2D is the two-dimensional elastic constant determined by fitting the energy-strain curve to quadratic polynomial (see Fig. S1), our calculated values for C 2D are consistent with previous reported values41. T represents the temperature, m * is the effective mass in the transport direction, and m d is calculated as m d  = \(\sqrt{{m}_{x}{m}_{y}}\). Here m x and m y are the effective masses along armchair and zigzag directions. E 1 is the deformation potential constant defined by E 1 = \(\partial {E}_{{edge}}/\partial \delta \), where E edge is the conduction band minima (CBM) and δ is uniaxial strain. The shift in CBM by applying uniaxial strain is shown in Supplementary Fig. S2. The relaxation time is evaluated from mobility using the following relation:

$$\tau =\frac{{m}^{\ast }\mu }{e}$$
(3)

The calculated E 1, C 2D , m * and temperature-dependent μ and τ are tabulated in Table 2. We found an anisotropic behaviour in mobility and relaxation time for these monolayer. The carrier mobility of SnSe is in a good agreement with ref. 38. The GeS has the highest carrier mobility and relaxation time along the armchair direction among these four monolayer due to the low deformation constant and the low effective mass.

Table 2 Deformation potential constant (E 1), two dimensional elastic constant (C 2D), effective mass (m *), carrier mobility (μ) and relaxation time (τ) at 300 K, 500 K, and 700 K, in the armchair and zigzag directions of the group IV–VI compounds.
Full size table

The thermoelectric properties of SnSe, SnS, GeSe, and GeS are calculated using the Boltzmann transport equation for electrons under a constant scattering time. Boltztrap code calculates-relaxation-time dependent electrical conductivity (σ/τ) and electronic thermal conductivity (κ e /τ). Since there is no experimental data available for the electrical conductivity to calculate the exact value of the relaxation time of these monolayers, deformation potential theory is used to predict the temperature-dependent relaxation time for each material as compiled in Table 2. The electrical (σ) and electronic thermal (κ e ) conductivities are plotted as a function of carrier concentration (n) in Fig. 3(a–d) and (e–g) respectively, at 300 K, 500 K, and 700 K along armchair and zigzag directions. The carrier concentration shows the doping (positive values for p-type doping and negative values for n-type doping). The electrical and thermal conductivities increase by increasing carrier concentration. When the Fermi level occurs in the middle band gap region, the conductivities are increased slowly with respect to the carrier concentration and when it moves down into the valence band (for p-type) or moves up into conduction band (for n-type), the conductivities are increased quickly. GeSe has the highest electrical conductivity of 69.85 × 106 S/m at n = −8.9 × 1014/cm2 in the n-type doping among these compounds. The band gap is higher for the monolayer than bulk (see Table 1). The monolayers have lower electrical conductivity than the bulk due to the increase in the band gaps.

Figure 3
figure 3

Electrical (σ) and electronic thermal (κ e ) conductivities for the (a,e) SnSe, (b,f) SnS, (c,g) GeSe, and (d,h) GeS along the armchair (solid lines) and zigzag (dashed lines) directions at 300 K, 500 K, and 700 K.

Full size image

The Seebeck coefficients are calculated as a function of carrier concentration at different temperatures along the armchair and zigzag direction, as shown in Fig. 4. As the temperature decreases, the Seebeck coefficient also increases because of bipolar conduction42. The Seebeck coefficients of the two-dimensional monochalcogenides are two times greater than those of the bulk material as shown in the Table 3. This results from, the increase in the band gaps and the density of states near the Fermi level. The GeS has the largest Seebeck coefficient of 2810 μ VK−1 at 300 K because of the large band gap and the flatness in the band structure. The Seebeck coefficient (S) is calculated with the expression,

$$S=\frac{{\int }_{-\infty }^{\infty }dE\,g(E)(E-\mu )(-\frac{\partial f(E,\mu ,T)}{\partial E})}{T{\int }_{-\infty }^{\infty }dE\,g(E)(\frac{\partial f(E,\mu ,T)}{\partial E})},$$
(4)

where E, g(E), f(E, μ, T), μ, and T are energy, transport function, Fermi function, chemical potential, and temperature, respectively43. The transport function is:

$$g(E)=N(E){v}^{2}(E)\tau (E),$$
(5)

where N(E) is the density of states, v(E) is the Fermi velocity and τ(E) is the scattering time43. The Seebeck coefficient changes dramatically near the Fermi level because of the term \(\frac{\partial f}{\partial E}\) in Eq. 4, which behaves like a Dirac delta function.

Figure 4
figure 4

Calculated Seebeck coefficients (S) as a function of the carrier concentration (n) along the armchair (solid lines) and zigzag (dashed lines) directions at 300 K, 500 K, and 700 K for the group IV–VI monolayers (a) SnSe, (b) SnS, (c) GeSe, and (d) GeS.

Full size image
Table 3 The largest values of Seebeck coefficients (S) of bulk and monolayer SnSe, SnS, GeSe, and GeS at 300 K (unit: μ VK−1).
Full size table

Phonon dispersions of SnSe, GeSe, SnS, and GeS were computed to examine the thermal stability using density functional perturbation theory44, as shown in Fig. 5. There is no imaginary line in the dispersion curves, which means that these materials are vibrationally stable. There are twelve modes of vibrations: three loweer modes that are acoustic (a transverse acoustic mode, a longitudinal acoustic mode, and a flexural acoustic mode), and the others are optical modes. The flexural acoustic mode is an out-of-plane transverse acoustic mode similar to other two-dimensional materials like graphene, phosphorene, and stanene, quadratic near Γ point45,46,47. The flexural mode vibrational direction is exactly perpendicular to the plane. It is an important mode in order to understand thermal and mechanical properties of two-dimensional materials.

Figure 5
figure 5

Phonon dispersions along high symmetry k-points for (a) SnSe, (b) SnS, (c) GeSe, and (d) GeS.

Full size image

The lattice conductivities are calculated by solving the Boltzmann transport equation for phonons (BTEP) using the iterative method and the relaxation time approximation (RTA). The iterative method exactly solves the BTEP, while RTA is a good approximation for low conductivity compounds. The results in Fig. 6 show good agreement with recently reported results30. All four of the materials have very low lattice thermal conductivity compared to other two-dimensional materials like graphene, phosphorene, and monolayers of MoSe2 and WSe2, and are comparable to bismuth monolayer as shown in Table 4  25, 45, 46, 49. We also found the different lattice thermal conductivity along the armchair and zigzag directions. Because of the heavy masses of Sn and Se, SnSe has the lowest lattice conductivity of 2.44 Wm−1 K−1 and 2.63 Wm−1 K−1 at room temperature along the armchair and zigzag directions, respectively.

Figure 6
figure 6

Lattice thermal conductivity (κ l ) for the group IV–VI monolayers are calculated as a function of the temperature using iterative (solid lines) and SMRTA (dashed lines) method.

Full size image
Table 4 Comparison of lattice thermal conductivities κ l of group IV–VI monolayers with other two-dimensional materials at room temperature.
Full size table

According to glass dynamical theory, the lattice thermal conductivity is calculated as \({\kappa }_{l}=1/3{C}_{v}l{v}_{s}\), where C v is the heat capacity, l is the mean free path, and v s is the sound velocity. As the temperature increases, the lattice softens and the stiffness decreases, which reduces the sound velocity and hence the lattice thermal conductivity50. This trend is shown in Fig. 6.

Finally, using the Seebeck coefficient and the electrical and thermal conductivities, we calculated ZT as a function of the carrier concentration along the armchair and zigzag directions at 300 K, 500 K, and 700 K, as shown in Fig. 7. These monolayers had very high ZT. SnSe had the highest ZT of 2.63 along the armchair direction at 700 K because of the high electrical conductivity and Seebeck coefficient and the low lattice thermal conductivity. In the case SnS, a high ZT of 1.88 is predicted along the zigzag direction.

Figure 7
figure 7

Calculated Figures of merit (ZT) as a function of the carrier concentration (n) for the monolayer of (a,b) SnSe, (c,d) SnS, (e,f) GeSe, and (g,h) GeS along armchair and zigzag directions at temperature 300 K, 500 K, and 700 K.

Full size image

Conclusion

We analyzed the structural, electronic, thermoelectric, and phonon-transport properties of the two-dimensional monochalcogenide compounds SnSe, SnS, GeSe, and GeS using density functional theory combined with Boltzmann transport theory for electrons and phonons. These compounds are energetically and vibrationally stable, and SnSe, SnS, and GeS have indirect band gaps while GeSe has a direct band gap. The Seebeck coefficients of these two-dimensional materials are two times larger than those of their bulk structures, and two-dimensional GeS has the largest Seebeck coefficient of 2810 μ V K−1 at room temperature. These monolayer materials have very low lattice thermal conductivities in comparison to other two-dimensional materials. ZT of SnSe, GeSe and GeS along the armchair direction was 2.63, 1.99, and 1.85, respectively, while that of ZT of SnS along the zigzag direction was 1.88. These ZT values are higher than those of their bulk analogs. Hence, the materials are very promising for thermoelectric applications.

Methods

Our calculations are based on density functional theory combined with Boltzmann transport theory and were performed using the Vienna Ab initio Simulation Package (VASP) and the Boltztrap code51,52,53. The generalized gradient approximation proposed by Perdew-Burke-Ernzerhof was chosen as an electronic exchange correlation functional54. The vdW-DF scheme is used to include the van der Waals interaction55. A Monkhorst mesh of 10 × 10 × 1 k-points is used for lattice optimization and 450 eV is used as a plane wave cutoff energy. Structures are optimized until the Hellmann–Feynman force on each atom is less than 0.001 eV/Å. A vacuum region of 15 Å in the z-direction is produced to avoid the interaction between periodic images.

Thermoelectric properties were computed by solving the Boltzmann transport equation under a constant relaxation time (τ) and a rigid band approximation performed in the Boltztrap code53. We used a very dense k-point mesh of 60 × 60 × 1 to obtain a convergent density of states. The Seebeck coefficient S(T, n), electrical (σ/τ) and electronic thermal (κ e /τ) conductivities devided by the relaxation time were calculated. Boltztrap code used Wiedemann-Franz law to calculate electronic thermal conductivty from electrical conductivity.

To calculate the lattice thermal conductivity (κ l ), we used the ShengBTE code56. The second-order (harmonic) and the third-order (anharmonic) interatomic force constants (IFCs) are required to calculate lattice thermal conductivity. In order to calculate the second-order IFCs, we used the Phonopy code with a supercell of 5 × 5 × 1 and a k-mesh of 10 × 10 × 1. For the third-order IFCs, a supercell of 4 × 4 × 1 was used with the of interactions up to the 15th nearest neighbors57,58,59.