Non-Monotonic Variation of Acoustic Spectrum with the Mass or Thickness of a Layered Structure

Abstract

We are examining the behavior of resonance frequencies and their response to variations of material parameters such as thicknesses, masses, and bulk velocities for certain Rayleigh–Lamb acoustic modes in a multilayered structure. The treatment is based on recent explicit analytic solutions that have allowed us to explore the entire parametric space using dimensionless ratios. This exploration has revealed a complex parametric dependence of the phase velocities and their mass loading response. Specifically, for the fundamental flexural modes in a bilayer, we have shown that both quantities change in a strongly non-monotonic way with thickness, density, or bulk velocity ratios. Even in the regime of thin coating, commonly encountered in acoustic sensing applications, we have found important differences from previously known results, e.g., that response to loading may switch its sign multiple times when the velocity of the deposited material is increased. We have also discovered that the fundamental dilatational modes can be highly effective in stabilizing resonant frequencies against even large variations of the thickness or mass of the exposed layer. This property is demonstrated in an explicit form by the derived expression for the mass coefficient of frequency for an arbitrary number of layers.

1. Introduction

Due to the ubiquitous nature of acoustic vibrations and a direct dependence on the properties of the propagation medium, ability to interact with other excitations and need to keep a low attenuation rate, they demonstrate a broad diversity of applications ranging from the macro- to nanoscale, showing classical or quantum behavior In the commonly encountered bounded layered structures like thin films, membranes, plates, etc., acoustic spectra are described by the waves of the Rayleigh–Lamb (RL) type [1,2,3,4,5] and consist of an involved arrangement of frequency–wavenumber dispersion curves ${\omega }_{n,p}\left(\mathbf{q}\right)$, where q is the in-plane wavenumber, $n=0,1,2,\dots$ is the band index and p is the polarization. At long wavelengths, the curves terminate at finite cutoff values except for a few fundamental branches, such as $n=0$ and the gapless shear horizontal, flexural (F) and dilatational (D) modes (we have dropped the p subscript). At shorter wavelengths, most of the curves converge towards bulk shear and longitudinal velocities and some evolve into evanescent surface or interface-guided modes, including the famous Rayleigh surface acoustic waves. However, knowledge of the band structure is insufficient for applications where one also needs to know the way it changes with various parameters. For example, consider semiconductor or metal-insulator wafers with a thickness of a few dozen nanometers that are the main components of electronic ultra-low power detectors or on-chip microcoolers [6,7]. One of the key processes is the transfer of heat carried by electrons and its dissipation by the RL phonon modes of the numerous dispersion branches, which in turn depend on a set $\left\{{\xi }_m\right\}$ of material and structural parameters, consisting of thicknesses $h_{\gamma }$ of the layers $\gamma =1,2,\dots$, bulk shear $s_{\gamma }$ and longitudinal ${\mathcal{l}}_{\gamma }$ velocities, elastic constants, mass densities ${\rho }_{\gamma }$, etc. [8,9]. The intensity of the electron–phonon heat transfer can be increased by an order of magnitude by reducing layer thicknesses; however, there exists a large variation between specimens of different materials [10]. Because any particular set of parametric values corresponds to just a single point in a multidimensional space, there exists a fundamental and practical interest in finding more explicit dependence of the band structure on parametric coordinates using analytic or semi-analytic results. For instance, in [11], it has been shown that by increasing the number of layers with an alternating ratio of elastic constants, the rearrangement of the band structure can result in the accumulation or depletion of the density of states in a specified spectral region. This property is relevant for, e.g., thermal transport and for tuning of electron–phonon and magneto-acoustic coupling in microelectronics and spintronics [12,13]. On the other hand, many applications (e.g., in acoustic sensing) spotlight the dependence of separate RL modes on external factors due to their effect on ${\xi }_m$ [4,14]. This principle is implemented in many types of modern microsensors on the basis of thin layer technologies measuring mass, temperature, pressure, humidity, gas concentration, etc. [15,16,17,18,19,20]. The relative simplicity of the layered structure gives advantages in fabrication and ensures low power consumption, excellent interfacing with electronics, compatibility with wireless data transmission and the ability of passive (battery-less) operation in both harsh and biological environments. Thus, such sensors using Rayleigh or Sezawa modes are successfully used in industry, biology, chemistry and medicine [21,22,23]. Analyzing the frequency response of RL cutoff modes are instrumental in the detection of defects or determination of the degree of corrosion in nondestructive testing [24]. The effect of mass loading on the frequency shift generated by metallic electrodes in film bulk acoustic resonators are used in the operation of modern 5G technology pass band filters [25]; the high-sensitivity resonators of this kind can detect masses down to femtograms [26].

