A differential thin film resistance thermometry method for peak thermal conductivity measurements of high thermal conductivity crystals

Thermal management has become a major challenge for microelectronic devices because the high-density electrical power is dissipated as Joule heating.^1–3 High operating temperatures at transistor hot spots can lead to performance degradation and thermomechanical failures in devices. Advanced thermal management approaches play an essential role in future progress of microelectronics. As an effort to find potential solutions to the mounting thermal management challenges, ultrahigh thermal conductivity materials, including diamond,^4,5 graphitic materials,^6–8 and semiconducting zinc blende boron arsenide (BAs),⁹ have been studied for heterogeneous integration with microelectronic devices.

Accurate thermal transport measurements are essential in the study of ultrahigh thermal conductivity materials. Several steady-state and transient measurement methods have been established for measuring the thermal conductivity of bulk crystals, thin films, and nanostructures. The standard guarded hot plate method¹⁰ is suitable for samples with a large cross section. The accuracy of this measurement method has been improved using a guard and numerical modeling to minimize edge heat losses and precisely control the temperature profile. For a high thermal conductivity sample, however, the thermal resistance across the thickness of the plate sample can be small compared to the contact thermal resistance between the thermometer and each end of the plate. In addition, it is a significant challenge to grow ultrahigh thermal conductivity single crystals with a sufficiently large size for guarded hot plate measurements.

Small-size, high thermal conductivity bulk crystals can be prepared into a rod sample for steady-state thermal conductivity measurement. In this approach,^11,12 one end of the sample rod is attached to a heater and a heat flux sensor made of a reference rod material with known thermal conductivity, while the other end is attached to a heat sink. The temperature drops along the heat flux sensor and the sample rod are measured with the use of either differential thermocouples or resistance thermometers attached to the sample with either thermal grease or pressure contact. Due to the presence of heat loss into these surface-mounted thermometers with electrical lead wires and the contact thermal resistance between the thermometers and sample, the thermometer temperature is different from the sample temperature. This difference can become appreciable compared to the small temperature drop along a high thermal conductivity crystal rod sample. The resultant measurement error is especially important at a low temperature where the peak thermal conductivity appears as a result of the interplay between intrinsic phonon–phonon scattering processes and extrinsic phonon-defect and phonon-boundary scattering processes.¹³ Although the peak thermal conductivity contains valuable insight into the competition between these intrinsic and extrinsic phonon scattering mechanisms, a prior steady-state measurement of small BAs rod samples was not able to obtain the peak thermal conductivity because of the uncertainty due to this measurement error as well as the decreased Seebeck coefficient of the differential thermocouples at low temperatures.¹⁴ 

Among the transient thermal transport measurement methods, laser flash measurement apparatus is commercially available for thermal diffusivity measurement of relatively large disk-shaped samples. In this apparatus, transient optical heating is applied to one side of the sample disk, while the temperature evolution of the other side is monitored using an optical non-contact measurement technique, such as infrared thermometry. It is often necessary to coat the two sample surfaces with a black graphite spray to confine the optical heating and temperature measurement on the sample surfaces. Although the transient measurement method is relatively insensitive to contact thermal resistance compared to the steady-state guarded hot plate measurement of a high thermal conductivity sample, the contact thermal resistance between the graphite spray and the sample can still become a source of error in laser flash measurement of a high thermal conductivity thin sample. In addition, the specific heat of the sample needs to be measured separately and used with the measured thermal diffusivity to obtain the thermal conductivity. Accurate measurement of the specific heat of small-size samples is a challenge itself at both high and low temperatures.

The time-domain thermal reflectance (TDTR) measurement technique has allowed for the extraction of both the sample effusivity and the thermal interface conductance of a thin film metal transducer coated on a bulk or thin film sample.^15,16 The specific heat of the sample needs to be determined from separate measurements or a theory for the conversion of the effusivity into the thermal conductivity. Because the thermal penetration depth of the ultrafast laser heating is usually smaller than tens of micrometers, the TDTR method probes thermal transport in the surface layer of a bulk sample. In high thermal conductivity crystals at low temperatures near the peak thermal conductivity, phonons with mean free paths longer than the penetration depth traverse the temperature gradient ballistically and do not contribute to the thermal conductivity measured by the TDTR.^17,18 Consequently, TDTR measurements have not obtained the peak thermal conductivity of high thermal conductivity crystals including diamond, BAs, and boron phosphide (BP).^14,19–22

