The “discretized” version of the two-dimensional “Generalized Skettrup Model” is presented and implemented at simulations on the effect of alterations in the aspect ratio(s) of rib lengths on temperature-dependent harmonic and anharmonic specific lattice thermal capacities of rectangular graphene nano-ribbons. The obtained simulation results are discussed in comparison with appropriate experimental data and their counterparts reported elsewhere for the square graphene flakes [V. Ligatchev, ECS J. Solid State Sci. Technol. 9, 093014 (2020); V. Ligatchev, “Generalized Skettrup Model” and Lattice Thermal Capacity of Spatially Non-Homogeneous and Low-Dimensional Semiconductors and Insulators. Lambert Academic Publishing, Saarbrücken, Germany, 2021].

Thermal properties of graphene and other two-dimensional (2D) solids as well as of 1D semiconductors and insulators with crystalline, amorphous, and amorphous–crystalline (e.g., “polymeric”) atomic structures emerged in recent decades as attractive branches of material and biomedical sciences due to unique size-dependent and temperature-dependent features of those low-dimensional materials.1–6 In particular, harmonic and anharmonic isobaric lattice heat capacity, Cp, thermal conductivity,7,8 and even melting temperature9,10 of many 2D nano-sized solids were found to be considerably size-dependent, while their low-temperature dependences deviate significantly from the famous Debye’s law: Cp(T) ∝ T3.11 Those unusual thermal properties make the low-dimensional solids to be attractive for implementation in physics,3 material science,7,8 biology, medicine, and so on.4,5

Different computational techniques (e.g., density-functional-perturbation-theory,12 Peierls-Boltzmann transport equation, and density functional theory7,8) are currently in use for fairly accurate “first-principles” evaluations on the key features of the phonon dispersion and spectra, isochoric lattice thermal capacity, and thermal conductivity of solids of different natures and spatial dimensionalities. Such techniques allow one to evaluate many important thermal characteristics of real solids with remarkable accuracy (often well comparable with accuracies of corresponding experimental techniques), and they are still computationally very demanding even for some elemental solids, like alpha- and beta-modifications of the boron, comprising 12 and 106 (!) atoms (respectively) in their unit cells.2,8,9,12 Therefore, quests for alternative (and less computationally-demanding) “first-principle” and semi-empirical approximations and techniques are still far from over.

Herein, the two-dimensional (2D) “discretized” version of the so-called “Generalized Skettrup Model” (2D GSM)1,2,6,13 is presented and implemented at quantitative evaluations on the temperature-dependent harmonic and anharmonic fractions of the lattice thermal capacity of the rectangular graphene nano-flakes (ribbons) with different aspect ratios of their orthogonal rib lengths.

Section II depicts basic equations of the “discretized” 2D GSM, whereas Sec. III of this article exhibits obtained simulation results. Section IV summarizes the outcomes of the entire work.

The basic equation(s) of 2D GSM1,2 links both single-particle (N ≡ 1) harmonic as well as many-particle (N = 2, 3, 4, …) anharmonic components of the specific lattice thermal capacity of a rectangular (in general) flake (of submicron- or micrometer-sized rib length) to features of the (idealized) quasi-continuous spectrum of its in-plane transverse and longitudinal thermal acoustic vibrations (isotropic TA and LA phonons); for more details, see Refs. 1 and 6 and sub-Sections 5.1 and 6.1 of Ref. 2. In particular, the aforementioned components of the isobaric specific lattice thermal capacity of 2D GSM1,2,6 read

(1a)

where Lx and Ly denote rib lengths of the rectangular (in general) flake (ribbon), C2 = [(kBθD)2/(LxLy)] * {LxLy/[2π2 (h ceff)2]}3/2 is a dimensional (of eV−1cm−2) coefficient, ET stands for the energy of confined multi-phonon state, Emin is its low limit, imposed by the flake’s morphology (and defined shortly below), and ceff denotes the effective sound velocity of the given 2D solid: 2/ceff = (1/cl + 1/ct); here cl and ct are its longitudinal and transverse sound velocities, respectively, whereas kB and θD symbolize the Boltzmann constant and Debye temperature, respectively.1,2

Furthermore, F2(Lx, Ly) stands for a dimensionless function, which is dependent entirely on (Lx/Ly) ratio,

(1b)

LxLy is presumed herein. The Emin energy in Eq. (1a) is defined as follows:1,2

