The effect of ultrasmall grain sizes on the thermal conductivity of nanocrystalline silicon thin films

Nanocrystallization has been an important approach for reducing thermal conductivity in thermoelectric materials due to limits on phonon mean-free path imposed by the characteristic structural size. We report on thermal conductivity as low as 0.3 Wm−1K−1 of nanocrystalline silicon thin films prepared by plasma-enhanced chemical-vapor deposition as grain size is reduced to 2.8 nm by controlling hydrogen dilution of silane gas during growth. A multilayered film composed by alternating growth conditions, with layer thicknesses of 3.6 nm, is measured to have a thermal conductivity 30% and 15% lower than its two constituents. Our quantitative analysis attributes the strong reduction of thermal conductivity with decreasing grain size to the magnifying effect of porosity which occurs concomitantly due to increased mass density fluctuations. Our results demonstrate that ultrasmall grain sizes, multilayering, and porosity, all at a similar nanometer-size scale, may be a promising way to engineer thermoelectric materials. Thermoelectric materials convert heat into electricity and their performance is determined by their figure of merit ZT, which is generally too small in many materials for practical applications. Here, the authors demonstrate that a reduction in grain size for nanocrystalline Si can reduce thermal conductivity and potentially be used as a method to engineer greater ZT in Si for thermoelectric applications.

S ilicon is one of the most abundant and environmentally friendly materials on earth. Despite its excellent charge carrier transport properties and dominance in the semiconductor industry, it is not yet clear if it can also become an important player in thermoelectric applications 1 . A thermoelectric material would need to have a high Seebeck coefficient S, a high electrical conductivity σ, and a low thermal conductivity κ, in order to be competitive as measured by the dimensionless thermoelectric figure-of-merit (ZT), given by ZT = S 2 σT/κ, where T is the absolute temperature. The power factor, S 2 σ, of optimally doped crystalline silicon 2 can be as good as materials with a much higher ZT (>1) 3 .
Despite its high thermal conductivity, doped silicon and its alloys have a ZT that increases almost linearly with T 4-6 . Silicon also has a high melting point (1414°C). These properties could make silicon an appealing thermoelectric candidate that is usable in the ultrahigh temperature region (~1000°C), where suitable thermoelectric materials are scarce 7,8 , provided that its thermal conductivity can be reduced significantly.
In the kinetic expression, thermal conductivity can be written as: where ω is the phonon angular frequency; C(ω, T), v(ω), and Λ eff (ω, T) are the spectral specific heat, the group velocity, and the effective mean-free path (MFP) of the phonons, respectively. Phonon dispersion determines C(ω, T) and v(ω), while Λ eff (ω, T) depends on the phonon scattering mechanisms. The high thermal conductivity of silicon is mostly due to its large group velocity and MFP; at room temperature κ ≈ 150 Wm −1 K −1 , the group velocity characterized by the polarization averaged speed of sound at the low-frequency limit, v s = 6084 ms −1 (ref. 9 ), and the estimated spectral and polarization averaged MFP, Λ ≈ 260 nm, depending on the details of dispersion and modeling 10,11 . With the advent of nanotechnology, a steady reduction of thermal conductivity in silicon up to two orders of magnitude has been achieved using thin films [12][13][14] , nanowires and nanotubes [15][16][17][18][19][20][21] , nanostructures with porosity [22][23][24][25][26] , and nanocrystals [27][28][29] . What all of these experimental studies have in common is that in each category at least one of the relevant characteristic structural sizes L c , either thickness, diameter, or material spacing (e.g., between adjacent voids), has approached 10-20 nm (~5 nm in refs. 14,20,25,30 ), thus limiting Λ eff (ω, T) close to L c by means of phonon boundary scattering. Although modifications, such as surface roughness 14,18 , defects 19,21 , and composition, including oxidation 4,17 , all affect thermal conduction, the L c accounts for most of the experimentally observed reduction in κ 22 .
A reduced group velocity has also been invoked to explain results, where κ is smaller than the amorphous limit of the corresponding material 20,25 . The so-called minimum thermal conductivity, κ min , which limits phonon MFP to no less than one-half of their wavelength, while keeping the group velocity constant at its low-frequency limit, gives a reference lower bound for a given material 31 . So far, values of κ smaller than κ min at room temperature have been observed in nanocrystalline Si 80 Ge 20 composites (0.75 Wm −1 K −1 ) prepared by ball milling with average grain sizes of 9-12 nm (ref. 32 ), in multilayered Si/Ge nanodot thin films (0.9 Wm −1 K −1 ) with layer spacings of 3 nm (ref. 33 ), in polycrystalline Si 80 Ge 20 nanomeshes (0.55 Wm −1 K −1 ) with pore diameters of 31 nm (ref. 24 ), in porous silicon nanowires (0.33 Wm −1 K −1 ) with pore boundary spacing <5 nm (ref. 25 ), and in a number of other nanostructured materials having a characteristic structural size L c ≤ 10 nm (refs. [34][35][36] ). These encouraging results demonstrate the potential to reduce κ further in silicon. A recent theoretical work proposes that most of the heat is transported diffusively in amorphous, as well as nanostructured solids by substituting v(ω)Λ eff (ω, T) in Eq. (1) with a diffusivity D(ω, T) 37 . This approach reduces κ min by 37% or more if v is reduced.
Among those nanostructuring efforts, nanocrystallization has been the most promising approach to make silicon and silicongermanium (SiGe) alloys into practical thermoelectric materials [4][5][6]32 . Nanocrystalline silicon (nc-Si) is highly scalable, and can be prepared in both bulk 4,32 and thin film forms 5,6,38,39 , without the practical challenges associated with nanowires and nanomeshes. Even nc-Si material harvested from thin film processing can be pressed into bulk form in applications 5,39 . However, only one published work of nc-SiGe has obtained a κ lower than κ min 32 due to an exhaustive ball milling effort that reduced its grain size to 9 nm. Obviously, the success of the nanocrystallization approach depends critically on how small grain sizes can be reduced.
Grain size-dependent thermal conductivity has been studied systematically in bulk nc-Si for grain sizes from 550 to 76 nm, where the room temperature κ of bulk nc-Si decreases from 81 to 24 Wm −1 K −1 (ref. 27 ). Another study showed that the room temperature in-plane κ of a 1-μm thick free-suspended nc-Si film with grain size of 190 nm is~22 Wm −1 K −1 (ref. 28 ). Recently, we found that the room temperature κ of a 590-nm thick nc-Si thin film with grain size of 9.7 nm is 4.5 Wm −1 K −1 (ref. 29 ). Ultrasmall grain sizes of 3 nm have been achieved by a recent epitaxial growth and oxidation approach with grains coated with a monolayer of thermal oxide, resulting in a κ of 1.1 Wm −1 K −1 (ref. 40 ). The strong grain size dependence of κ in nc-Si is a convincing piece of evidence that phonon scattering at grain boundaries dominates heat conduction, demonstrating the potential of this approach to effectively reduce thermal conductivity. Alternative and more practical nc-Si production methods that may reach smaller grain sizes need to be explored.
In this work, we introduce a plasma-enhanced chemical vapor deposition (PECVD) technique to produce thin films of nc-Si, in which the grain sizes are controlled effectively by deposition parameters, specifically by hydrogen dilution of silane. We are able to reach an average grain size of 3 nm and to generate a layered structure with layer thickness of 3.6 nm by varying the grain sizes periodically during growth. The lowest κ we reached is 0.3 Wm −1 K −1 , about a factor of three below κ min .

