Molecular dynamics study of thermal transport in amorphous silicon carbide thin film

Abstract

The emergence of amorphous silicon carbide (a-SiC) thin film based photovoltaic applications has provoked great interest in its physical properties. In this work, we report the first comprehensive study of thermal transport in the a-SiC thin film from 10 nm to 50 nm under various conditions using empirical molecular dynamic (MD) simulations. The thermal conductivity increases from 1.38 to 1.75 W m⁻¹K⁻¹ as temperature increases from 100 K to 1100 K. A similar increase in the thermal conductivity from 1.4 to 2.09 W m⁻¹K⁻¹ is obtained with densities from 2.7 to 3.24 g cm⁻³. Besides, a slight increase in the thermal conductivity (15%) with calculation domain from 10 nm to 50 nm is observed, indicating that the size dependence of thermal transport also exists in nanoscale amorphous structures. For the physical interpretation of simulation results, the phonon mean free path (MFP) and specific heat are calculated, which are responsible for the temperature dependence of the thermal conductivity. The phonon group velocity is the key factor for the change in thermal conductivity with density. The results also show that the phonon MFP decreases rapidly with temperature and is subject to the Matthiessen's rule.

---

1. Introduction

As an alternative material to silicon, silicon carbide (SiC) has been studied extensively because of its great potential in industrial applications.^1–8 For example, crystalline SiC (c-SiC) can be used as a surface coating material in micro-sensors in extreme harsh environments.⁹ In addition, amorphous silicon carbide (a-SiC) has great potential as a photovoltaic material.¹⁰ The c-SiC has different forms, such as 3C–SiC, 6H–SiC, etc., with excellent physical properties. The a-SiC thin film is also attractive for industrial applications due to its less crucial deposition conditions while preserving most of the extraordinary properties of c-SiC, which is much more expensive to synthesize.¹¹ One drawback for amorphous materials is their thermal properties. Crystalline materials exhibit higher thermal characteristics due to their ordered lattice vibration (phonon), whereas amorphous structures introduce phonon scattering inside the material. Besides the crystalline/amorphous structures, the multi-scale dimensions (from bulk to nanoscale) of the material exhibit different thermophysical characteristics.^12–15 In addition, it is believed that the size effect in nanoscale thermal transport is only applicable in crystalline structures rather than amorphous materials due to their phonon boundary scattering. However, the dominance of boundary scattering of phonon transport in amorphous materials, in other words, whether the size effect still exists in nanoscale amorphous structures, still remains a question.

To date, not much information is available for the thermal transport of nanoscale a-SiC thin film, in spite of its importance for the thermal design of nanoscale devices. In an experimental study, T. Jeong et al. used a thermoreflectance technique to characterize the thermal conductivity of a-SiC films with a thickness from 100 to 2500 nm. They pointed out that the thickness dependence of the thermal conductivity was attributed to thermal boundary resistance between Au (substrate) and a-SiC during the measurements. However, the size of 100 nm used in their work is still too large to demonstrate the size effect.¹⁶ M. Mazumder et al. studied the temperature dependence of the thermal conductivity in 22 nm period Si/SiC amorphous multilayer films using the 3ω method. They assumed that the discrepancy between Debye's model and the dispersion model of phonon mean free path (MFP) of a-SiC was the same as that of the crystalline Si, and the phonon MFP of a-SiC was calculated according to others' experimental results because the intrinsic thermal conductivity of a-SiC could hardly be extracted.¹⁷ Besides experimental work, modeling tools have been used to study the thermal transport in SiC. Molecular dynamics (MD) simulation is very convenient and cost-effective to study thermal energy transport. Especially with the rapid development of computation capacity, MD simulation is frequently used to study the phonon transport behavior of materials at nanoscale.¹⁸ To date, not many simulation results are available for the thermal transport study of a-SiC compared with other nanostructured materials. L. J. Porter et al. studied the thermal properties of bulk β-SiC, which is one of the polymorphs of SiC with a zincblende crystal structure. In their work, both lattice dynamical calculations and MD simulations were performed. However, the size effect was not discussed because they only focused on the bulk materials based on equilibrium MD simulation.^6,7