(1c)

while the N-particle “statistical sum,” Z2DN, reads1,2

(1d)

It is noteworthy that Eqs. (1a)(1d) (as well as many other variants of the basic equations of the 2D GSM and especially those of the 3D GSM13) are based essentially on (Debye’s11) assumption on linear phonon dispersion(s) for the in-plane TA and LA phonons with the static plane-wave basis and quasi-continuous phonon spectrum, although the former assumption apparently not valid even for the in-plane acoustic phonons as their wavevectors approaching edges of corresponding Brillouin Zone (BZ), and especially for the “outplane” (or “flexural,” ZA) mode(s).1–3 

As it is just mentioned, Eq. (1a) ultimately implies the existence of the (quasi)-continuous vibrational (single-phonon) spectrum of the (2D) solid. Thus, it is rightfully valid for monolayers of large (as compared to inter-atomic distances) spatial extents only. However, spatial extents of real graphene monolayers and especially those of the nano-flakes and nano-ribbons are always limited by definition, which requires some modification in the (anisotropic in general) basic equation of 2D GSM.

In particular, the continuous ET quantity (integration variable) in Eq. (1a) has to be replaced with an appropriate set of “discretized” energy levels of pEmin, where p stands for an integer index (variable), ranging from 1 to pmax = floor(kBθD/Emin); the floor function routinely yields the highest integer quantity, which equals to or just below the (NkBθD/Emin) ratio.

In addition, a combination of “normalizing” C2 and Z2DN terms in Eq. (1a) yields the flake-size-independent harmonic fraction of the lattice heat capacity for the square flakes (i.e., for those of Lx = Ly = L), while their anharmonic fractions become size-dependent based on Eqs. (1b)(1d), and its dominant term (at N = 2) is inversely proportional to the amendment in the rib length of the square flake(s): CAH2(T) ∝ (1/L)1,2 at elevated temperatures.

However, it could be verified readily that for the rectangular flakes (ribbons) with a single variable rib length (e.g., at the variable Lx extent and fixed Ly one), Eq. (1a) does not yield the size-independent specific lattice thermal capacity anymore. Indeed, even its harmonic fraction varies at the ribbon size alterations, which is quite unusual for the specific thermal capacity of 2D solids of macroscopic spatial extents (sizes).

Thus, the basic equation of the 2D GSM for the lattice heat capacity of the (rectangular in general) flakes (ribbons) originating from the (confined) in-plane LA and TA phonons has to be modified appropriately in order to address all mentioned above drawbacks of Eq. (1a).

First of all, let us replace the quasi-continuous spectra of the in-plane LA and TA phonons with the more appropriate (for nano-ribbons) “discretized” one

(2a)

Thereafter, a combination of the Emin term [expressed by Eq. (1c) above] with Eq. (2a) would yield

(2b)

It is noteworthy that the latter equation does not contain the obsolete “normalizing” C2 and Z2DN terms anymore. These terms have to be omitted as well even for the square flakes but only within the framework of the “discretized” version of the 2D GSM, while both those C2 and Z2DN items have to be retained for the basic equation of the previous version(s) of the 2D GSM, presented, for instance, elsewhere in Refs. 1 and 2. As it will be shown shortly below, Eqs. (1b), (1c), and (2b) emerge in a size-independent “in-plane” harmonic fraction of the thermal lattice capacity of square flakes, while a similar fraction still remains (slightly) affected by the sizes of the rectangular flakes (ribbons).

Let us introduce the following dimensionless function:

(2c)

This would allow one to simplify (formally) the expression for the denominator in the last term of Eq. (2b). The expression for the latter function has some similarity with the expression for the F2(Lx, Ly) one, defined by Eq. (1b) above, although right-hand-side of Eq. (2c) is much simpler as compared to that of Eq. (1b) one.

Figure 1 shows the behavior of the F2(Lx, Ly) and Fxy(Lx, Ly) functions, as defined by Eqs. (1b) and (2c), respectively. As it is seen readily from the figure, both these functions depend considerably on the (Lx/Ly) ratio at a relatively low (Lx/Ly) fraction [1 ≤ (Lx/Ly) ≤ 5] only.

FIG. 1.

Dependence of the F2(Lx, Ly) and Fxy(Lx, Ly) functions on the (Lx/Ly) ratio evaluated by the numerical integration(s) of Eqs. (1b) and (2c), respectively.

