Anisotropic thermal expansion of silicon monolayer in biphenylene network

Abstract

Materials with a negative thermal expansion property are of great importance in the emerging family of two-dimensional materials. For example, mixing two materials with negative and positive coefficients of thermal expansion avoids volume changing with temperature. In this work, based on first-principles calculations and Grüneisen's theory, we investigated the thermal expansion properties of a silicon monolayer in biphenylene networks. Our results show that the thermal expansion is greatly negative and anisotropic, as the linear thermal expansion coefficient along the a-direction is significantly smaller than the one along the b-direction, even at high temperatures. At 300 K, the thermal expansion coefficients along the two lattice directions are −17.010 × 10⁻⁶ K⁻¹ and −2.907 × 10⁻⁶ K⁻¹, respectively. By analyzing the Grüneisen parameters and the elastic compliance, we obtained an understanding of the giant negative thermal expansion of the material. Rigid unit modes are also responsible for the negative thermal expansion behavior. Our work provides fundamental insights into the thermal expansion of silicon monolayer in biphenylene networks and should stimulate the further exploration of the possible thermoelectric and thermal management applications of the material.

---

Introduction

Since the successful exfoliation of graphene, two-dimensional (2D) materials have gained widespread attention and become one of the most active areas of nanomaterials research.^1–3 Among them, the 2D group-IV materials in the same family as C are also developing rapidly due to their numerous remarkable physical properties and their great application potential.^4–8 For example, the quantum Hall effect, valley polarization, tunable band gap and fast endowment mobility have been discovered in silicene and germanene with buckled hexagonal lattices.^1,9–15 This indicates that they have great promise for applications in electronics, optics and fundamental physics research.

With the advancement of 2D materials research, quasi-one-dimensional (1D) carbon nanoribbons were synthesized in a biphenylene network on the surface of Au(111) by Fan et al.¹⁶ consisting of adjacent octagonal (o), hexagonal (h), and square (s) rings, called the ohs structure. It is a great achievement to realize a metallic, non-hexagonal, nanoscale monolayer of carbon proposed/predicted in various previous theoretical studies,^17–20 paving the way for the exploration of new two-dimensional allotropes and their properties. The available studies show that the C ohs monolayer has excellent stability, electronic properties and thermal transport properties,^21,22 and its thermal transport properties can be suppressed by hydrogenation.²³ The C ohs monolayer is predicted to have an ultrahigh melting point up to 4500 K,²⁴ negative thermal expansion (NTE)²⁵ and a negative differential resistance effect.²⁶ Moreover, C ohs monolayer is a good catalyst for oxygen reduction reaction,²⁷ an effective gas sensor for NO and NO₂,²⁶ and a promising anode material for high-performance sodium-ion batteries.²⁸ The superconducting critical temperature of the C ohs monolayer can be increased from 0.03 K to 11.22 K by adsorption of Li atoms, and the adsorption of alkaline metal atoms increases the transport anisotropy by more than two times.²⁹ These findings suggest that C ohs have great potential applications, which prompted us to turn our research attention to the counterpart of other 2D group-IV materials in the same family.

As a cognate of carbon, the Si ohs monolayer has also attracted the attention of researchers. Salih et al. showed that the Si ohs monolayer is metallic and can be stacked to form vertical heterostructures, which can form metallic junctions and thus be applied to diodes.³⁰ Tylan et al. showed that by dehydrogenation or deoxygenation stripes on Si ohs monolayers that are completely covered by H or O, various lateral composite structures can be constructed. And these structures consist of commensurately bare and atom-covered regions of ohs monolayers. Moreover, depending on the pattern and size of these regions, novel materials with unusual electronic and magnetic properties can be constructed, such as alloys, heterostructures, quantum wells, antipoints, and anticyclic rings.³¹

