nanomaterials Article Anisotropic Thermal Conductivity in Few-Layer and Bulk Titanium Trisulphide from First Principles Fernan Saiz 1 , Jesus Carrete 2 and Riccardo Rurali 1, * 1 Institut de Ciència de Materials de Barcelona (ICMAB–CSIC) Campus de Bellaterra, Bellaterra, 08193 Barcelona, Spain; fsaiz@icmab.es 2 Institute of Materials Chemistry, TU Wien, A-1060 Vienna, Austria; jesus.carrete.montana@tuwien.ac.at * Correspondence: rrurali@icmab.es  Received: 6 March 2020; Accepted: 2 April 2020; Published: 8 April 2020  Abstract: We study the thermal conductivity of monolayer, bilayer, and bulk titanium trisulphide (TiS3 ) by means of an iterative solution of the Boltzmann transport equation based on ab-initio force constants. Our results show that the thermal conductivity of these layers is anisotropic and highlight the importance of enforcing the fundamental symmetries in order to accurately describe the quadratic dispersion of the flexural phonon branch near the center of the Brillouin zone. Keywords: thermal transport; first-principles calculations; lattice dynamics; 2D materials 1. Introduction The study of two-dimensional materials is a very active research field, as these compounds exhibit enhanced properties with respect to those shown by their bulk counterparts. For example, the thermal conductivity of graphene [1,2] has been reported to be as high as 5300 W m−1 K−1 , which is an order of magnitude higher than graphite’s [3]. Beyond graphene, the properties of new two-dimensional materials are being investigated as new methods to exfoliate flakes with thicknesses of several nanometers are being developed. Among these materials, transition metal trichalcogenides (MX3 with M = Ti, Zr, or Hf and X = S, Se, or Te) have attracted much attention due to their narrow band gap of around 1.0 eV [4,5], which makes them not only suitable for optical and electronic properties [6], but also for thermoelectric applications [7]. Furthermore, the structure of these chalcogenides is anisotropic, leading to different transport coefficients in the flake plane, from which undesirable backscattering from hot electrons can be minimized [8]. Within this trichalcogenide family, recent studies have shown that titanium trisulphide (TiS3 ) is an excellent candidate to build electronic and optic nanodevices [6], electrodes for Li and Na ion batteries [9], anodes for hydrogen photogeneration [10], and photodetectors [11–13] due to a direct band gap of 1.0 ± 0.1 eV for thin films with a thickness of around 300 nm [14]. TiS3 is an n-type semiconducting material formed by abundant and low toxicity elements that is also attractive for thermoelectric applications because of its large Seebeck coefficient S of −650 µV/K at 300 K and a low thermal conductivity, κ, of 3.7 W m−1 K−1 at 325 K in bulk form [15]. Similar measurements have reported a thermopower of nearly −700 µV/K [16] at 300 K with an electrical conductivity σ at around 200 K between 0.5 and 5.0 (Ω cm)−1 . The efficiency of a thermoelectric material is characterized by the figure of merit, which is calculated for a given temperature T as S2 σ zT = T, (1) κ where σ and κ are the electrical and thermal conductivities, and S is the Seebeck coefficient. Therefore, any degree of anisotropy in all transport properties present in Equation (1) can be exploited Nanomaterials 2020, 10, 704; doi:10.3390/nano10040704 www.mdpi.com/journal/nanomaterials Nanomaterials 2020, 10, 704 2 of 10 to optimize the conversion of heat into electricity in these compounds. In particular, we ask ourselves if anisotropy can be a degree of freedom to increase the efficiency of few-layer and bulk TiS3 , particularly where the reduction of the thermal conductivity is concerned. Estimates of the thermal properties of monolayer TiS3 from first-principles have been previously reported by Zhang et al. [17]. In this work, however, we also compute the thermal conductivity of bilayers and of the bulk form of this material to help fill the gap in the experimental data between the well-established results for the single layer and for (bulk) multilayers. The goal of this work is then to calculate the thermal conductivity tensor of monolayers and bilayers of TiS3 , as well as that of its bulk phase. This calculation is performed employing a first-principles approach whereby density functional theory (DFT) is used to relax the geometry of these TiS3 structures to then predict their thermal properties with a lattice dynamics scheme. Our intention here is to track the evolution of the thermal conductivity tensor as the thickness of a nanosheet is doubled and compare it with that obtained in the bulk, as limit of a thick multilayer. In this work we pay special attention to the correct enforcement of the translation and rotation symmetries of free space for the two nanosheets and to their effect on the thermal conductivity. Conserving these symmetries requires correcting the second-order force constants by defining a set of internal coordinates for the monolayer and bilayer that span the same space as the original Cartesian coordinates of the crystalline system except for rigid translations and rotations. This has been shown to lead to a correct description of out-of-plane vibrations [18]. Therefore, we apply the correction formulated in Ref. [18], where it was found that by enforcing these symmetries in the case of borophene, its computed thermal conductivity decreases by 50% and reverses its anisotropy as well as that the lowest-lying phonon branch becomes quadratic in the neighborhood of the Γ point in reciprocal space. The manuscript is organized as follows: in Section 2 we describe the recipe used to run our first-principles calculations, starting from the geometry optimisations of the TiS3 systems to then compute the thermal properties with lattice dynamics; next, in Section 3, we compare the phonon dispersions of these systems before and after correcting the interatomic force constants and we analyse the impact of this correction on the anisotropy and values of the thermal conductivity tensor; finally, Section 4 concludes with our main findings and the implications of this work. 2. Theoretical Methods Our methodology starts by employing the Vienna Ab Initio Simulation Package (VASP) [19–22] to relax the positions and cell parameters of these nanosheets. The monolayers are first represented as a unit cell belonging to the P21/m space group as shown in Figure 1. Reciprocal space integrations are performed using a mesh of 10 × 14 × 1 k-points centered at Γ in the Brillouin zone. The geometry relaxation is carried out setting thresholds of 0.01 eV/nm for the forces and 1 × 10−6 eV for the self-consistent solution of the wavefunction. We use the generalized gradient approximation to the exchange-correlation potential in the Perdew-Burke-Ernzerhof (PBE) flavour [23] and the projector augmented wave method (PAW) [24,25]. We expand the valence orbitals on a plane wave basis with an energy cutoff of 350 eV. Once the geometry of the TiS3 monolayer is optimized, we obtain the following parameters for the monolayer: a = 0.496 nm, b = 0.339 nm; for the bilayer: a = 0.495 nm, b = 0.340 nm; for the bulk: a = 0.494 nm, b = 0.339 nm, c = 0.887, and β = 96.10 degrees. As a reference, bulk values from Ref. [26] are a = 0.458 nm, b = 0.340 nm, c = 0.878 nm, and β = 97.32 degrees. Long-range van der Waals forces are included using the zero damping Grimme DFT-D3 scheme [27]. We choose this method because it yields an interlayer separation of 0.313 nm, in better agreement with the bulk value than the predictions of other van der Waals forces descriptions available in VASP. Nanomaterials 2020, 10, 704 3 of 10 z x y Figure 1. Representation of a TiS3 unit cell with atoms painted in blue for titanium and in orange for sulfur. Starting from the optimized primitive cell, we compute the interatomic force constants (IFCs) in 4 × 5 supercells for the monolayer and the bilayer and 4 × 5 × 3 supercells for the bulk using a finite differences method. For the harmonic displacements we use the P HONOPY code [28] considering all neighboring interactions. THIRDORDER . PY [29,30] is used to characterize the anharmonic interactions, neglecting those beyond seventh neighbors. The IFCs are then used as an input to solve the Boltzmann Transport Equation (BTE) iteratively with the almaBTE code [31] and the lattice thermal conductivity is obtained as κij = ∑ κij,λ = C ∑ f λ ( f λ + 1)(hνλ )2 vi,λ Fj,λ , (2) λ λ where i and j are the spatial directions x, y, and z, C −1 = k B T 2 ΩN and k B is Boltzmann’s constant, h the Planck’s constant, Ω the unit cell’s volume, and N the number of q-points. The summation in Equation (2) runs over all phonon modes λ; each mode has a frequency νλ and a group velocity vλ , and at thermal equilibrium at temperature T its occupancy follows the Bose-Einstein distribution f λ . The mean free displacement is initially calculated as Fj,λ = τλ v j,λ , where τλ is the lifetime of mode λ within the relaxation time approximation (RTA). Starting from this guess, the solution is then obtained iteratively and Fj,λ takes the general form τλ (v j,λ + ∆ j,λ ), where the correction ∆ j,λ captures the changes in the heat current deriving from the deviations in the phonon populations with respect to the solution at the RTA level and that result from the iterative process [32,33]. Scattering from isotopic disorder is also included considering the natural distributions of Ti and S isotopes within Tamura’s model [34]. 3. Results and Discussion Our first step is to evaluate the dependence of the thermal conductivity tensor on the number of q-points included in the solution of the BTE. This is a preliminary analysis necessary to assess the accuracy of our calculations. Figure 2 shows that, while the values of κ xx at 300 K are virtually insensitive to the number of q-points, we obtain variations of κyy of less than 2.0% with a grid of least 10 × 10 q-points, if compared to the results obtained with a much thicker and computationally expensive 28 × 28 mesh. This accuracy is comparable with the dispersions typically reported by experimental measurements. We observe a similar behavior at 100 K, where convergence might be more difficult to achieve because of the larger relative weight of long-wavelength phonons. At this Nanomaterials 2020, 10, 704 4 of 10 temperature we observe that κ xx oscillates at around an average value of 33.31 W m−1 K−1 for all grid resolutions, whereas κyy requires at least a mesh of 16 × 16 q-points to stabilize with fluctuations smaller than 7% for finer meshes. Based on these convergence tests, hereinafter we present all our results in this manuscripts using a mesh of 16 × 16 q-points for the monolayer and bilayer and 16 × 16 × 16 q-points for the bulk phase, which both give an optimal compromise between accuracy and computational load. 100 κxx,100 K 90 κyy,100 K 80 κxx,300 K κyy,300 K 70 κ (W m-1 K-1) 60 50 40 30 20 10 0 4 8 12 16 20 24 28 q x, q y Figure 2. Thermal conductivity as a function of the number of q-points in the x (qx ) and y (qy ) directions at 100 K and 300 K. We next illustrate in Figure 3 the dispersion relations calculated using the second-order force constants for these three structures. The most striking finding is that the lowest acoustic branches near the Γ point exhibit a pocket of negative frequencies along the Γ-Y direction for both nanosheets. For the monolayer, this pocket yields imaginary frequencies as low as −0.19 THz at q = (0.06, 0.0, 0.0) and for the bilayer, this singularity is extended farther away from Γ with lower frequencies that decrease down to −0.91 THz at q = (0.13, 0.0, 0.0). The increased number of low-frequency bands in close proximity to one another in the case of the bilayer makes the artifact worse by increasing the sensitivity of the results to spurious contamination between polarizations, but the problem completely disappears in the 3D solid once the monolayer bands are hybridized to give the familiar bulk acoustic branches. These imaginary frequencies appear because the flexural phonon band near Γ is inaccurately described, as the periodic boundary conditions imposed along all three axes in DFT calculations of layered systems cause a breakdown of the continuous rotation symmetry at the root of that quadratic character of the ZA branch. [18,35] Instead, this band should have a quadratic dispersion with a vanishing derivative in q, i.e., ∂ω/∂q → 0. Nanomaterials 2020, 10, 704 5 of 10 4.0 4.0 4.0 3.5 Corrected 3.5 3.5 3.0 Uncorrected 3.0 3.0 a b c 2.5 2.5 2.5 Frequency (THz) 2.0 2.0 2.0 1.5 1.5 1.5 1.0 1.0 1.0 0.5 0.5 0.5 0.0 0.0 0.0 -0.5 -0.5 -0.5 -1.0 -1.0 -1.0 Γ Y/4 Γ Y/4 Γ Y/4 18 18 18 16 16 16 14 14 14 d e f Frequency (THz) 12 12 12 10 10 10 8 8 8 6 6 6 4 4 4 2 2 2 0 0 0 Γ Y C Z Γ Γ Y C Z Γ Γ Y C Z Γ Figure 3. Zoomed phonon dispersions along the Γ-Y direction before and after correcting the harmonic IFCs according to the procedure given in Ref. [18] for the monolayer (a) and bilayer (b) as well as the uncorrected spectra for the bulk (c). Full corrected dispersion relations for the monolayer (d) and the bilayer (e) and uncorrected for the bulk (f). The following step is then to enforce rotational symmetry by applying the correction proposed in Ref. [18]. This correction adds a post-processing step to the second-order force-constant calculation, where it proceeds by systematically building a basis of independent internal coordinates (distances, angles, and dihedrals) that are explicitly scalar, i.e., rotationally invariant. The Cartesian force constants are then projected on that basis and back. Since those internal coordinates can express any atomic movement except for the three continuous translations and the three continuous rotations, the result is compatible with the homogeneity and isotropy of free space by construction, leading to the expected quadratic behavior of the flexural branch. As shown in Figure 3a,b, the pocket is also exhibited by the other acoustic bands close to Γ, due to cross-contamination between polarizations; the correction solves the issues for those bands too, which is especially noticeable in the case of the bilayer. After correcting the second-order force constants, the pocket of negative frequencies disappears and the flexural band exhibits the expected quadratic dependence in the full dispersion relations illustrated in Figure 3c,d. With the corrected harmonic force constants, we now calculate the thermal conductivities for both nanosheets as a function of temperature between 100 K and 600 K under the RTA and by solving the full BTE. Figure 4 indicates that this quantity decreases with decreasing temperature, indicating that Umklapp processes are the dominant scattering mechanism for the RTA and the full iterative solution. We also find that solving the full BTE produces higher conductivities than those determined under the RTA. For instance, this increase for the monolayer is 29.40% for κ xx and 13.34% for κyy at 300 K. This enhancement gives an estimate of the role of momentum conserving Normal processes as a function of temperature and the underestimation of RTA derives from the fact that they are erroneously considered as resistive collisions. Furthermore, we evaluate the influence of symmetrizing Nanomaterials 2020, 10, 704 6 of 10 the IFCs on the thermal conductivity. Figure 4b,d show that only the x-component is significantly affected by this correction since the pocket of imaginary frequencies appears in the Γ-Y direction with variations of 6.24 % for the monolayer and −6.44 % for the bilayer at 300 K. 1.15 70 κxx,RTA 60 κxx,BTE 1.10 a b κyy,RTA 50 κxx, κyy (W m-1 K-1) κyy,BTE 1.05 κcorr / κuncorr 40 1.00 30 20 0.95 10 0.90 0 0.85 100 200 300 400 500 600 100 200 300 400 500 600 1.15 70 60 1.10 c d 50 κxx, κyy (W m-1 K-1) 1.05 κcorr / κuncorr 40 1.00 30 20 0.95 10 0.90 0 0.85 100 200 300 400 500 600 100 200 300 400 500 600 Figure 4. Corrected values of the thermal conductivity as a function of temperature after correcting the IFCs of the TiS3 monolayer (a) and bilayer (b). Ratios of the corrected thermal conductivity with respect to their uncorrected values for the monolayer (c) and bilayer (d). Our results evidence a significant anisotropy of thermal transport in layered TiS3 , as opposed to other 2D materials such as graphene [36] or transition metal dichalcogenides [35,37]. For instance, our BTE values at 300 K are κ xx = 8.24 W m−1 K−1 and κyy = 16.28 W m−1 K−1 for the monolayer and κ xx = 6.36 W m−1 K−1 and κyy = 11.57 W m−1 K−1 for the bilayer, yielding anisotropy ratios κyy /κ xx of 1.98 for the former and 1.82 for the latter. For the bulk, we find BTE values at 300 K are κ xx = 8.28 W m−1 K−1 and κyy = 14.34 W m−1 K−1 (see Figure 5). We believe that this anisotropy in all three systems is a consequence of the atomic structure of TiS3 that, to some extent, can be considered as consisting of covalently bonded chains along the y-axis, with weaker interchain couplings along the x-axis. Furthermore, our calculations suggest that the thermal conductivity reduction when doubling the monolayer thickness is caused by the van der Waals forces. Although these forces are weak, they create small distortions in the atomic structure of the two joined layers that perturb the original degeneracy of their phonon branches. As a result, this perturbation increases the phase space for anharmonic scattering and hence its rates. This trend is inverted when moving from the bilayer to the bulk, where the thermal conductivity increases and the values of the monolayer are almost recovered. In this case the higher symmetry of an infinitely thick multilayer (i.e., the bulk), narrows again the phase space and thus the thermal conductivity can increase. As expected, we predict a much more Nanomaterials 2020, 10, 704 7 of 10 inefficient thermal transport in the out-of-plane direction, where phonon propagation is mediated by the low van der Waals interlayer interaction, with values of κzz that are 7 times smaller than those of κ xx and 12 times than those of κyy . 70 κxx,RTA κxx,BTE 60 κyy,RTA κyy,BTE 50 κxx, κyy (W m-1 K-1) 40 30 20 10 0 100 200 300 400 500 600 Temperature (K) Figure 5. Evolution of the components of the thermal conductivity tensor as a function of temperature of bulk TiS3 . Last but not least, we evaluate the contributions of all phonon frequencies to the thermal conductivity. Figure 6a illustrates that the component κ xx plateaus for all three systems at around 5 THz. At this frequency, the individual contributions drop to a value below 0.1 W m−1 K−1 THz−1 given that most of the heat is transferred by the acoustic branches drawn in Figure 3a,b. In contrast, the component κyy for both layers plateaus near 7 THz as the second-lowest lying acoustic branch is slightly steeper near Γ in the direction Γ-Z than along the Γ-Y. This variation produces a slightly higher group velocity, which in turn leads to higher contributions to the thermal conductivity. The plots of Figure 6 also reveal that the higher thermal conductivity of the monolayer essentially derives from the larger contribution of low-frequency phonons, particularly in the case of κ xx where its cumulative value dramatically rises at around 1 THz. 1.0 1.0 1.0 a d g κcum/κcum,max κcum/κcum,max κcum/κcum,max 0.8 0.8 0.8 0.6 κxx 0.6 κxx 0.6 κxx 0.4 κyy 0.4 κyy 0.4 κyy 0.2 0.2 0.2 0.0 0.0 0.0 100 b 100 e 100 h 10-1 10-1 10-1 κxx κxx κxx 10-2 10-2 10-2 10-3 10-3 10-3 100 c 100 f 100 i 10-1 10-1 10-1 κyy κyy κyy 10-2 10-2 10-2 10-3 10-3 10-3 0 5 10 15 20 0 5 10 15 20 0 5 10 15 20 Frequency (THz) Frequency (THz) Frequency (THz) Figure 6. Modal decomposition of the components of the thermal conductivity tensor κ xx and κyy (both in W m−1 K−1 THz−1 ) for the monolayer (a–c), bilayer (d–f), and bulk (g–i) with respect to the phonon frequency. Panels (a,d,g) show the normalised cumulative function of both components. Nanomaterials 2020, 10, 704 8 of 10 4. Conclusions In conclusion, we have investigated the thermal transport in monolayers and bilayers of titanium trisulphide as well as in its bulk phase with first-principles calculations. The second-order force constants are corrected to enforce rotational symmetry, which has been recently shown to be lost by imposing periodic boundary conditions in systems with vacuum gaps. Our results indicate that the thermal conductivity is highly anisotropic with values at 300 K of 8.24 W m−1 K−1 in the x direction and 16.28 W m−1 K−1 along the y-axis for the monolayer. For the bilayer, our model predicts significant reductions of 22.72 % with 6.36 W m−1 K−1 in the x direction and 28.95 % with 11.57 W m−1 K−1 along the y-axis. This decrease is attributed to the van der Waals forces that introduce small coupling between the modes of both monolayers, slightly break the degeneracy of the phonon bands, and therefore increase the phase space for anharmonic scattering. In addition, these in-plane conductivities show the highly anisotropic nature of thermal transport in layered titanium trisulphide with ratios of 1.98 for monolayer and 1.82 for the bilayer. We believe that our results will help the community to recognize the importance of taking into account the loss of rotational symmetry in first-principle calculations of 2D systems with a vacuum gap, as well as to compare more carefully numerical data from monolayers and the experimental measurements taken samples with thicknesses of several layers. These results are also important for emerging nanodevices based on 2D materials and many layered van der Waals crystals, which are expected to play an important role in several applications, ranging from optoelectronics to thermoelectricity. Author Contributions: Conceptualization, F.S., J.C. and R.R.; methodology, F.S., J.C. and R.R.; software, F.S. and J.C.; validation, F.S. and J.C. and R.R.; formal analysis, F.S., J.C. and R.R.; investigation, F.S., J.C. and R.R.; resources, F.S.; data curation, F.S., J.C. and R.R.; writing—original draft preparation, F.S. and R.R.; writing—review and editing, F.S., J.C. and R.R.; visualization, F.S.; supervision, R.R.; project administration, R.R.; funding acquisition, F.S. and R.R. All authors have read and agree to the published version of the manuscript. Funding: This research received the financial support by the Ministerio de Economía, Industria y Competitividad (MINECO) under grant FEDER-MAT2017-90024-P and the Severo Ochoa Centres of Excellence Program under Grant SEV-2015-0496 and by the Generalitat de Catalunya under grants no. 2017 SGR 1506. The authors acknowledge the funding received from the European Union’s Horizon 2020 research and innovation programme under grant agreement No. 793726 (TELIOTES—Thermal and electronic transport in inorganic-organic thermoelectric superlattices). Acknowledgments: We acknowledge the support of The Supercomputing Center of Galicia (CESGA) where the calculations have been carried out. The authors also thank fruitful discussions with Andrés Castellanos-Gómez. Conflicts of Interest: The authors declare no conflict of interest. References 1. Novoselov, K.S.; Geim, A.K.; Morozov, S.V.; Jiang, D.; Zhang, Y.; Dubonos, S.V.; Grigorieva, I.V.; Firsov, A.A. Electric Field Effect in Atomically Thin Carbon Films. Science 2004, 306, 666–669. [CrossRef] 2. Novoselov, K.S.; Jiang, D.; Schedin, F.; Booth, T.J.; Khotkevich, V.V.; Morozov, S.V.; Geim, A.K. Two-dimensional atomic crystals. Proc. Natl. Acad. Sci. USA 2005, 102, 10451–10453. [CrossRef] [PubMed] 3. Balandin, A.A.; Ghosh, S.; Bao, W.; Calizo, I.; Teweldebrhan, D.; Miao, F.; Lau, C.N. Superior thermal conductivity of single-layer graphene. Nano Lett. 2008, 8, 902–907. [CrossRef] [PubMed] 4. Flores, E.; Ares, J.R.; Ferrer, I.J.; Sánchez, C. Synthesis and characterization of a family of layered trichalcogenides for assisted hydrogen photogeneration. Phys. Status Solidi-Rapid Res. Lett. 2016, 10, 802–806. [CrossRef] 5. Biele, R.; Flores, E.; Ares, J.R.; Sanchez, C.; Ferrer, I.J.; Rubio-Bollinger, G.; Castellanos-Gomez, A.; D’Agosta, R. Strain-induced band gap engineering in layered TiS3. Nano Res. 2018, 11, 225–232. [CrossRef] 6. Jin, Y.; Li, X.; Yang, J. Single layer of MX3 (M = Ti, Zr; X = S, Se, Te): A new platform for nano-electronics and optics. Phys. Chem. Chem. Phys. 2015, 17, 18665–18669. [CrossRef] 7. Saeed, Y.; Kachmar, A.; Carignano, M.A. First-principles study of the transport properties in bulk and monolayer MX3 (M = Ti, Zr, Hf and X = S, Se) compounds. J. Phys. Chem. C 2017, 121, 1399–1403. [CrossRef] Nanomaterials 2020, 10, 704 9 of 10 8. Abudukelimu, A.; Kakushima, K.; Ahmet, P.; Geni, M.; Tsutsui, K.; Nishiyama, A.; Sugii, N.; Natori, K.; Hattori, T.; Iwai, H. The effect of isotropic and anisotropic scattering in drain region of ballistic channel diode. In Proceedings of the 2010 10th IEEE International Conference on Solid-State and Integrated Circuit Technology, Shanghai, China , 1–4 November 2010; pp. 1247–1249. [CrossRef] 9. Wu, J.; Wang, D.; Liu, H.; Lau, W.M.; Liu, L.M. An ab initio study of TiS3: A promising electrode material for rechargeable Li and Na ion batteries. RSC Adv. 2015, 5, 21455–21463. [CrossRef] 10. Barawi, M.; Flores, E.; Ferrer, I.J.; Ares, J.R.; Sánchez, C. Titanium trisulphide (TiS3) nanoribbons for easy hydrogen photogeneration under visible light. J. Mater. Chem. A 2015, 3, 7959–7965. [CrossRef] 11. Island, J.O.; Buscema, M.; Barawi, M.; Clamagirand, J.M.; Ares, J.R.; Sánchez, C.; Ferrer, I.J.; Steele, G.A.; van der Zant, H.S.J.; Castellanos-Gomez, A. Ultrahigh Photoresponse of Few-Layer TiS3 Nanoribbon Transistors. Adv. Opt. Mater. 2014, 2, 641–645. [CrossRef] 12. Island, J.O.; Barawi, M.; Biele, R.; Almazán, A.; Clamagirand, J.M.; Ares, J.R.; Sánchez, C.; van der Zant, H.S.J.; Álvarez, J.V.; D’Agosta, R.; et al. TiS3 Transistors with Tailored Morphology and Electrical Properties. Adv. Mater. 2015, 27, 2595–2601. [CrossRef] [PubMed] 13. Molina-Mendoza, A.J.; Barawi, M.; Biele, R.; Flores, E.; Ares, J.R.; Sánchez, C.; Rubio-Bollinger, G.; Agraït, N.; D’Agosta, R.; Ferrer, I.J.; et al. Electronic Bandgap and Exciton Binding Energy of Layered Semiconductor TiS3. Adv. Electron. Mater. 2015, 1, 1500126. [CrossRef] 14. Ferrer, I.; Ares, J.; Clamagirand, J.; Barawi, M.; Sánchez, C. Optical properties of titanium trisulphide (TiS3) thin films. Thin Sol. Films 2013, 535, 398–401. [CrossRef] 15. Guilmeau, E.; Berthebaud, D.; Misse, P.R.; Hebert, S.; Lebedev, O.I.; Chateigner, D.; Martin, C.; Maignan, A. ZrSe3-Type variant of TiS3: Structure and thermoelectric properties. Chem. Mat. 2014, 26, 5585–5591. [CrossRef] 16. Hsieh, P.L.; Jackson, C.; Grüner, G. Disorder effects in the linear chain compound TiS3. Sol. Stat. Commun. 1983, 46, 505–507. [CrossRef] 17. Zhang, J.; Liu, X.; Wen, Y.; Shi, L.; Chen, R.; Liu, H.; Shan, B. Titanium Trisulfide Monolayer as a Potential Thermoelectric Material: A First-Principles-Based Boltzmann Transport Study. ACS Appl. Mater. Interfaces 2017, 9, 2509–2515. [CrossRef] 18. Carrete, J.; Li, W.; Lindsay, L.; Broido, D.A.; Gallego, L.J.; Mingo, N. Physically founded phonon dispersions of few-layer materials and the case of borophene. Mater. Res. Lett. 2016, 4, 204–211. [CrossRef] 19. Kresse, G.; Hafner, J. Ab initio molecular dynamics for liquid metals. Phys. Rev. B 1993, 47, 558–561. [CrossRef] 20. Kresse, G.; Hafner, J. Ab initio molecular-dynamics simulation of the liquid-metal–amorphous-semiconductor transition in germanium. Phys. Rev. B 1994, 49, 14251–14269. [CrossRef] 21. Kresse, G.; Furthmüller, J. Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set. Comp. Mat. Sci. 1996, 6, 15–50. [CrossRef] 22. Kresse, G.; Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B 1996, 54, 11169–11186. [CrossRef] 23. Perdew, J.P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865–3868. [PubMed] 24. Blöchl, P.E. Projector augmented-wave method. Phys. Rev. B 1994, 50, 17953. 25. Kresse, G.; Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B 1999, 59, 1758–1775. [CrossRef] 26. Furuseth, S. L. Brattas and A. Kjekshus. Acta Chem. Scand. A 1975, 29, 623. 27. Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 2010, 132, 154104. doi:10.1063/1.3382344. [CrossRef] 28. Togo, A.; Tanaka, I. First principles phonon calculations in materials science. Scr. Mater. 2015, 108, 1–5. 29. Li, W.; Mingo, N.; Lindsay, L.; Broido, D.A.; Stewart, D.A.; Katcho, N.A. Thermal conductivity of diamond nanowires from first principles. Phys. Rev. B 2012, 85, 195436. 30. Li, W.; Carrete, J.; Katcho, N.A.; Mingo, N. ShengBTE: A solver of the Boltzmann transport equation for phonons. Comp. Phys. Commun. 2014, 185, 1747–1758. Nanomaterials 2020, 10, 704 10 of 10 31. Carrete, J.; Vermeersch, B.; Katre, A.; van Roekeghem, A.; Wang, T.; Madsen, G.K.; Mingo, N. almaBTE: A solver of the space-time dependent Boltzmann transport equation for phonons in structured materials. Comput. Phys. Commun. 2017, 220, 351–362. 32. Li, W.; Lindsay, L.; Broido, D.A.; Stewart, D.A.; Mingo, N. Thermal conductivity of bulk and nanowire Mg2 Six Sn1− x alloys from first principles. Phys. Rev. B 2012, 86, 174307. 33. Torres, P.; Torelló, A.; Bafaluy, J.; Camacho, J.; Cartoixà, X.; Alvarez, F.X. First principles kinetic-collective thermal conductivity of semiconductors. Phys. Rev. B 2017, 95, 165407. 34. Tamura, S. Isotope scattering of dispersive phonons in Ge. Phys. Rev. B 1983, 27, 858–866. 35. Torres, P.; Alvarez, F.X.; Cartoixà, X.; Rurali, R. Thermal conductivity and phonon hydrodynamics in transition metal dichalcogenides from first-principles. 2D Mater. 2019, 6, 035002. 36. Kong, B.D.; Paul, S.; Nardelli, M.B.; Kim, K.W. First-principles analysis of lattice thermal conductivity in monolayer and bilayer graphene. Phys. Rev. B 2009, 80, 033406. 37. Li, W.; Carrete, J.; Mingo, N. Thermal conductivity and phonon linewidths of monolayer MoS2 from first principles. Appl. Phys. Lett. 2013, 103, 253103. c 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).