FIG. 1.

Dependence of the F2(Lx, Ly) and Fxy(Lx, Ly) functions on the (Lx/Ly) ratio evaluated by the numerical integration(s) of Eqs. (1b) and (2c), respectively.

Close modal

In particular, the F2(Lx, Ly) function descending from (π/2) ≅ 1.5708 at Lx = Ly until F2(Lx, Ly) ≅ 1.0505 at (Lx/Ly) = 5, while the Fxy(Lx, Ly) one diminishes from 1 at Lx = Ly until ∼0.72111 within the same argument range (see Fig. 1). At higher (Lx/Ly) ratios, both the F2(Lx, Ly) and Fxy(Lx, Ly) functions diminish slowly, approaching unity in the (Lx/Ly) → ∞ limit in the former case and the 0.5√2 ≅ 0.707 quantity in such a limit for the latter case. As a result, the F2(Lx, Ly)/Fxy(Lx, Ly) ratio becomes weakly-size-dependent, descending from the ∼1.5708 at Lx = Ly until ∼1.4194 at (Lx/Ly) = 20: i.e., it drops just by ∼11% at the variation of a single rib length of the nano-ribbon in the range from the (Lx/Ly) = 1 until (Lx/Ly) = 20. Further “elongation” of a nano-ribbon would not affect the F2(Lx, Ly)/Fxy(Lx, Ly) ratio considerably since both these functions are “saturating” gradually at (Lx/Ly) → ∞ (see Fig. 1). However, for the square flakes (of Lx = Ly), both harmonic and anharmonic fractions of the “in-plane” lattice thermal capacity obviously become size-independent. In contrast, substantial size-dependent effects are expected to emerge for the anharmonic fractions of the lattice thermal capacity of the 2D GSM of the elongated flakes (ribbons) since—based on Eq. (2b)—the enlargement in the (Lx/Ly) ratio would yield an enhancement in [F2(Lx, Ly)/Fxy(Lx, Ly)]N one, “powered” by the integer indices N = 2, 3, 4….

Since the F2(Lx, Ly)/Fxy(Lx, Ly) ratio in the last term of Eq. (2b) is affected by the (Lx/Ly) fraction, both conventional (at N = 1) dimensionless vibrational density-of-states (VDOS) function, expressed by that last term, and its many-particle counterpart(s) (at N = 2, 3, 4, …),2,13 as well as entire Eq. (2b) become anisotropic in general. Thus, both harmonic and anharmonic contributions to the lattice thermal capacity expressed by Eq. (2b) turn out to be anisotropic for elongated nano-ribbons but remain isotropic for square nano-flakes.

Furthermore, the 2D and 3D versions of the GSM might take into account possible interactions among (at least partially occupied) “in-plane” multi-phonon states, via implementation of the two “model parameters,” r1 and r2,2,13,14

(2d)

In the latter case, those model parameters are “embedded” into the basic equation of the 2D GSM, initially expressed by Eq. (2a) above in this section.

As it can be seen readily from the very last term in the square brackets of Eq. (2d), an amendment in the r1 quantity would directly affect the effective number of quasi-particle in the in-plane ensemble. Subsequently, based on the Fock space formalism,15 which is expected to be relevant both for ensembles of interacting and/or non-interacting particles,14,15 the total number of the (quantum) states available for a multi-phonon system would be amended as well when the phonon–phonon interactions are “switched on”: see, for instance, a brief discussion on this issue in Appendix D of Ref. 13 and the third section of Ref. 2. This implies that the existence of (different kinds of) interactions between the LA and/or TA phonons might be (in principle) reflected in an enlargement (attained at (r1 < 1), for the case of the 2D GSM,1,2 and at (r1 < 2), for the 3D one13) in the effective number(s) of the TA and LA phonons confined within the 3D crystallite, and/or square or rectangular 2D flake, and associated with enlarged numbers of multi-phonon states (of the given aggregate energy) available for the appropriate phonon ensembles2,13,14 For example, any alteration in the quantity of the first model parameter, r1, in Eq. (2d) would rather amend a number of the many-phonon states (corresponding to N = 2, 3, 4, …), potentially involved in the phonon–phonon interactions, but not the number of the fundamental states (corresponding to the N ≡ 1) of the “discretized” version of the 2D GSM.2 