Besides these four commonly used measurement techniques, several other thermal transport measurement techniques have also been reported. As one example, a thermal mirage measurement technique was used for earlier measurements of the thermal conductivity of diamond samples.²³ When the sample is heated by a laser beam, a temporarily varying gradient in the index of refraction in both the sample and the surrounding air is produced and can be detected from the measured deflection of a probe beam at a grazing incidence. This method yields the thermal diffusivity. The sensitivity and accuracy of this measurement method decrease at low temperatures because of the long thermal penetration depth in the sample, decreasing diffusivity of the air, diffusivity mismatch between the sample and the metal coating on top, and decreasing interface thermal conductance.²⁴ Due to the small temperature rise in the ultrahigh thermal conductivity sample and the limited temperature sensitivity, this and other measurement methods have not been able to obtain the peak thermal conductivity of high-quality diamond and other ultrahigh thermal conductivity materials.

Here, we report a steady-state method that is designed for peak thermal conductivity measurements of high thermal conductivity rod samples. This method is based on the deposition and patterning of thin film resistive heaters and thermometers directly on the sample surface to control the interface thermal resistance between the sample and thermometers. The electrical lead wires for the heater and thermometers are arranged to be outside of the sample region of interest to eliminate the heat loss through the thermometers and the error caused by the interface thermal resistance. Heat loss from the heater into the lead wires is obtained directly from two measurements with different numbers of lead wires attached to the heater. Moreover, a differential Wheatstone bridge circuit is employed to measure the small temperature drop between two resistive thermometers. Meanwhile, the steady-state heating is stepped according to a sinusoidal pattern to allow Fast Fourier Transform (FFT) analysis of the data to improve the signal to noise ratio. A silicon sample is used as the calibration standard to investigate the accuracy of this method. A comparison between the obtained measurement results and literature reports and theoretical calculations is used to evaluate the effectiveness of this method for peak thermal conductivity measurements.

Instead of attaching thermocouples to the sample surface with thermal grease or glue, thin film Au/Cr serpentine resistive thermometers are deposited and patterned on the sample surface that is covered by a thin silicon nitride (SiN_x) dielectric film, as shown in Figs. 1(a)–1(f). A thin film Au/Cr serpentine resistive heater is patterned on the suspended end of the sample. The other end of the sample is joined with silver epoxy to a copper heat sink. Because of the use of microfabrication to pattern the heater and thermometers directly and precisely on the sample bar, this approach is particularly useful for mesoscopic samples with a lateral dimension in the range between hundreds of micrometers and several millimeters.

The thermal resistance of the 300-nm-thick SiN_x dielectric layer under the thermometer line at environmental temperature T₀ = 300 K is calculated to be 91 K W⁻¹ as

where A_c is the 300 × 2 μm² cross-sectional area for vertical heat conduction across the dielectric between the thermometer line and the sample and $κSiNx$ is the thermal conductivity of the SiN_x film and reported to be about 5.5 W m⁻¹ K⁻¹.²⁵ In addition, the thermal interface resistance (r_i) for an evaporated or sputtered thin film is typically on the order of 3 × 10⁻⁸ K m² W⁻¹.²⁶ Thus, the interface resistance of contact area A_c is

The total contact thermal resistance (R_c) between the thermometer and the sample is ∼141 K W⁻¹ based on

Similarly, the thermal resistance across the metal layer is 0.7 K W⁻¹ calculated based on

where κ_m and t_m are the thermal conductivity and thickness of the metal film. In comparison, the resistance of radiation (R_rad) from the surface of the metal thermometer line is about 1.6 × 10¹⁶ K W⁻¹ based on

where ϵ_m and σ are the emissivity of gold and Stefan–Boltzmann constant, respectively. Hence, the (R_c + R_m) value is about fourteen orders of magnitude smaller than R_rad. Consequently, the metal film temperature is expected to be essentially the same as the sample temperature based on the thermal resistance circuit representation of the sample [Fig. 1(g)].

There are four contacts to the serpentine heater and six contacts to the two thermometers. The sensing current flows between two contacts through the two thermometers that are connected in series. The contacts are connected to a printed circuit board (PCB) with 25.4 μm-diameter bonding wires that are made of 99% Al and 1% Si. The bonding wires for the serpentine heater and two thermometer lines are arranged to be at the two ends of the sample to minimize surface heat loss from the region between the heater and thermometer 2.