Sensitivity, the key characteristic of such devices, is quantified by the shift of resonance frequency in response to an external factor (e.g., temperature, T) and can be defined in the form of the temperature coefficient of frequency (TCF):

$$TCF={\displaystyle \frac{dln{\omega }_n}{dT}}.$$

(1)

Then, the relative variation of ${\omega }_n$ or of the phase velocity $C_n$ is driven by the temperature dependence of ${\xi }_m$:

$$TCF\left(\xi \right)={\displaystyle \frac{\partial ln{\omega }_n\left({\xi }_m\right)}{\partial {\xi }_m}}{\displaystyle \frac{d{\xi }_m}{dT}}={\displaystyle \frac{\partial C}{C\partial {\xi }_m}}{\displaystyle \frac{d{\xi }_m}{dT}}-\alpha .$$

(2)

Here, $\alpha$ is the linear thermal expansion coefficient due to the elongation of the device in the direction of wave propagation, i.e., it accounts for the change in the excited wavelength. The derivatives of $\xi$ are related to bulk thermal coefficients of, e.g., thickness, elastic constants, sound velocities, etc. [27,28], and similar considerations apply to pressure, humidity or other types of sensors. Therefore, to maximize the sensitivity, one needs to know the rate $\partial {\omega }_n/\partial {\xi }_m$ at which the acoustic band structure changes with parameters. For some applications, on the contrary, a reduced sensitivity is required to ensure the frequency stability of the resonance mode against a change in, e.g., temperature, pressure or mass [29]. Respectively, the condition $\partial {\omega }_n/\partial {\xi }_m=0$ identifies surfaces of relative stability towards variations in parameters. On the other hand, these surfaces are closely related to the monotonicity property of the $\omega -\xi$ dependence, i.e., the existence of a one-to-one correspondence between parameters and resonance frequencies, a basic characteristic of measurement and sensing devices. The non-monotonic behavior is signaled by the presence of maxima or minima of the ${\omega }_n\left({\xi }_m\right)$ with the condition that the second derivative does not switch its sign on the critical surface. These critical surfaces are therefore of significant interest both for applications exploiting the stability of the resonance frequencies and acoustic sensing, as higher values of $\partial {\omega }_n/\partial {\xi }_m$ and higher sensitivities are achieved in the regions further away from the surfaces. For instance, a nearly linear behavior of resonance frequencies with temperature, pressure or gas concentration is reported for sufficiently large working intervals in high-performance acoustic sensors based on Rayleigh, Sezawa and some Lamb waves [16,30,31,32]. The monotonic dependence on parameters has also been found in an analytical form by the perturbative calculation of the mass loading sensitivity in a two- and three-layer isotropic composite plate [15,33,34]. In particular, for sensors based on $SH$, fundamental dilatational and flexural waves, the sensitivity increases for thinner and lighter structures or with a higher contrast between the velocities of the layers.For a single layer, monotonic thickness dependence of the whole spectrum is well established, e.g., the phase velocities of the symmetric modes are decreasing functions of H, as illustrated in the Figure 1 for an iron plate.

For a bulk medium, parametric monotonicity can be understood by considering a transverse (or longitudinal) wave traversing a sequence of layers with respective velocities $s_{\gamma }$ (or ${\mathcal{l}}_{\gamma }$ ). For the normal incidence, the effective velocity is determined by considering the total propagation time [35]:

$$C_{eff}={\displaystyle \frac{H}{{\sum }_{\gamma =1}{h}_{\gamma }/s_{\gamma }}},$$

(3)