Although the approach depicted just above allows one to incorporate (formally) interactions among the in-plane LA and TA phonons and makes the whole framework of the 2D GSM to be more appropriate for quantitative description of the thermal characteristics of real 2D solids, it is actually comprised of the (mentioned above) model parameters r1 and r2, which are difficult to define based entirely on the “first principles.” Furthermore, in general, the r1 and r2 quantities could be treated even as temperature-dependent ones,14,16 although their actual temperature dependence(s) are hard to evaluate merely based on phonon statistics;14,16 still such dependences could rather be “tailored” in a more or less reasonable way(s).14,16

On the other hand, the appropriate equation(s) of the Fock space formalism15 comprises all suitable “permutations” of many-particle (Fock) states, which (formally) correspond to contribution(s) from mixed many-particle states of a fermionic and/or bosonic ensemble into its given total energy, and apparently incorporates interaction effects among (quasi-) particles of different eigenenergies within the ensemble.

Within the framework of the “discretized” 2D GSM, such “permutations” could be taken into account by carrying out an additional summation over all allowed p quantities within the last (square) brackets of Eqs. (1a) and (2b) for N = 2, 3, 4,…. On the other hand, a number of the fundamental states of the “in-plane” acoustic phonons (corresponding to N ≡ 1) should not be affected by the aforementioned permutations (of the fixed many-phonon—or Fock—energy) since (by their definition) those permutations contribute to the many-particle mixed states only.15 This means that this additional summation has to be fulfilled only for external summation terms in Eqs. (1a) and (2b) with N ≥ 2. Appropriate sum over quantized energies of the Fock states of the “discretized” version of the 2D GSM could be evaluated readily based on the well-known rule for the sum of an “arithmetic progression,”17 which would eventually yield the following factor, κ = (pmax/2) [1 + (pmax/4)], placed within the last (square) brackets of Eqs. (1a) and (2b). Then, an appropriate number of the mixed multi-phonon states of the given Fock energy, pEmin (indexed with the integer p), will be taken into account aptly. The same κ term has to be incorporated into Eq. (1d) at N = 2, 3, 4 … as well but it will eventually amend Eqs. (1a) and (2a) only. Therefore, the right-hand side of Eq. (2b) has to be re-normalized with the p! term, which yields total number of permutations for the Fock level of pEmin energy.17 

Consequently, a more appropriate form of the basic equation of the “discretized” 2D GSM for the isobaric anharmonic lattice capacity reads

(3)

The latter equation does not contain any model parameters (e.g., r1 and r2) anymore. It is worth to be mentioned as well that the harmonic fraction of the lattice thermal capacity of the “discretized” 2D GSM still could be evaluated based on the first term (corresponding to N ≡ 1) of external summation(s) in Eq. (2b).

The “in-plane” LA and TA phonon modes dominate the lattice heat capacity of 2D materials at their relatively high (typically above 100 K) temperatures, while their low-temperature specific lattice heat capacity is usually “controlled” by the “outplane” (or “flexural,” ZA) mode(s).1–3 Furthermore, similar to the spectra of the in-plane LA and TA phonons, confined within the nano-sized flakes, the spectrum of their “outplane” acoustic ZA mode is also expected to be “discretized” for the flakes (ribbons) of finite spatial extents. Moreover, the phonon energy “quanta,” Emin, is expected to be defined as well by Eq. (1c) above. This eventually yields the following counterpart of the quasi-continuous version of Eq. (7) in Ref. 1 (see also Eq. (29) in the sub-section 6.1 of Ref. 2):

(4)

with the integer index q varying in the range from 1 to qmax = floor(hωZAM/Emin). Here, hωZAM stands for the cut-off energy of the “flexural” ZA mode of the (nano-sized) flake, and ZZA = [LxLy (hωZAM)2/(2 π h ζZA)], while ζZA denotes a dimensional (of m2/s) coefficient, which is characterizing quadratic dispersion relation for the “flexural” ZA mode.1–3 It is noteworthy that Eq. (4) takes into account the single-particle states of ZA phonons only.1,2 Typical simulation results obtained based on the basic equations of the “discretized” version of the 2D GSM will be revealed and discussed in Sec. III.