In this work, we employed the MD simulation to conduct a comprehensive study of thermal transport in a-SiC thin film. The application of the material might be subjected to different conditions, for example, working temperature can be as high as 1000 K. Therefore, information on the thermal properties of nanoscale a-SiC thin film at various temperatures is important. In addition, different synthesis methods/processes to produce a-SiC with various states are important. Whether the inner pressure of the materials impacts the thermal transport of a-SiC thin film is also of great interest. Meanwhile, the application of a nanoscale a-SiC thin film is usually in the regime of above 10 nm. Whether sample dimension has any effect on the thermal transport is another point of interest to be focused on. The following section describes the methods including the principle and details about the MD simulation, followed by the description of the high temperature annealing process used for the construction of the amorphous structures. In the results and discussion section, density, temperature and size dependence of the thermal transport of a nanoscale a-SiC thin film are presented with their corresponding physical explanations. Finally, the phonon scattering/transport theory is used to comprehend the simulation results.

2. Molecular dynamics (MD) simulation details

2.1 MD simulation principle

In MD simulation, the individual atoms of silicon and carbon are treated as ideal mass points with no volume, interacting with each other through a given set of interatomic potentials.¹⁹ The macroscopic variables can be deduced from the coordination and velocities of atoms integrated numerically under the rule of classic mechanics.²⁰ The system is first stabilized at a constant temperature. The time integration is performed on the Nose–Hoover style non-Hamiltonian equations of motion, which is designed to generate the positions and velocities sampled from the canonical (nvt), isothermal–isobaric (npt) and isenthalpic (nph) ensembles.^20,21 The interatomic potential used in this work is a many-body description proposed by Tersoff, emphasizing the significance of the local coordination of atoms in determining bond strength between two atoms. The Tersoff potential has been adopted in a wide range of silicon and carbon-based materials, including graphene, carbon nanotubes and different polymorphs of SiC.^6,7,22–27

The direct method is employed to calculate the thermal conductivity of a-SiC based on its Fourier's law: , where J_u is one component of thermal current, ∂T/∂x_v is the temperature gradient along the x_v direction and k_uv is the thermal conductivity tensor. Each component of the thermal conductivity tensor is identical due to its isotropic characteristics. In the simulation, a cube is divided into 30 slices of identical size, and the slices 1/4 L away from the front and back boundary are set as the heat source and heat sink, respectively, as shown in Fig. 1(a).^28–31 Periodic boundary conditions are used in every direction.

Fig. 1 (a) Schematic description of our system model. The boundaries of all directions are periodic, and the slices of 1/4 L from each boundary act as the heat source and heat sink. (b) Illustration image of the atoms and bonds of SiC before amorphization. (c) Illustration image of atoms and bonds of SiC after amorphization.

2.2 High temperature annealing of a-SiC structures

The amorphous structure features the long-range random distribution of atoms rather than an orderly distribution of atoms with a crystalline structure. The commonly used methods for building amorphous structures are ion implantation and high temperature annealing. In addition, unlike crystalline structures, the density of the amorphous structure varies depending on the preparation pathways. V. Heera et al. studied the structural and density changes of a-SiC prepared by ion implantation experiments and found that the density of a-SiC decreases about 10% (as 2.75 ± 0.09 g cm⁻³) compared with c-SiC.³² P. Vashishta et al. showed that the density of a-SiC without internal pressure prepared by the annealing method is 3.076 g cm⁻³.³³ MD simulation can be used in the annealing method to prepare amorphous structures, and the zincblende structure of crystalline 3C–SiC is used for the sample preparation. First, the crystalline 3C–SiCs are maintained with different lattice constants (from 4.35 Å to 4.62 Å, corresponding to different densities from 3.229 g cm⁻³ to 2.695 g cm⁻³). Then, the system is heated to 8000 K at a rate of 2 × 10¹³ K s⁻¹ with a short time-step of 0.2 fs. After thermalization and relaxation, the system is carefully controlled and cooled to 1100 K at a rate of 1 × 10¹³ K s⁻¹. After the second equilibration lasting 400 ps, the amorphous structure of SiC is formed, as shown in Fig. 1(c). Our results show that the density of a-SiC at zero internal pressure is 3.016 g cm⁻³, which is consistent with the density calculated by the ab initio method by Finocchi et al.³⁴

