Fernan
Saiz
*a,
Jesús
Carrete
b and
Riccardo
Rurali
c
aInstitut de Ciència de Materials de Barcelona (ICMAB-CSIC), Campus de la UAB, Bellaterra, 08193, Spain. E-mail: fsaiz@icmab.es
bInstitute of Materials Chemistry, TU Wien, Getreidemarkt 9, 1060 Vienna, Austria. E-mail: jesus.carrete.montana@tuwien.ac.at
cInstitut de Ciència de Materials de Barcelona (ICMAB-CSIC), Campus de la UAB, Bellaterra, 08193, Spain. E-mail: rrurali@icmab.es
First published on 14th October 2020
The goal of this work is to investigate the influence of mechanical deformation on the electronic and thermoelectric properties of ZrS3 monolayers. We employ density functional theory (DFT) calculations at the hybrid HSE06 level to evaluate the response of the electronic band gap and mobilities, as well as the thermopower, the electrical conductivity, the phononic and electronic contributions to the thermal conductivity, and the heat capacity. Direct examination of the electronic band structures reveals that the band gap can be increased by up to 17% under uniaxial strain, reaching a maximum value of 2.32 eV. We also detect large variations in the electrical conductivity, which is multiplied by 3.40 under a 4% compression, but much smaller changes in the Seebeck coefficient. The effects of mechanical deformation on thermal transport are even more significant, with a nearly five-fold reduction of the lattice thermal conductivity under a biaxial strain of −4%. By harnessing a combination of these effects, the thermoelectric figure of merit of strained ZrS3 could be doubled with respect to the unstrained material.
(1) |
Another feature to be considered with a view to maximising the thermoelectric efficiency of layered trichalcogenides is the anisotropy of their thermal and electronic properties. For instance, in a recent work17 we have shown that the thermal conductivity of TiS3 monolayers and bilayers is highly dependent on the in-plane direction, with values at 300 K of 8.85 W m−1 K−1 in the x-direction and 13.56 W m−1 K−1 along the y-axis (the reference system is the same used here for ZrS3 and illustrated in Fig. 1). Moreover, under a 7% compression the electrical conductivity of TiS3 monolayers, and therefore their power factor, can nearly double in value with respect to the undeformed case.18 Apart from TiS3, little is known about the thermoelectric potential of other trichalcogenides. A recent study has compared the thermoelectric performances of the bulk phase and a monolayer of zirconium trisulphide (ZrS3),19 showing that the thermoelectric properties of the latter are anisotropic, a feature also exhibited by TiS3 in a number of experimental and theoretical studies.17,20,21 However, it is not clear if and how the interplay of anisotropy and strain can be harnessed to optimise the conversion of heat into electricity in ZrS3 monolayers.
Fig. 1 Ball-and-stick representation of the ZrS3 monolayer, with sulphur atoms depicted in blue and zirconium atoms in green. |
The goal of this work is thus to evaluate how the thermoelectric properties of zirconium trisulphide monolayers behave under strain. We use total energy first-principles calculations to obtain electronic and thermal properties such as electronic mobilities and band gaps, Seebeck coefficients, and electrical and thermal conductivities for a number of strains applied either uniaxially or biaxially to ZrS3 nanosheets. Those nanosheets are represented as a vacuum gap in contact with a unit cell whose atomic structure is optimised within density functional theory (DFT). We employ plane waves to represent the electronic wave functions, and the generalised gradient approximation (GGA) to the exchange and correlation energies as a starting point. We then recalculate the electronic structure using the hybrid exchange-correlation functional HSE06,22 which corrects the underestimated electronic band gap of the GGA calculation and substantially improves the agreement with the experimental values reported for ZrS3 nanostructures. Moreover, we compute the thermal properties of ZrS3 monolayers with a similarly fully ab initio approach starting from the relaxed atomic structure, paying attention to dynamical stability. ZrS3 is an n-type23 semiconductor with a direct optical band gap of (2.02 ± 0.02) eV, as measured in multilayered samples with an average thickness of (37 ± 2) nm.24 Previous measurements,23 performed on single-crystal nanostructures 10 mm in length, 1 mm in width, and a few microns in thickness, showed that this material has an outstanding thermopower of −850 μV K−1 measured at 300 K, with an electrical conductivity of 15 (Ω cm)−1 and an electron mobility of 25.98 cm2 V−1 s−1 for a charge carrier concentration of 4.5 × 1018 electrons per cm3. Those measurements correspond to transport along the atomic chains, i.e., the direction represented as the y axis in Fig. 1. A few years later, the electrical properties of ZrS3 crystals on a WSe2 substrate at 289 K were reported as σ = 0.038 (Ω cm)−1 and μ = 10.3 cm2 V−1 s−1 at a concentration of 2.33 × 1016 e− per cm3,25 along with a maximum Hall mobility of 29.1 cm2 V−1 s−1 at 144 K. More recently, suspended plates of ZrS3 yielded an electron mobility of 1132.9 cm2 V−1 s−1 and an electrical conductivity of 0.12 (Ω cm)−1 at a concentration of 0.66 × 1015 electrons per cm3.26 In this case, the suspended samples were fixed with silver paste on four contacts, which highlights to the key role of the substrate employed in the mobility measurements. We use those values as a reference to validate our calculations.
After this introduction, this manuscript is organised as follows. In Sec. 2, we detail the methodology employed to compute the charge mobilities, effective masses, Seebeck coefficients and other relevant thermoelectric quantities of ZrS3 monolayers under uniaxial and biaxial deformation. Next, in Sec. 3 we study the influence of stretching and compression on the electronic and thermal transport properties of the monolayer and evaluate the resulting figure of merit in order to detect the best strategies to increase the thermoelectric efficiency of those nanostructures. Finally, Sec. 4 presents our main findings and their implications for future work.
Once the geometry of the ZrS3 monolayer is relaxed, we obtain the following values for the moduli of the lattice vectors: a = 0.512 nm and b = 0.363 nm. We then apply uniaxial and biaxial strains from −4.0% to 7.0% along the x and/or y directions in increases of 1%, keeping the lattice vector along the z axis, i.e. the vacuum layer thickness, constant. In this range of strains we always find exclusively real phonon frequencies, which points to the dynamical stability of the strained systems investigated. For uniaxial strains exerted along the x (or y) direction, the cell is allowed to adjust along the y (or x) axis, thus enabling a straightforward calculation of the in-plane components of the Poisson's ratio. We then employ the HSE06 functional22 to obtain the electronic structure of the optimised structures, from which we can obtain more reliable thermoelectric transport properties using the BoltzTraP code.35 This code uses semiclassical Boltzmann theory under the uniform single-mode relaxation time approximation, and therefore the results for the electrical conductivity and the electronic contribution to the thermal conductivity are proportional to that parameter. The carriers' relaxation time is defined as
(2) |
(3) |
(4) |
(5) |
We now describe the workflow employed to calculate the thermal transport properties of ZrS3 monolayers. We use the almaBTE37 package, which needs matrices of second- and third-order interatomic force constants (IFCs) from supercells. These force constants are obtained with VASP using the same computational parameters described above with the exception of the reciprocal-space integration grid, which now consists only of the Γ point for each one of the displaced supercell configurations. Those configurations are generated by replicating the unit cell 4 × 5 times with Phonopy38 and thirdorder.py.39 A cutoff of 0.44 nm (i.e., up to the seventh nearest neighbors) is chosen for the third-order calculation. Due to the periodic boundary conditions and the corresponding breakdown of continuous rotation symmetry, the phonon frequencies of 2D systems calculated using this method are known to include numerical errors that can lead to the appearance of artifactual imaginary acoustic frequencies near Γ.18,40,41 Therefore, we need to apply a correction42 to enforce the rotational symmetry in our calculations. After this correction is applied, analysis of the phonon dispersion curves for all strains applied reveals that the structure becomes unstable for uniaxial and biaxial compressive strains higher than 4%, with non-negligible negative acoustic frequencies near Γ.
We also computed the Born effective charges in order to include the non-analytic corrections that, in 3D systems, result in a splitting of the transverse and longitudinal optical phonon branches at the Γ point. However, recent theoretical and experimental studies have shown that such a splitting is suppressed in 2D materials because of the fundamentally different physics of polar 2D vs. 3D systems regarding the effect of the dipole–dipole interactions.43,44
Once all force constants are computed, we use them as inputs to solve the Boltzmann Transport Equation (BTE) self-consistently45 and obtain the lattice thermal conductivity as
(6) |
Strain | Carrier | C | E 1 | m * | μ | τ |
---|---|---|---|---|---|---|
ε x | Hole | 5.48 | −1.03 | −0.54 | 4038.27 | 1235.28 |
ε y | Hole | 8.09 | −7.16 | −0.97 | 37.90 | 20.97 |
ε x | Electron | 5.48 | 0.39 | 1.65 | 3004.82 | 2826.12 |
ε y | Electron | 8.09 | 3.43 | 0.29 | 1810.53 | 302.45 |
Strain | Carrier | C | E 1 | m * | μ | τ |
---|---|---|---|---|---|---|
ε b | Hole | 15.81 | −8.04 | −0.54 | 192.31 | 58.83 |
ε b | Hole | 15.81 | 4.36 | −0.97 | 199.18 | 110.24 |
ε b | Electron | 15.81 | −8.04 | 1.65 | 20.34 | 19.13 |
ε b | Electron | 15.81 | 4.36 | 0.29 | 2186.61 | 365.27 |
Regarding the calculation of mobility from DFT using eqn (3), the disagreement may be caused by an incorrect calculation of the effective mass m* or of the constant E1 from deformation potential theory (or both) or by limitations of the formalism itself versus a more detailed treatment of electron-phonon scattering. However, it seems reasonable to assume that the HSE06 functional likely reproduces the shape of the conduction band in reciprocal space fairly well for two reasons. First, we find a very good agreement of the experimental optical and electronic band gaps; second, our hole effective masses are in excellent agreement with the values of −0.87 (±0.09) me along the x axis and −0.49 (±0.20) me along the y axis recently measured on ZrS3 nano-whiskers using nanospot angle resolved photoemission spectroscopy (nanoARPES).49
If we analyse the conditions at which simulations and experiments are carried out, we find the following factors that can be responsible for the disagreement: first, the ZrS3 samples in the experiments are multilayers rather than the monolayer we model with DFT; second, these samples can contain crystallographic defects or surface impurities as well as slightly different compositions that deviate from their ideal stoichiometry; and third, the samples were measured with different substrates. The first cause can be expected to have only a minor impact because the band gap of the experimental multilayer24 is very similar to that reported in the bulk19 and computed for the monolayer in this work. The second factor is likely to play a more significant role since the presence of crystallographic defects cannot be completely ruled out. Moreover, the stoichiometry of exfoliated ZrS3 layers has been reported to differ from the ideal ratio S/Zr = 3/1, which might a priori produce relevant changes in the electronic structure. Additionally, a strong influence is also exerted by the substrate employed to encapsulate the sample in the laboratory setup. For instance, the data provided in ref. 23 and 25 differs from that published in ref. 26. The first study used samples suspended on four gold contacts, while the second work employed WSe2 as a substrate, and the third fixed the samples with four silver contacts. Since silver is the best conductor among the three materials, it is to be expected that the highest mobility be recorded with this material.
Fig. 2 Electronic band gaps Egap, calculated using the HSE06 functional, as a function of the uniaxial strains, (εx) and (εy), and the biaxial strain, (εb). |
Another crucial electronic property in n-type semiconductors, such as ZrS3, is the electron mobility, which dictates how fast an electron can move through the lattice, in effective terms, under the effect of an electric field. In this study, we are interested in quantifying the values of the electron mobility and its anisotropy, but also how much it varies with mechanical deformation. Fig. 3 shows a strong growth of μy with increasing εx and a slower increase of μx with increasing εy. In contrast, μx decreases with increasing εx and μy decreases with increasing εy. Therefore, the mobility increases under perpendicular strains but is reduced under parallel uniaxial deformations. The enhancements upon compressions of 4% are 51.73% along the x axis and 45.00% along the y axis, while the correspondent percentual reductions upon 7% stretching are 49.70% and 51.49%. For the biaxial cases, both components of the mobility decrease with increasing strain, which corresponds to the widening the band gap as shown in Fig. 2.
Fig. 3 Electronic mobilities vs. uniaxial strains using the HSE06 functional along the x (a) and y (b) directions and biaxial strains (c). |
Fig. 4 Effect of the electron concentration on the directional Seebeck coefficient Sy when strain is applied in the y direction of ZrS3 monolayers. |
We now discuss how we choose the relaxation times necessary to obtain numerical values of the electrical conductivity. Our DFT workflow, based on band deformation potentials, predicts electron relaxation times τx,DFT = 2825.12 fs and τy,DFT = 302.45 fs, whereas if we accept the experimental mobility reported in ref. 26, μexp = 1132.9 cm2 V−1 s−1, and plug it into eqn (2), we obtain the corresponding relaxation times τx,DFT-exp = 2826.12 fs and τy,DFT-exp = 302.45 fs, which in turn yield electrical conductivities at 0.66 × 1015 e− per cm3 of σx,DFT-exp = 0.46 (Ω cm)−1 and σy,DFT-exp = 0.17 (Ω cm)−1. Note that the latter is in excellent agreement with the 0.12 (Ω cm)−1 reported in ref. 26. We therefore opt for this approach, which is justified in view of the inability of the theoretical calculations to capture the complex experimental conditions in terms of substrate, defects, grain boundaries and other kind of lattice imperfection.
Note though that experimental data can also be taken to calculate electron mobilities following the method proposed in ref. 51. This method starts by setting the experimental Seebeck coefficient S(T;μ) as an input to determine the reduced chemical potential ηeff from its definition in BoltzTraP
(7) |
(8) |
In our case, we use the Seebeck coefficient component Sy = −850 μV K−1 reported in ref. 23; in eqn (7) to find a reduced chemical potential of ηeff = −9.74. We then plug this value in the following expression
(9) |
Having chosen appropriate relaxation times, we now illustrate the characteristic response of the electrical conductivity with respect to electron concentration and strain. The response of the in-plane component σyvs. εy is depicted in Fig. 5, and all combinations of in-plane components σx and σyvs. εx, εy and εb are shown in Fig. S2.† For uniaxial strains, each component experiences higher changes when the compression is parallel to the mechanical deformation. For instance the component σx increases by a factor of 2.87 with an εx of −4% and by 3.40 with an εb of −4%, and the component σy by 2.27 under an εy of −4%, at 0.66 × 1015 electron per cm3. Therefore, these enhancements evidence not only the possibility of improving the electronic transport by narrowing the band gap but also by increasing electron mobility. Specifically, the most remarkable improvements in electrical conductivity are detected when the electron mobility grows by nearly 50% upon 4% compression of the monolayer.
Fig. 5 Effect of electron concentration on the directional electrical conductivity σy when strain is applied in the y direction of ZrS3 monolayers. |
κel = LσT, | (10) |
Fig. 6 and S3† show that strains applied along the x and y axes induce a similar response of the components κxxel and κyyel to those exhibited by their electrical conductivity counterparts, i.e., shrinking the ZrS3 nanosheet results in an enhancement of thermal transport by electrons, whereas stretching it reduces the amount of heat carried by these particles. Note that for concentrations higher than 1017 electrons per cm3, the heat transferred by electrons is comparable to that carried by phonons, as we show next. For instance, for the undeformed case, κxxel requires a concentration of 1.20 × 1017 electron per cm3 to reach a value of 1.15 W m−1 K−1, while κyyel requires a concentration of 5.75 × 1017 electron per cm3.
Fig. 6 Effect of the electron concentration on the directional electronic contribution to the thermal conductivity κel,yy when strain is applied in the y direction. |
We now analyse the effects of deforming ZrS3 on phonon transport. We first show in Fig. 7 how the phonon dispersion curves change under biaxial strains of −4% and 7%. In particular, we draw the lowest-lying branches between Γ = [0.0, 0.0, 0.0] and Y = [0.125, 0.0, 0.0] to illustrate how these curves change near the zone centre when the rotational symmetry is restored, as described above. The results underscore the importance of this correction since, going by the raw results of the phononic calculations, a pocket of negative frequencies forms which would point to the presence of mechanical instabilities. Those disappear when the proper correction is applied and would also look unlikely in view of the fact (see Fig. S4†) that we are working in the elastic regime, with perpendicular deformations in x and y created in the monolayer as a response to the correspondent strain exerted in y and x that follow a linear relationship. As a result, we can extract the Poisson's ratios νxy = 0.12 and νyx = 0.18. The correction is also important to obtain reasonable values of the lattice thermal conductivity.18,42
Fig. 7 Details of the phonon dispersions from Γ = [0, 0, 0] to Y = [1/2, 0, 0] before (in black) and after (in red) correcting the harmonic IFCs according to the procedure given in ref. 40 for the ZrS3 monolayer for biaxial strains of −4% (a), 0% (b), and +7% (c) and their corresponding full views in (d), (e), and (f) along the path Γ–Y–C–Z–Γ, where C = [1/2, 1/2, 0] and Z = [0, 1/2, 0]. |
We then plot the components κph,xx and κph,yy of the thermal conductivity tensor (as calculated from the corrected force constants) in Fig. 8 to illustrate their different responses to uniaxial and biaxial stresses. When the monolayer is compressed by −4% along the x axis, both components decrease, down to 56.75% for κph,xx and 82.25% for κph,yy with respect to their zero-strain values of 18.65 W m−1 K−1 along the x-axis and 35.12 W m−1 K−1 along the y-axis. Meanwhile, a tensile εx has no significant effect on either in-plane conductivity. Along the y-axis, both components decrease with increasing strain, down to 7.49% of the original value for κph,xx and 19.20% of the no-strain value for κph,yy with εy = −4%. Nonetheless, the most significant effects are seen under biaxial strains εb: a compression of 7% causes both conductivities to decrease by 49.62% for κph,xx and 74.93% for κph,yy, and a compression εb = 4% decreases κph,xx by 30.05% and κph,yy by 34.28%. We note that the thermal transport along the x direction is enhanced for intermediate stretching as κph,xx increases by a factor of 1.29 with εx = 1% and 1.70 for with εb = 4%.
Since our study is intended to find those cases in which the figure of merit can be enhanced due to a reduction in the thermal transport, we now perform a deeper analysis of how strain reduces the thermal conductivity at 4% compression and at 7% stretching. According to eqn (6), changes in the thermal conductivity must be correlated to variations in the frequencies, group velocities or lifetimes of phonons. We show in Tables S1–S4† that the group velocities of the six lowest-lying branches (as a sample of the modes that carry most of the heat) vary significantly near Γ for uniaxial and biaxial deformations. However, not all of these group velocities are reduced upon compression. Therefore, we make a modal decomposition of the thermal conductivities (see Fig. S6 and S7† for details) and find that the lower frequencies contribute less when the monolayer is shrunk by 4%. Additionally, we show in Fig. S6† that heat capacities can be enhanced by a factor of up to 1.07 with εb = −4%, decreasing with increasing strain in both uniaxial and biaxial scenarios.
A similar picture is also exhibited by the figure of merit in Fig. 10, with the addition now that the component (zT)y is also enhanced under εx = εy = −4% as a result of depressed thermal transport caused by a shrinking of the monolayer. Although the absolute values of S2σ and zT depend on the choice of the relaxation time, our calculations predict that compression is successful at increasing the maximum value of zT with respect to the case with no stress applied on the atomic structure. In particular, the maximum (zT)x increases from 0.55 with 1.10 × 1019 electron per cm3 under no strain to 0.69 with 6.72 × 1018 electron per cm3 under a compressive strain εx = −4% and to 0.69 with 5.06 × 1018 electron per cm3 and εb = −4%, representing improvements of nearly 25% in both cases. A more striking result is obtained for the maximum (zT)y, which grows from 0.25 with 1.66 × 1019 electron per cm3 under no deformation to 0.47 with 1.05 × 1019 electron per cm3 applying εx = −4% and to 0.52 at 9.96 × 1019 electron per cm3, representing a rise of 86.06% and 107.00%, respectively. Based on these relative results, our DFT methodology is able to suggest a recipe to produce substantial improvements in the thermoelectric efficiency of ZrS3 due mostly to the reduction of the phonon thermal transport but also to a higher electrical conductivity.
Overall, our calculations suggest that using strain to engineer the band gap and, more generally, the thermoelectric behaviour of a monolayer requires careful consideration of the direction along which deformation is applied. We note that possible buckling under compressive strain is outside the scope of this study, although the lack of imaginary phonon frequencies in our calculations suggests that the shrinkages that we propose should be achievable without inducing buckling.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/d0na00518e |
This journal is © The Royal Society of Chemistry 2020 |