where $H={\sum }_{\gamma =1}{h}_{\gamma }$ is the total thickness. Such elementary consideration turns out to be sufficiently reliable for estimating the crossover temperature of a thermal phonon gas from bulk to the confined (or quasi-2D) regime, even for nanometrically thin multilayers [6]. Note that the number of relevant independent parameters can be reduced by switching to scaled quantities, e.g., taking $\gamma =1,2$ in Equation (3) gives $C_{eff}/s_1=\left(1+\delta \right)/\left[1+\delta /\left(s_2/s_1\right)\right]$, with just two independent parameters, $s_2/s_1$ and $\delta =h_2/h_1$. This illustrates the physically transparent relation: $C_{eff}$ unconditionally increases with $s_{\gamma }$ as a faster material of any layer means a faster overall propagation. Also, $C_{eff}$ either increases or decreases with a layer thickness depending on the given set of values for the other parameters (see below).

Nevertheless, the above monotonicity property can not be generally true, even for the spectrum of an isotropic traction-free plate, as recently demonstrated in [36]. One notable implication of the aforementioned exceptions is that acoustic plate waves can propagate with a higher speed in a slower medium (i.e., with a lower bulk velocity). The $S_1$ modes of the calculated band structure in the Figure 2 are an example of such an unconventional parametric dependence that can be compared with the totally regular one in Figure 1. It should also be mentioned that due to the large dimensionality of the parametric space, the frequency–wavenumber $\omega \left(q\right)$ spectrum is known to have a number of different anomalies in this extended space, e.g., “diabolic” and exceptional crossing points, avoided crossing (osculation or veering) points related to degeneracy of quasi-degeneracy of the spectral surfaces [37,38,39,40]. The latter crossings have a conical character in the vicinity of such special points that may well involve a non-monotonic variation of frequency with $\xi$ [36,37]. The non-monotonicity of $\omega$ versus q dependence is also associated with anomalous zero and negative group velocity (ZNGV) modes [41,42,43,44], which have recently attracted special interest due to negative refraction and its application in the design of metamaterials or ultra-resolution acoustic lenses [45,46].

However, non-monotonicity is not necessarily related to singularities or degeneracy of spectral lines and may take place in a large area of the parametric space, as has been known from calculations and experimental studies for particular materials and parameters of multilayer structures [41,47]. The analytical solutions considered below allow us to obtain and analyze explicit expressions for the response functions in the multidimensional parametric space and reveal some general properties of the RL modes relevant to sensing applications.

2. Parametric Non-Monotonicity of the RL Spectrum

In the short wavelength limit, $Q=qH\gg 1$, the spectrum becomes non-dispersive, group and phase velocities coincide, dependence on total thickness H disappears and, by analogy with bulk media, the monotonic parametric behavior can be expected. Indeed, this is easily verified by the asymptotic approximation in [48,49] for the velocity of the Rayleigh SAW in a single-layer plate. Its parametric derivatives are all strictly positive and only vanish at $\nu \simeq -1.3$ or $\nu \simeq 0.86$ outside the physical domain (where $\nu$ is Poisson’s ratio). The sensitivity (1) is then found as:

$$TC{F}_R={\displaystyle \frac{\partial s}{s\partial T}}+f\left(\nu \right){\displaystyle \frac{\partial \nu }{\partial T}}-\alpha ,$$

(4)

with $f\left(\nu \right)\simeq 0.1-0.2$ (the lengthy explicit expression is not shown for brevity). Theoretical arguments show that the temperature dependence of bulk parameters in many solids is dominantly linear in sufficiently large temperature or pressure intervals [27,28] so that, according to (4), one can expect that values of the $TCF$ or $PCF$ (pressure coefficient of frequency defined by analogy with (1)) for the Rayleigh mode remain nearly constant at sufficiently short wavelengths. However, there are limited data on the measured temperature or pressure dependence of parameters such as s and $\nu$ in the literature, and it therefore could be possible to use similar explicit expressions for different modes in obtaining such complementary information.

In the opposite, long wavelength, limiting the fundamental dilatational mode appears to be the physically closest generalization of the bulk waves in multilayers (3) to finite H because of its linear frequency–wavenumber dispersion. Then, the general expression for the velocity of the $D-$ mode in a traction-free multilayer [11] can be expressed in a form reduced to dimensionless scaled quantities:

$$C_D\left(N\right)/V_1={\left({\displaystyle \frac{{\sum }_{\gamma =1}^N{\varrho }_{\gamma }{\delta }_{\gamma }{\sigma }_{\gamma }^2}{{\sum }_{\gamma =1}^N{\varrho }_{\gamma }{\delta }_{\gamma }}}\right)}^{1/2},$$

(5)

where

$${\varrho }_{\gamma }={\rho }_{\gamma }/{\rho }_1,{\delta }_{\gamma }=h_{\gamma }/h_1,{\sigma }_{\gamma }=V_{\gamma }/V_1,$$

(6)

N is the total number of layers, $V_{\gamma }=2s_{\gamma }\sqrt{1-J_{\gamma }}$ is the velocity of the fundamental symmetric $S_0$ Lamb mode in a single free plate $\gamma$ and $J_{\gamma }={\left(s_{\gamma }/{\mathcal{l}}_{\gamma }\right)}^2$. Alternatively, related parameters can be used, $s=\sqrt{\mu /\rho },$$\mathcal{l}=2\sqrt{\left(\lambda +2\mu \right)/\rho },m=\rho h,{\nu }^{-1}=2\left(\lambda +\mu \right)/\lambda =\left(1-J\right)/\left(1/2-J\right)$, where $\lambda$ and $\mu$ are the Lamé elastic constants (the $\gamma -$ index is omitted for brevity wherever confusion can be avoided). By following the arguments in the above discussion of bulk waves, we represent (5) in a suggestive form, revealing its similarity to the parametric dependence in (3).

$$C_D^2\left(N\right)={\displaystyle \frac{{\sum }_{\gamma =1}^N{m}_{\gamma }{V}_{\gamma }^2}{M}},$$

(7)

where $M={\sum }_{\gamma =1}^N{m}_{\gamma }$ is the total mass per unit area ($m_N$ exposed as top layer). Indeed, here we have an additive or linear superposition of mass and “kinetic energy” contributions that demonstrates a completely monotonic behavior of the $D-$ mode in the multidimensional parametric space. This form is convenient for the analysis of mass loading response, and we readily obtain a generalization of the earlier result [15] to the sensitivity function of a structure with an arbitrary number of layers:

$$MC{F}_D\left(N,\xi \right)={\displaystyle \frac{\partial ln{C}_D}{\partial m_N}}={\displaystyle \frac{{\sum }_{\gamma =1}^{N-1}{m}_{\gamma }\left(V_N^2-V_{\gamma }^2\right)}{2M{\sum }_{\gamma =1}^N{m}_{\gamma }{V}_{\gamma }^2}}.$$

(8)

Importantly, the condition $MCF=0$ does not involve the mass of the exposed layer, and the equation explicitly demonstrates how one can achieve a perfect mass loading stabilization effect for the dilatational wave in a multilayer. For a bilayer, this possibility is irrelevant because it imposes an unrealistic condition that the two materials are of the same velocity $V_1=V_2$. However, for a larger number of layers, this restriction is not necessary, and mass stabilization of the resonance frequency can be achieved by an appropriate adjustment of the thickness or mass ratios for the layers other than $\gamma =N$. Similar to bulk waves, the physical cause of such a cancellation is that in this way, effective velocities of the sub-sets of N layers can become equal to each other, as follows from (7) and (8). For example, when $N=3$, the condition $MCF\left(N=3\right)=0$ can be satisfied by choosing the masses or thicknesses of the layers 1 and 2 in such a way that:

$$m_2/m_1=\left(V_1^2-V_3^2\right)/\left(V_3^2-V_2^2\right).$$

Note that when such a condition is fulfilled, $\partial C_D\left(N\right)/\partial m_N=0$, all of the higher order derivatives in $m_N$ vanish as well. The above similarity with the bulk waves can be proven further: if the effective velocity of a double layer $C_D\left(2\right)$ is the same as $C_D\left(N-2\right)$, then it also equals the velocity of the composite structure $C_D\left(N\right)$.