It is well known that the thermal expansion properties of materials are crucial. While the temperature changes, the accumulated thermal stresses and strains remarkably affect the performance and life of almost all devices. Mismatched thermal expansion of multiple materials may lead to severe device damage such as wire breakage and interface spalling. For example, Si is recognized as the next-generation anode material for lithium-ion batteries due to its extremely high theoretical capacity, low working potential and abundant natural abundance.³² Whereas, particle fragmentation due to volume expansion and crushing of Si particles are the major drawbacks to its commercialization.^32–34 Therefore, pre-growing void space³⁵ or adopting a finer particle structure³⁶ is one of the ways to buffer volume expansion. However, it is even more crucial to select materials with appropriate expansion coefficients.

Although most materials exhibit positive thermal expansion when heated due to atomic vibrations,^37–40 such as conventional single-crystal silicon with a linear coefficient of thermal expansion of about 2.6 × 10⁶ K⁻¹ at room temperature,^41,42 NTE has been reported for many 2D materials.^43,44 For example, theoretical calculations predicted negative linear thermal expansion coefficients (TECs) for graphene, silicene, germanene, and phosphene,^45–47 and several theories have been proposed to reveal possible physical mechanisms. The special temperature-responsive behavior of NTE, which violates the common sense of the “thermal expansion and cold contraction” effect, is of great scientific and technological significance. Some amount of NTE can facilitate the scaffolding of different materials into multicomponent structures, as NTE tends to counteract the usually harmful positive expansion and helps reduce thermal strain.⁴⁸ Control of thermal expansion also expands the possibilities of designing nanoscale devices that can be expanded/contracted across a required temperature range. Therefore, materials with NTE properties are of great importance in the emerging family of 2D materials, and accurate measurement and prediction of TECs will also provide an important reference for the design of related functional devices. However, studies on the TECs of Si ohs monolayer structures are still lacking.

In this paper, we investigate the TECs of the Si ohs monolayer by first-principles calculations based on the Grüneisen theory, which can save a significant amount of computational time compared to routine quasi-harmonic approach (QHA), solving thermal expansion by the direct minimization of the total free energy. The results show that, similar to C ohs,⁴⁹ the Si ohs monolayer contracts up to 800 K and has a higher thermal expansion anisotropy than C ohs. In addition, temperature-dependent lattice constants, generalized mode and macroscopic Grüneisen parameters, and the effect of elastic compliance tensor on TECs were investigated. It is found that the large and anisotropic NTE of the Si ohs monolayer over the entire temperature range is attributed to its large and anisotropic value of its macroscopic Grüneisen parameter (G). Rigid unit modes are also responsible for thermal expansion behavior. The revealed thermal expansion properties contribute to a further understanding of the study of thermal expansion and transport properties of 2D materials.

Computational and theoretical methods

In this work, we used Grüneisen's theory^7,49 to study the linear TECs of Si ohs monolayer. This method can save a lot of computational time, showing the efficiency and validity for anisotropic materials, compared to the routine QHA of direct minimization of the total free energy.^50–53 Specifically, in the routine QHA method, a series of free energies calculations need to be performed on a grid of lattice parameter points, the dimension of which is determined by the number of independent lattice parameters.^4,54 This implies the calculations of dozens or even hundreds of volumes of phonon spectra, requiring huge computational resources and time consumption. Moreover, it is well known that ZA modes are soft near the Γ point and their frequencies may become negative at larger strains. However, in order to accurately fit the energy data points to the equation of state, these data points should span a considerable wide energy range.⁵⁵ This contradiction then affects the accuracy and validity of the method for two-dimensional anisotropic materials.⁵⁵ Even with a much smaller computational cost, Grüneisen's theory can obtain good results of thermal expansion properties and has a wide range of applications.^43,51,56,57 Therefore, we have investigated linear TECs for Si ohs monolayers using Grüneisen's theory.

The linear TECs of 2D materials with a rectangular lattice can be expressed as^49,57:

where α₁ and α₂ denote the linear TECs along the a- and b-directions, respectively. C_V is the constant volume specific heat. A₀ is the area of the primitive cell and S_ij is the matrix element of the elastic compliance tensor, which is the inverse of the elastic stiffness tensor:where C_ij is the matrix element of the elastic stiffness tensor. In addition, G(a) denotes the macroscopic Grüneisen parameter along a-direction, which is the average of the generalized mode Grüneisen parameter γ_λ(a) along the same direction, weighted by c_λ (c_λ is the contribution of phonons with angular frequency ω_λ to the constant volume specific heat):where N is the total number of wave vectors in Brillouin zone (BZ).

