Robust Metallic Nanolaminates Having Phonon-Glass Thermal Conductivity

Heat transfer phenomena in multilayer structures have gained interest due to their promising use in thermal insulation and thermoelectricity applications. In such systems, nanostructuring has been used to introduce moderate interfacial density, and it has been demonstrated that interfacial thermal resistance plays a crucial role in reducing thermal conductivity κ. Nevertheless, the main constraint for actual applications is related to their tiny size because they are extremely thin to establish appreciable temperature gradients. In this work, by severe plastic deformation process of accumulative roll bonding (ARB), a 250 µm thick Cu-Nb multilayer containing more than 8000 interfaces with periods below 40 nm was obtained, enabling the production of bulk metallic nanolaminates with ultralow κ. Multilayers present an ultralow κ of ∼0.81 W/mK at 300 K, which is 100 times smaller than its Cu-Nb bulk counterpart, and even lower than the amorphous lattice limit for the Cu-Nb thin film system. By using electron diffusive mismatch model (EDMM), we argue that both electrons diffusively scattering at interface and those ballistically crossing the constituents are responsible for heat conduction in the Cu-Nb multilayers at nanoscale. Hence, ARB Cu-Nb multilayers are intriguing candidate materials which can prove avenues to achieve stable ultralow κ thermal barriers for robust applications.

Determination of the interface physical properties for actual surfaces such as real area in contact, the number of contact points per unit area, and the average thickness of the void gaps via image analysis, as well as calculation of the thermal contact resistance RTCR. Determinations of the materials properties (Cu and Nb) such as electron mean free path, density of states, Fermi velocity, and transmission coefficient, as well as calculation of the interfacial thermal resistance RITR. Detailed discussion of the EDMM and TCR approaches, as well as calculations of the minimum thermal conductivity.

Experimental Procedure
The phase diagram of binary Cu-Nb shows that this system is immiscible, and no intermixing should occur in the equilibrium state [1]. This is due to the difference in crystal structure between Cu and Nb, and the large atomic radii mismatch between the constituent atoms. Therefore, such materials are good candidates for the development of bulk nanolamellar Cu-Nb composites for robust thermal barrier applications. The starting materials consist of single-phase laminates of commercial purity Cu and Nb which are degreased, wire-brushed, and stacked. Iterating these steps increases the number of layers exponentially while decreasing layer thickness. A schematic representation of the general ARB process is shown in Figure S1a. ARB processing was carried out in a rolling mill with a maximum separating force of 20 metric tons, equipped with D2 steel rollers of 10 cm in diameter, following the procedure as described in reference [2]. Before stacking, surfaces of Cu-Nb-Cu sheets of 500, 1000 and 500 microns thick, respectively, were brushed and cleaned using acetone in an ultrasonic bath. Then, stacks were preheated in an inert atmosphere furnace at 700 °C for 4 min, and the rolling linear speed was set at 5 × 10 −3 m/s. After each ARB pass, the bonded laminate was cut transversally in half using a mechanical shear. The two halves were then subjected to a repetitive process of cleaning, stacking and rolling until obtaining the multilayers with ultrahigh content of interfaces. In Figure S1b, optical images of the multilayers after eight sequential ABR process are shown. Evidently, after each ABR stage, the number of interfaces is increased following the power law = 2 , where stands for the process number. Based on this, the optical image with = 8 presents a transversal section of the sample containing 256 interfaces, which is a considerable number of interfaces. The inset shows the EDS chemical analysis of the multilayers, Cu and Nb were identified as alternating layers. Although ARB was applied to several different bimetal systems [3][4][5][6], the Cu-Nb system offers low solid solubility between the two phases. These characteristics result in excellent microstructural stability during ARB processing, allowing production of bulk metallic nanolaminates. Likewise, Figure S2 shows SEM images for samples with = 11, = 12 and = 13 iterations, which theoretically correspond to = 2048 , = 4096 and = 8192 interfaces, respectively. Well-defined interfaces can be seen formed by continuous layers of some nanometers in thickness, e.g., in sample with = 13, a period below 40 nm is identified. Insets present HRTEM images at the interfaces.