To confirm our amorphous structure, radial distribution function (RDF) analysis was conducted on the c-SiC, and the structure obtained after the annealing process is shown in Fig. 2. RDF is an effective method to distinguish the difference between amorphous and crystalline structures by studying the peaks shown at various positions inside the materials (representing orderly distributed atoms).³⁵ From Fig. 2, we can clearly observe the orderly distributed peaks at regular positions for c-SiCs, while there are no regular peaks in the a-SiC, which implies that the amorphous structure was constructed successfully. Meanwhile, the C–C curve, C–Si curve and Si–Si curve have their first peaks at 0.157 nm, 0.19 nm and 0.243 nm, respectively, consistent with the work of Finocchi et al. and Tersoff et al.^34,36 In addition, the areas under the first peaks correspond to the number of neighboring atoms and is defined as , where g(r) is the RDF, N the number of atoms, V the volume of system, r₁ the lower boundary of the first peak, and r₂ the upper boundary of the first peak.^36–38 Despite the un-defined partial coordination of Si atoms, as explained by Tersoff, the number of nearest neighboring atoms of carbon atoms are computed for different densities and are listed in Table 1. The coordination of the carbon atoms is always four no matter whether the density is changed, indicating that threefold sp²-bonded carbons are not present inside the amorphous structure.

Fig. 2 Calculated radial distribution function for crystalline SiC (red line) and a-SiC (black line). The density is 3.229 g cm⁻³. The partial pair correlation functions of Si–Si, C–C and C–Si are presented.

Table 1 The partial coordination of carbon atoms calculated from the partial radial distribution function and degree of chemical disorder of a-SiC

Density (g cm^−3)	Pair	Partial coordination	Total coordination of carbon atoms	Peak position Angstrom
3.229	C–C	1.42	4.00	1.57^+0.02−0.02
	C–Si	2.58		1.91^+0.02−0.02
3.016	C–C	1.44	4.01	1.57^+0.02−0.02
	C–Si	2.57		1.91^+0.02−0.02
2.695	C–C	1.55	3.98	1.59^+0.02−0.02
	C–Si	2.43		1.93^+0.02−0.02

---

2.3 Simulation details for thermal transport study

The temperature gradient and heat flux are required to calculate the thermal conductivity of materials. The temperature of each slice is monitored by averaging the data from the last 250 ps. The heat flux is controlled by changing the non-translational kinetic energy added to the heat source, and the same amount of energy is subtracted at the heat sink piece by rescaling the atoms' velocities using the method described by Jund and Jullien.³⁹ The whole system is controlled at constant N (number of atoms), V (volume) and E (energy) as a microcanonical ensemble with the heat sink and heat source for 1 ns with a time-step of 0.5 fs. In this work, systems with six different sizes ranging from 2 × 2 × 10 to 2 × 2 × 50 nm³ and nine different densities ranging from 2.695 to 3.229 g cm⁻³ are simulated with temperatures from 100 K to 1100 K. Fig. 1 shows an example of one system with L = 15 nm and a = 2 nm. In the simulation, heat flux is first applied and after the temperature profiles along the z axis become stable (Fig. 3) by linear fitting of the temperature profiles between the heat source and heat sink, the thermal conductivity can be calculated. To ensure that the transverse area is large enough that our results are insensitive to it, the k values are estimated using different transverse areas as displayed in Table 2, indicating that a transverse area of 5 × 5 cells is sufficient for this work.

Fig. 3 Temperature profiles of a-SiC with the size of L = 15 nm and a = 2 nm along the direction of the heat current under different ambient temperatures after the stable temperature field forms. The heat flux of the system is 16.23 GW m⁻² for temperature profiles of 100 K and 200 K, and 81.15 GW m⁻² for all the other temperature profiles.

Table 2 The thermal conductivities of the system with different transverse areas. The heat flux per area through cross-plane with constant length along the heat current direction

Cell size (1100 K, 2.695 g cm^−3)	q (W m^−2)	k (W m^−1K^−1)
2 × 2 × 35	81.15	1.760^+0.084−0.076
4 × 4 × 35	81.15	1.347^+0.019−0.018
5 × 5 × 35	81.15	1.400^+0.013−0.013
10 × 10 × 35	81.15	1.405^+0.008−0.011

---

3. Simulation results and discussions

3.1 Temperature dependence of thermal transport in a-SiC thin film

