Robust Modeling of Acoustic Phonon Transmission in Nanomechanical Structures

TL;DR

This paper compares two finite-element numerical methods for modeling acoustic phonon transmission in nanomechanical structures, extending understanding beyond low-frequency limits and enabling analysis of complex mesoscopic device geometries.

Contribution

It introduces and validates two robust numerical techniques for calculating phonon transmission in arbitrary structures, including mode conversion effects, beyond traditional analytic solutions.

Findings

Methods agree with thin-plate theory at low frequencies

Simulations accurately model abrupt junction transmission

Tools applicable for heat transfer analysis in nano-resonators

Abstract

The transmission of acoustic phonons is an important element in the design and performance of nano-mechanical devices operating in the mesoscopic limit. Analytic expressions for the power transmission coefficient, T, exist only in the low-frequency (quasi-static) limit described by thin-plate elastic theory, and for well-defined elastic wave-guiding geometries. We compare two numerical techniques based on finite-element computations to determine the frequency dependence of T for arbitrary phonon scattering structures. Both methods take into account acoustic mode conversion to acoustic and optical modes. In one case, phase and amplitude of complex-valued reflected waves are determined and related to transmission through a Fresnel equation, while in the other the magnitude of the transmitted mechanical power is directly calculated. The numerical robustness of these methods is demonstrated by considering the transmission across an abrupt junction in a rectangular elastic beam, a well-known problem of considerable importance in mesoscopic device physics. The simulations presented extend the standard results for acoustic phonon transmission at an abrupt junction, and are in good agreement with analytic predictions from thin-plate elastic theory in the long-wavelength limit. More generally, the numerical methods developed provide an effective tool for calculating acoustic mode energy loss in nano-mechanical resonators through mode conversion and heat transfer in arbitrary mesoscopic structures.

Equations

$$-\rho {\omega }^2\mathbf{u}=\nabla \cdot \mathrm{T}+\mathbf{F}$$

$${\mathbf{u}}_{\mathbf{s}}={\mathbf{u}}_{\mathbf{t}}(1+\gamma e^{i\varphi }),$$

$${\mathbf{u}}_{\mathbf{r}}={\mathbf{u}}_{\mathbf{t}}(1+e^{i\varphi }).$$

$$\gamma =\frac{\widetilde{u_s}-1}{\widetilde{u_r}-1},$$

$$\mathbf{I}=-\frac{1}{2}Re[i\omega {\mathbf{u}}^{\ast }\cdot \mathrm{T}].$$

$$P={\int }_{\mathcal{S}}\mathbf{I}\cdot \hat{\mathbf{n}}d\mathcal{S},$$

$$\mathcal{T}=\frac{P_{\mathrm{J}}}{P_{\mathrm{T}}},$$

$$G(T)=\frac{k_B^{2}T}{2\pi \hslash }\sum _{\alpha }{\int }_{x_c^{\alpha }}^{\infty }{\mathcal{T}}_{\alpha }(x)\frac{x^2\text{exp}(x)}{[\text{exp}(x)-1{]}^2}\text{d}x$$

Full text

I Introduction

The rate at which mechanical (phonon) energy dissipates via nano- and micro-mechanical support structures plays a critical role in determining the quality of an optomechanical deviceCole et al. (2011); Bochmann et al. (2013); Rostem et al. (2016). At low temperatures where thermal wavelengths are large compared to the minimal dimension of an elastic waveguide structure, phonon energy transfer is dominated by acoustic phonon modesRego and Kirczenow (1998); Schwab et al. (2000). The physical principles that govern phonon transport in mesoscopic structures are well understood; however, development of full mathematical models for the propagation of phonon modes is only tractable in limited casesCross and Lifshitz (2001); Photiadis and Judge (2004, 2007); Wilson-Rae (2008). Even the relatively simple geometry of a narrow rectangular beam joining a wide rectangular cavity cannot be modelled analytically for the full wave vector. The nearest mathematical analysis is that performed by Cross and LifschitzCross and Lifshitz (2001), which evaluates phonon transmission across an abrupt geometric boundary and produces analytic results only for the quasi-static long-wavelength limit. Photiadis and JudgePhotiadis and Judge (2004, 2007) expand on this work, analytically evaluating 3-dimensional cantilevered resonators, finding that neglecting thickness does in fact have a significant effect on phonon transmission. This suggests that even for a low-quality factor support structure Schwab et al. (2000), assuming 2-dimensional geometries may be an oversimplification. Furthermore, the analysis of Photiadis and JudgePhotiadis and Judge (2004, 2007) is only valid for long wavelengths. Thus, a method of accurately analyzing geometrically realistic structures beyond the quasi-static limit remains highly desirable.