Figures 2(a) and 2(b) reveal the results of the simulation on the temperature dependencies of the “in-plane” [Fig. 2(a)] and “outplane” [Fig. 2(b)] harmonic fractions of the specific lattice heat capacity of the rectangular graphene flakes (nano-ribbons) with the different (from 1:1 to 10:1) aspect ratios of the (orthogonal) rib lengths Lx and Ly, while their temperature-dependent [CAH(T)/CH(T)] ratios are shown in Fig. 2(c).

FIG. 2.

(a) Temperature dependencies of the “in-plane” harmonic specific lattice thermal capacity of rectangular graphene nano-ribbons of Ly = 10 nm and with the different aspect ratios (indicated directly in the figure), evaluated based on Eqs. (1b), (1c), (2c), and (3). Some other simulation parameters are given elsewhere in Refs. 1 and 2 and shown directly in the figure. The open star represents room-temperature experimental data from Ref. 18. The log–log slope for the CH(T) dependence simulated for the aspect ratio of 10:1 is evaluated via a standard least-square-fit procedure17 within the temperature range 80 K ≤ T ≤ 200 K, and it is indicated directly in the figure with the blue symbols. See the main text for more details. (b) Temperature dependencies of the “out-plane” harmonic specific thermal capacity of rectangular graphene nano-ribbons of Ly = 10 nm and with the different aspect ratios (indicated directly in the figure), evaluated based on Eqs. (1c) and (4). Log–log slope for the CHZA(T) dependence simulated for the aspect ratio of 10:1 is evaluated via standard least-square-fit procedure17 and indicated directly in the figure with the blue symbols. (c) Temperature dependencies of the anharmonic fraction of the specific lattice heat of rectangular graphene nano-ribbons with the different aspect ratios, (Lx/Ly), evaluated based on Eqs. (1b), (1c), (2c), and (3); see main text for details. High-temperature log–log slope for the curve obtained at (Lx/Ly) = 10 is indicated directly in the figure with the blue symbols.

FIG. 2.

(a) Temperature dependencies of the “in-plane” harmonic specific lattice thermal capacity of rectangular graphene nano-ribbons of Ly = 10 nm and with the different aspect ratios (indicated directly in the figure), evaluated based on Eqs. (1b), (1c), (2c), and (3). Some other simulation parameters are given elsewhere in Refs. 1 and 2 and shown directly in the figure. The open star represents room-temperature experimental data from Ref. 18. The log–log slope for the CH(T) dependence simulated for the aspect ratio of 10:1 is evaluated via a standard least-square-fit procedure17 within the temperature range 80 K ≤ T ≤ 200 K, and it is indicated directly in the figure with the blue symbols. See the main text for more details. (b) Temperature dependencies of the “out-plane” harmonic specific thermal capacity of rectangular graphene nano-ribbons of Ly = 10 nm and with the different aspect ratios (indicated directly in the figure), evaluated based on Eqs. (1c) and (4). Log–log slope for the CHZA(T) dependence simulated for the aspect ratio of 10:1 is evaluated via standard least-square-fit procedure17 and indicated directly in the figure with the blue symbols. (c) Temperature dependencies of the anharmonic fraction of the specific lattice heat of rectangular graphene nano-ribbons with the different aspect ratios, (Lx/Ly), evaluated based on Eqs. (1b), (1c), (2c), and (3); see main text for details. High-temperature log–log slope for the curve obtained at (Lx/Ly) = 10 is indicated directly in the figure with the blue symbols.

Close modal

As one can see from Figs. 2(a)2(c), both harmonic (including the “outplane” one) and anharmonic fractions of the thermal lattice capacity of those nano-ribbons are essentially size-dependent. Indeed, the low-temperature (T < 50 K) parts of the CH(T) curve in Fig. 2(a) and those of CHZA(T) curves in Fig. 2(b) are apparently affected by the morphological “cut-off” effects, which yield a gap in the spectra of the “in-plane” and “out-plane” phonons, confined within those nano-ribbons. The gap width and “cut-off” temperature(s) are affected [based on Eq. (1c) from Sec. II] by the sizes of the rectangular nano-ribbons. Therefore, nano-ribbons with lower (Lx/Ly) ratios are expected to exhibit the “cut-off” effects (or so-called “low-temperature anomaly”1,2) at higher temperatures due to smaller (average) sizes of such nano-ribbons and a wider low-energy gap in their LA, TA, and ZA spectra.1,2