Moreover, the generalized mode Grüneisen parameter γ_λ(a) can be calculated as:

where a₀ is the equilibrium lattice constant.

Our calculations were performed using the Vienna Ab initio Simulation Package (VASP)^58,59 based on density generalized function theory. The projector-augmented-wave (PAW)^60,61 pseudopotentials and Perdew–Burke–Ernzerhof (PBE)⁶² exchange–correlation functionals were chosen. Based on the existing structural model of the C ohs monolayer, the structural model of the Si ohs monolayer is determined by replacing C atoms with Si atoms. Unlike the planar structure of the C ohs monolayer, the Si ohs monolayer has a buckled structure. Therefore, a random perturbation was added to the z-axis of the lattice using VASPKIT⁶³ before the structure optimization. To ensure the convergence, the plane-wave energy cutoff was set to 500 eV, much higher than the recommended value. A vacuum space of 18 Å was taken in the direction perpendicular to the crystal plane to prevent the interaction of two adjacent layers. In all calculations, the energy convergence was 10⁻⁶ eV and the maximal residual Hellmann–Feynman forces were reduced to 10⁻⁴ eV Å⁻¹. And a 4 × 6 × 1 k-mesh was used during the structural relaxation. The phonon dispersions of Si ohs monolayer is obtained using the finite displacement method implemented in the PHONOPY code.⁶⁴ To obtain the convergent phonon properties in the calculation of the harmonic interatomic force constants, a 5 × 5 × 1 supercell is used. After convergence tests, a 51 × 51 × 1 q-mesh is used for the calculations of harmonic properties.

Results and discussions

The optimized Si ohs monolayer is shown in Fig. 1(a) and (b), consisting of octagons, hexagons and squares. The primitive cell of Si ohs structure contains six silicon atoms. The calculated equilibrium lattice constants are a = 7.13 Å, b = 5.73 Å, and the total buckling height is Δh = 1.04 Å, which are in good agreement with the previous work.³⁰ The Si–Si bond lengths are d₁ = 2.323 Å, d₂ = 2.285 Å, d₃ = 2.292 Å and d₄ = 2.287 Å as shown in Fig. 1. Based on the optimized structure, the phonon dispersion of Si ohs calculated using the PHONOPY code is shown in Fig. 1(c). As there are six silicon atoms in the primitive cell, eighteen phonon branches in the phonon spectra, including three phonon branches and fifteen photon branches. The absence of negative frequencies in the phonon spectrum indicates that the structure of the material is dynamically stable.

Fig. 1 (a) Top view and (b) side view of Si ohs monolayer, where the rectangle represents the unit cell with a and b being the lattice constants. (c) Phonon dispersion of Si ohs monolayer. The coordinates of the high-symmetry points in the Brillouin zone are Γ (0,0,0), X (0.5,0,0), L (0.5,0.5,0), and Y (0,0.5,0).

Essentially, the Grüneisen parameters are directly accountable for the TECs of the material. For 2D materials, the ZA mode is so soft that a slight change in the lattice constant may lead to a negative phonon frequency near the Γ point. This implies that the applied strain should be small enough. Therefor we expanded the lattice by 1% along the a- and b-directions, respectively, and used the same parameter settings to calculate the phonon dispersion of the Si ohs monolayer after stretching the lattice. Moreover, there are no negative frequencies in any of the three phonon spectra, ensuring the availability of the data.

The calculated frequency-dependent generalized mode Grüneisen parameters (γ) are shown in Fig. 2. The γ of the Si ohs monolayer shows a clear anisotropy. It can be seen from the figure that the overall trends of γ in both directions are very similar, but the distribution range of γ in the a-direction is much wider than that in the b-direction. Moreover, in the low frequency region, the γ values of phonons in the a-direction are much lower than those in the b-direction. And not only that, in the range of 3–6 THz, the γ values in the a-direction and b-direction also show obvious differences. Along the b-direction, it shows greater positive values of γ. In the high frequency region, γ changes from negative to positive values. Generally speaking, the negative value of γ indicates NTE. It means that the negative value of γ is the reason for the NTE of Si ohs monolayer, while its anisotropic thermal expansion is affected by the anisotropy of γ.