Numerical modeling is one of the most powerful tools available for analysis of phonon transport in structures beyond the simplest cases. This paper explores numerical methods for determining acoustic phonon transmission probabilities across an abrupt junction in a rectangular beam (see Fig. 1) as a function of frequency using a finite element (FE) approach based on elastic wave theory. The numerical results are compared to the analysis performed by Cross and LifshitzCross and Lifshitz (2001), and explored in terms of phononic heat transfer in a self-suspended mesoscopic deviceSchwab et al. (2000). More importantly, the numerical methods established can be readily applied to any geometry and to a wide variety of elastic wave transmission calculations where acoustic modes are scattered into acoustic and optical modes.

II Methodology

The governing equation for elastic wave propagation expressed in the frequency domain isAuld (1973)

[TABLE]

where $\rho$ is density, $\omega$ is angular frequency, and $\mathbf{u}$, $\bm{{\rm T}}$, $\mathbf{F}$ are the displacement vector, stress tensor, and forcing vector. Numerical solutions to Eq. 1 are readily determined using commercial finite element software. However, using solutions of this formulation to calculate power transmission poses two challenges. First, the phase and amplitude of a propagating displacement wave must be extracted from the complex-valued displacement field solution. The phase-counting method discussed below outlines a procedure for achieving this, while the total power method avoids this issue by using the stress and velocity fields to calculate power transmission directly. Second, in most if not all available finite-element tools, excitation and detection of pure acoustic modes is not possible. Instead, finite geometric loads are typically used, which excite the appropriate acoustic mode, as well as evanescent modes. This implies that phase of the propagating mode cannot simply be assumed to be the product of the propagation constant of the mode and distance away from the source load. Quantification of the error caused by the interaction of reflected waves with the source load is discussed below.

Finally, we note that relative to time-domain solvers, the frequency-domain methods described here require significantly less computational resources for transmission calculations.

II.1 Complex Amplitude Reflection (Phase-Counting)

Figure 1a illustrates the reflection of elastic modes from an abrupt junction. To obtain the amplitude of the complex displacement field $\bf u_{s}$ at the reference plane $\Gamma$, an infinite beam of uniform cross-section (Fig. 1b) is modelled and the response $\bf u_{t}$ is obtained for a given excitation as a function of frequency. $\bf u_{t}$ provides a normalization for the finite-element model, which is always in the elastic limit but otherwise scaled in response depending on the magnitude of the source load, and the reactance of the elastic waveguide at the excitation frequency.

The displacement field at $\Gamma$ in the presence of the abrupt junction is related to $\bf u_{t}$,

[TABLE]

where $\gamma$ is the amplitude reflection coefficient at the cavity boundary, and $\phi$ is an arbitrary phase offset. To solve for $\gamma$, a second semi-infinite beam with one fixed end is also modeled (Fig. 1c), producing the response $\bf u_{r}$,

[TABLE]

The final expression for $\gamma$ is,

[TABLE]

where $\tilde{u_{s}}=\bf\|u_{s}\|/\|u_{t}\|$ and $\tilde{u_{r}}=\bf\|u_{r}\|/\|u_{t}\|$ are normalized displacement amplitudes evaluated point-wise at the reference plane $\Gamma$. In the absence of loss, the Fresnel relationship $\mathcal{T}=1-\left|{\gamma}\right|^{2}$ can be used to obtain the power transmission across the beam-cavity junction.