Results
Material characterizations. It has been generally accepted that hydrogen dilution plays an important role in removing disorder and reconfiguring the Si-Si network with improved structural ordering 41,42 , although the exact growth mechanism of nc-Si is still debated 43 . The ratio of flowing hydrogen (H 2 ) to silane (SiH 4 ) gases, expressed as R = [H 2 ]/[SiH 4 ], is the most important film deposition parameter that determines the grain size, the crystalline content, and the porosity of the nc-Si in this work. For R < 10, the structure is still almost completely amorphous containing bonded hydrogen~10 at.% or more. As R increases, but remains <30, nanocrystallites start to appear in the amorphous matrix with a dramatic increase in both bonded and unbonded hydrogen to as high as 20 at.%, which concentrates at the interfaces between amorphous matrix and crystalline grains. As R continues to increase >30, the crystalline content increases rapidly and the nanograins begin to contact each other, forming heat conduction percolation paths 44 . The structure transforms from mostly amorphous to mostly nanocrystalline followed by a large reduction of hydrogen concentration to~6 at.% or less 41,42 . This is the point where thermal conductivity becomes interesting and relevant to thermoelectric applications; we thus start our thermal conductivity study from R = 25 and up.
We list R, other relevant film parameters, and room temperature thermal conductivities in Table 1. The sample of R = 80/100 is a multilayered film prepared by alternating the growth condition between R = 80 and R = 100 for 51 times with each layer thickness of 3.6 nm. The sample with R > 100 was made with a reduced SiH 4 flow below its control limit and its exact R is unknown. Despite their varying amorphous content, we call all of our PECVD-grown films nc-Si for simplicity. Technically, these films are structurally different from those studied in refs. 27,28 , and should be called hydrogenated nc-Si or nc-Si:H, due to their incorporation of hydrogen. For convenience, we drop "hydrogenated" in this article. Figure 1 shows the grazing-incidence X-ray diffraction (XRD) results of all seven PECVD films after background substraction. The five diffraction peaks are associated with the (111), (220), (311), (400), and (331) planes of silicon. Their positions and relative peak intensities are in good accordance with the diamond cubic crystal structure of silicon (ICDD no. 03-065-0565; plotted at the bottom for reference), indicating no preferred orientation for all films. The data for R > 100 is much noisier than the others due to its smaller film thickness and lower diffracting volume. The R dependence of the average grain sizes, calculated from the three main peaks at (111), (220), and (311), using the Scherrer equation, is shown in the inset where we assign R = 90 for the film with R = 80/100 and represent the R > 100 film with a right arrow. The same conventions are used for the rest of the paper. The average grain size increases with R starting when the film is predominantly amorphous (R = 25) and reaches a maximum at R = 80, when the structure is mostly nanocrystalline. This trend is similar to what we have observed earlier in nc-Si films prepared by hot-wire chemical vapor deposition (HWCVD) 44 . However, the grain size of the PECVD films increases only by~50% at most and remains within the range of 3-4 nm. This contrasts with the R dependence in HWCVD films, where the average grain size increases by a factor of 4 and reaches 9 nm at its maximum. Depending on deposition conditions, a wide range of grain sizes have been produced with PECVD 42 . The ability to produce ultrasmall grains makes PECVD-grown nc-Si films an appealing material for thermal conductivity studies and for potential thermoelectric applications.
The results of Raman spectra are shown in Fig. 2a for all seven PECVD films. A peak at~510 cm −1 associated with transverse optical modes is visible for all films, even for R = 25, indicating the appearance of the crystalline phase, consistent with the XRD results shown in Fig. 1. A faint peak at~615 cm −1 caused by vibrations of Si-H bonds is also visible. In order to estimate their crystalline content, we use three Gaussian peaks located at 480, 505, and 520 cm −1 to fit the Raman spectra, which originate from the amorphous, grain boundary, and crystalline phases, respectively. An example of the fitting process for R = 100 is shown in Fig. 2b. We calculate the integrated areas of the three peaks weighted by the scattering cross-section ratio for amorphous and crystalline silicon, which depends on grain sizes listed in Table 1. The crystalline content is approximated by the ratio of the combined weighted areas of both grain boundary and crystalline peaks to the total weighted areas of the three peaks; the results are shown in Fig. 2c. The crystalline content increases rapidly with R reaching a maximum of 73% at R = 100 and drops again for R > 100. Similar R dependence has previously been observed in PECVD-and HWCVD-grown Si films 42,44 . Upon closer inspection of the crystalline peak, one finds it is shifted to the right  Porosity for R = 100 is estimated from R = 80 and R = 80/100. The exact R value for R > 100 could not be determined. See "Results" section for details. Fig. 1 Grazing-incidence X-ray diffraction patterns. For all seven plasmaenhanced chemical vapor-deposited nc-Si films, the backgrounds are subtracted, integrated peak areas are normalized, and y-axes are shifted for clarity. The positions and the relative peak intensities of polycrystalline silicon (ICDD no. 03-065-0565) are at the bottom for reference. The Rdependent grain size is in the inset. The right arrow indicates the uncertainty in R for the R > 100 sample explained in the "Results" section.
progressively from R = 25 up to R = 80-100. For R > 100, this trend is reversed. The redshift of Raman crystalline peak with decreasing grain size has been associated with the phonon confinement effect due to phonon localization in nanocrystals 45 . Thus, the concomitant R dependence on grain size and crystalline content in our PECVD nc-Si films demonstrate the effect of controlling hydrogen dilution in both promoting crystalline growth and limiting the grain size. Figure 3 shows the high-angle annular dark-field images (HAADF) of cross-sectional scanning transmission electron microscopy (STEM) for three films: R = 80 in Fig. 3a-c grown on silicate glass, R = 80/100 in Fig. 3d-f and R > 100 in Fig. 3g-i, both grown on single-crystal silicon. A columnar growth structure can be observed in the lower magnification images in Fig. 3a, d, which is not obvious in Fig. 3g. Note that the white dots in upper part in Fig. 3g are artifacts due to the STEM crosssection preparation. Unlike the HWCVD-grown nc-Si films studied earlier, PECVD films have no amorphous incubation layer and silicon nanocrystallites grow almost instantaneously on substrates. This is advantageous in terms of thermoelectric applications. The higher magnification images in Fig. 3c, f, i reveal lattice planes within individual 2-5 nm crystallites, demonstrating that the crystallites are randomly oriented. This is in agreement with the XRD results, despite the columnar growth structure that is related to density variation. A particularly interesting point is that in Fig. 3e, we can clearly see the layered structure due to the alternation of deposition conditions between R = 80 and 100. The visible thickness is also consistent with measurements from the N&K Spectrophotometer. The HAADF intensity variation also indicates that the films are not fully dense. Pores between grains are visible as dark regions without visible lattice planes in the higher magnification images. The image intensity across the film-substrate interface provides a basis for estimation of the film porosity, at 2.5%, 7%, and 22% for the R = 80, R = 80/100, and R > 100 samples, respectively, also listed in Table 1.
Thermal conductivities. The temperature-dependent thermal conductivity results κ(T), measured by differential 3ω method (see "Methods" section), of all seven films are shown in Fig. 4. Additional time-domain thermoreflectance (TDTR) experiments were performed on R = 80 at University of Illinois (Urbana-Champaign) to measure its longitudinal speed of sound and to verify 3ω results, which are within 5% (smaller than the thickness uncertainty mentioned later) 46 . For comparison, the results of five HWCVD films we studied earlier are included 29,44 . Also plotted are κ of a 0.52-μm thick sputtered a-Si film grown at 230°C (ref. 47 ). It represents the thermal conductivity of a typical a-Si film, as the κ of a variety of a-Si films fall within a narrow range of it, such as a 0.6-μm thick HWCVD a-Si 29 and a 0.13-μm thick electron-beam evaporated a-Si 48 , although some special type of HWCVD a-Si can exhibit much higher κ 49,50 . The κ D min and κ BvK min curves are the conventional minimum thermal conductivity calculated, using the Debye dispersion and an alternative one that takes a reduction of the group velocity of high-frequency phonons into account (to be introduced later). On the top x-axis, we show the dominant phonon wavelength λ dom = v s /ν dom , responsible for heat transport at T, where ν dom = 90[GHz K −1 ]T is the frequency of dominant heat carrying phonons.
In the temperature range investigated, grain boundary scattering dominates the thermal transport as in the HWCVD films 29,44 . This is evident from the small values of κ and their overall positive temperature dependence. However, the R dependence of PECVD films is dramatically different from that of HWCVD films. While most of HWCVD films have a κ higher than the a-Si, all PECVD films have a κ lower than the a-Si. To demonstrate the difference, a comparison of the R-dependent κ at~294 K between HWCVD and PECVD is shown in Fig. 5, on the top and bottom x-axes, respectively, because different ratios of R are required to reach the onset of crystallinity 51 , depicted by a vertical blue bar. The room temperature κ of a typical a-Si film, κ D min and κ BvK min shown in Fig. 4, and κ Diff min from ref. 37 are presented as solid or patterned horizontal red lines. Starting from the small R, where the film is mostly amorphous, κ increases with R for HWCVD, but decreases for PECVD. In HWCVD, κ reaches a maximum at R = 8 and goes down for R = 10. For PECVD, κ decreases monotonically with R; passes κ D min at R = 80, and drops below κ BvK min for R > 100. An interesting observation is that the κ of R = 80/100 is 30% and 15% lower than R = 80 and R = 100, respectively, although both its average grain size and crystalline content are in-between. This suggests that alternating the grain sizes at an interval of 3.5 nm further reduces heat conduction despite an otherwise polydispersed grain distribution is not expected to affect grain boundary scattering much 52 . The most remarkable observation is that the κ of R > 100 is more than a factor of five smaller than that of a-Si and a factor of three smaller than the amorphous limit of silicon, κ D min . Compared to all previous results, this is achieved with a very modest porosity.
As the crystalline content in both HWCVD and PECVD increases with R in a similar way (see Fig. 3 in ref. 44 for HWCVD), it is obvious that the different R dependence of κ can be explained by the greatly reduced grain sizes in PECVD. In HWCVD, the increase of κ with R has been interpreted with an enhancement in nanocrystalline percolation paths that outweighs the increased interface scattering 44 . In PECVD, the opposite would be true. The reduction in grain sizes greatly increases the grain boundary content and hence their scattering of heat carriers. It affects boundary scattering transmissivity and specularity, as well as the speed of sound, rendering crystalline percolation paths devitalized. Studies show that high-frequency phonons or diffusive modes, which are dominant heat carriers at room temperature and above, are particularly sensitive to a reduction of L c 33, 36 .
We note that although the nc-Si films studied in this work are hydrogenated with varying hydrogen content, we don't think the inclusion of hydrogen or hydrogen bonds plays a meaningful role in thermal conduction. In our previous HWCVD nc-Si films, removal of hydrogen by a 600°C anneal increases κ by 35%, mostly attributable to an increase in grain size 29 . In addition, no hydrogen content dependence has been found in thermal conductivity of sputtered a-Si with H varying from 1 to 20 at.% (ref. 47 ).
Analytical modeling. In order to preserve their ultrasmall grain sizes for the interest of the thermal conduction study, none of the PECVD films are annealed to remove their amorphous content. At R > 80, their crystalline content reaches~60% and grain boundary scattering is expected to dominate heat conduction. We extend the model of grain boundary scattering with the Born-von Karman (BvK) dispersion to examine whether this approach is still adequate to explain our results. In our previous work on HWCVD nc-Si with a grain size of 9.7 nm, we also used the Debye and the Morelli-Callaway dispersions 29 . As the Debye model fails to give a satisfactory fit and the Morelli-Callaway model gives a result similar to BvK's, we only focus on the BvK dispersion in this work, as it is also considered to be most accurate among analytical models 11 . We previously considered both frequency-dependent and -independent grain boundary scattering 29 . However, as pointed out by Wang et al. 27 , the frequency-dependent model only affects the heat conduction of low-frequency phonons at T < 100 K. We only consider the frequency-independent (gray) grain boundary scattering in this work, which gives phonon MFP Λ B (ω) as, here α and D denote a fitting parameter and the average grain size, respectively. In addition, we also include impurity/defect scattering and Umklapp scattering as detalied in Methods section. For the demonstration purposes, we apply the model to R = 80 (D = 4.2 nm) and R > 100 (D = 2.8 nm), representing the largest and smallest D in this work, respectively. Both films also have relatively high crystalline content. Their BvK thermal conductivity results are shown in Fig. 6 together those of D = 76 (ref. 27 ) and 9.7 nm (ref. 29 ). We thus demonstrate that it is possible to describe the temperature-dependent κ of nc-Si with grain sizes varying by two orders of magnitude from 550 to 2.8 nm with a single set of parameters, where the only grain size-dependent parameter is α. As we noticed earlier in HWCVD nc-Si, the measured κ rises with increasing temperature stronger than a typical a-Si and BvK modeling. Similar behaviors have been observed in epitaxially grown nc-Si thin films with 3-40 nm grains coated with a monolayer of thermal oxide 40 . Higher temperature measurements would help to show how far this trend continues. Figure 7 summarizes the room temperature κ vs. grain size D for nc-Si studied so far, including the HWCVD nc-Si thin films 29,44 , the bulk nc-Si 27 , the polycrystalline Si thin films 28 , the epitaxial nc-Si 40 , and the PECVD nc-Si thin films presented in this work. The overall trend resembles those modeled in refs. 25,30,53,54 ; approaching saturation for D > 100 nm and decreasing steadily for D < 10 nm. However, we didn't observe the expected κ ∝ D at the small D limit 36 . The reason may be that the α we obtained by least squares fitting decreases with D, e.g., we have α = 0.26 and 0.15 for D = 4.2 and 2.8 nm, respectively. This trend has been observed at larger grain sizes as well 27,29 . So it is likely that the presumed linear dependence is not sufficient to describe the enhanced grain boundary scattering with decreasing D and a reduction in α compensates the deficiency. To improve our modeling, we set Λ B (ω) = Λ B as a fitting parameter in Eq. (10) and plot it vs. D in the inset of Fig. 7 covering three orders of magnitude. Here, we only use those data from the bulk, HWCVD and PECVD for which the values of α have been determined. We find that a superlinear dependence, as ; ð4Þ where κ rest = 145 Wm −1 K −1 represents the large grain asymptotic limit of κ, which is close to the value of crystalline silicon and κ ≈ 0.12D 1.3 (here κ is in unit of [Wm −1 K −1 ] and D in [nm]) is the small grain size limit at D < 10 nm. As α is determined by the the transmissivity t and specularity p of phonons traveling perpendicular and parallel to an interface respectively, and has been shown in Si/Ge superlattice nanowire that 9,27 α À1 ¼ 3t=4 it is reasonable to understand that a decrease of both t and p with D may be the source of grain size dependence of α.
Minimum thermal conductivity. The minimum thermal conductivity is a reference lower bound for a given material, signaling that a minimum in phonon MFP has been reached at a given frequency ω due to localization. Thermal conductivity smaller than κ D min would suggest that mechanisms other than a reduction of Λ eff take effect, such as a reduction in phonon density of states and/or speeds of sound.
In our previous study, we found that κ of nc-Si approaches that of a-Si with decreasing grain size and concluded that κ of nc-Si does not differ much from that of a-Si 29,44 . This observation is supported by molecular dynamic simulations, showing that κ  The R values for hot-wire chemical vapor deposited and plasma-enhanced chemical vapor deposited are on the top and bottom axes, respectively. The right arrow indicates the uncertainty of R for R > 100 explained in the "Results" section. The vertical blue bar represents the R range beyond which the films become mostly crystalline. The room temperature κ of a sputtered a-Si film, 0.52 μm thick, grown at 230°C (ref. 47 ) and the two κ min shown in Fig. 4 and an additional diffusive κ Diff min (ref. 37 ) are presented as patterned horizontal red lines. The error bars for plasma-enhanced chemical vapor-deposited nc-Si films are given (unless they are smaller than the data points), as discussed in detail in Supplementary Note 2.
decreases almost linearly with D and approaches that of a-Si for D = 2-3 nm (refs. 55,56 ). In this work, we demonstrate that with a further reduction of grain sizes the κ of many PECVD nc-Si films becomes smaller than that of a-Si and the amorphous limit, thus questioning such similarity. Indeed, Wingert et al. found that the κ of crystalline Si nanotubes with wall thickness of~5 nm is 30% smaller than that of a-Si nanotubes with the same dimensions 20 . Studies of nanostructured a-Si confirm that the amorphous limit is approached but not broken by a-Si with decreasing L c to~5 nm (ref. 57,58 ). It appears that the characteristic structural size, which dominates nanostructured crystalline materials, is not as important for a-Si.
In a-Si, heat is conducted approximately equally by propagating modes with MFP up to 1 μm, called propagons, and interactive nonpropagating modes with a MFP on the order of interatomic spacings, called diffusons 59 . In principle, a similar heat conduction mechanism applies to nc-Si with small grain sizes 55 . However, Zhao et al. studied κ of polycrystalline nanowires with grain sizes from~1 to 3.5 nm, using molecular dynamics simulations and found κ can be reduced to one-third of a-Si nanowires having the same structure dimensions 60 . Their detailed spectroscopic analysis shows a severe reduction of heat conduction in the effective phonon MFP range of~0.2-3 nm, pointing to a strong reduction of propagon contribution to heat conduction in nc-Si. They suggested that the presence of grain boundaries disrupts the continuity of atomic bonds causing an increased scattering and a decreased velocity of propagons, consistent with our observations. What is still unclear and beyond the scope of molecular dynamics simulations is the contribution of propagons with wavelength larger than the simulated structure and grain sizes. The strong grain size depenence observed in this work indicates their contributions to be minimal. We suggest that the combination of a polydispersed grain distribution and the occurrence of density fluctuation or porosity (to be discussed below) may have reduced their contribution. A detailed frequency-dependent thermal conductivity study may help to provide more insights 49 .
A comparison of Debye and BvK models in analyzing κ of nc-Si 27,29 demonstrate that the Debye model overestimates the contribution of high-frequency phonons by about a factor of two at room temperature. Thermal conductivity smaller than κ D min has been observed in an increasing body of nanostructured materials in recent years 24,[32][33][34][35][36] . While changes other than phonon MFP may account for some of the observations, overestimation of high-frequency phonon group velocity may have contributed to some of them 20,36 .
In this work, we insert the same Λ eff (ω, T) = πv s /ω as in κ D min , the smallest allowable phonon MFP at a given frequency ω, into Eq. (9) in "Methods" section to obtain the minimum thermal conductivity by the BvK dispersion, κ BvK min . The result is shown in Figs. 4 and 5 together with the standard κ D min with the Debye dispersion. This means that if a more realistic dispersion is used, no violation of Ioffe-Regel criterion nor other mechanisms would be needed for most PECVD nc-Si films except for R > 100. Our κ BvK min is also lower than the diffuson-dominated thermal conductivity limit κ Diff min 37 , also shown in Fig. 5. So far, we have not taken a reduction of sound speeds due to nanostructuring into account as has been observed in some recent studies, as a result of reduced structure size 20 or increased porosity 25 . The TDTR measurement on R = 80 shows that its longitudinal speed of sound is only~50% of the bulk silicon. If we Fig. 6 Fitting thermal conductivity vs. temperature. Two plasmaenhanced chemical vapor-deposited nc-Si films, R = 80 and R > 100, with the average grain sizes of 4.2 and 2.8 nm, respectively, are shown together with a hot-wire chemical vapor-deposited and annealed nc-Si film with grain size of 9.7 nm (ref. 29 ) and a bulk nc-Si sample with grain size of 76 nm (ref. 27 ). The fitting curves using the Born-von Karman dispersion and frequency-independent (gray) grain boundary scattering explained in the "Analytical modeling" section are shown in dashed blue lines. The error bars for plasma-enhanced chemical vapor-deposited nc-Si films are given, as discussed in detail in Supplementary Note 2. Fig. 7 The room temperature thermal conductivity vs. the average grain size D. The hot-wire chemical vapor-deposited nc-Si thin films 29, 44 , the bulk nc-Si 27 , the polycrystalline silicon 28 , the epitaxial nc-Si 40 , and the plasma-enhanced chemical vapor-deposited nc-Si thin films presented in this work are compared. The solid curve is the result of the Born-von Karman model calculation using the empirical superlinear D-dependent phonon mean-free path shown in the inset and explained in the "Analytical modeling" section. It is in contrast to the reference line labeled T 1 at the small D end. The error bars for plasma-enhanced chemical vapor-deposited nc-Si films are given (unless they are smaller than the data points), as discussed in detail in Supplementary Note 2.
reduce v s by half (assuming both v l and v t reduce proportionally), κ BvK min and κ Diff min would be reduced by 42% and 50%, respectively. This is lower than the measured κ of all our PECVD films and we expect sound speeds in R = 100 and R > 100 would be even lower.
The effect of porosity. At R = 80, porosity starts to appear in HAADF and it reaches 22% for R > 100. According to molecular dynamics simulations, porosity reduces κ by a reduction in phonon conduction channels and an increase in phonon scattering at the pore surfaces 61,62 . Crystalline silicon with high porosity (64-89%) has pushed down κ to as low as 0.04 Wm −1 K −1 (ref. 63 ). Single-crystalline silicon nanowires with 43% porosity have reached κ of 0.33 Wm −1 K −1 (ref. 25 ). In contrast to early theories where porosity is the only parameter, many experiments suggest a strong size effect on κ, such as pore sizes and their spacings 23,24,26,30,[64][65][66] .
For a quantitative analysis, we first apply the phononhydrodynamics approach which considers both the porosity φ and strucutre size L c (grain or pore) and κ eff is given by refs. 65,67 where C 0 ¼ 2C 1 þ C 2 expðÀγ=K n Þ and the Knudsen number K n = Λ/L c is the ratio of phonon MFP of silicon and its structure size. So the first and the second terms in the parenthesis reduce κ eff by an increase of φ or K n , respectively. This quantitative evaluation depends on the choice of parameters C 1 , C 2 , γ, and even Λ 65,67 . We use C 1 = 0.864, C 2 = 0.290, γ = 1.25, and Λ = 40 nm for our R > 100 film as in ref. 65 , where L c = 1-2.8 nm as the size range of porous regions are generally smaller than D (Fig. 3h-i). We reach κ eff = 2.1-5.8 Wm −1 K −1 . A Monte Carlo simulation work that considers an array of parallel pores along the heat conduction direction shows κ decreases not only with increasing porosity, but also with decreasing pore size as L c < 10 nm (ref. 66 ). However, κ remains above~1 Wm −1 K −1 even when φ = 30% and L c = 0.5 nm. Therefore both approaches, although completely different, severely underestimate the observed reduction of κ in our PECVD nc-Si. This suggests that a combination of the above discussed effects (grain size, phonon softening, and porosity) may have enhanced the reduction in κ.