Fig. 2 Frequency-dependent generalized mode Grüneisen parameters for the strain along the a- and b-directions for Si ohs monolayer. The insets show a large range of mode Grüneisen parameters along the corresponding directions.

The variation of G with temperature is shown in Fig. 3, and the inset shows the whole range of G. It is worth noting that the G are remarkably negative along both the a- and b-directions in the range 0–800 K. The reason is that a large number of phonons with negative γ are distributed over a rather wide range of frequencies, and the absolute values of them are large, up to 10² orders of magnitude, while the positive values are below 10. At low temperatures, the predominantly excited phonon modes are low-frequency phonon modes with a large negative γ. As the temperature increases, high-frequency phonon modes, most of which have positive γ, are gradually excited. Thus, G(a) and G(b) increase rapidly with increasing temperature and reach their saturation values around 200 K. The saturation value of G along the a-direction is −1.161, which is lower than that of −0.594 along the b-direction. Moreover, in general, the range of G in the a-direction is much wider than that in the b-direction, consistent with the case of the γ.

Fig. 3 Variation of the macroscopic Grüneisen parameters with increasing temperature, along the a- and b-directions. The inset shows the results for the overall range.

In order to clearly quantify the distribution of the Grüneisen parameters, we introduced the mode Grüneisen parameters density for uniaxial deformation along different directions.⁵⁰ The density of the γ is expressed as ⁵⁰ and displayed in Fig. 4(a). Note that g₁ and g₂ are the average densities of the γ along the a- and b-directions, respectively. As can be seen from the figure, the main distribution of γ ranges from −5 to 6. In the a-direction, the density reaches a maximum value of 0.57 when γ is 0.4, and γ reaches a maximum density value of 0.46 when γ is 0.5 in the b-direction. In addition, the percentages of phonons with positive γ values are more than 50% in both directions. In fact, we find that phonons with negative γ values occupy only 34.8% and 34.4% along the a-direction and b-direction, respectively. Combining with Fig. 2, we can notice that although the percentage of phonons with negative γ values is relatively low, the phonons with negative γ are dominant in thermal expansion.

Fig. 4 (a) The mode Grüneisen parameters density for Si ohs monolayer. (b) The frequency dependence of the mode Grüneisen parameters.

To show the frequency-dependent information of the γ, we introduced Γ, the frequency-resolved phonon density of states weighted by the Grüneisen parameters, defined as ⁵⁰ Here, Γ₁ and Γ₂ also correspond to the strains applied along the a- and b-directions. And D denotes the width of the interval in which the frequency is equally divided into multiple intervals. The frequency dependence information of the Γ is shown in Fig. 4(b). It can be seen that in the low frequency region (below 3 THz), the Γ values are all negative. Above 3 THz, the vast majority of Γ values are positive. It is also found the absolute value of negative Γ is larger than the one of positive Γ significantly, resulting in the negative values of the macroscopic Grüneisen parameters. And in general, the Γ values in the b-direction are generally larger than those in the a-direction, which is consistent with the same trend of the macroscopic Grüneisen parameters.

In Fig. 5(a), we give the calculated ratios of the lattice constants a and b to the values of 0 K at different temperatures. As the temperature increases, both a and b are compressed, and a compressed much more than b. The ratios of the lattice constants are anisotropic like γ and G, and the anisotropy increases with the increase of temperature. From 0 to 800 K, it can be seen that the lattice parameters of the Si ohs monolayer shrink up to about 1.5 percent, which is greater than the 0.8 percent of C ohs⁴⁹ and exhibits stronger anisotropy. And more precisely, at 300 K, the lattice constants a and b contract by 0.54% and 0.16%, respectively, and the percentage contraction increases to 1.36% and 0.28% at 800 K, respectively.