To quantify the interaction between the source load and reflected fields, a uniform beam with one free end is simulated Rostem et al. (2016); Gaiewski, Dunleavy, and A Castro (1997). The simulation accuracy can be quantified by measuring the deviation of $|\gamma|^{2}$ from unity. This test provides a robust procedure for refining the finite-element mesh until the desired accuracy is achieved. For the models presented here, an accuracy of $<\,$-30 dB for out-of-plane shear excitation, and less that -50 dB for the other three acoustic modes is readily achieved. This error can be further reduced with the use of a finer mesh, at the expense of increased computation time.

II.2 Total Power Transmission

An inherent assumption in the phase-counting method is in the symmetry of the elastic scattering region, which excites a single mode both in the beam and in the cavity. As such, the method as illustrated in Fig. 1 fails for geometries producing multi-moded reflections (e.g. an asymmetric junction). Therefore, a second model is developed to directly evaluate the total transmitted power. The average power (in W m*-2*) is given by the Poynting vectorAuld (1973)

[TABLE]

The total transmitted power is determined by evaluating the following integral over a surface $\mathcal{S}$ (see Fig. 1) in the cavity enclosing the geometric junction of interest,

[TABLE]

where ${\rm\hat{\bf n}}$ is the outward-pointing normal vector along $\mathcal{S}$. Power transmission is determined as

[TABLE]

where subscripts J and T refer to the Junction and Thru mechanical power transmission across the reference plane $\mathcal{S}$. This method makes no assumption about the type and number of modes transmitted across the elastic scattering region.

It is important to note that the normalization of the vector field $\bf u$ in Eq. 4 and scalar field $P$ in Eq. 7 is necessary to obtain accurate scattering coefficients from finite-element computations. At best, neglecting this step may produce severely skewed transmission profiles, and in most cases greater than unity transmission is observed Yu, Chen, and Wang (2013); Wang, Wang, and Zhang (2014); Lim et al. (2015).

II.3 Model Setup & Scaling

The model of a beam connected to a wider rectangular cavity has received much interest in optomechanical applications Schwab et al. (2000); Cross and Lifshitz (2001); Glavin (2001); Blanc et al. (2014). In the example pursued here, the beam is 300 nm wide with a square cross-section. The elastic material properties $E=280$ GPa, $\rho=3140$ kg m*-3* and $\nu=0.28$ are used to approximate silicon nitride Osman et al. (2014). In developing an FE model, the port geometry and mesh should be scaled to capture the propagating modes with the longest wavelengths of interest. In Fig. 3a, the longest wavelength that can be simulated is set by the beam and cavity lengths that are terminated by a perfectly matched layer. Given the dispersion relation for each section, the cavity length is typically much larger than the beam length in a FE model. Here a symmetrically loaded abrupt junction model is considered and only one mode is excited by the beam source and reflected by the abrupt junction.

Properly establishing a numerical model of the cavity is challenging due to the nontrivial relationship between the cavity modes and excitation frequency of the beam. Moreover, there is the possibility that multiple modes are simultaneously excited within the cavity Royer and Dieulesaint (1999).

To determine which modes are present in the cavity, an in-situ spectroscopic study is performed by exciting the beam at a range of frequencies for each of the four acoustic modes. Figure 2 shows the dispersion relation of the cavity modes and the cutoff imposed by the FE mesh and cavity length. Modes lying between these two limits are effectively captured by the FE model.

For every simulation represented in Fig. 2, a Fourier transform of the displacement field is taken along a longitudinal edge of the cavity. The power spectral density output by these transforms are plotted as a function of the beam source frequency. The results clearly show which cavity modes are excited by which beam modes. As is predicted by mode coupling theoryRoyer and Dieulesaint (1999), most cases show efficient coupling of a beam mode to a single cavity mode. In general, for any multi-moded elastic waveguiding structure, this type of study is needed to accurately model and capture the response of acoustic and optical modes.