Based on the temperature profiles in Fig. 3, the thermal conductivities are calculated at different temperatures, as shown in Fig. 4. The simulation of systems with three different densities is conducted. These three densities represented the positive internal pressure, zero internal pressure and negative internal pressure, respectively. Fig. 4 shows that the thermal conductivity of a-SiC increases with temperature. In addition, the results for the different densities have the same tendency, and there is a limit for the thermal conductivity before reaching 1100 K. The saturation of the thermal conductivity at high temperatures can be explained by our simulation of specific heat, which will be discussed later. The thermal conductivity calculation results can be checked with the modified Slack's minimum thermal conductivity proposed by Cahill, assuming that the minimum phonon mean free path (MFP) equals the half of the average interatomic spacing.⁴⁰ The tendency of our results is similar to the modified Slack's minimum thermal conductivity model.

Fig. 4 Relationship between the thermal conductivity of a-SiC and temperature. The results of three different densities are presented with the purple line representing negative pressure, the orange line representing zero pressure and the red line representing positive pressure. The blue line (lowest) is the minimum thermal conductivity calculated by the modified Slack's model.

The phonon density of states (DOS) from the Fourier transform of the velocity autocorrelation function can be obtained by:³³g(ω) = ∫e^iωt〈v⃑(t)v⃑(0)〉/〈v⃑(0)v⃑(0)〉dt, where v⃑(t) is the velocity vector of an atom and v⃑(0) is the initial velocity vector of the same atom. The angular brackets denote an average of the canonical ensemble. Fig. 5(a) shows the velocity autocorrelation function of the first 1 ps converging to zero when the duration is more than 0.3 ps. Fig. 5(b) shows the phonon DOS of a-SiC with a density of 3.016 g cm⁻³ at 1100 K (the black line). The blue and red line indicate the vibrational DOS of the silicon atoms and carbon atoms, respectively. Due to the molecular weight, the peak frequency of silicon is almost two times lesser than the carbon atoms. In Debye's model, the specific heat at a constant volume C_v can be determined as: ,³³ where k_B is the Boltzmann's constant, ℏ is the Plank's constant divided by 2π and g(ω) is the phonon DOS. From Fig. 5(c), for phonon DOS at different temperatures, we find that the vibration of carbon atoms (higher phonon frequency) is enhanced a little as the temperature increases, while the vibration of the silicon atoms (lower phonon frequency) remains constant. In above mentioned equation, two factors that determine specific heat are phonon DOS and temperature. Due to the higher proportion of low frequency phonons, a small change in high frequency phonon would not significantly change on the integration of the phonon spectrum (less than 10%), explaining that the phonon DOS is not the main factor for specific heat. In fact, temperature is the key to specific heat, as shown in Fig. 5(d), it increases from 112.6 J kg⁻¹K⁻¹ to around 1032.1 J kg⁻¹K⁻¹ with an increase in temperature from 100 K to 1100 K.

Fig. 5 (a) The velocity autocorrelation function (VACF) converges to zero; (b) phonon density of states (DOS) of Si atoms (blue), C atoms (red) and the summation (black); (c) phonon DOS of a-SiC at different temperatures (density: 3.016 g cm⁻³); and (d) specific heat per mass calculated by the phonon density of states. A saturated value of 1.2461 × 10³ J K⁻¹ kg⁻¹ at infinite temperature is listed for comparison.

Thermal conductivity can be derived from kinetic theory as: k = C_vvΛ/3,⁴¹ where C_v is the specific heat per volume, v is the phonon group velocity and Λ is the phonon mean free path (MFP). Thermal conductivity is a combined effect of the specific heat, phonon group velocity and phonon MFP. Treated as an elastic mechanical problem, the average phonon group speed can be calculated by ,⁴² where P is pressure and ρ is the mass density. The results reveal that the average phonon group velocity is independent of temperature with little variance from 8555 m s⁻¹ to 8794 m s⁻¹ (around 2.8%). Therefore, phonon MFP is another main factor, despite the specific heat, which is responsible for the temperature dependence of the thermal conductivity. Detailed analysis of the phonon MFP is presented in the following sections.

3.2 Density (internal pressure) dependence of a-SiC thermal transport

