Investigation of phonon coherence and backscattering using silicon nanomeshes

Phonons can display both wave-like and particle-like behaviour during thermal transport. While thermal transport in silicon nanomeshes has been previously interpreted by phonon wave effects due to interference with periodic structures, as well as phonon particle effects including backscattering, the dominant mechanism responsible for thermal conductivity reductions below classical predictions still remains unclear. Here we isolate the wave-related coherence effects by comparing periodic and aperiodic nanomeshes, and quantify the backscattering effect by comparing variable-pitch nanomeshes. We measure identical (within 6% uncertainty) thermal conductivities for periodic and aperiodic nanomeshes of the same average pitch, and reduced thermal conductivities for nanomeshes with smaller pitches. Ray tracing simulations support the measurement results. We conclude phonon coherence is unimportant for thermal transport in silicon nanomeshes with periodicities of 100 nm and higher and temperatures above 14 K, and phonon backscattering, as manifested in the classical size effect, is responsible for the thermal conductivity reduction.

Supplementary Figure 2 | Validating Matthiessen's rule using ray tracing. Points are ray tracing simulations explicitly including both Umklapp scattering and diffuse boundary scattering, while the dashed line is the analytical Matthiessen's rule prediction for the given boundary scattering B and Umklapp scattering U . The very good agreement (7% maximum error for these dimensions of =30 nm, =90 nm, =100 nm and =0.8 m for pitch aspect ratio AR=1 and AR=8 nanomeshes) indicates that Matthiessen's rule accurately describes the combined scattering processes for the silicon nanomeshes.  Figure 2 of the main text and L=1.6 μm. For P  0.85, the periodic and aperiodic nanomeshes have the same within the error bars of the ray tracing (standard deviation 1 %, calculated as described in Supplementary Note 3). For larger specularities, the periodic conductivity is larger than the aperiodic conductivity, which we attribute to a "critical angle" phenomena as described in Supplementary Note 7. The experiments (see main text) are modeled using the P=0 results. is surprisingly accurate for all aspect ratios (<3% error even for 1:1 aspect ratio), while the commonly used approximation B,RNW = 1.12√ℎ is accurate within 5% for all < 2.5. Figure 7 | Multiple backscattering. The average number of collisions per transmitted phonon 〈 c 〉 increases with decreasing aspect ratio AR for the mesh with rectangular holes, indicating that multiple backscattering effects are important. Here =45 nm, = 60 nm, and = 100 nm. 〈 c 〉 is normalized by the AR=10 value 〈 c, =10 〉; the inset shows that this ratio is independent of for all aspect ratios.

Supplementary Figure 8 | Particle model comparison with previous experiment.
Particle BTE modeling results are compared to previously reported measurements 3 of thermal conductivity mesh for structures with different interpenetrating small hole arrays, normalized to the value of the original supercell SC . Ray tracing simulations predict that adding the interpenetrating small holes decreases B and mesh , in contrast to the previous constant-B assumption 3 which predicts a constant mesh (dotted line). The experimental trend is well explained by the BTE model without any modified dispersion or scattering due to coherence effects.