III Results

The power transmission $\mathcal{T}$ of the four acoustic modes across an abrupt junction characterized by $B/b=5$ is shown Fig. 3. Additionally, transmission into a cavity of infinite width is also shown. All four modes show rejection at low frequencies and rapidly approach unity as frequency is increased and more modes propagate within the cavity. Results of both the phase-counting and total power methods are in good agreement.

Cross and Lifshitz Cross and Lifshitz (2001) derive analytic solutions for the transmission of the four acoustic modes across an abrupt junction using elastic beam theory in the long-wavelength limit. As shown in Fig. 4, the theoretical prediction and the results derived from the FE methods described here are in excellent agreement in the long-wavelength limit. In the case of the infinitely wide cavity, the analytic solutions quickly diverge from simulated results as frequency increases, which illustrates the effects of mode coupling to acoustic and optical modes. These effects cannot be captured using elastic beam theory in the long-wavelength limit, where the equations have simple plane wave-like solutions. It is notable that in Fig. 4, the total power method shows slightly closer agreement with theoretical predictions than the phase-counting method.

In general, the phase counting method demonstrated quicker convergence than the total power method and produces smoother curves at a given mesh size, as can be seen in Fig. 3. The main advantage of the total power method lies in its robustness as discussed in Section II.2. It also tends to show closer agreement with elastic theory in the long-wavelength limit (see Fig. 4).

IV Thermal Conductance

The abrupt junction is of particular interest in the context of heat transfer in mesoscopic structures Cross and Lifshitz (2001); Cole et al. (2011); Blanc et al. (2014); Puurtinen, Neittaanmaki, and Baskin (2012). Thermal conductance across a junction can be estimated by summing over the transmission of the propagating modes in the beam:

[TABLE]

where $x\equiv\hbar\omega/k_{B}T$, $\omega_{c}^{\alpha}=x_{c}^{\alpha}k_{B}T/\hbar$ is the angular cutoff frequency of mode $\alpha$, and $\mathcal{T}_{\alpha}$ is the power transmission through the junctionRostem et al. (2016). Evaluating this expression produces the relationship between thermal conductance and temperature presented in Fig. 5 for various values of $B/b$. The conductance of a uniform beam ($B/b=1$) is included as a point of reference, and at very low temperatures, is defined by the quantum of conductance Schwab et al. (2000). In evaluating Equation 8 we compute transmission for $\alpha=[1,4]$ up to $\omega_{c}^{5}$ in the beam. Past $\omega_{c}^{5}$, transmission is assumed to be unity. Optical modes of the beam up to 20 GHz are included in all calculations. As $B/b$ is increased, the thermal conductance quickly converges to the value obtained for a cavity of infinite width. A notable result from the inset in Fig. 5 is that the minimum in conductance in the 0.05 – 0.15 K range occurs at $B/b\approx 5$ due to the strong rejection of low-frequency waves at this cavity width, particularly for the in-plane shear mode.

V Conclusions

Two methods have been explored for obtaining robust simulations of acoustic phonon transmission in elastic waveguiding structures using FE computations. The accuracy of each method is evaluated using the model of an abrupt geometric junction. In the long-wavelength limit, the behavior of the power transmission is in close agreement with the available analytic results for a square beam connected to a wider cavity. The two methods can be successfully applied to compute the transmission profiles of acoustic modes in a beam at frequencies below the first optical mode. More generally and in the context of high sensitivity applications, the numerical methods described here provide versatile and reliable tools in the analysis of heat transport in phononic crystal structures, and energy storage in mesoscopic resonators.

Acknowledgements.

We gratefully acknowledge financial support from the NASA Astrophysics Research and Analysis (NNX17AH83G) grant. We thank Tuomas Puurtinen for useful discussions.