Thermal conductivity measurements of samples
Thermal response measurement of the samples was carried out via the hot-plate method. In such technique, the sample is installed between the hot plate (heater element) and cold plate (Aluminum heat sink), as it is shown in Figure S3. The heater element is placed on Silica aerogel tiles to reduce the heat leakage by conduction. Additionally, the whole system is put inside of a conventional vacuum system to avoid heat losses by convection. A radiation shield is also used in order to protect the equipment from thermal fluctuations coming from the environment. Furthermore, silver paste is used to enhance thermal contact between sample surfaces and heater element as well as cold plate. Then, once steady state conditions have been established at room temperature, a DC electrical current with regular increments is applied into a 1 kΩ resistive heater by using the Keithley 6221 AC/DC power supply (company, city, state, country). Measurements are performed until steady state is accomplished for each DC electrical current level applied. This action will cause several temperature rises ΔT on the sample which are sensed via micro-thermocouples embedded into the hot and cold plates, respectively. Temperatures from thermocouples are logged by using a data acquisition module from national instruments model NI-9213. In this way, a set of ΔT and their corresponding heat fluxes Q are obtained, then finally the transfer curve Q−ΔT can be plotted. To verify the correct performance of the system, before each ARB sample measurement, a reference sample of SiO2 was used for thermal conductivity system calibration. System calibration is performed via thermal conductivity measurement of a reference sample consisting of 1.1 cm x 1.1 cm x 1 mm of a corning glass slide. Figure S4 shows the ΔT-Q transfer curve for SiO2 corning slide. By using the relation = ⁄ , the thermal conductivity of the sample is calculated. Here L is the sample thickness, A the area, and R the thermal resistance which is estimated from the slope of the ΔT-Q curve. A value of κ ~ 1.37W/mK is obtained for SiO2 calibration sample. Such result is in good agreement with previously reported values for SiO2 [7]; hence, the system is reliable for measuring thermal conductivity. Moreover, to ensure measurement repeatability, each measurement was performed at least 5 times, results shown in Figure S4 are the average of the whole set of measurements for the SiO2 sample, whereas Figure S5 shows the results for ARB sample with n = 13. By a similar procedure, the thermal conductivity of the multilayer samples from = 1 to = 13 was calculated. Figure S5 shows the result for sample with = 13 iterations. As you can see, a value of κ~0.81W/mK is obtained.

Determination of parameters for modeling and calculations
In order to determine , , and , we analyze the interfaces of the two contacting surfaces using SEM images, as well as the specialized Digital Micrograph GATAN software (Gatan Microscopy Suite Software, version 3, Gatan, Inc., Pleasanton, CA, USA). With the aim to show the general procedure employed, an arbitrary sample with = 2 iterations is used as example. However, it is worth mentioning that all samples with period length in the mm range, i.e., 1 ≲ ≲ 6 were analyzed following the next procedure. Figure S6 shows an illustrative scheme of the parameters estimated so as to determine which is the squared root of the ratio of the real area of contact Ac to the total contact area A.
is a parameter which expresses the equivalent idealized gap thickness in terms of the average heights , and of the voids of actual surfaces in the x and y direction, respectively. Additionally, the number of contact points per unit area n is obtained by dividing the number of contacts on each pair of profiles in the x and y directions, i.e., nxi and nyi over an area defined by lxly. Clearly, Figure S6 shows an illustrative and idealized representation of actual interfaces. For most contacts, the height of the void is small compared with is width. Under this condition, it is reasonable to neglect both radial and convection heat conduction in the voids, thus, the thermal contact resistance associated with these type of interfaces can be estimated using the method proposed by W.M. Rohsenow [8] as This expression for is the sum of two terms. The first term represents the heat flow through the metallic contact, and the second term with the square brackets in the numerator represents the heat flow across the voids.
represents the thermal conductivity of the metals forming the interface, and the thermal conductivity of the fluid in the void. Likewise, is the squared root of the ratio of the real area of contact Ac to the total contact area A which is approximately = (2) Figure S7a shows the period length as function of the iteration number n for the ARB Cu-Nb composites. It can be seen that period length approximately spans over three scales; mm range for 1 ≲ ≲ 6, μm range for 7 ≲ ≲ 10, and nm range for 11 ≲ ≲ 13. Thus, only samples for ≳ 11 must experience significant ITR effects. Therefore, the drastic reduction observed in κ for the ARB Cu-Nb composites with ≲ 10 is linked to the TCR. Inset in Figure 7a shows an optical image for sample with = 1 iteration, contact points and voids are evident at the interfaces. Moreover, Figures 7b, c and d show SEM images for samples with = 2, = 4 and = 6 iterations, respectively. Clearly, when pressure on the interface is increased via ARB process, the points in contact are deformed and they increase both in size and number, as highlighted by the red arrows in the images for samples with = 2 and = 4. Hence, the interfaces between metals become gradually in total contact because of the reduction of the surface's roughness, as shown for sample with = 6.
To gain a quantitative understanding of these results, we have used the method proposed by W.M. Rohsenow [8]. In such an approach, the height of the void must be small compared with its width so as to neglect both radial and convection heat conduction in the voids.