Different densities represent the different internal pressures applied onto the system during the annealing process. It provides the pathway for studying the relationship between thermal transport and internal pressure. Fig. 6(a) is the calculation result at 1100 K (as an example). It indicates that the thermal conductivity is strongly affected by density, which also agrees well with the lower-limit value and the trend calculated by the modified Slack's model. Fig. 6(b) shows the phonon spectra for different densities. Obviously, the effect of density on the phonon DOS is less significant than that of temperature by comparing Fig. 6(b) with Fig. 5(c). Through calculation, we find that the specific heat per mass is less affected (about 5%) when the density increases from 2.695 g cm⁻³ to 3.229 g cm⁻³, which is due to the relatively increased portion of carbon atoms' vibration in the phonon DOS, displayed as a red line in Fig. 6(c). However, the specific heat per volume increases ∼15% at a higher density, shown as a black line in Fig. 6(c). Meanwhile, the sound velocity increases rapidly with the density, following the same tendency as the pressure, as shown in Fig. 6(d). From the tendency and magnitude of the variance of k (∼50%) and sound velocity (∼30%), we can infer that the sound speed is the main reason for the density dependence of k, and that the specific heat per volume is the secondary source.

Fig. 6 (a) Thermal conductivity of amorphous SiC at different densities from 2.695 g cm⁻³ to 3.229 g cm⁻³; (b) phonon DOS of a-SiC for different densities; (c) specific heat per mass (red) and per volume (black) of a-SiC; (d) the internal pressure of a-SiC after annealing and the sound velocity as a function of density.

3.3 Size effect in thermal transport of nanoscale a-SiC thin film

Size effect is a common phenomenon in the thermal transport of crystalline materials on a nanoscale.^29,43–45 It is attributed to the phonon scattering effect at the boundary when the characteristic length is comparable to the phonon MFP. Whether size effect appears in the amorphous solid is described in this section. Fig. 7(a) shows that the phonon DOS varies with the sample size ranging from 10 nm to 50 nm in the z axis. It does not show much difference for the phonon DOS spectra. It is also confirmed by our calculation that the sound speed and specific heat of different sizes remain constant as the size changes. However, the thermal conductivity results reveal that the size effect does exist, as shown in Fig. 7(b). The thermal conductivity increases from 1.38 W m⁻¹K⁻¹ to 1.57 W m⁻¹K⁻¹ as the characteristic length increases from 10 nm to 50 nm. This phenomenon is similar with the size effect encountered in crystalline solids explained by Matthiessen's rule but is not significant.⁴⁶ Thermal transport inside the a-SiC thin film is dominated by phonon scattering. Phonon scattering in a limited space consists of intrinsic and boundary scattering, which can be considered as independent events, as expressed in this relation: 1/L_eff = 4/L + 1/L_∞,²⁸ where L_eff is the effective phonon MFP in our system, L is the characteristic length of our system and L_∞ is the phonon MFP of the bulk materials. The factor of 4 stems from the half distance between the heat source and heat sink, and the average position where the phonons collide at the center between the heat source and heat sink. Based on the thermal conductivity derivation, the reciprocal of thermal conductivity can be derived as: 1/k = a/L + b, where a = 12/C_vv, b = 3/C_vvL_∞ are constants, which can be determined from the linear relationship between 1/k and 1/L. The result is displayed as the embedded figure in Fig. 7(b). The phonon MFP of a-SiC bulk can be estimated through extrapolation with the computed average specific heat of 1051.1 J (kg⁻¹K⁻¹) and sound speed of 7499 m s⁻¹ as 0.225^+0.014_−0.018 nm, which is very close to the interatomic spacing of 0.231 nm. The corresponding thermal conductivity of a-SiC bulk is calculated as 1.59^+0.10_−0.12 W m⁻¹K⁻¹, which is very close to the value of 1.49 W m⁻¹K⁻¹ measured by Dongsik Kim et al. at 500 K.¹⁵ The deviation may stem from our higher temperature and different densities from each other as well as the uncertainty during measurement. These results reveal that the phonon MFP of amorphous solid is close to the interatomic spacing.

Fig. 7 (a) Phonon DOS of a-SiC with characteristic length from 10 nm to 50 nm; (b) thermal conductivity of a-SiC changes with characteristic length (density: 2.695 g cm⁻³), the embedded figure is the linear fitting of 1/k as a function of 1/L.

Fig. 8 (a) Phonon MFP decreases with temperature; (b) phonon MFP increases with characteristic length; (c) phonon MFPs of a-SiC thin film (purple) and bulk materials (orange) changes with density, and the calculated average interatomic distance (shown as the green line) is listed for comparison.

3.4 Phonon mean free path analysis