Fig. 5 (a) The ratios of lattice constants for Si ohs monolayer. (b) Linear TECs for Si ohs monolayer.

Fig. 5(b) shows the linear TECs with increasing temperature. As can be seen from the figure, the linear TECs of Si ohs monolayer decrease at first, then increases and finally stabilizes with increasing temperature. The linear TECs of many other 2D materials have the similar trend, such as graphene,^45,65,66 hexagonal silicene,^45,47 and C ohs monolayer.⁴⁹ It is notable that the Si ohs monolayer exhibits giant NTE and anisotropy. At 300 K, the linear TECs reach −17.0 × 10⁻⁶ K⁻¹ and −2.9 × 10⁻⁶ K⁻¹ along a- and b-directions. For comparison, the linear TECs of C ohs are −11 × 10⁻⁶ K⁻¹ and −9 × 10⁻⁶ K⁻¹ along a- and b-directions.⁴⁹ It is found the linear TECs of Si ohs monolayer are much lower in a-direction but higher in b-direction compared with C ohs monolayer. Thus, it shows more significant anisotropy. This significant in-plane anisotropy can be used for the thermal management of nanoelectronic devices and may attract more attention to orientation-dependent thermal devices.^67,68

Here, we present the thermal expansion data of several related materials in Table 1. It is found that unlike the isotropic thermal contraction of graphene and silicene, the Si ohs monolayer and the C ohs monolayer are anisotropic. As shown in eqn (1) and (2), this is attributed to the anisotropy of their Grüneisen parameters, as the elastic compliance tensor S of Si and C ohs monolayers are nearly isotropic. The highest values of G for Si ohs monolayer in the a- and b-directions are −1.2 and −0.6 respectively, leading to the high anisotropy of the monolayer. And the absolute values of Grüneisen parameters for the Si ohs monolayer are much smaller than those of C ohs monolayer (−3 and −2.5 in the a and b directions, respectively). The much more significant anisotropy of G also results in more anisotropy of thermal expansion for Si ohs monolayer. Although the absolute values of G for the Si ohs monolayer are smaller, the S_ij are about 6.5–8.5 times larger than those of C ohs monolayer, finally leading to more remarkable NTE. Thus, it can be included that the elastic properties have an important effect on the thermal expansion, and soft materials may possess more remarkable positive/negative thermal expansion. Our calculations show that the S₁₁ and S₂₂ are positive, while S₁₂ is negative. And the S₁₁ and S₂₂ are 2.38 and 2.47 times the absolute value of S₁₂, respectively. The negative linear TECs are mainly determined by G. Both G(a) and G(b) are negative throughout the temperature region, which give rise to the negative linear TECs.

Table 1 Elastic compliance tensor and thermal expansion data for graphene, hexagonal silicene, C ohs monolayer, and Si ohs monolayer. Note that the superscripts of a and b here indicate the data along a-direction and b-direction, respectively

Material	Elastic constants (N m^−1)		Sij (m N^−1)		α (×10^−6 K^−1)	Reference
Graphene	C11 = C22 = C12	312	—		−3.7	55 and 69
Hexagonal silicene	C11 = C12	68.9	S11	−0.0074	−5.3	55
	C22	23.3	S22 = S12	−0.0219
C ohs monolayer	C11	283.2	S11	0.0040	−11(a), −9(b)	49
	C22	236.0	S22	0.0048
	C12	88.5	S12	−0.0015
Si ohs monolayer	C11	39.9	S11	0.0302	−17(a), −3(b)	This work
	C22	38.4	S22	0.0314
	C12	16.1	S12	−0.0127

---

In rigid unit mode theory, the NTE is also related to the vibrational modes of the rigid unit.^49,70 So, we examined the real-space vibrations of all phonon modes at the Γ point. Fig. 6 shows the rigid cell modes (square and hexagonal rings) observed at ω = 3.09 THz and ω = 5.19 THz. We can find that the rigid unit modes are significantly rotated. During the rotation of rigid units such as the square/hexagonal ring in Fig. 4(a) and (b), other rings are deformed to fill the empty space, causing the strong shrinkage of the material. This phenomenon also occurs in C ohs monolayer.⁴⁹