Determination of the real area of contact Ac
The real area of contact Ac is directly obtained by adding the projected width on the plane of contact of all contact points for each directional profile and by multiplying the total contact width in one direction by the total contact width in the perpendicular direction. Thus, calling the width of the contacts in the x direction a1, a2, . . ., an, and in the y direction b1, b2, . . ., bn (see Figure S6), the real area of contact and the total contact area A are given by where lx and ly are the total lengths of the recorded profiles in each direction. Figure S8 shows the length estimations for the and for an interface of the sample with = 2 iterations. Figure S8. Size of the contact points in the x and y directions, respectively.
By using Equations (3) and (4), as well as the data extracted from Figure S8, a value of ~0.021 was calculated along with Equation (2).

Determination of the parameter
Furthermore, in Equation (1), is a parameter which expresses the equivalent idealized gap thickness in terms of the average thicknesses , and of the voids of actual surfaces in the x and y directions, respectively, as shown in Figure 6. The subscript i refers to metal 1 or 2, and j refers to the number of the void. Thus, the following expression for is defined by where is then calculated using the relations: = + (1 − ) and = + (1 − ). Figure S9 shows the measurements of the heights , and .

Sample with n=2 iterations
Height of the gaps estimations (μm) The summation extending over all the void segments xi and yi of the lx and ly profiles, respectively. Figure S10 shows the measurements of the width of the voids in each direction. By using Equations (6) and (7), as well as the data extracted from Figure 9 and Figure 10

Determination of the number of contact points per unit area n
Moreover, the number of contact points per unit area n is obtained by counting the number of contacts on each pair of profiles in the x and y directions, i.e., nxi and nyi. Since the two directions are perpendicular to one another, one deduces that over an area lxly, an approximation for n is ≅ (8) Figure S11. Number of contact points per unit area n.
By using Equation (8), as well as the data extracted from Figure 11, a value of ~6.9 10 11 −2 was calculated.
The term C is a factor which is only used to compact Equation (1), and it is given by By applying the values obtained for , n and , a value of ~0.573 is calculated. Here, ( ) and ( ) stand for the bulk thermal conductivity of the Nb and Cu, respectively.
Finally, with the parameters , n, and , and the factor , the magnitude of is calculated by using Equation (1). A summary of the parameters used for calculations of the thermal contact resistance of the sample with = 2 iterations is shown in Table S1. An effective thermal conductivity value of ~13.5 / was estimated using the value of = 4.7 × 10 / for the sample with = 2 iterations.
The above procedure was performed for samples with 1 ≲ ≲ 6 iterations; here, only the procedure for a sample with = 2 is illustrated.
The problem with Equation (3) Yi and Ji are the Bessel functions of the first and second kind. For practical purposes, the ratio of the squared root of the real area of contact to the total contact area must be ≲ 0.1. Under this condition, ( ) ≈ 1 , otherwise, ( ) drops drastically well below 1. The physical meaning of ( ) ≪ 1 implies that the real area of contact to the total contact are very similar, thus ≈ 1, i.e., no more voids or a very high contact points density; hence, the model fails in such a situation. To overcome such condition, the first term in Equation (1), with the square brackets in the numerator, can be neglected because it represents the heat flow across the voids. By using such mathematical artifice, Equation (1) becomes independent of ( ), and the model can be applied for samples with periods at microscale where the interfaces have much better contact, i.e., 7 ≲ ≲ 10.

Determination of parameters for modeling and calculations
In the EDMM approach, the interfacial thermal resistance is given by where = , and is the electronic heat capacity of the metal in side i, as well as being the electron Fermi velocity. Besides, → is the transmission coefficient which is given by Here, electronic heat capacity is and denotes the density of states of the metal in side i, which is given by Furthermore, the electron velocity on Fermi surface is Moreover, the mean speed of all the electrons can be calculated from the Fermi velocity by Hence, the mean free path in turn is By using Equations (13) to (17), the material properties for Cu and Nb were estimated. Data are shown in Table S2. Table S2. Calculated parameters for Cu and Nb using the free electron model. Additionally, the interfacial thermal resistance and the transmission coefficient were estimated via Equations (11) and (12), respectively, values are shown in Table S3. Table S3. Calculated data for the interface thermal resistance between copper and niobium interfaces.