The sample assembled on the heat sink and the PCB is loaded into a high vacuum cryostat with two radiation shields to reduce the surface radiation and heat loss to gas molecules. The heat sink is attached directly to the cryostat stage with screws to improve the temperature uniformity between the heat sink and the sample stage. The sample stage temperature is in the range from 350 K to as low as 12 K for these measurements. The measured DC heating current (I) and voltage drop (V) across the serpentine are used to obtain the heating power as Q = IV.

Instead of measuring the electrical resistances (R₁ and R₂) of the two thermometer lines to obtain their temperature rises and the temperature drop as the small difference between two large values, we use the differential Wheatstone bridge circuit shown in Fig. 2 to measure the temperature drop between them directly. In this circuit, the reference resistor (R_r) is significantly larger than the sum of all other resistors in series so as to limit the sensing current (i) in the two resistance thermometers as

where R_l is the resistance of the lead wires including the Al:Si bonding wires, R_L is the 50 Ω built-in output impedance of a SRS 830 lock-in amplifier, V_L is the voltage output of the lock-in amplifier, and V_r is the voltage drop measured by a lock-in amplifier across R_r. Here, the leakage current through the lock-in amplifier is ignored because its 10 MΩ input impedance is much larger than R_r = 100 k Ω.

In Fig. 2, the bridge circuit voltage signal is related to the electric resistances of the two thermometers according to

where g₁ and g₂ are the gains of the two SR560 voltage preamplifiers, which are used to measure the voltage drops (V₁ and V₂) across the two thermometers, R₀₁ and R₀₂ are the built-in output impendence of the two preamplifiers, R_1a, R_1b, and R_2a are three precision resistors, and R_2b is a potentiometer to be adjusted to minimize the following zero-point offset of the Wheatstone bridge output prior to the application of a heating current,

Here, V_s,0 and i₀ are the V_s and i values and R_1,0 and R_2,0 are the resistances of the two thermometer lines at Q = 0.

Under steady-state Joule heating in the serpentine heater, the temperature rise (ΔT₁ and ΔT₂) of the two resistance thermometers leads to the following voltage signal:

where Δi ≡ i − i₀, ΔR₁ ≡ R₁ − R_1,0, and ΔR₂ ≡ R₂ − R_2,0. In addition, $α1≡dR1,0R1,0dT$ and $α2≡dR2,0R2,0dT$ are the measured temperature coefficients of resistances (TCR) of the two thermometers. Thus, the temperature difference between the two thermometers can be obtained as

Due to the much larger R_r than ΔR₁ and ΔR₂, $Δii0Vs,0$ is much smaller than ΔV_s. The use of the Wheatstone bridge to nullify V_s,0 helps to make $Vs,0ii0α2θ2$ much smaller than ΔV_s. When α₁ and α₂ are close to each other for the two thermometer lines made of the same metal, the second term proportional to $α1−α2$ is small so that ΔV_s becomes the dominant signal source in the δθ expression. The resistance of the second thermometer adjacent to the heat sink is measured with the use of a lock-in amplifier to obtain θ₂. Meanwhile, the heater temperature rise $θh$ is obtained from the change in the four-probe resistance obtained from the measured DC current–voltage (I–V) curve.^27,28 Joule heating in the Al:Si bonding wires of the heater can be determined based on the voltage drop (V_w) measured across the Al:Si bonding wires with a voltage preamplifier in a four-probe configuration.

Based on the measured temperatures, the thermal conductivity of the sample can be obtained from a radiative fin model, in which the heat loss to the bonding wires of the heater and from the sample surface are accounted for. In this model, the total heat flow in the sample under the heater toward thermometer 1 is

where Q_end is calculated as the radiation heat loss from the suspended end surface of the heater segment and Q_loss,i is the heat loss into each of the four bonding wires. To evaluate the heat loss into the bonding wires, we treat the bonding wire as a fin with known temperatures T_h and T_l at two ends (x = 0 and x = l). The steady-state heat diffusion equation for such a fin is²⁹ 