Similar “low-temperature anomalies” are, however, almost “suppressed” in the temperature dependencies of the CAH(T)/CH(T) ratios in Fig. 2(c) because of (at least partial) compensation(s) of the “cut-off” effects in the low-temperature parts of the CAH(T) and CH(T) dependencies for the given nano-ribbon.

On the other hand, entire CAH(T)/CH(T) dependences [within the whole temperature range revealed in Fig. 2(c)] are affected significantly by the aspect ratios of the nano-ribbons rib lengths: larger CAH(T)/CH(T) ratios are routinely obtained for the nano-ribbons with the larger aspect ratios within the whole temperature range, represented in this figure.

Such (temperature-dependent in general) enlargement in the anharmonic fraction of the lattice heat capacity of the elongated rectangular graphene nano-ribbons originate from interplay among (more than seven-fold!) increment in κ the factor [apparently caused by lower Emin energy for the nano-ribbon with bigger aspect ratios, see Eqs. (1b), (3) in the Sec. II], and diminishment in the F2(Lx, Ly)/Fxy(Lx, Ly) ratio for elongated nano-ribbons (see Fig. 1 in the previous section and comments to that figure therein). Their CAH(T)/CH(T) fraction also reveal significant (and size-dependent) “low-temperature anomalies” at (T < 100 K), see Fig. 2(c).

In general, anharmonic lattice thermal capacity is related intimately to interactions among (different kinds of) the “in-plane” and “outplane” acoustic phonons1,2,19 as well as to so-called “quasi-harmonic” effects,20 see also references therein. Such interactions define the mean free path (MFP) of acoustic phonons and also dominate the thermal conductivity of solids.20 Thus, in general, both lattice thermal capacity and thermal conductivity of rectangular flakes (nano-ribbons) are expected to be size- and temperature-dependent due to the interplay among temperature-dependent phonon MFP quantity and (variable) lengths of their orthogonal ribs. Such interplay is the main origin for the size-dependent harmonic and anharmonic CH(T) and CAH(T) functions for nano-ribbon(s), evaluated based on Eqs. (1b), (1c), (2c), and (3) [see Figs. 2(a)2(c) above]. In addition, scattering of the “in-plane” acoustic phonons on the edges of the rectangular nano-flakes (known as “Casimir scattering”21) is expected to affect both CH(T) and CAH(T) functions of the “discretized” version of the 2D GSM.

However, its basic equations yield size-independent harmonic and anharmonic fractions of the lattice thermal capacity of square graphene flakes at their elevated temperatures (T > 100 K), while features of their “low-temperature anomalies” on appropriated CH(T) dependencies would be essentially size-dependent even for the square graphene flakes of L < 100 nm: see, for instance, Figs. 1(a) and 1(b) in Ref. 1 and Figs. 20(a) and 20(b) in Ref. 2. This independence has to be attributed predominantly to the fact that Eqs. (1b), (1c), (2b), (2c), and (3) are actually valid only for the nano-sized (graphene) flakes and ribbons with the typical lengths of their orthogonal ribs significantly shorter than the experimental room-temperature MFP of the acoustic phonons in graphene (of ∼0.8 µm).1,2,22 Casimir’s scattering21 becomes isotropic for the square graphene flakes: its rate is expected to be the same along two orthogonal dimensions of such flakes.

In contrast, the effect of Casimir scattering on the edges of the rectangular nano-ribbons is anisotropic at LxLy. This eventually causes a weak (as long as both Lx and Ly lengths remain much shorter than the MFP) effect of changes in the rib length of (elongated) nano-ribbons on their CH(T) and CAH(T) functions [see Figs. 2(a)2(c)]. From a more formal viewpoint, the anisotropic harmonic and anharmonic lattice thermal capacities of the nano-ribbons at LxLy could be attributed as well to the anisotropic conventional and many-particle VDOS functions embedded into Eqs. (2b) and (3).

On the other hand, temperature-dependent harmonic and anharmonic fractions of the lattice thermal capacity of the (rectangular) nano-flakes with micrometer-sized rib lengths have to be evaluated using a different version(s) of the basic equations of the 2D GSM [see Refs. 1 and 6 and sub-section 6.1 of Ref. 2; see also Eq. (1a) in Sec. II]. The harmonic fraction of the lattice heat capacity evaluated based on those equations [“re-normalized” with the C2 and Z2DN terms—see comments to Eq. (1a)] becomes size-independent at elevated temperatures with increment, while anharmonic counterpart become size-dependent and diminishes as ∝ L−1 with increments in the rib length of those square nano-flake,1,2,6 reflecting result(s) of the interplay among the rib length of such flakes and the MFP.