Discussion
It would be desirable if all the factors encountered in this work could be analyzed in a comprehensive model. Recent theoretical work suggests hierarchical nanostructuring involving nanograins, multilayering and voids can all be beneficial to the power factors 68 . It is worth mentioning that the grain size reduction from R = 80 to R > 100 is~50%, but the reduction of κ is more than a factor of 3. The multilayered structure in R = 80/ 100 shows that a polydispersed grain distribution is effective, but this has so far not been understood theoretically. The low κ is further aided by an emergence of porosity and a reduction of sound speeds. However, neither a higher porosity 25 nor a larger reduction of sound speeds 20 yields a smaller κ in silicon than shown in this work. Based on our STEM observation, we find the porosity in PECVD nc-Si manifests more as isolated low-density regions than actual pores with sizes no larger than the crystalline grains. It is likely that such a density fluctuation is more effective than pores in disrupting heat conduction, via either propagons or diffusons. We thus attribute the strong reduction of κ, as well as v s to the magnifying effect of increased porosity that originates from the density fluctuation in our high R PECVD films to the effect of ultrasmall grain sizes, in agreement with recent Monte Carlo simulations 69 .
To summarize, we show nanocrystallization by hydrogen dilution using PECVD and HWCVD has quite different effects in terms of the resulting grain sizes and thermal conductivity. The smaller grain sizes and the lower thermal conductivities achieved by PECVD make this material a compelling choice for thermoelectric applications. We find that a superlinear dependence of phonon MFP on the average grain size best describes grain boundary scattering for grain sizes varying more than two orders of magnitude from 2.8 to 550 nm. With that, we provide a consistent description of both temperature and grain size dependences with the BvK dispersion. The observation that κ of a multilayered nc-Si is smaller than its two constituents not only emphasizes the importance of boundary scattering, but also reminds us of the limitation of our current analytical modeling in dealing with a broad or nonuniform distribution of grain sizes. Our results support a different heat conduction mechanism between nanocrystalline and amorphous materials due to the strong characteristic structural size dependence observed in the former. The origin of such difference is not known and more studies are needed. Applying the BvK dispersion to the minimum thermal conductivity effectively lowers its value by more than a factor of two. This evaluation method allows a reduction of the group velocity of high-frequency phonons, while keeping phonon MFP not smaller than half of their wavelength, a key criteria of minimum thermal conductivity 31 . Compared with the similarly low values of κ achieved in silicon 25 , our nc-Si films have a rather low level of porosity that is the natural consequence of our film production process. We conclude that the low κ achieved in this work is due to a combination of ultrasmall grain size, and the appearance of a density fluctuation at the same length scale of their grain sizes. We therefore propose that ultrasmall grain sizes, multilayering, and a moderate amount of porosity, all at comparable size scales, may be the ideal ingredients to control heat conduction in thermoelectric materials.