where κ_w, A_w, P, h, and $q̇$ are the thermal conductivity, cross-sectional area, perimeter, surface heat transfer coefficient, and Joule heating density of the Al:Si bonding wire. For the two bonding wires with Joule heating, $q̇=VwIlAw$⁠, whereas $q̇=0$ for the other two wires without Joule heating. The solution of this fin equation yields the heat loss from the sample to the bonding wire as

where $θi≡Ti−T0−q̇AhP$ for i being h or l, fin parameter $m=hP/κA$ with $h=4ϵσT∞3$⁠, ϵ being the emissivity, and T₀ being the sample stage temperature that is assumed to be the same as the temperature of the inner radiation shield mounted on the sample stage.³⁰ The heat loss difference between two identical bonding wires with and without internal Joule heating in the wire is thus

In this model, the thermal conductivity of the bonding wire can be evaluated according to Wiedemann–Franz law from the measured electrical resistance (R_w) as

where L_o is the Lorenz number. The surface emissivity of the aluminum wire can be calculated from its resistivity (ρ_w) as³¹ 

where the wavelength λ is chosen as $λ=bT0$ with b ≈ 2898 µm K based on Wien’s displacement law.

The heat flow across the SiN_x film under the surface heater is much larger than that under each thermometer line. At low temperatures, the thermal resistance for heat flow from the surface serpentine heater across the SiN_x film into the high thermal conductivity sample is appreciable compared to the thermal resistance of the sample. Thus, the temperature rise (θ_b) in the sample right underneath the serpentine heater is lower than the measured heater temperature rise (θ_h). In particular, there exists a three-dimensional (3D) temperature distribution in the sample in the vicinity of the surface heater. Because the surface radiation thermal resistance is many orders of magnitudes larger than the thermal resistance of the SiN_x film, in comparison, the temperature distribution in the bar sample is expected to become one-dimensional (1D) at a region sufficiently far away from both the surface heater and the copper heat sink. For our samples, the distance between each of the two thermometer lines and the surface heater or the copper heat sink is about ten times the sample thickness. As such, 1D temperature distribution is established in the segment between the two thermometers to yield the same temperature of each thermometer line and the sample right beneath the thermometer line.

By treating the segment of the sample from the heater to the second thermometer as a fin, we use the measured temperature rise (θ₁ = θ₂ + δθ and θ₂) at the first and second thermometer lines to calculate²⁹ 

where L₁ is the distance between the heater and the first thermometer line, L₂ is the distance between the two thermometer lines, and m_s is the fin parameter of the sample and calculated with emissivity values for Si³² and SiN_x. For the semi-transparent thin film SiN_x layer, the emissivity can be calculated based on the real part of the refractive index and the extinction index of the SiN_x thin film at each temperature.^33,34

The fin heat transfer rate from the heater toward the two thermometers is

which can be used to obtain the sample thermal conductivity κ_s. The obtained κ_s is the weighted average between the silicon bar, the SiN_x film, and the metal thermometer layer. Due to the much smaller cross section of the SiN_x film and the thin metal film, the difference between κ_s and the silicon thermal conductivity is well within 0.01%.

As a comparison, a pure heat conduction model is also employed to analyze the measurement data to obtain the sample thermal conductivity as

where heat loss to the bonding wires is similarly accounted for. To ensure that the heat loss through the bonding wires is small compared to Q_in as suggested by the model calculation, in addition, we removed the two bonding wires 3 and 4 that were used to measure the voltage drop of the heater and repeated the measurements to obtain the temperature drop (δθ′) between the two thermometers. Based on the heat conduction model, removal of these two bonding wires reduces the heat loss by

The measured $δθ′−δθδθ$ is within 2.5% at above 20 K, below which the uncertainty in the measured δθ′ and δθ increases due to reduced TCR as a result of the Kondo effect associated with electron scattering by magnetic impurity in Au.³⁵ 

The accuracy of this method is investigated with the use of the silicon bar sample shown in Fig. 1 as the calibration standard. The silicon bar sample was diced from an undoped silicon wafer. A 320-nm-thick SiN_x thin film was sputtered on top of the sample for electrical insulation of the resistive thermometer from the semiconducting silicon sample. Both the heater and thermometers were then patterned from a 10-nm-thick chromium adhesion layer and a 100-nm-thick gold film deposited on the sample surface. Scanning electron microscopy (SEM) measurements show that the silicon sample is 624.4 μm wide, 514.4 μm thick, and the length between the two thermometer lines is 5.1 mm.