Supplementary Note 1: Porosity correction factor
Even in the absence of subcontinuum size effects, introducing porosity decreases thermal conductance G. To account for this bulk effect in interpreting measurements, several different conventions of a bulk porosity correction factor f have been taken in previous works when converting G to k . Yu et al. 4 assumed the effective cross-sectional area of the nanomesh to be the same as a nanowire array, resulting in f=n/py. This method yields an effective upper bound on k due to the neglected area of the bridging necks. Tang et al. 5 used the experimental porosity ϕ to calculate a volumetric reduction f= (1 -ϕ), which yields a lower bound on k. El Kady et al. 3 and Jain et al. 6 used finite element method (FEM) simulations to obtain the porosity correction factor f. Among these conventions, using the FEM f factor is most meaningful because it allows direct comparison between diffusive theory predictions and experimental measurements of k.
To quantify the porosity correction factor f, we performed FEM simulations using COMSOL software for arrays of square holes and circular holes. The porosity correction factor can be readily obtained by the ratio of the thermal conductance G of mesh structures to G of completely filled structures, i.e. a homogenous plate with the same overall external dimensions but no holes. As shown in Supplementary Figure 1(a), the porosity correction factor for arrays of square holes and circle holes decreases as porosity increases and the porosity correction factor difference between arrays of square holes and circle holes is minimal (less than 3% for all porosities below 0.57). Next, we compared the FEM simulation results with the 2D Maxwell-Garnett effective medium model (EMM) 7 f = (1 -ϕ)/(1 + ϕ) and find that results from FEM simulation and EMM are in excellent agreement (Supplementary Figure 1(a)).
We also computed the porosity correction factor f for rectangular holes with various aspect ratios using COMSOL, as shown in Supplementary Figure 1(b). It can be seen that f depends on aspect ratio and larger aspect ratios result in larger f for a given porosity. These FEM simulation results for rectangular holes agree well with a 2D EMM, suitably generalized for aligned elliptic cylinders 8 , f = (1 -ϕ)/(1 + ϕ * (py/px)). In this work, we calculated f for our experiments using the measured porosity and 3 rd order polynomial fitting between the FEM results for f as a function of porosity shown in Supplementary Figure 1 Lastly, we illustrate the effect of different conventions for the porosity correction factor f on the reported conductivity k in Supplementary Figure 1(c). For the same experimentally measured quantities G, L and A, the choice of f can dramatically change the reported k value by a factor of 2 for typical nanomesh dimensions (here, the 100 nm by 100 nm mesh). Therefore, different f conventions must be taken into account when comparing data sets from multiple sources.

Supplementary Note 2: Details of BTE model
The inputs to the particle BTE model are the phonon dispersion and the scattering parameters. We chose an isotropic Born-van Karman (BvK) sine-type dispersion approximation = 0 sin ( 2 0 ) , where 0 is the cutoff wavevector set by the atomic number density and 0 is determined by the weighted speeds of sound such that ( ) is exact in the low boundary scattering limit 9 . Though silicon's dispersion is not isotropic, of bulk and nanostructured silicon is often successfully modeled with an assumed isotropic dispersion 10 , and the isotropic dispersion approximation is essential in interpreting the ray tracing results. The LA and two TA branches are lumped into one effective branch. Our previous studies showed that either an isotropic full dispersion model that accounts for six independent phonon branches 11 or the Born von Karman (BvK) model with lumped branches 9,10 adequately describes phonon transport in bulk and nanostructured silicon. Different choices of dispersion can lead to different predictions of reduction due to boundary scattering. For example, for the mesh with =30 nm, = 90 nm, and = =100 nm, the BvK dispersion model for the pure silicon predicts around a 20% larger room temperature conductivity reduction than predicted by one first-principles approach 10 . The BvK model was selected because its predictions agree well with previous measurements of silicon nanostructures 10 from 20 to 300 K, and reduces to give the correct ( ) in the boundary scattering regime by construction.
The Umklapp scattering was determined by a least-squares fit 10 of experimental ( ). We use the Umklapp scattering rate U −1 = 2 exp (− U ) with = 1.53 * 10 −19 s K -1 and U = 144 K. The impurity scattering rate for the literature silicon was I −1 = I 4 with I = 2.54 * 10 −45 s 3 . We fit the impurity scattering for our devices to the periodic and aperiodic data shown in Figure 2c and obtain I = 1.96 * 10 −44 s 3 . This factor of 10 larger phonon impurity scattering rate appears reasonable given the literature silicon room temperature resistivity was 2000 Ω-cm (p-type) 12 and our boron-doped silicon resistivity is 14-22 Ω-cm, indicating larger dopant concentrations and scattering rates for our silicon. The impurity scattering modifies k by at most 15% at T=300 K, and for all T below 60 K the impurity correction to k is less than 5%. We use Matthiessen's rule to combine ray tracing boundary scattering with frequency dependent intrinsic scattering when predicting ( ). To test the validity of Matthiessen's rule for nanomesh structures, we also performed several ray tracing simulations that incorporate Umklapp scattering as well as boundary scattering. After a phonon is input to the domain or undergoes a scattering event, a free path is determined by drawing a random number from a uniform distribution between 0 and 1. The phonon free path is then -U ln ( ), where U is the Umklapp scattering mean free path. If the free path is smaller than the distance to the next surface, an Umklapp scattering event occurs at the corresponding point within the domain, with the phonon's outgoing direction uniformly randomized among 4 sr. Supplementary Figure 2 compares the calculated from the simulations explicitly including both Umklapp and boundary scattering with the predictions from Matthiessen's rule, for different input values of B / U . Both the AR=1 and AR=8 structures have = 30 nm, =90 nm, =100 nm, and =0.8 m. The discrepancy between the points that include the coupled scattering and the dashed line showing the Matthiessen's rule calculation is never larger than 7%, indicating that using Matthiessen's rule is not a major source of error in the high temperature analysis. Similarly good agreement between Matthiessen's rule and rigorously coupled scattering in nanostructures has also been observed in previous works 10,13 .