Similarly, $TCF$ and other coefficients can be easily obtained from (7):

$$TC{F}_D\left(N,\xi \right)={\displaystyle \frac{1}{M{C}_D^2}}\sum _{\gamma =1}^N{m}_{\gamma }{V}_{\gamma }^2\left({\displaystyle \frac{\partial s_{\gamma }}{s_{\gamma }\partial T}}+{\displaystyle \frac{1}{2\left(1-\nu \right)}}{\displaystyle \frac{\partial \nu }{\partial T}}\right)-\alpha .$$

(9)

In particular, the $TC{F}_D$ for a single layer (when $D=S_0$, the symmetric Lamb mode) may provide information complementary to (4) in finding the temperature coefficients of bulk materials:

$$TC{F}_{S_0}={\displaystyle \frac{\partial s}{s\partial T}}+{\displaystyle \frac{\partial \nu }{2\left(1-\nu \right)\partial T}}-\alpha ,$$

(10)

as all the derivatives remain finite. It can be shown that the rest of the long wavelength modes follow a qualitatively different, non-monotonic behavior in the parametric space. Here, we focus on the example where the explicit solution has been found recently.

Let us consider the fundamental flexural wave (F), the lowest velocity branch of the RL spectrum in a bilayered traction-free plate, $\gamma =1,2$, of total thickness H. The analytical expression obtained in [50] allows us to consider this issue in the full parametric range and extend the earlier perturbative results for the response functions to, e.g., finite masses and thicknesses:

$$C_F/C_{A_0}=\sqrt{{\displaystyle \frac{{\left(1-\varrho {\left(\delta \sigma \right)}^2\right)}^2+4\varrho \delta {\sigma }^2{\left(1+\delta \right)}^2}{\left(1+\varrho \delta {\sigma }^2\right)\left(1+\varrho \delta \right)}}}.$$

(11)

Here, $C_{A_0}=q{h}_1{s}_1\sqrt{\left(1-J_1\right)/3}$ is the reference velocity of the fundamental $A_0$ antisymmetric mode corresponding to a single free plate, $\gamma =1$ of thickness $h_1$, so that $H=h_1\left(1+\delta \right)$ (the layer index at ${\delta }_2$ has been dropped; henceforth, ${\delta }_2=\delta$). In Equation (11), the dimensionality of the whole set of parameters $\xi$ has been reduced to only three independent variables, and the parametric response functions can be calculated in an explicit form. The result for the mass loading sensitivity $MCF=\partial ln{C}_F/\partial m_2$ is given in Appendix A. We also note that due to the known restriction on Poisson’s ratio, $0\le \nu \le /2$, one of the contributing quantities has a limited variation interval, $1/2\le \left(1-J_2\right)/\left(1-J_1\right)\le 2$. Moreover, for most of commonly used solids, their values are close to unity, $0.8÷1.2$, and we therefore assume $J_1=J_2$. Then, examination of the expression in Appendix A of the first derivative with respect to the $\sigma$ ratio (or, essentially, with respect to shear velocity $s_2$ if $\gamma =1$ is fixed as a reference layer) shows that it never crosses zero and proves that $C_F$ is strictly monotonic in this parameter. However, dependence of the flexural mode on the thickness or mass ratios shows a rather different behavior. Note that if the mass densities remain constant, then mass and thickness variations are equivalent $d{m}_{\gamma }={\rho }_{\gamma }d{h}_{\gamma }$. As discussed above, the essential properties of response functions can then be understood by considering the surfaces of the vanishing parametric derivative $MCF=0$, shown in Figure 3. To represent the full interval of thickness ratios $\delta$, it is convenient to switch to the scaled variable, $\epsilon =h_2/H=\delta /\left(1+\delta \right)$.

First, it can be easily shown that expression (11) confirms the known result for the thin coating limit $\delta \to 0$ ($D_1$ in Reference [15]):

$$\underset{h_2\to 0}{lim}MC{F}_F\left(N=2,\xi \right)={\displaystyle \frac{1}{m_1}}\left(3V_2^2/V_1^2-1\right).$$

(12)