Fig. 6 Real-space visualization of typical rigid unit modes observed in Si ohs monolayer with (a) ω = 3.09 THz, and (b) ω = 5.19 THz. Red circles indicate locally rigid units. The equilibrium structure is shown in the left box, while the right box depicts a snapshot of rigid units rotated clockwise.

Conclusions

In summary, we have investigated the thermal expansion properties of Si ohs monolayer through first-principles calculations. Our calculations show that the Si ohs monolayer exhibits giant NTE behavior in the temperature range of 0 to 800 K. Moreover, the absolute value of the linear TECs in the a-direction of the Si ohs monolayer is much larger than that of the C ohs monolayer over the entire temperature range, while the one in the b-direction is smaller than that of the C ohs monolayer. The linear TECs of Si ohs monolayer at 300 K are −17.0 × 10⁻⁶ K⁻¹ and −2.9 × 10⁻⁶ K⁻¹ along the a- and b-directions, respectively. The linear TECs along the a-direction is approximately 6 times greater than that along the b-direction. The in-plane thermal expansion exhibits a much stronger anisotropy than the C ohs monolayer. Our investigations suggest that the large negative value of the G and the rigid unit modes are responsible for its strong NTE, and the anisotropy of G mainly leads to anisotropic thermal expansion in the Si ohs monolayer. Furthermore, the elastic compliance is studied, confirming that soft materials usually possess significant thermal expansion. Our results on the NTE of Si ohs monolayer enrich the range of NTE 2D materials and will have a significant impact on the fundamental understanding and potential applications of Si ohs monolayer in electronic devices.

Data availability

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

Author contributions

Aiqing Guo contributed to the writing of one original draft. Gang Liu, Fengli Cao, and Xiaodong Qiu contributed to the writing, review, and editing of this article. Weiwei Ju and Zhibin Gao contributed to the conceptualization and supervision. All authors contributed to this article and approved the submitted version.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (Grant No. 11974100 and 61874160).

References