Supplementary Note 3: Description of ray tracing computational implementation
We use a ray tracing method with Monte Carlo integration to find the average transmission coefficient 〈 〉 for the 3D nanomesh structures considering only boundary scattering. A phonon is launched at the device-lead contact with a randomly selected initial position. The azimuthal angle of injection is also uniformly random, while the polar angle is sampled from a cos sin weighted distribution 14 . The phonon is propagated along a straight line to the nearest intersecting surface, and a uniform random number between 0 and 1 is drawn to determine the reflection angle. If the surface specularity parameter is less than (greater than) the random number, then the particle is reflected diffusely (specularly). The ray tracing continues until the phonon exits the computational domain by reaching either thermalizing contact. For all results presented, at least 3*10 5 phonons were launched, and 〈 〉 is calculated as the number of transmitted phonons divided by number of total phonons injected.
As a measure of uncertainty, the standard deviations of 〈 〉 and thus k are estimated following a common subsampling method 15 . The phonons are uniformly distributed into J subsamples (here, J = 40). Each of these subsamples is now regarded as an independent measurement, and the difference between these subsamples provides an estimate of the variance of the total transmission coefficient 〈 〉. The mean transmission coefficients of the n subsamples 〈 〉 are calculated, and the variance 2 of 〈 〉 is found as The computational domain is determined based on the symmetry of the structure. We take the global heat flow direction to be in-plane in ̂. Since the experimental sample widths are much larger than the pitch, we consider only one unit cell in the in-plane direction ̂ and apply fully specular boundary conditions on the periodic unit cell in ̂ due to the mirror symmetry. The top and bottom surfaces of the film form two additional walls of the domain in the planes. The neck n between pores along ̂ and ̂ is always identical. When comparing the periodic and aperiodic mesh, we explicitly place distinct mesh pores along the length of the device. For the aperiodic structure, and are held constant while is varied randomly and uniformly by up to ±20% of the nominal spacing by drawing a different random number for each pore. We perform multiple trials with different aperiodic structures to achieve ensemble statistics.
For simulations where all structures of interest have a discrete translational symmetry along ̂, we populate only one unit cell of surfaces to save computational time calculating the nearest collision. The mirror symmetry plane at the midpoint of the unit cell in ̂ allows us to convert the periodic boundary condition to specular boundary conditions at = 0 and = . We track the total distance the phonon has travelled along and terminate the simulation when the phonon reaches either of the device-lead contact planes.
A characteristic behavior of diffusive transport is a mean free path that saturates for long . The mean free path calculated from the transmission coefficient 〈 〉 describes the combined effects of diffusive nanomesh boundary scattering and a length-independent ballistic resistance (which can be represented through a transmission coefficient ball ) that depends on the configuration of the device-contact connection. To find the intrinsic B , defined as lim →∞ ( ), from finite-L simulations, we follow the technique of previous works 16

Supplementary Note 4: Ray Tracing Validation for Nanomesh Thermal Conductivity
As an additional validation of our modeling technique, we compare the ray tracing simulation predictions against published computational results for a silicon nanomesh 1 . The literature results were obtained using a variance-reduced Monte Carlo (VRMC) algorithm to solve the BTE while including both boundary and intrinsic scattering in the 3D diffuse nanomesh. We performed ray tracing simulations using the same dimensions and specularity to obtain B for a mesh with square holes ( = 34 nm, = 23 nm, = 22 nm) and a mesh with circular holes ( = 34 nm, diameter =11 nm, = 22 nm). For the purpose of this comparison only, we use the same Si [100] dispersion and intrinsic scattering times as in the original work 1 to make a fair comparison between numerical predictions. We note that the plotted quantity for the mesh in the reference 1 is eff = = / , and we use the EMM result that matched well with our FEM analysis, = (1 − )/(1 + ), for both circular and square holes. As shown in Supplementary Figure 4, the ray tracing and the VRMC predictions for eff ( ) show excellent agreement within 3% for all temperatures for the circular and square holes, validating our ray tracing code. The small remaining discrepancy between techniques may be due to our use of Matthiessen's rule or slight differences in the dispersion fits. A benefit of using a ray tracing approach as compared to solving the BTE with a known mean free path spectrum as an input is that in the ray tracing, the material-independent B is determined once and for all for a given geometry and surface specularity. Once B is found from the ray tracing, different isotropic dispersions, intrinsic scattering rates, or temperatures are then easily incorporated using Matthiessen's rule without requiring the additional computationally expensive simulations of a BTE solver. We note that B could also be found by solving the BTE assuming that all phonons have the same mean free path. In addition, detailed BTE solvers can also obtain the full temperature profile in addition to the heat flux.