To understand the temperature and size dependence of k and to investigate the thermal transport mechanism of amorphous solids, the phonon MFP of a-SiC bulk and thin film at all temperatures and densities is calculated. Fig. 8(a) shows the relationship between phonon MFP and temperature for a-SiC thin film simulated in this work and the bulk materials. It shows that there is a rapid decrease of phonon MFP with temperature for both a-SiC bulk and film. Considering the aforementioned result that the thermal conductivity increases with temperature, this figure reveals that at low temperatures the phonon MFP is dominant and specific heat is low, while at high temperatures specific heat is the key factor for higher thermal conductivity. Fig. 8(a) shows that the phonon MFP decreases significantly with temperature, and the size effect is apparent in the low-temperature region in both the thin film and bulk amorphous structures. It is because the phonon MFP becomes closer to the characteristic length of the computation domain. Fig. 8(b) shows that the phonon MFP changes with size for the a-SiC bulk and thin film. Despite the uncertainty in the calculation of the bulk material, the phonon MFP of a-SiC thin film increases about 10% as size increases from 10 nm to 50 nm, suggesting that size effect exists but is very weak. The reduced phonon MFP is responsible for the reduced thermal conductivity at smaller sizes.

The calculation results show that the phonon MFP is comparable with the interatomic spacing in our amorphous structures. We are questioning if the slightly changed interatomic distance impacts the phonon MFP. Fig. 8(c) shows the relation between phonon MFP and density. The green line represents the average interatomic distance, which is calculated based on the density. There is no obvious trend for any relations between the MFP and average interatomic distance. Einstein's original concept of a random walk of localized oscillations can explain that the phonon MFP in amorphous solids is not determined by the global parameters but the localized interaction and geometry, including its disordered structure and unchanged minimum distance between the atoms.⁴⁷

4. Conclusion

In this work, we have employed the nonequilibrium molecular dynamic simulation to perform a comprehensive study regarding the thermal transport properties of a nanoscale a-SiC thin film. Temperature, density and the size dependence of the thermal conductivity are analyzed. The results show that the thermal conductivity increases with temperature until it approaches a certain value, which is attributed to the variance of the specific heat and phonon mean free path. The sensitivity of the thermal conductivity to density is explained by the change of phonon group speed and specific heat. A weak but existent size effect in the thermal transport of nanoscale a-SiC thin film is observed with the thermal conductivity increasing by about 15% with characteristic length from 10.0 nm to 50 nm. The extremely small phonon MFP is determined by localized oscillation rather than the global geometrical parameters in amorphous solids. Our study of thermal transport in nanoscale a-SiC thin film could be beneficial to the thermal design of all a-SiC based nanodevices.

Acknowledgements

The financial support from the Start-up of Wuhan University and National Science Foundation of China (no. 51206124) are gratefully acknowledged.

References