Figure 3 shows the measurement results of the silicon bar sample. The measured electrical resistance at the sample stage temperature without electrical heating is shown in Fig. 3(a) and can be used to obtain the TCR, which is then used to convert the measured electrical resistance rise at each heating current to the temperature rise. Figures 3(b) and 3(c) show the measured voltage drop of thermometer 2 and the bridge voltage signal as a function of the heating current. The results can be used to obtain the electrical resistance change of thermometer 2 as shown in Fig. 3(b) and temperature rise (θ₂) and temperature difference (δθ) shown in Fig. 3(d) at the stage temperature T₀ = 35 K, respectively.

At each temperature, the heating current from a DC current source is stepped to take 16 discrete steady values according to a sinusoidal pattern as shown in Fig. 4(a), which is repeated for 20 cycles. Each discrete and steady current is held for 12 s to allow the thermal response to reach the steady state and be measured. Each cycle takes 192 s to yield an equivalent modulation frequency (f_m) of about 0.005 Hz. FFT analysis is used to obtain the modulation amplitudes of the measured I, V, ΔR_e,2, and δθ, as shown in Fig. 4. The FFT amplitudes are used to calculate the second harmonic amplitudes of Joule heating (Q), temperature rise of thermometer 2 (θ₂), and temperature difference (δθ) with Eq. (12). Figures 4(d)–4(f) show clear modulation at the second harmonic component corresponding to the modulation frequency (f_m) of the heating current. These second harmonic amplitudes are used in the model calculation of the thermal conductivity. Compared to averaging the measurement results of the 20 cycles of data, this FFT approach can improve the signal to noise ratio when the temperature drop along the sample becomes small near the peak thermal conductivity.²⁸ 

However, both the FFT spectra in Fig. 4(f) and the average temperature rise data in Fig. 3(d) show higher random noise in δθ measured by the Wheatstone bridge with unitary gain set for the two differential amplifiers than in θ₂ measured by a lock-in amplifier. For this measurement where the sample stage temperature was highly stable at 35 K, the noise floor in the measured θ₂ is as small as about 5 × 10⁻⁴ K. With either an increased stage temperature or reduced time for stabilizing the stage temperature, the noise floor in the measured θ₂ can exceed 3 × 10⁻³ K, as shown in Fig. 5. Consequently, the noise floor can exceed 4 × 10⁻³ K in the calculated (θ₁ −θ₂), where both θ₁ and θ₂ are obtained from the voltage drops (V₁ and V₂) measured directly across the two thermometer lines. With the gain for the two differential amplifiers increased from 1 to 100, in comparison, the noise floor in δθ measured by the bridge circuit can be reduced to below 2 × 10⁻³ K. This comparison indicates that the resolution of the differential measurement method can be enhanced to the level of 1 × 10⁻³ K for measuring the small temperature drop along a high thermal conductivity material, even when the temperature fluctuation in the sample stage makes it impossible to obtain this resolution by measuring θ₁ and θ₂ directly to obtain (θ₁ − θ₂).

Figure 6 shows the thermal conductivity measurement results based on the radiative fin model and the conduction model. The uncertainty in the thermal conductivity was calculated by a Monte Carlo simulation based on the uncertainties in the various measurements.²⁸ For each measured modulated parameter such as the ΔV_s/i term in Eq. (12), the local noise of the FFT amplitude spectrum is treated as a folded normal distribution to calculate the standard deviation of the distribution. The obtained FFT peak amplitudes are perturbed with normally distributed random generator with the corresponding standard deviations in the Monte Carlo calculation of δθ from Eq. (12). While only random uncertainty is shown for δθ and θ₂ in Fig. 3(d), additional systematic errors are accounted for by the Monte Carlo simulation. Specifically, the preamplifier gain (g₁) is similarly perturbed with a normally distributed random generator to account for the 1% gain error. In addition, the dominant systematic uncertainty source for δθ is the uncertainty in the measured $α1−α2$ that shows a slight dependence on the actual temperatures (T₁ and T₂) of the two thermometer lines. At each stage temperature T₀, T₁, and T₂ differ from T₀ by a small amount depending on the heating level. Thus, in the Monte Carlo calculation, the thermometer temperature is perturbed with a uniformly distributed random generator within this amount to determine α₁ and α₂ corresponding to the perturbed thermometer temperature. The dimensions are measured several times to obtain the standard deviations of the dimensions. The perturbated values of δθ, Q, and the dimensions are then used to calculate the thermal conductivity based on Eqs. (19) and (20). This process is repeated for 5000 times to achieve independency on the number of Monte Carlo simulations. The uncertainty in the thermal conductivity is calculated based on the standard deviation of the 5000 calculated values of thermal conductivity, at 95% confidence interval.