Supplementary Note 5: Porosity Calculation
The hole shapes are neither perfect squares nor perfect circles due to the complexity of microfabrication processes. As a result, quantification of porosity is crucial for obtaining the porosity correction factor, which is required to finally convert experimentally measured thermal conductance to thermal conductivity. Thus, we performed image processing on SEM images of Si nanomesh structures using MATLAB to obtain the porosity. First, we convert SEM images to binary images using a grayscale intensity threshold. Then, the pixels corresponding to the empty and the filled state are each counted on the binary images. The porosity ϕ can be readily obtained from ϕ = (number of pixels of empty state)/(total number of pixels). The obtained porosity is used to calculate the porosity correction factor, as described in Supplementary Note 1.

Supplementary Note 6: Experimental Uncertainty
We start with the thermal conductivity expression = (Eq. 1) to estimate the uncertainty of the experimental data using propagation of error. The relative uncertainty of k can be written as

(Supplementary Equation 3)
The nanomesh conductance G is determined 18 by measuring the temperature rise of the heating membrane ΔTh and the temperature rise of the sensing membrane ΔTs to an input heating power Q. Since ΔTh >> ΔTs in our experiments, the relative uncertainty of G is

(Supplementary Equation 4)
The dominant term of δG/G is , which is estimated to be 1~3 %. The length L, width w, and thickness t are measured by SEM. While δL/L and δw/w are negligible, δt/t is ~4 % considering the uncertainty (~3 nm) of the SEM measurement. f is estimated as described previously in Supplementary Notes 1 and 5, and due to the uncertainty in the SEM measurement, δf/f is ~4 %. Consequently, δk/k is estimated to be ~6 % by error propagation.

Supplementary Note 7: Periodic and Aperiodic Thermal Conductivity
A main conclusion of this work is that the particle BTE model predicts very similar k for the periodic and aperiodic nanomeshes that are experimentally measured. We presented ray tracing results in the main text (Fig. 2d) showing the periodic and aperiodic as a function of L for two surface specularities, P=0 and P=0.8. Supplementary Figure 5 shows the periodic/aperiodic thermal conductivity ratio p / ap in the boundary scattering regime as a function of specularity P for nanomeshes with the same dimensions in Figure 2 and fixed L=1.6 μm. For very large specularities (P>0.9), the periodic conductivity is significantly larger than the aperiodic conductivity. We attribute this result to a "critical angle" effect: for the periodic structures, phonons incident at critical angles can be specularly reflected forward through the nanomesh for many periods before being backscattered, as illustrated in the specular cartoon in Supplementary Figure 5. Introducing aperiodicity can backscatter these critical angle phonons due to the shifted positions of the holes. Of course, adding the aperiodicity can also cause some phonons who were previously not transmitted in the periodic structure to transmit; however, the increased for periodic structures at P=1 indicates that the net effect of the aperiodicity is to increase backscattering in this scenario.
These ray tracing results reveal that at very low temperature (<5 K) where the dominant phonon wavelengths are very large compared to surface roughness so that surfaces are well approximated as fully specular, the particle model also predicts that p and ap can differ. But this clearly has nothing to do with coherence effects on the phonon dispersion relation. Thus, in this regime great care must be taken to properly identify the coherent contributions to thermal transport, since differences cannot be simply attributed to wave effects. Lastly, we note that our experimental trends and comparisons thoughout the main text are well described by the P=0 ray tracing predictions. Such fully diffuse behavior in our temperature range is expected from specularity model predictions such as Ziman's commonly used expression 19 assuming surface roughness >1 nm, as well as previous measurements of silicon nanostructures 20,21 .