Other important factors affecting the thermal conductivity of spatially non-homogeneous 1D, 2D, and 3D semiconductors and insulators, flakes, nano-wires (etc.) are phonon scattering on point and linear defects23,24 as well as on their edges,21 grains boundaries, roughness on the edges of the flakes,23,24 and so on. The latter mechanisms were discussed in detail elsewhere in Refs. 23 and 24 for the case of graphene, while the former approach incorporates not only phonon scattering on atomic-size impurities but also their “isotopic” scattering in graphene flakes and monolayers fabricated using different techniques.23 The basic idea of the latter model resembles that of one of the well-known Casimir’s model,21 although the isotropic scattering model incorporates a specific “specularity” parameter, which characterizes phonon scattering rate on the “non-ideal” edges of the graphene flakes.23,24

It is noteworthy as well that the heat capacity of bilayers of the graphene (BLG) and twisted bilayers of the graphene, T-BLG, could be affected as well by interactions among phonons located within the adjacent graphene layers; such effects were studied and reviewed thoroughly elsewhere in Refs. 24–26. In contrast to predictions of the (“discretized” version) of 2D GSM [see Fig. 2(c) above], additional phonon scattering caused by interactions among the neighboring graphene layers (of formally infinite in-plane sizes) of the T-BLGs rather reduce significantly their lattice heat capacity23 and affect considerably features of its temperature dependence. Indeed, the low-temperature log–log slopes of the temperature dependencies plotted separately for the “in-plane” and “out-plane” components of the harmonic lattice thermal capacity become ∼2.3 for its “in-plane” fraction for the T-BLG with 13.2° rotation among the adjacent layers [instead of ∼2.0 for the single nano-ribbon, see Fig. 2(a) herein] and ∼1.3 for the “outplane” fraction for such T-BLG [instead of ∼1.0 for the single flake; see also Fig. 5(b) in Ref. 24 and Fig. 2(b) herein]. Based on simulation results presented elsewhere in Ref. 25, the Cv(T) ∝ T1.3 emerges for low-temperature (T < 15 K) volumetric specific lattice heat capacity of the BLG, while Cv(T) ∝ T1.6 occurs for the T-BLG with the interlayer rotation angle of 21.8°.25 This indicates the unambiguous significance of the interaction and anisotropic effects,25,26 which are affecting not only the magnitudes of the low-temperature volumetric heat capacity of the BLG and T-BLG structures but also their “vibrational density-of-states” (VDOS) functions.23 The anisotropic conventional and many-particle VDOS functions emerge for a single rectangular graphene flake within the framework of the 2D GSM1,2,13 [see also comments to Eq. (1a) in Sec. II], although this anisotropy arises only for those flakes with non-equal rib lengths. Effects of features of the conventional phonon spectrum on thermal properties of the bi-layered and three-layered graphene structures (with a comparison to those of graphene monolayers) were reviewed elsewhere in Ref. 26.

Simulated temperature-dependent harmonic and anharmonic components of the specific lattice thermal capacity of rectangular graphene nano-ribbons of different (from 1:1 to 10:1) aspect ratios have been obtained within the framework of the “discretized” two-dimensional “Generalized Skettrup Model” (2D GSM). Its anisotropic basic equations incorporate “discretized” conventional and many-particle spectra on confined (within the nano-ribbons) LA, TA, and ZA phonons and ensure weakly ribbon-size-dependent “in-pane” harmonic and anharmonic fractions of the lattice thermal capacity of the elongated (graphene) nano-ribbons at their elevated temperatures. The anharmonic fraction(s) of the lattice thermal capacity of the “discretized” 2D GSM are affected more significantly (as compared to alterations in the harmonic fraction) by alterations in the ratio of the lengths of orthogonal ribs of those rectangular nano-ribbons. However, both harmonic and anharmonic “in-plane” lattice thermal capacities become size-independent at elevated temperatures for the square flakes. All of these prompt implementation of the “discretized” version of 2D GSM for realistic simulations on temperature-dependent harmonic and anharmonic fractions of the specific lattice thermal capacity of rectangular flakes of various chemical compositions and sizes.