In Figure 3, the term in parentheses corresponds to the straight line in the base plane $h_2=0$. The regions separated by the colored surface in this three-dimensional graph correspond to a monotonic behavior (constant sign of the derivative $\partial C_F/\partial m_2$). If velocity of the F-mode would change monotonically with the thickness of the added layer [33], then the critical surface $\partial C_F/\partial m_2=0$ in the figure should be strictly vertical, because the solution of this equation would not depend on $h_2$ or $m_2$ (Table 1 in [15]), similar to the dilatational mode discussed earlier. However, our expression leads to a different behavior and, instead, we obtain a strongly curved and folded $3D$ surface when $\epsilon >0$. In Figure 4 and Figure 5, a few cross-sections of this surface are shown. These results demonstrate that the response function $\partial C_F/\partial h_2$ can also have multiple extrema and change its sign several times as a function of thickness, density and bulk velocity ratios. Surprisingly, these new extrema are found even for the thin coating regime ($\epsilon <0.1$) in a large area of $\left(\varrho ,\sigma \right)$ values, and they are found away from those described by the slanted surface that crosses the $\epsilon =0$ plane in Figure 3 (also seen as the lowest curve in the Figure 5). This crossing corresponds to the known result in (12), $V_2/V_1=1/\sqrt{3}$. Note that the obtained non-trivial behavior of the critical surface does not reveal itself in a perturbative treatment. The slanted part of the surface can be approximated by $3{\sigma }^2\simeq 1/\left(1+4\delta \left(1-\varrho /6\right)\right)$ for $\varrho <6$ and $\epsilon <0.1$ (as shown in Appendix A). Moreover, to clarify whether these extrema correspond to maxima (minima) or to the inflection points of the $C_F\left(h_2\right)$, we also show in Figure 6 the critical surface of the second derivatives ${\partial }^2{C}_F/\partial h_2^2=0$. From these figures, it can be concluded that the mass and thickness response of the F-mode shows a strongly non-monotonic behavior for a broad range of material parameters at small values of $\epsilon =h_2/H$, while at some larger values, there exists an area of enhanced stability of the resonance frequencies (Figure 7). The reason becomes clear if one notices that when both first and second parametric derivatives vanish at the intersection of gray and colored surfaces, Figure 6, the frequency of the F-wave shows the weakest sensitivity to mass or thickness variation. As soon as the thickness of the second layer becomes sufficiently large, $\epsilon \gtrsim 0.3$, the frequency varies monotonically with $h_2$, and crossing the gray surface only signals a change in curvature in the parametric dependence. The exception belongs to the case when the material of the second layer has a smaller bulk velocity, $\sigma <1/\sqrt{3}.$ Then, the vanishing of the mass or thickness response can take place at any ratio $\delta =h_2/h_1$. However, as the second derivative in $\delta$ does not cross zero in this area (see Figure 6), the phase velocity $C_F$ behaves in a strictly non-monotonic way, i.e., it has parametric maxima and minima, Figure 8. For the domain limited by the conditions $\sigma >1/\sqrt{3}$ and $\varrho <3$, the resonance frequency of the F-wave varies in a monotonic way in the full range of all of its parameters and, respectively, shows a non-vanishing response to their variation.

3. Discussion and Conclusions

We have analyzed the behavior of acoustic resonance frequencies and their response to the variation of material parameters for certain Rayleigh–Lamb modes of multilayered structures. Based on explicit analytic expressions, we have considered the entire parametric space in terms of dimensionless ratios, allowing us to reveal some qualitatively new features relevant to thin-membrane acoustic sensing applications. Thus, we have found that the flexural wave’s mass loading sensitivity depends non-monotonically on thickness, density, or bulk velocity ratios, even in the thin coating regime. We have identified the regions of parametric space where one can expect enhanced or reduced sensitivity to mass loading, with purely monotonic and predictable dependence on material parameters. We have demonstrated that the dilatational mode in multilayers offers unique possibilities in applications that require the stabilization of acoustic resonance frequencies with respect to thickness or mass variation and have explained the physical nature of this property. We have also discussed the possibility of using the analytic form of the obtained expressions in exploring the dependence of material parameters on external factors (temperature, pressure, humidity, etc.), which can stimulate further research.