Supplementary Note 8: Backscattering definition
Our definition of the backscattering fraction as the fraction of the total boundary scattering events that change the ̂ component of the phonon velocity from positive to negative for the phonons emitted from the hot reservoir differs slightly from previous discussions of backscattering in the literature. In the previous work by Ravichandran and Minnich 1 , backscattering is determined by whether an energy packet, which is emitted at the entrance of the repeated unit cell and which has one or more boundary scattering events within the cell, exits the repeated unit cell at the far end or returns to the entrance. Moore et al 22 use the term backscattering to describe the heat flux diffusely reflected from a single sawtooth on a nanowire surface. Though our definition is similar to these two scenarios, our method counts each boundary scattering event individually rather than summarizing the effect of one periodic unit cell or one sawtooth on the flux. In our definition, the choice of periodic unit cell has no influence on the backscattering fraction calculation, which is desirable because there are an infinite number of valid choices of primitive unit cell.
An alternate "both directions" definition of backscattering would consider both +̂ to −̂ and −̂ to +̂ boundary scattering events for phonons emitted from the hot reservoir as backscattering, since the phonon leaves the boundary in a different direction than it entered. The backscattering fraction would then be defined as the number of backscattering events divided by the number of total boundary scattering events for phonons emitted from the hot reservoir. Due to the high symmetry of the nanomesh, the fraction of boundary scattering events in which phonons heading in +̂ are scattered into the −̂ direction is identical to the fraction of boundary scattering events in which phonons heading in −̂ are scattered into the +̂ direction for long nanomeshes. For short nanomeshes, edge effects can result in backscattering fraction asymmetries; for all results shown here, we have verified that the backscattering fraction changes by less than 5% (e.g. from 0.5 to 0.5025) with an increase of a factor of 2 in length. Therefore all numerical results in the main text are unchanged whether the backscattering fraction definition of the main text is used or the "both directions" definition discussed here is used.
Lastly, we note that the transmittance 〈 〉 calculated in the ray tracing fundamentally quantifies the total backscattering of the device including multiple backscattering effects over the entire length, since 1-〈 〉 is the probability that a phonon injected at the entrance returns to the entrance due to any combination of backscattering events. Therefore conductance G in the boundary scattering regime is itself a direct measure of the total device backscattering, and the backscattering capability of different devices can be compared using the product GL/A. Figure 4 in the main text showed that the ray tracing for diffuse surfaces predicts that conductance G is independent of pitch aspect ratio. Since G is proportional to B in the boundary scattering regime, B of a mesh with square or rectangular holes and B of a nanowire array (large aspect ratio structure) with the same and are similar for the nanomesh dimensions considered in this work. This indicates that the boundary scattering mean free path of a rectangular nanowire B,RNW can be used to estimate the mean free path of the nanomesh  [23][24][25] is to take B,RNW ≈ 1.12√ℎ in analogy with the square nanowire of diameter result B,SNW = 1.12 . This error in the 1.12√ℎ approximation is less than 5% for < 2. 5.

Supplementary Note 9: Rectangular Nanowire Equivalent Conductance
Surprisingly, we see that the large limit, ) is within 3% of the exact result for all > 1 and more accurate than the small limit for > 1.5. This large limit is reminiscent of the logarithmic scaling 17 Table 1. For all results shown here, the neck n between holes is the same in both in-plane directions ̂ and ̂. We also calculate the percentage error in the equivalent nanowire approximation (discussed in Supplementary Note 9) B = B,RNW ( / ) , where B,RNW is calculated from the known analytical solution 2 for a rectangular nanowire with the same n and t as the mesh. Supplementary Table 1 shows results for a mesh with square holes arranged in a square array, Supplementary Table 2 shows results for a mesh with rectangular holes in a rectangular array, and Supplementary Table 3 shows results from a square array of circular holes. We note that the rectangular nanowire assumption is generally accurate at the 10% level for square and rectangular holes, but that using the smallest dimension between circular pores as an effective nanowire is a poor predictor of the circular mesh conductance, as shown in Supplementary Table 3. Figure 4(d) in the main text showed that the increase in backscattering fraction for the nanomesh did not affect the conductance G, which is directly proportional to B in the boundary scattering regime. There are several possible hypotheses we considered to describe the apparent lack of correlation between backscattering fraction and conductance for diffuse surfaces. For example, in small AR meshes, the lateral necks (n) between pores in ̂ are relatively large compared to the hole length (px-n) along ̂, which could feasibly provide additional channels for phonons to transmit forward instead of hitting a hole in a large AR nanomesh. This argument would predict that transmitted phonons hit fewer surfaces in small AR meshes: even though more of the boundary scattering events result in backscattering and are efficient at reducing heat transfer, there may be fewer scattering events overall. Another hypothesis is that multiple backscattering effects dominate. Phonons which backscatter once inside the mesh necks can scatter again off of the opposite wall before exiting the mesh in the original direction, or indeed may backscatter many times within the necks before reemerging. The multiple backscattering argument predicts that phonons transmitted through short AR meshes will hit more surfaces on average than phonons transmitted through the long AR meshes.