1. S. Cahangirov, M. Topsakal, E. Aktürk, H. Şahin and S. Ciraci, Phys. Rev. Lett., 2009, 102, 236804.
2. K. S. Novoselov, A. K. Geim, S. V. Morozov, D.-e. Jiang, Y. Zhang, S. V. Dubonos, I. V. Grigorieva and A. A. Firsov, Science, 2004, 306, 666–669.
3. A. K. Geim and K. S. Novoselov, Nat. Mater., 2007, 6, 183–191.
4. N. Mounet and N. Marzari, Phys. Rev. B: Condens. Matter Mater. Phys., 2005, 71, 205214.
5. N. Liu, G. Bo, Y. Liu, X. Xu, Y. Du and S. X. Dou, Small, 2019, 15, 1805147.
6. A. Acun, L. Zhang, P. Bampoulis, M. v. Farmanbar, A. van Houselt, A. Rudenko, M. Lingenfelder, G. Brocks, B. Poelsema and M. Katsnelson, J. Phys.: Condens.Matter, 2015, 27, 443002.
7. P. Schelling and P. Keblinski, Phys. Rev. B: Condens. Matter Mater. Phys., 2003, 68, 035425.
8. K. Takeda and K. Shiraishi, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 14916.
9. X. Li, J. T. Mullen, Z. Jin, K. M. Borysenko, M. B. Nardelli and K. W. Kim, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 87, 115418.
10. Z. Ni, Q. Liu, K. Tang, J. Zheng, J. Zhou, R. Qin, Z. Gao, D. Yu and J. Lu, Nano Lett., 2012, 12, 113–118.
11. N. Drummond, V. Zolyomi and V. Fal'Ko, Phys. Rev. B: Condens. Matter Mater. Phys., 2012, 85, 075423.
12. W.-F. Tsai, C.-Y. Huang, T.-R. Chang, H. Lin, H.-T. Jeng and A. Bansil, Nat. Commun., 2013, 4, 1500.
13. M. Ezawa, Phys. Rev. Lett., 2012, 109, 055502.
14. C.-C. Liu, W. Feng and Y. Yao, Phys. Rev. Lett., 2011, 107, 076802.
15. H. Şahin, S. Cahangirov, M. Topsakal, E. Bekaroglu, E. Akturk, R. T. Senger and S. Ciraci, Phys. Rev. B: Condens. Matter Mater. Phys., 2009, 80, 155453.
16. Q. Fan, L. Yan, M. W. Tripp, O. Krejčí, S. Dimosthenous, S. R. Kachel, M. Chen, A. S. Foster, U. Koert and P. Liljeroth, Science, 2021, 372, 852–856.
17. A. Balaban, in Symmetry 2, Elsevier, 1989, pp. 397–416.
18. X.-Q. Wang, H.-D. Li and J.-T. Wang, Phys. Chem. Chem. Phys., 2013, 15, 2024–2030.
19. O. Rahaman, B. Mortazavi, A. Dianat, G. Cuniberti and T. Rabczuk, FlatChem, 2017, 1, 65–73.
20. M. A. Hudspeth, B. W. Whitman, V. Barone and J. E. Peralta, ACS Nano, 2010, 4, 4565–4570.
21. A. Bafekry, M. Faraji, M. Fadlallah, H. Jappor, S. Karbasizadeh, M. Ghergherehchi and D. Gogova, J. Phys.: Condens. Matter, 2021, 34, 015001.
22. H. P. Veeravenkata and A. Jain, Carbon, 2021, 183, 893–898.
23. P. Zhang, T. Ouyang, C. Tang, C. He, J. Li, C. Zhang, M. Hu and J. Zhong, J. Mater. Chem. C, 2021, 9, 16945–16951.
24. Y. Luo, C. Ren, Y. Xu, J. Yu, S. Wang and M. Sun, Sci. Rep., 2021, 11, 19008.
25. B. Mortazavi and A. V. Shapeev, FlatChem, 2022, 32, 100347.
26. G. Liu, T. Chen, X. Li, Z. Xu and X. Xiao, Appl. Surf. Sci., 2022, 599, 153993.
27. T. Liu, Y. Jing and Y. Li, J. Phys. Chem. Lett., 2021, 12, 12230–12234.
28. T. Han, Y. Liu, X. Lv and F. Li, Phys. Chem. Chem. Phys., 2022, 24, 10712–10716.
29. Y. Mu and S.-D. Li, Appl. Surf. Sci., 2023, 628, 157255.
30. S. Demirci, Ş. Çallıoğlu, T. Görkan, E. Aktürk and S. Ciraci, Phys. Rev. B, 2022, 105, 035408.
31. T. Gorkan, S. a. Çallıoğlu, S. Demirci, E. Aktürk and S. Ciraci, ACS Appl. Electron. Mater., 2022, 4, 3056–3070.
32. Z. Chen, A. Soltani, Y. Chen, Q. Zhang, A. Davoodi, S. Hosseinpour, W. Peukert and W. Liu, Adv. Energy Mater., 2022, 12, 2200924.
33. S. Chae, M. Ko, K. Kim, K. Ahn and J. Cho, Joule, 2017, 1, 47–60.
34. K. Feng, M. Li, W. Liu, A. G. Kashkooli, X. Xiao, M. Cai and Z. Chen, Small, 2018, 14, 1702737.
35. W. An, P. He, Z. Che, C. Xiao, E. Guo, C. Pang, X. He, J. Ren, G. Yuan and N. Du, ACS Appl. Mater. Interfaces, 2022, 14, 10308–10318.
36. H. Li, Z. Chen, Z. Kang, W. Liu and Y. Chen, Energy Storage Mater., 2023, 56, 40–49.
37. L. F. Huang, P. L. Gong and Z. Zeng, Phys. Rev. B: Condens. Matter Mater. Phys., 2014, 90, 045409.
38. Y. Aierken, D. Çakır, C. Sevik and F. M. Peeters, Phys. Rev. B: Condens. Matter Mater. Phys., 2015, 92, 081408.
39. L.-F. Huang and Z. Zeng, J. Phys. Chem. C, 2015, 119, 18779–18789.
40. H. Sun, G. Liu, Q. Li and X. Wan, Phys. Lett. A, 2016, 380, 2098–2104.
41. H. Watanabe, N. Yamada and M. Okaji, Int. J. Thermophys., 2004, 25, 221–236.
42. D. S. Kim, O. Hellman, J. Herriman, H. Smith, J. Lin, N. Shulumba, J. Niedziela, C. Li, D. Abernathy and B. Fultz, Proc. Natl. Acad. Sci. U. S. A., 2018, 115, 1992–1997.
43. P. Anees, M. Valsakumar and B. Panigrahi, Phys. Chem. Chem. Phys., 2017, 19, 10518–10526.
44. C. Sevik, Phys. Rev. B: Condens. Matter Mater. Phys., 2014, 89, 035422.
45. S. Mann, R. Kumar and V. Jindal, RSC Adv., 2017, 7, 22378–22387.
46. L. Wang, C. Wang and Y. Chen, J. Phys.: Condens. Matter, 2019, 31, 465003.
47. L.-F. Huang, P.-L. Gong and Z. Zeng, Phys. Rev. B: Condens. Matter Mater. Phys., 2015, 91, 205433.
48. L.-F. Huang, X.-Z. Lu and J. M. Rondinelli, Phys. Rev. Lett., 2016, 117, 115901.
49. Q. Li, J. Zhou, G. Liu and X. Wan, Carbon, 2022, 187, 349–353.
50. G. Liu, J. Zhou and H. Wang, Phys. Chem. Chem. Phys., 2017, 19, 15187–15193.
51. B. Arnaud, S. Lebègue and G. Raffy, Phys. Rev. B, 2016, 93, 094106.
52. C. K. Gan and Y. Y. F. Liu, Phys. Rev. B, 2016, 94, 134303.
53. C. K. Gan, J. R. Soh and Y. Liu, Phys. Rev. B: Condens. Matter Mater. Phys., 2015, 92, 235202.
54. L.-F. Huang, X.-Z. Lu, E. Tennessen and J. M. Rondinelli, Comput. Mater. Sci., 2016, 120, 84–93.
55. X.-J. Ge, K.-L. Yao and J.-T. Lü, Phys. Rev. B, 2016, 94, 165433.
56. G. Liu and J. Zhou, J. Phys.: Condens. Matter, 2018, 31, 065302.
57. G. Liu, Z. Gao and J. Ren, Phys. Rev. B, 2019, 99, 195436.
58. G. Kresse and J. Furthmüller, Comput. Mater. Sci., 1996, 6, 15–50.
59. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169.
60. G. Kresse and D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758.
61. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953.
62. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865.
63. V. Wang, N. Xu, J.-C. Liu, G. Tang and W.-T. Geng, Comput. Phys. Commun., 2021, 267, 108033.
64. A. Togo and I. Tanaka, Scr. Mater., 2015, 108, 1–5.
65. W. Gao and R. Huang, J. Mech. Phys. Solids, 2014, 66, 42–58.
66. S. Sarikurt, Y. Z. Abdullahi, E. Durgun and F. Ersan, J. Phys. D: Appl. Phys., 2022, 55, 315303.
67. M. Zeraati, S. M. V. Allaei, I. A. Sarsari, M. Pourfath and D. Donadio, Phys. Rev. B, 2016, 93, 085424.
68. Z. Gao, Z. Zhang, G. Liu and J.-S. Wang, Phys. Chem. Chem. Phys., 2019, 21, 26033–26040.
69. E. Cadelano, P. L. Palla, S. Giordano and L. Colombo, Phys. Rev. Lett., 2009, 102, 235502.
70. M. G. Tucker, A. L. Goodwin, M. T. Dove, D. A. Keen, S. A. Wells and J. S. Evans, Phys. Rev. Lett., 2005, 95, 255501.

---