At temperature above 20 K, the difference between the conduction model and the radiation fin model is within 6.2%, which includes less than 3.3% radiation loss from the sample surface and less than 3.1% due to the difference between the fin model calculation of the heat loss into the bonding wires and the measurement with different numbers of bonding wires. At low temperatures, the radiation loss becomes negligible. At temperatures from 400 to 100 K, the differences between the current measurement results and literature results^36,37 for larger size silicon samples are within 5% and 10%, respectively. At low temperatures, the increased surface-phonon scattering processes in the smaller size sample measured here lower the peak thermal conductivity and increase the peak temperature compared to the literature reports.

To verify that the lower peak thermal conductivity measured for the mesoscale silicon sample than the larger bulk sample is due to boundary scattering instead of measurement errors, we incorporate appropriate boundary scattering and defect scattering into three-phonon scattering in a theoretical calculation of the lattice thermal conductivity. This first-principles calculation is based on the density functional theory as implemented in Quantum Espresso,^38,39 which yields results in good agreement with measurement results of a number of pure crystals without the use of fitting parameters.^40,41 Norm-conserving pseudopotentials and the Perdew and Zunger⁴² exchange-correlation functional were used. An energy cutoff of 70 Ry was taken, which was sufficient to ensure that forces were accurately described. After structural relaxation, the calculated lattice constant was 5.40 Å, close to the measured value of 5.43 Å. Harmonic interatomic force constants (IFCs) were calculated using density functional perturbation theory, from which phonon dispersions were obtained, which showed good agreement with the experiment. Third-order anharmonic IFCs were calculated using a finite displacement approach with a 250 atom supercell. The harmonic and anharmonic IFCs were used to calculate three-phonon scattering rates. Here, a mass disorder model with mass variance parameter, g₂,⁴³ was used to describe phonon scattering by the naturally occurring silicon isotope mixture as well as by other point defects. Phonon scattering from sample boundaries was described by an empirical scattering rate formula, $τb,λ−1=vλ/Lb$⁠, where v_λ is the magnitude of the phonon velocity in mode λ, while L_b is an averaged sample size. The linearized Peierls–Boltzmann transport equation was solved using our home-grown code to obtain the silicon thermal conductivity. Both g₂ and L_b are taken as adjustable parameters and used to give a best fit to the measured thermal conductivity data over the range of temperatures studied. The measurement data can be fitted with the use of L_b = 584 μm and g₂ = 2.07 × 10⁻⁴. As the L_b value is close to the (wt)^1/2 = 567 μm of the silicon sample with width w = 624.4 μm and thickness t = 514.4 μm, diffuse surface scattering can explain the measured lower peak thermal conductivity of the mesoscale silicon sample than literature reports^36,37 of larger silicon crystals.⁴⁴ 

These results demonstrate the accuracy and effectiveness of this differential thin film resistance thermometry approach for peak thermal conductivity measurement of high thermal conductivity materials. Direct deposition and patterning of the thin film thermometers on the sample surface allow for judicious arrangement of the electrical lead wires to essentially eliminate both the heat loss from the sample into the thermometers and the corresponding error associated with the interface thermal resistance. Measurements with different numbers of electrical connection wires to the heater directly yield the heat loss into the wires. The small temperature drop along a high thermal conductivity sample is resolved with a combination of Wheatstone bridge measurement of the responses of two resistive thermometers and FFT analysis of the steady-state data obtained at different heating rates. As this method is validated here with a silicon calibration standard, we expect that this method will find increasing use in the ongoing search of high thermal conductivity materials for thermal management, including diamond, BAs, BP, boron nitride (BN),⁴⁵ BC₂N,⁴⁶ and theta-phase tantalum nitride (θ-TaN).⁴⁷ 

This work was supported by Office of Naval Research MURI Award No. N00014-16-1-2436.

The data that support the findings of this study are available from the corresponding author upon reasonable request.