Supplementary Note 11: Phonon backscattering
Supplementary Figure 7 shows the average number of boundary scattering collisions per transmitted phonon 〈 c 〉, normalized by the AR=10 value 〈 c, =10 〉. The inset shows 〈 c 〉 as a function of length for different aspect ratios, indicating that the ratio 〈 c 〉/〈 c, =10 〉 is independent of . 〈 c 〉 is larger for smaller aspect ratios, meaning that phonons transmitted through small AR meshes on average hit more surfaces than phonons transmitted through large AR meshes. This result supports the multiple backscattering interpretation.

Supplementary Note 12: Particle model explanation for previous experimental result using ray tracing
A recent experimental work 3 measured of silicon nanomeshes with circular holes ( =1100 nm, =366 nm, =250 nm, and various diameters) and suggested that coherence effects arising from the mesh periodicity are important for heat transfer at 300 K. Here, we use ray tracing simulations to compute rigorous incoherent model predictions for these geometries assuming diffuse scattering, and find that the particle model also predicts the experimental trend previously attributed to phonon wave effects.
The experiment 3 was designed to manipulate the phononic band structure by interpenetrating an additional periodic array of small circular holes (radius s =102.5 nm) within the larger periodic mesh. The work calculates a modified dispersion relation for sub-10 GHz phonons due to changes in this secondary periodicity for five meshes with different interpenetrated arrays. We note that 10 GHz is relatively low energy compared to the characteristic thermal frequency at 300 K, estimated as th = b ℎ ≈ 6,250 GHz, where h is Planck's constant and b is Boltzmann's constant. Two trends in the experimental results were interpreted as suggestive of room temperature coherence effects. First, the measured mesh normalized by a measured tf of a film of the same thickness was smaller than a "lower-bound" prediction of a BTE model using an assumed B = . Second, the experimental mesh trend showed a stronger decrease with increasing porosity than predicted by a BTE model with constant B between structures. A crucial approximation in the reported BTE analysis is that B was essentially unchanged between the different meshes even though additional holes have been introduced.
We performed ray tracing simulations using the geometries and dimensions of the five reported mesh structures. Contrary to the previous assumption 3 , we find that B decreases monotonically with the number of additional small holes, and that the decrease in B from the original supercell mesh and the "1x1" mesh with the highest density of small holes is around 20% (from B = 320 nm to B = 256 nm). We use these B values as inputs to the BvK model using Matthiessen's rule with the same silicon scattering parameters used in the main text here, to find mesh at 300 K using our definition of from Eq. 1 in the main text. All surfaces are taken fully diffuse. To focus on the BTE model predictions for the change in mesh with additional holes, we normalize the simulation and experimental values by their respective supercell conductivity SC (see leftmost schematic of Supplementary Figure 8). We present results using the BvK model for clarity since the six different models analyzed in Yang and Dames 10 yield mesh / SC that differ by no more than 4%.
Supplementary Figure 8 shows mesh / SC as measured in the experiment 3 and as predicted here from the ray tracing and BTE model described in the preceding paragraph. The labeling notation follows the original work, and the experimental fractional percentage uncertainty in mesh / SC is taken to be the same as the reported percentage uncertainty 3 in mesh / tf . In the BTE model, the reduced mesh with increasing small hole density is due exclusively to increased boundary scattering because the BTE model uses bulk dispersion and bulk intrinsic scattering. This result shows that the BTE model predicts the observed experimental trend in mesh / SC once boundary scattering is rigorously considered, without any appeal to coherent effects.