1. R. Madar, Nature, 2004, 430, 974–975.
2. D. Nakamura, I. Gunjishima, S. Yamaguchi, T. Ito, A. Okamoto, H. Kondo, S. Onda and K. Takatori, Nature, 2004, 430, 1009–1012.
3. A. Elasser and T. P. Chow, Proc. IEEE, 2002, 90, 969–986.
4. Y. T. Yang, K. L. Ekinci, X. M. H. Huang, L. M. Schiavone, M. L. Roukes, C. A. Zorman and M. Mehregany, Appl. Phys. Lett., 2001, 78, 162–164.
5. M. Mehregany, C. A. Zorman, N. Rajan and W. Chien-Hung, Proc. IEEE, 1998, 86, 1594–1609.
6. L. J. Porter, J. Li and S. Yip, J. Nucl. Mater., 1997, 246, 53–59.
7. J. Li, L. Porter and S. Yip, J. Nucl. Mater., 1998, 255, 139–152.
8. R. Maboudian, C. Carraro, D. G. Senesky and C. S. Roper, J. Vac. Sci. Technol., A, 2013, 31, 050805.
9. Á. Barna, S. Gurban, L. Kotis, J. Lábár, A. Sulyok, A. L. Tóth, M. Menyhárd, J. Kovac and P. Panjan, Appl. Surf. Sci., 2012, 263, 367–372.
10. K. S. Lim and O. Shevaleevskiy, Pure Appl. Chem., 2008, 80, 2141–2150.
11. W. Wang, R. K. Kalia, A. Nakano and P. Vashishta, 2007 MRS Spring Meeting, San Francisco, 2007.
12. G. A. Slack, J. Appl. Phys., 1964, 35, 3460–3466.
13. Properties of Silicon Carbide, ed. G. L. Harris, Institute of Electrical Engineers, London, 1995.
14. N. Papanikolaou, J. Phys.: Condens. Matter, 2008, 20, 135201.
15. S. R. Choi, D. Kim, S.-H. Choa, S.-H. Lee and J.-K. Kim, Int. J. Thermophys., 2006, 27, 896–905.
16. T. Jeong, J.-G. Zhu, S. Mao, T. Pan and Y. J. Tang, Int. J. Thermophys., 2012, 33, 1000–1012.
17. M. Mazumder, T. Borca-Tasciuc, S. C. Teehan, E. Stinzianni, H. Efstathiadis and S. Solovyov, Appl. Phys. Lett., 2010, 96, 093103.
18. C.-L. Tien, R. L. Jennifer and F.-C. Chou, Microscale Thermophys. Eng., 1998, 2, 133–137.
19. Z.-Y. Ong and E. Pop, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 81, 155408.
20. S. Plimpton, J. Comput. Phys., 1995, 117, 1–19.
21. D. Frenkel and B. Smit, Understanding Molecular Simulation: From Algorithms to Applications, Elsevier, 2001.
22. S. Berber, Y.-K. Kwon and D. Tománek, Phys. Rev. Lett., 2000, 84, 4613–4616.
23. W. J. Evans, L. Hu and P. Keblinski, Appl. Phys. Lett., 2010, 96, 203112.
24. J. Tersoff, Phys. Rev. B: Condens. Matter Mater. Phys., 1989, 39, 5566–5568.
25. J. Tersoff, Phys. Rev. Lett., 1988, 61, 2879–2882.
26. J. Tersoff, Phys. Rev. Lett., 1986, 56, 632–635.
27. M. Ishimaru, I.-T. Bae, Y. Hirotsu, S. Matsumura and K. E. Sickafus, Phys. Rev. Lett., 2002, 89, 055502.
28. P. K. Schelling, S. R. Phillpot and P. Keblinski, Phys. Rev. B: Condens. Matter Mater. Phys., 2002, 65, 144306.
29. D. Sellan, E. Landry, J. Turney, A. McGaughey and C. Amon, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 81, 214305.
30. X.-L. Feng, Microscale Thermophys. Eng., 2003, 7, 153–161.
31. Y. W. Yang, X. J. Liu and J. P. Yang, Mol. Simul., 2008, 34, 51–56.
32. V. Heera, F. Prokert, N. Schell, H. Seifarth, W. Fukarek, M. Voelskow and W. Skorupa, Appl. Phys. Lett., 1997, 70, 3531–3533.
33. P. Vashishta, R. K. Kalia, A. Nakano and J. P. Rino, J. Appl. Phys., 2007, 101, 103515.
34. F. Finocchi, G. Galli, M. Parrinello and C. M. Bertoni, Phys. Rev. Lett., 1992, 68, 3044–3047.
35. E. Kaxiras, Atomic and Electronic Structure of Solids, Cambridge University Press, 2003.
36. J. Tersoff, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 49, 16349–16352.
37. P. C. Kelires, Europhys. Lett., 1991, 14, 43.
38. V. I. Ivashchenko, P. E. A. Turchi, V. I. Shevchenko, L. A. Ivashchenko and G. V. Rusakov, Phys. Rev. B: Condens. Matter Mater. Phys., 2002, 66, 195201.
39. P. Jund and R. Jullien, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 13707–13711.
40. D. G. Cahill and R. O. Pohl, Phys. Rev. B: Condens. Matter Mater. Phys., 1987, 35, 4067–4073.
41. G. Chen, Int. J. Therm. Sci., 2000, 39, 471–480.
42. S. A. Bludman and M. A. Ruderman, Phys. Rev., 1968, 170, 1176–1184.
43. C. Jianwei, Ç. Tahir and A. G. William, III, Nanotechnology, 2000, 11, 65.
44. M.-H. Bae, Z. Li, Z. Aksamija, P. N. Martin, F. Xiong, Z.-Y. Ong, I. Knezevic and E. Pop, Nat. Commun., 2013, 4, 1734.
45. W. Jang, Z. Chen, W. Bao, C. N. Lau and C. Dames, Nano Lett., 2010, 10, 3909–3913.
46. J. M. Ziman, Electrons and Phonons: The Theory of Transport Phenomena in Solids, Oxford University Press, 2001.
47. D. G. Cahill and R. O. Pohl, Annu. Rev. Phys. Chem., 1988, 39, 93–121.

---