All simulation results exhibited in this article are obtained with utilities created using Octave 3.6.4 free software: Copyright © 1996–2011 John W. Eaton [email protected]; see also http://www.octave.org.

The authors have no conflicts to disclose.

Valeri Ligatchev: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Methodology (lead); Software (lead); Writing—original draft (lead); Writing—review & editing (lead).

The data that support the findings of this study are available within the article.

1.
V.
Ligatchev
,
ECS J. Solid State Sci. Technol.
9
,
093014
(
2020
).
2.
V.
Ligatchev
,
‘Generalized Skettrup Model’ and Lattice Thermal Capacity of Spatially Non-Homogeneous and Low-Dimensional Semiconductors and Insulators
(
LAP – Lambert Academic Publishing
,
Saarbrücken, Germany
,
2021
), pp.
1
152
.
4.
M.
Pyda
and
B.
Wunderlich
,
J. Therm. Anal.
49
,
685
(
1997
).
5.
M. L.
Di Lorenzo
,
G.
Zhang
,
M.
Pyda
,
B. V.
Lebedev
, and
B.
Wunderlich
,
J. Polym. Sci., Part B: Polym. Phys.
37
,
2093
(
1999
).
6.
V.
Ligatchev
,
ECS J. Solid State Sci. Technol.
10
,
119001
(
2021
).
7.
L.
Lindsay
,
C.
Hua
,
X. L.
Ruan
, and
S.
Lee
,
Mater. Today Phys.
7
,
106
(
2018
).
8.
A. J. H.
McGaughey
,
A.
Jain
,
H.-Y.
Kim
, and
J.
Bo Fu
,
Appl. Phys.
125
,
011101
(
2019
).
9.
S. K.
Singh
,
M.
Neek-Amal
, and
F. M.
Peeters
,
Phys. Rev. B
87
,
134103
(
2013
).
10.
J. H.
Los
,
K. V.
Zakharchenko
,
M. I.
Katsnelson
, and
A.
Fasolino
,
Phys. Rev. B
91
,
045415
(
2015
).
12.
S.
Baroni
,
S.
de Gironcoli
,
A.
Dal Corso
, and
P.
Giannozzi
,
Rev. Mod. Phys.
73
,
515
(
2001
).
13.
V.
Ligatchev
,
Polycrystalline and Spatially Non-homogeneous Amorphous Semiconductors and Insulators
(
Nova Science Publishers
,
Hauppauge, NY
,
2017
), pp.
1
289
.
14.
V.
Ligatchev
, “
Anharmonic lattice thermal capacity and size-dependent melting temperature of square graphene nano-flakes
,” in
Invited Presentation on the NANOMACH 2022
,
Liberty Hotels Lykia, Oludeniz, Mugla, Turkey
,
22–28 April 2022
.
16.
V.
Ligatchev
, “
Anharmonic lattice thermal capacity and size-dependent melting temperature of square graphene nano-flakes
,”
AIP Proc.
(submitted).
17.
Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables
, edited by
M.
Abramowitz
and
I. A.
Stegun
(
Dover
,
New York
,
1965
), pp.
1
1046
.
18.
E.
Pop
,
V.
Varshney
, and
A. K.
Roy
,
MRS Bull.
37
,
1273
(
2012
).
19.
P. G.
Klemens
and
D. F.
Pedraza
,
Carbon
32
,
735
(
1994
).
20.
22.
S.
Ghosh
,
I.
Calizo
,
D.
Teweldebrhan
,
E. P.
Pokatilov
,
D. L.
Nika
,
A. A.
Balandin
, and
W.
Bao
,
Appl. Phys. Lett.
92
,
151911
(
2008
).
23.
D. L.
Nika
,
E. P.
Pokatilov
,
A. S.
Askerov
, and
A. A.
Balandin
,
Phys. Rev. B
79
,
155413
(
2009
).
24.
A. I.
Cosemasov
,
D. L.
Nika
, and
A. A.
Balandin
,
Nanoscale
7
,
12851
(
2015
).
25.
D. L.
Nika
,
A. I.
Cocemasov
, and
A. A.
Balandin
,
Appl. Phys. Lett.
105
,
031904
(
2014
).
26.
D. L.
Nika
and
A. A.
Balandin
,
Rep. Prog. Phys.
80
,
036502
(
2017
).