Material
0.42 9.7 × 10 −11 m 2 ·K/W 9.7 × 10 −11 m 2 ·K/W 4.07 × 10 −11 m 2 ·K/W By analyzing the results, evidently, electrons move ballistically from the copper side to the niobium side, where then they scatter diffusively depending on the surface roughness. The effective thermal interface resistance is thus the sum of the contribution of both the electrons ballistically crossing the Cu-Nb interface ( ( → ) = ), and those diffusively scattering in Niobium ( ). Therefore, the effective ITR is given by is a term derived from ballistic diffusive equations [10] and gives the exponential decay of heat flux across interface, d is half the period thickness and λe is an effective mean free path of the electrons in the Cu-Nb system seen as an effective medium. Besides, the ballistic component is related to the diffusive component through

Minimum thermal conductivity calculations
Here, following Cahill and co-workers, we use "ultra-low" to describe a κ value lower than what is predicted by the minimum κ model developed by Cahill et al. Therefore, the existence of a lower limit to the thermal conductivity of disordered crystals based on the idea that lattice vibrations in those solids are essentially the same as those of an amorphous solid can be estimated by Equation (20) can be reduced to where and are the transversal and longitudinal speeds of the sound, respectively, na is the number of density of atoms, and the constant of Boltzmann. For comparison, Table 1 shows the minimum thermal conductivity values predicted by Equation (21), as well as the reported experimental thermal conductivity for amorphous bulk materials. Clearly, thermal conductivity value of 0.81W/mK at room temperature which is even lower than the amorphous lattice limit for the Cu-Nb thin film system and is determined to be well below that of the amorphous dielectrics like Al2O3 or SiO2.

Experimental interfacial thermal resistance estimations at nanoscale
The interfacial thermal resistance at nanoscale was estimated from the experimental ΔT-Q transfer curves. The sample thermal resistance can be expressed as where terms in brackets represent the thermal resistance due to the Cu and Nb films, whereas the last term represents the thermal resistance due to the total number N of interfaces. Evidently, the interfacial thermal resistance due to a single Cu-Nb interface is given by It is clear that Cu and Nb films are not monocrystalline; in this sense, it is well known that energy carriers experience severe size confinement and scattering at grain boundaries as the mean free path accomplish or even overcome the sample size [11].
For the sake of simplicity, by using a previously reported model [12] we estimated the sizedependent electron thermal conductivity in Cu(kCu-FILM) and Nb(kCu-FILM) films. Results are shown in Figure 12, and to confirm the validity of the modeling results, experimental thermal conductivity values for Cu were obtained from the existing literature [13,14], whereas size-dependent thermal conductivity data for Nb are completely absent in literature. However, due to the modeling results being consistent with the existing data for Cu, we assumed that the experimental data for Nb can be predicted with reasonable uncertainty. Based on this, the thermal conductivity data obtained for polycrystalline films of 20nm thick of Cu, and 30nm thick of Nb are kCu-FILM=13.7W/mK, and kNb- ·K/W; such value corresponds to the value of the interfacial thermal resistance for sample n = 13. By a similar procedure, the interfacial thermal resistance values were obtained for sample n = 11, and n = 12 from the experimental data of the corresponding slope of the Q vs ΔT plot. Figure S12. Size-dependent thermal conductivity for Cu and Nb polycrystalline films. Solid lines represent the results obtained from model reported in reference [12], whereas solid points represent the existing experimental data for Cu.

Uncertainty assessment
In connection to the uncertainty sources, a common method developed by Moffat [15] was used in the present research for estimating the overall uncertainty in the measurements. This method is based on the specification of the uncertainties in the various primary experimentally measured parameters. For instance, If the uncertainty present in Q−ΔT plots is a function of several measured parameters (P1, P2, …, Pn) each with a corresponding uncertainty (ΔP1, ΔP2, …, ΔPn), an overall uncertainty on the slope of the (ΔS) of the Q−ΔT plots using the method of Moffat can be estimated: Considering that, the parameters Pi experience small perturbations ΔPi, the partial derivatives can be approximated as: The fitting error ΔSi was performed assuming a t-distribution due to the small number of Q-ΔT points (Q,ΔT) taken at each measurement, below 10 points. This type of error including noise source errors ΔSnoise only accounts for the random error in the experiment and was simply added to the systematic error found from a sensitivity analysis. Table S6 summarizes the uncertainty values for each parameter that was varied in the analysis, the uncertainty in the slope of a typical Q−ΔT measurement can be determined based on the parameters' outline. Cu Experimental data [8] Cu experimental data [7] Cu k(d) modeling Nb k(d) modeling Cu-Bulk value Using the formulation ∆ , the uncertainty in the slope can be determined based on the overall system parameters outlined in the above table. The uncertainty was found to be ∆ ~ 6.56% for a determined Q−ΔT plot. Hence, the uncertainty in the measurements performed in the present research does not go beyond 10%. Based on this, uncertainty bars were added to the measurements present in the main manuscript.