Methods
Thin film preparation. Thin films of nc-Si were prepared by PECVD with an Oxford PlasmaPro100 system. For all depositions, the chamber pressure was maintained at 700 mTorr and the substrate temperature was held at 250°C with a radio-frequency (13.56 MHz) plasma power of 90 W. Films were grown on sapphire, single-crystal silicon, and borosilicate substrates for 3ω thermal conductivity measurements and structural characterizations. Besides R, other deposition parameters, such as chamber pressure, plasma power, and substrate temperature, have effects as well 43,70 . We optimized and kept them constant in this work. The desired R was achieved by throttling [H 2 ] in the range from 25 to 100 s.c.c.m., while keeping [SiH 4 ] = 1 s.c.c.m.. As the flow rate for both H 2 and SiH 4 can only be adjusted between 1 and 100 s.c.c.m. in this system, the R value is limited to a maximum of 100. Most of the films were prepared with a constant R between 25 and 100. One film was prepared by alternating the growth condition between R = 80 and R = 100, for a growth time of 160 and 190 s, respectively, in each cycle based on a predetermined growth rate to reach 3.6 nm thickness for each layer. The process was repeated 51 times, resulting in a multilayered structure with varying grain sizes over a total thickness of 366 nm. Another film was targeted to grow at R = 100. However, an unintentional mechanical problem caused a reduced SiH 4 flow by an unspecified amount, resulting in an effective R that is larger than 100. We include this film for its significant reduction of κ. Systematic investigation of hydrogen dilution for R > 100 is, however, a subject of future study.
Film thickness was measured with an N&K Spectrophotometer at locations where 3ω thermal conductivity measurements took place. The results were verified by stylus profilometry on selected films after a photolithography and reactive-ion etching process to create a depth profile. The results agree to within 8%, which constitutes one of the main sources of error in the thermal conductivity evaluation.
Experimental. Although R-dependent studies of PECVD-grown a-Si and nc-Si abound [41][42][43]70 , we start with a detailed structural characterizations of our films. XRD measurements were performed with a Rigaku SmartLab X-ray diffractometer using fixed CuKα radiation in parallel-beam mode. The Raman spectra were acquired with a home-built micro-Raman setup using a single-mode 488 nm laser with typical spot size <1 μm and power~0.5 mW at the sample. Cross-sectional STEM studies of selected films were performed with a Nion UltraSTEM 200-X operated at 200 kV. HAADF images were collected with a nominal probe size of 120 pm, and probe current of 100 pA. Estimates of the film porosities were obtained from quantitative analysis of the HAADF image intensity, using the growth substrate as an in situ intensity calibration standard. The estimates included an assumption that variations in sample thickness across the interface are negligible, and the crystallite density is equivalent to bulk silicon. Two of the films (R = 80 and R = 80/100) were prepared by focused ion beam lift-out with a Helios G3 FIB, while the third (R > 100) was prepared by conventional mechanical polishing and low-angle Ar ion milling with a Fishcione Nanomill.
To measure κ, we employed the differential 3ω method which was established by Cahill et al. 47 . More details are given in Supplementary Note 1. A planar metal line, made of 5 nm Ti/120 nm Au, was deposited directly on the surface of the nc-Si thin film by a photolithography and lift-off process, see Supplementary Fig. S1. This patterned metal film served as both a heater and a thermometer, which supplied an AC current with a frequency ω in the range 3-1500 Hz and probed temperature oscillations in the film generated by an AC current at 3ω. This temperature oscillations is proportional to the sum of the thermal resistances of the film, the sapphire substrate, and their interface. The contribution from sapphire substrate can be subtracted as it is measured separately in the absence of a film, the contribution of which to the total varies from 10% at 85 K and 1500 Hz to 70% at 294 K and 3 Hz, see Supplementary Fig. S2. The interface thermal resistance is not known and is assumed to be negligible. This is based on the lack of significant thickness dependence of similarly prepared films and the excellent agreement between results on the same R = 80 sample obtained with 3ω and TDTR measurements, where thermal waves are confined to the film and no interfacial effect should exist. Four 3ω devices were fabricated together on each sample, with nominal metal film length of 1.2 mm and width of 10, 15, 20 and 30 μm, The exact length and width were determined after fabrication under a calibrated microscope. Thicknesses of our films were much smaller than the heater width, satisfying the one-dimensional heat flow condition of the differential 3ω technique. For confirmation purposes, most of the films were measured with more than one device and the differences in measured κ are typically within error bars given in Supplementary Note 2, where more details of the error estimation are given.

Data availability
The data that support the findings of this study are available from the authors on reasonable request, see author contributions for specific data sets.