Analysis of membrane phononic crystals with wide bandgaps and low-mass defects

TL;DR

This paper develops models and fabrication techniques for membrane phononic crystals with low-mass defects, optimizing for force sensing by balancing acoustic bandgap size and dissipation, and provides design principles for low-dissipation optomechanical systems.

Contribution

It introduces new modeling and fabrication methods for low-mass defect phononic crystals, analyzing the tradeoffs between mass contrast, bandgap size, and dissipation.

Findings

High mass contrast enhances acoustic isolation.

Low dissipation defect modes are achievable at low temperatures.

Design principles for integrating low-mass mechanical modes into optomechanics.

Abstract

We present techniques to model and design membrane phononic crystals with low-mass defects, optimized for force sensing. Further, we identify the importance of the phononic crystal mass contrast as it pertains to the size of acoustic bandgaps and to the dissipation properties of defect modes. In particular, we quantify the tradeoff between high mass contrast phononic crystals with their associated robust acoustic isolation, and a reduction of soft clamping of the defect mode. We fabricate a set of phononic crystals with a variety of defect geometries out of high stress stoichiometric silicon nitride membranes, and measured at both room temperature and 4 K in order to characterize the dissipative pathways across a variety of geometries. Analysis of these devices highlights a number of design principles integral to the implementation of low-mass, low-dissipation mechanical modes into optomechanical systems.

Equations

$$S_{FF}=\frac{4k_{\mathrm{B}}Tm{\omega }_{\mathrm{m}}}{Q}=\frac{2\hslash k_{\mathrm{B}}T}{Q{x}_{\mathrm{zp}}^2}$$

$$\frac{1}{Q}=\frac{1}{Q_{\mathrm{bend}}}+\frac{1}{Q_{\mathrm{rad}}}$$

$$\dfrac{d^{2}}{dx^{2}}\big{[}I(x)E\dfrac{d^{2}u(x)}{dx^{2}}\big{]}-\mathcal{T}\dfrac{d^{2}u(x)}{dx^{2}}-\rho_{\mathrm{1D}}\dfrac{d^{2}u(x)}{dt^{2}}=0$$

$$v_{p,t}=\sqrt{\frac{\mathcal{T}}{\rho h{\tilde{\mathrm{w}}}_{p,t}}}$$

$$\arccos S=ak$$

$$\begin{split}S=\frac{(V+1)^{2}}{4V}\cos{\big{[}\omega(t_{p}+t_{t})\big{]}}\\ -\frac{(V-1)^{2}}{4V}\cos{\big{[}\omega(t_{p}-t_{t})\big{]}}\end{split}$$

$$S_{\mathrm{max}}=\frac{(V+V^{\scalebox{0.65}[1.0]{$-$}1}+2)}{4}$$

$$\Delta=\frac{\Delta\omega}{\omega_{\pi}}=\ \frac{2}{\pi}\arccos{\bigg{(}\scalebox{0.65}[1.0]{$-$}\frac{1-6V+V^{2}}{1+2V+V^{2}}\bigg{)}}$$

$$Q_{\mathrm{bend}}=Q_{\mathrm{int}}(h)\frac{24(1-\nu^{2})}{Eh^{3}}L_{\mathrm{s}}^{\scalebox{0.65}[1.0]{$-$}1}=\frac{24(1-\nu^{2})}{E_{2}(h)h^{3}}L_{s}^{\scalebox{0.65}[1.0]{$-$}1}$$

$$L_{\mathrm{s}}=\int \kappa (x,y)dxdy$$

$$\kappa =\frac{1}{U}({\partial }_x^{2}u(x,y)+{\partial }_y^{2}u(x,y){)}^2$$

$$g_{1D}(x_{\mathrm{l}})={\int }_{G}g(x_{\mathrm{l}},x_{\mathrm{t}})d{x}_{\mathrm{t}}$$

$${\rho }_p=\frac{\sqrt{3}}{2}{w}_{p}h\rho +\frac{l_t{t}_{w}h}{w_p}$$

$${\rho }_t=\frac{4}{\sqrt{3}}{w}_{p}h\rho$$

$${\tilde{l}}_p=l_p$$

$${\tilde{l}}_t=\frac{\sqrt{3}}{2}{l}_t$$

$$\tilde{\mathcal{T}}=\frac{4}{\sqrt{3}}\sigma t_{w}h$$

$$\tilde{v}_{p}=\big{[}\frac{3w_{p}\rho}{8t_{w}\sigma}+\frac{\sqrt{3}l\rho}{4w_{p}\sigma}\big{]}^{\scalebox{0.65}[1.0]{$-$}\frac{1}{2}}$$

$${\tilde{v}}_t=\sqrt{\frac{\sigma }{\rho }}$$

$$y_k=A_k{e}^{ikx-i\omega t}+B_k{e}^{-ik-i\omega t}$$

$$r_{pt}=\frac{v_t-v_p}{v_t+v_p}$$

$$t_{pt}=\frac{2v_t}{v_t+v_p}$$

$$r_{tp}=\frac{v_p-v_t}{v_t+v_p}$$

$$t_{tp}=\frac{2v_p}{v_t+v_p}$$

$$M=M_{hp}{M}_{tp}{M}_t{M}_{pt}{M}_{hp}$$

$$M_{hp}=\left(\begin{array}{cc}{e}^{i{k}_p{l}_p/2}& 0\\ 0& e^{-i{k}_p{l}_p/2}\end{array}\right)$$

$$M_{pt}=\begin{pmatrix}t_{pt}-\frac{r_{tp}r_{pt}}{t_{tp}}&\frac{r_{tp}}{t_{tp}}\\ \scalebox{0.65}[1.0]{$-$}\frac{r_{pt}}{t_{tp}}&\frac{1}{t_{tp}}\end{pmatrix}$$

$$M_t=\left(\begin{array}{cc}{e}^{i{k}_t{l}_t}& 0\\ 0& e^{-i{k}_t{l}_t}\end{array}\right)$$

$$M_{tp}=\begin{pmatrix}t_{tp}-\frac{r_{pt}r_{tp}}{t_{pt}}&\frac{r_{pt}}{t_{pt}}\\ \scalebox{0.65}[1.0]{$-$}\frac{r_{tp}}{t_{pt}}&\frac{1}{t_{pt}}\end{pmatrix}$$

$$\left(\begin{array}{c}{A}_1\\ B_1\end{array}\right)=e^{iKa}\left(\begin{array}{c}{A}_0\\ B_0\end{array}\right)$$

Full text

I Introduction

Many nanoscale sensing and transduction protocols are based on detection of the minuscule forces between a sample and a highly sensitive probe. Examples include microscopy methods such as atomic force microscopy (AFM) and magnetic resonance force microscopy (MRFM) rugar2004single ; poggio2010force-detected . Additionally, an increasing number of quantum protocols make use of a mechanical intermediary to transduce between two types of otherwise non-interacting quantum systems safavi-naeini2010proposal ; higginbotham2018harnessing ; rabl2010quantum . Readout of the mechanical motion can be achieved through a range of optical, electrical, or magnetic means aspelmeyer2014cavity ; ekinci2005nanoelectromechanical . Fundamentally, all these techniques are often limited by environmental, i.e. thermal or Brownian, noise coupled to the mechanical resonator. This noise leaks into the mechanical mode at the mechanical damping rate as dictated by the fluctuation-dissipation theorem, defining a force sensitivity noise floor.

The factors that define mechanical damping to the environment have been a long-standing question. Generally, loss mechanisms can be divided into two categories: Acoustic radiation into the external structure wilson2008intrinsic , and internal loss (dissipation) of acoustic motion to heat through bending – either at the clamping point, or within the structure itself. Historically, control of materials has been a main driver at reducing loss. Recently, progress is increasingly propelled by engineering of geometry and structural parameters to mitigate both mechanisms. Silicon nitride tensioned thin films are a platform that offer multiple control parameters verbridge2008megahertz ; zwickl2008high ; unterreithmeier2010damping ; yu2012control ; yuan2015silicon ; norte2016mechanical ; reinhardt2016ultralow-noise ; barasheed2016optically ; tsaturyan2017ultracoherent ; ghadimi2018elastic . Tension of the film results in dissipation dilution, which maintains low dissipation while increasing the oscillator energy. Acoustic isolation in the form of bandgap engineering of the surrounding substrate can control acoustic radiation alegre2010quasi-two-dimensional ; yu2014phononic ; tsaturyan2014demonstration , and in recent realizations nanopatterning of the tensioned SiN film itself can provide both (1) a phononic bandgap and (2) the phenomenon called soft clamping, which describes the gradual decay of the mechanical mode into the phononic crystal (PnC) structure reducing bending loss tsaturyan2017ultracoherent ; ghadimi2018elastic .

In this work we explore PnCs with varying mass contrast and associated defect designs in two-dimensional periodic structures. We observe that increasing mass contrast of the PnC widens the bandgap in the frequency domain [Fig. 1 (a,b)]; at the outset this improves isolation of the mechanical motion and avoids contamination of the spectrum due to thermal noise of the integrated motion of other modes. Further, designs with high mass ratio are amenable to analytic analysis, as has been carried out recently in the context of 1D PnC strings ghadimi2018elastic . This design methodology simplifies the design of a range of defects, in particular, smaller effective mass defects and hence better force sensitivity. Small-mass defects can also be incorporated into low-contrast PnCs, through insight from FEM simulations, and we study a range of defect designs and measurements in low-contrast PnCs tsaturyan2017ultracoherent .

When increasing mass contrast of PnCs mode localization is enhanced, but a consequence is a reduction in soft clamping, which increases susceptibility to internal (bending) loss. Here we consider this decrease in soft clamping quantitatively by examining the bending of multiple combinations of intertwined PnC and defect designs. We note that the choice of tradeoff between localization and soft-clamping will depend upon the application considered; when operating at room temperature soft-clamping is critical, but when at deeply cryogenic temperatures yuan2015silicon or in crystalline materials cole2014tensile ; buckle2018stress the associated small internal loss rates may allow for sufficiently low dissipation even when the mode is strongly localized.

For designs in this work, we specifically have in mind applications where high force sensitivity is desired for MHz-scale membrane resonators simultaneously coupled to an optical mode and a spin or electrical degree of freedom scozzaro2016magnetic ; fischer2019spin ; andrews2014bidirectional . For these functionalization purposes, we considered a few design constraints. This include designing defects large enough to couple to Fabry-Perot cavities purdy2012cavity , fabricating wide bandgaps that can incorporate a wide variety of defects robust to perturbations due to deposition of metallic or magnetic components, as well as limiting device size to avoid instabilities resulting from low frequency mechanical modes.

We begin the article with an analytic and numerical analysis of design considerations for optimizing both bending and radiative loss. In particular, we correlate the PnC mass contrast to bandgap width, phononic isolation, and variation in mechanical quality ($Q$) due to corresponding changes in soft clamping. In this way, we link the geometrical design considerations needed to achieve soft clamping to a relatively narrow bandgap. We then present measurements of a set of five devices with varying bandgap and defects; we study $Q$ at both room temperature and 4 K where internal loss is the dominant contributor, such that we can compare the experimental results to loss prediction from bending. As a primary metric, we examine the mechanical force sensitivity of all devices, developing design principles for combinations of both low mass and low loss mechanical defects. We also identify structures in which nuances of defect design has deleterious consequences, such as considerable defect bending, or large effective mass. In this study we utilize 100 nm thick SiN PnCs and demonstrate force sensitivities down to 12 $\mathrm{aN}/\sqrt{\mathrm{Hz}}$ at temperatures of 4 K and a frequency of 1.6 MHz. We expect that thinner films and colder temperatures will further enable optomechanics that combines high frequency membranes with excellent force sensitivity nichol2012nanomechanical ; yuan2015silicon ; tsaturyan2017ultracoherent ; reinhardt2016ultralow-noise ; norte2016mechanical .

II Concept and PnC Modeling

Here we present the design and analysis of two dimensional PnC structures and their defects. The defect is designed such that it has one or more mechanical mode within the bandgap of the surrounding PnC. Such a design supports a spatially localized mechanical mode within the larger patterned membrane structure. The design of all devices in this work is guided by an optimization of the force sensitivity saulson1990thermal :

[TABLE]

Here, $x_{\mathrm{zp}}$ is the resonator’s fluctuation at zero absolute temperature, $\omega_{\mathrm{m}}$ is the angular frequency of the mechanical mode, $T$ is temperature, $k_{B}$ is the Boltzmann constant, and $m$ is the effective mass of the mode. We emphasize that the appearance of $x_{\mathrm{zp}}$ is not an indication of quantum effects, but rather a convenient parameterization as $x_{\mathrm{zp}}$ is a fundamental parameter in both the force sensitivity and the optomechanical coupling aspelmeyer2014cavity . We see that force sensitivity can be enhanced by maximizing both $x_{\mathrm{zp}}$ and $Q$. In pursuit of the latter, we note that dissipation can be expressed as a sum of two contributions

[TABLE]

where $Q_{\mathrm{bend}}$ arises from the internal, or bending losses of the mechanical mode, while $Q_{\mathrm{rad}}$ encapsulates acoustic radiation that is lost to the substrate and wider environment. A SiN PnC can enhance both $Q_{\mathrm{rad}}$ and $Q_{\mathrm{bend}}$ via acoustic isolation and soft clamping, respectively.

II.1 1D Model for Phononic Crystals

In this work we are ultimately interested in the design and fabrication of 2D structures. However, the main design methodologies and concepts can be elucidated by considering a 1D analysis where we convert 2D unit cell (Fig. 2) to a 1D unit cell. In this work, we separate each unit cell into regions of thin tethers and wide pads as depicted in Fig. 1a. In this geometry, the pads (tethers) will have slower (faster) wave velocities than the unpatterned membrane. Periodic modulation of the wave velocity results in gaps in the acoustic spectrum. This procedure, in principle, integrates over the mass of the transverse dimension, and recalculates the appropriate wave velocities. Details are outlined in Appendix A.

We present analysis of both finite and infinite PnC structures. A complete description of finite structures involves solving the 1D Euler-Bernoulli equation for the out of plane displacement $u(x,t)$, which describes the bending of the PnC beam. Such an analysis captures the behavior of the mode shape near clamping points of the finite structure fedorov2019generalized

[TABLE]

where $I(x)$ is the geometric moment of inertia, $E$ is the Young’s Modulus, $\mathcal{T}$ is the tension, and $\rho_{\mathrm{1D}}$ is the linear mass density. It can be seen that for high aspect ratio devices, the 4th order term is negligible away from the clamps  fedorov2019generalized .

For infinite structures, we can assume periodic boundary conditions. Thus, any expected mode shape should have minimal contribution from the bending term in the Euler-Bernoulli equation. Under these conditions the Euler-Bernoulli equation reduces to the 1D wave equation, and thus the motion of the string is completely parameterized by its spatially dependent wave velocity:

[TABLE]

where $\rho$ is the bulk density of the membrane, $h$ is the membrane thickness, and $\tilde{\mathrm{w}}_{\mathrm{p}}$ ($\tilde{\mathrm{w}}_{\mathrm{t}}$) is the converted pad (tether) width. Definitions of the converted pad and tether widths appear in Appendix A.

We outline the solution for infinite structures in Appendix B. As is typical in these types of problems, the band structure is parameterized by a transcendental equation:

[TABLE]

[TABLE]

where $t_{\mathrm{p,t}}=l_{\mathrm{p,t}}/v_{\mathrm{p,t}}$ are the transit times of acoustic waves through the pad and tether respectively, $a$ is the unit cell length, $k$ is the Bloch wavenumber, and $V\equiv v_{t}/v_{p}$ is what we define as the contrast. As can be seen in Eq. 4, the contrast $V$ is geometrically defined by the mass contrast, or the relative converted width of the pads and tethers. Traveling wave solutions are prohibited when $|S|>1$; such a condition defines the range of $\omega$ that determines all bandgaps.

One goal of PnC design is to maximize the bandgap size in the frequency spectrum. For simplicity, we investigate the case where we wish to maximize the width of the first bandgap. Inspection of Eq. 6 reveals that this occurs when the second term is identically 0. In this case, $t_{p}=t_{t}$, and $|S|$ reaches its maximal value $S_{\mathrm{max}}$ when $\omega=\omega_{\pi}\equiv\pi/2t_{t}$:

[TABLE]

With this in mind, we can derive an expression for the normalized bandgap width:

[TABLE]

A strong agreement between the 1D model prediction of Eq. 8 (blue line) and FEM simulations of equivalent 2D unit cell (points) across a range of contrasts can be seen in Fig. 3. The main source of error in the 2D to 1D conversion arises from the inability of the 1D model to correctly account for the transition between the low mass tethers and high mass pads of the 2D PnC. In the high mass ratio limit, these regions account for a smaller fraction of the geometry, and thus the 1D model better captures the resulting band structure for large $V$. Figure 3 also displays results for fabricated devices. It is readily seen that the device illustrated in Fig 1(a) (red star) has a higher mass contrast compared to the device in Fig 1(b) (yellow square). This difference in mass contrast is also reflected in the corresponding bandgap width of both devices.

A technical point arises in the fabrication of devices in the limit of a very large bandgap. As said before, $V$ is ultimately is determined by the relative sizes of the pads and tethers, with wide pads and narrow tethers giving rise to large values of $V$. Therefore, the magnitude of $V$ is ultimately determined bt the minimum tether width that can be fabricated and the desired bandgap position. As an example, a PnC with $V=4$ would have pad widths of $40\text{\,}\mathrm{\SIUnitSymbolMicro m}$ assuming a $\mathrm{w}_{\mathrm{t}}=$ $1\text{\,}\mathrm{\SIUnitSymbolMicro m}$. The principle of equal transit times then sets the unit cell length to be $220\text{\,}\mathrm{\SIUnitSymbolMicro m}$, which in turn sets \omega_{\pi}=2\pi\times$$1\text{\,}\mathrm{MHz}. In general, higher contrast (higher $V$) phononic crystals will have lower bandgap center frequencies, because the larger mass ratio induces lower average wave velocities.

II.2 Effects of Contrast on Soft Clamping

The main advantages afforded by SiN PnCs include both soft clamping and phononic isolation of defect modes. We assert that $V$ (and thus also $\Delta$) is a geometrically defined indicator of the extent to which a PnC achieves these phenomena. To understand the effects of soft clamping on dissipation, we first examine the origin of bending loss for a given mechanical mode unterreithmeier2010damping ; schmid2011damping ; yu2012control :

[TABLE]

[TABLE]

[TABLE]

where $Q_{\mathrm{int}}(h)$ is the intrinsic $Q$ of silicon nitride, $E$ is the Youngs’ mondulus, $h$ is the nitride thickness, $\nu$ is the Poisson ratio, $U$ is the kinetic energy of the mode and $\mathrm{u(x,y)}$ is the mode shape. $E_{2}(h)$ is the imaginary part of the Youngs’ modulus, defined to be $E_{2}=E/Q_{\mathrm{int}}$. Here we call attention to the dependence of $Q_{\mathrm{int}}$, which has been shown to have a scaling $Q_{\mathrm{int}}\propto h$ at membrane thicknesses around 100 nm Villanueva2014evidence . We will call $L_{s}$ the loss factor, which quantifies how much bending losses a given mode has. We also define $\kappa$ as the loss density, which describes the spatial distribution of lossy motion of a mode. Therefore bending loss in membrane devices is most pronounced at clamping points where the mode shape curves strongly to adhere to the restrictive clamping boundary conditions.

Soft clamping and phononic isolation occur when the mechanical defect mode has a frequency within the acoustic bandgap of the PnC. Within the bandgap, $k$ become complex. The magnitude of $\operatorname{Im}(k)$ at the defect mode frequency determines the decay length $l_{0}=1/\operatorname{Im}(k)$ of the defect mode into the PnC.

To directly illustrate how PnC design affects the soft clamping of defect modes, we perform a series of 1D simulations where a single defect is placed in a series of PnCs with ascending contrast. From the calculated defect mode shapes, one can predict a value for $Q_{\mathrm{bend}}$ assuming that $Q_{\mathrm{int}}$ has a value of 6600 Villanueva2014evidence . In this 1D simulation, the defect mode’s frequency has a nearly constant position relative to the bandgap center [Fig. 4(a)]. Two effects as a function of increasing $\Delta$ are apparent. First, the normalized mechanical decay length $n_{0}=l_{0}/a$ decreases as a function of $\Delta$ [Fig. 4(b)], as predicted from our analytic analysis in the 1D model. Second, a decrease of $Q_{\mathrm{bend}}$ is apparent. This is indicative of the bending resulting from the strong decay into the PnC.

II.3 Effect of Contrast on Phononic Isolation

The numerical analysis above showed that high contrast PnCs exhibit short defect mode decay lengths, and we also see numerically that this short decay length provides enhanced suppression of radiative loss. In this analysis, we investigate the concept of phononic isolation in a two dimensional structure by FEM simulation using Comsol. For details, see Appendix D. Here, we calculate the ratio of the tensile energy stored in the entire structure to the energy stored in a $10\text{\,}\mathrm{\SIUnitSymbolMicro m}$ overhang around the PnC perimeter. We define this ratio $\Delta U$. As seen in Fig. 5, the higher contrast PnC has a faster energy decay compared to the low contrast crystal. We note that in the 1D case, $S_{\mathrm{max}}$ is reached when the defect mode frequency is at the center of the bandgap. However, even if this condition is not obtained, one can still achieve robust isolation in the high contrast limit.

II.4 Beyond the 1D Model: Incorporating Defects in 2D

Up to this point, we focused on PnC characteristics, and have shown that our 1D model reproduces the the basic behavior the equivalent 2D crystals. However, the same cannot be done for designing defects in two dimensions, where we generally find that the 1D results have at most a qualitative relationship to a 2D defect. This is mostly due to the complicated stress redistribution that occurs in intricate 2D structures around the defect (Fig. 6), which a 1D model cannot fully capture. This internal structure of the defect leads us to further divide the loss pathways in our 2D devices into losses that can encapsulated in 1D simulations ($Q_{\mathrm{1D-bend}}$) and losses that are inherently 2D ($Q_{\mathrm{2D-bend}}$). In this next section, for example, we designed 5 different devices with a various combinations of defect and PnC designs (Fig. 6). FEM simulations of the defect mode shape help illustrate the difference between the aforementioned separation between “1D loss” and “2D loss”. For example in the context of high-contrast PnCs, a linecut of the defect mode of device A resembles closely a string mode emanating radially from the PnC center. Device E has a large defect pad, and therefore there is considerable internal motion of the pad. This internal motion gives rise to a larger participation ratio of the “2D loss” as compared to device A. Therefore, one expects that 1D analysis is much more applicable to device A than to device E.

The stress redistribution of the 2D structure also strongly affects the position of the defect mode frequency within the bandgap, leading to discrepancies between the 1D and 2D simulations. While the 1D model allows positioning the defect modes directly in the bandgap center [Fig. 4(a)], the equivalent mode has a lower frequency in both 2D FEM simulations and the experimental results [Fig. 1(e)].

III Experimental Results

To test examples of the above ideas, we fabricated and characterized multiple devices (Fig. 6) utilizing the two different PnC designs illustrated in Fig. 1 and defined in Fig. 2, and multiple defect designs. The chosen devices illustrate the complex interplay of bending loss (both in the PnC and defect) and effective mass $m$ (characterized by $x_{\text{zp}}$) of the defect mode. The phononic devices were patterned into a 100 nm thick LPCVD silicon nitride layer on a silicon substrate, and suspended using a KOH wet etch. All Si substrates were 375 $\mu$m thick 5$\times$5 mm2 in the transverse dimensions. Full fabrication information appears in Appendix C.

The mechanical mode spectrum of the devices was characterized at both room temperature and 4 K by affixing the device to a mirror mounted on a piezoelectric transducer. The spectrum was imprinted on the amplitude fluctuations of light reflected from the etalon formed by the mirror and the membrane device. For all measurements, around 5 mW of 1064 nm laser power was incident on the sample. An absolute calibration of the mechanical displacement of our devices can be performed by assuming that the observed mechanical amplitude is entirely due to the expected thermal signal. However, many of our spectra were taken via white noise driving of the entire PnC chip via the piezo electric transducer. Therefore, all mechanical amplitudes are displayed relative to the shot noise limited noise floor the detection chain. The room temperature measurements were performed at a pressure less than $\mathrm{10^{-6}\,mBar}$. 4 K measurements were performed in a closed cycle cryostat with free space optical access.

To clearly understand the mode spectrum we characterize the motion in multiple regions of the PnC. Outlined in Fig. 1, the different colored spectra were obtained by probing near the defect (orange point) and in the crystal bulk (blue point). Modes that only appear when probing near the defect are confined by the PnC and considered to be defect modes [Fig 1(e,f)]. When probing in the crystal bulk [Fig 1(c,e)], one sees a wide range of frequencies with no discernible mechanical resonance (grey regions in Fig. 1(c-f). This frequency range is found to be broader for the high contrast PnC pictured in Fig. 1(a), and is accordance with the results of FEM simulations. We note that in order to adequately resolve the upper bandgap edge, we chose to measure on a tether as it accommodates more high frequency motion above the bandgap than in the pads of the PnC. Occasionally, spurious modes were found inside the bandgap (Fig. 1). These modes had quality factors on the order of 100, and therefore are inconsistent with pure silicon nitride membrane modes. This is consistent with the fact that these are hybridized modes between the membrane and its substrate or the chip-mirror assembly yu2014phononic .

Quality factors were determined via ringdown measurements of the phononic structures driven by the piezoelectric transducer. Figure 7(a) displays the quality factors for the symmetric modes of all devices at both room temperature and 4 K. Loss factors were calculated with Eq. 9 by FEM simulations of each device. The fits are based on the assumption that the dissipation is bending loss limited $Q_{\mathrm{rad}}\ll Q_{\mathrm{bend}}$. For the room temperature data, the fit value of E_{2}=\,$$80\text{\,}\mathrm{MPa} combined with the assumption of E=\,$$250\text{\,}\mathrm{GPa} infers $Q_{\mathrm{int}}=3125$. This has agreement tabulated values of $Q_{\mathrm{int}}$ found in the literature Villanueva2014evidence . The fit value of E_{2}=\,$$17\text{\,}\mathrm{MPa} at 4 K infers $Q_{\mathrm{int}}=14700$. We note that scatter of points, particularly at room temperature, indicates that this assumption may not hold for all devices. Notably, device C which has a rather large value of $\Delta U$, which may be indicative of the existence of measurable external loss for this device.

In Fig. 7(b) we see the force sensing performance across all devices at both room temperature and 4 K. As discussed previously, we use force sensing performance of these high-frequency resonators as a comparison metric. First, let us comparatively examine devices A and B. Device A displayed among the best performance among all devices, with a force sensitivity of 12 $\mathrm{aN}/\sqrt{\mathrm{Hz}}$ at 1.6 MHz at 4 K. Both devices A and B contain a well confined defect mode, leading to a relatively large $x_{\mathrm{zp}}>1\,\mathrm{fm}$ at MHz frequencies. This confinement is accomplished simultaneously with little dissipation. For device A (Fig. 6) the bending loss of the fundamental trampoline mode is most pronounced in the pads surrounding the defect itself. This type of bending profile is similar to that seen in SiN trampolines reinhardt2016ultralow-noise ; norte2016mechanical , while relaxing some of the clamping effects at the edge of the defect mode. Device B maintains low dissipation via the soft clamping of its surrounding low contrast PnC. Although both devices demonstrate relatively high force sensitivity, they differ considerably in design, and therefore include clear design trade-offs. Device A has low effective mass due to its low physical mass. Additionally, the low effective mass of Device A can also be attributed to its first order symmetric mode, which is generally lower in mass since it has a node of motion at the defect boundary. However, its frequency is close to the bandgap edge (Fig. 1). This frequency is most strongly influenced by the length of the defect tethers norte2016mechanical , and thus frequency cannot be easily pushed higher. Device B has a second order symmetric defect mode. In general, these second order modes are easier to place in the bandgap center, but are higher mass because the edge of the defect coincides with an antinode of motion.

To compare the effects of the PnC on the defect mode, we can compare devices A and C, which both have defects with trampoline geometries. Both devices also support a first asymmetric defect mode with frequencies between 1.5 and 2 MHz. We find overall that these devices exhibit a similar force sensitivity both experimentally [Fig. 7(b)] and in FEM simulations. There are however slight differences in how these devices achieve their sensitivity. Finite element simulations predict that device C will have more soft clamping and hence higher $Q$ than device A (although in experiment device A has higher Q, perhaps due to its larger number of unit cells). However, device A has lower mass, and hence larger $x_{\text{zp}}$, than device C. Note, even when the defect pad width of device C is set to be identical to that of A, the increase in pad width is associated with an increase in resonant frequency and the same conclusions hold.

The remaining devices (D and E) provide examples of how nuances in device design can cause adverse effects. For instance, device D exhibits a relatively low $x_{\text{zp}}$ [Fig. 7(b)]. Here, the defect accommodates the second symmetric mode, causing the PnC to carry a large portion of the defect motion (Fig. 6). This differs from the fundamental trampoline-like mode shape in that the second symmetric mode has an anti-node of motion at the PnC boundary. Thus the motion of the defect strongly drives the motion of the surrounding low contrast crystal adding considerable mass to the mode. This effect is not seen in device C (same PnC, different defect), because this defect was designed to accommodate the fundamental mode within the bandgap. Finally, device E is composed of a large defect pad intended to isolate the second symmetric mode through the PnC. which leads to a large amount of defect loss. The large pad is largely stress released, and thus the edges of the pad are allowed to freely bend (Fig. 6). We can compare this to device A which has a large released pad, but it is surrounded by a highly stressed boundary. The relieved stress regions accommodate extreme bending which leads to substantial internal loss, thus limiting the $Q$. This effect has ramifications for future PnC design, where a large pad is often desired for efficient coupling between an optical mode and the mechanical motion.

IV Conclusion

To summarize, we studied SiN membrane designs composed of a defect surrounded by a PnC structure. We demonstrated the effects that high contrast PnCs can have on producing wide bandgaps. Additionally, we presented numerical analysis that connects high contrast PnCs to robust energy isolation of defect modes. Conversely, high contrast crystal, compress the mode spatial distribution, induce more bending and increase the internal loss. However, the increased internal loss is not severe for certain designs, and may for some applications be worth the trade-off for phononic isolation and requirements on the number of unit cells. We have shown that a 1D model with proper conversion, based on the Euler-Bernoulli equation or its reduction to the wave equation, is suitable for analyzing the behavior of a 2D PnC structure. However, incorporating a defect to a 2D crystal requires a complex redistribution of the stress, which could be achieved only by a suitable FEM simulation.

Utilizing these tools, we designed, fabricated and measured 5 devices differing in PnC contrast and defect designs. These devices exhibit force sensitivity on the order 10 $\mathrm{aN}/\sqrt{\mathrm{Hz}}$ at 4K. We draw attention to both devices A and C, which we predict to have comparable force sensitivities when bending-loss limited. This similar performance is achieved despite the vast difference between the PnCs of the devices and we attribute this similarity in force sensitivity to the low-mass trampoline defect common to both devices. However, under conditions where both these devices are radiation-loss dominated, device A has the potential to outperform device C in force sensitivity due to its robust acoustic isolation, or equivalently achieve equal radiation isolation with a smaller number of unit cells. The low mass trampoline defect archetype also incorporates high stress tethers, which can be readily utilized to deposit metallic or magnetic materials for further device functionalization fischer2019spin .

For force sensing functionalization, it is generally required to add additional components to the device to induce a force between the resonator and the sample in order to couple the mechanical resonator to the force in question. These additional elements add mass and alter the designed mechanical mode shapes and properties. One particularly challenging aspect is how to position the defect mode’s frequency within the bandgap upon loading of additional mass. High contrast PnCs offer a larger design space in which to position the defect within the bandgap, while allowing a relatively large tolerance for mass deposition. However, the low mass defects studied here displayed a first symmetric mode at the lower end of the bandgap. This problem could be adverted by exploring more exotic defects beyond the trampoline archetype, such as delocalized double-pad defects, or by utilizing higher order modes.

As mentioned earlier, the low frequency modes, or notably fundamental mode frequencies ($\omega_{0}$ in Fig. 8) of a mechanical structure can induce instabilities when placed inside an optical cavity. High contrast PnCs can achieve mode localization and acoustic isolation with fewer unit cells (Fig. 2), and thus have higher fundamental frequencies (Fig. 8). However, the fundamental modes of these high contrast devices have low effective masses when compared to their low-contrast counterparts, and thus have more Brownian motion. Ultimately, the effect on the cavity stability will be determined by the complete mechanical spectrum, which requires knowledge of the dissipation properties of all mechanical modes, and the most relevant effects will be dependent on the optical cavity in question.

All devices reported in this work were fabricated with a thickness of 100 nm. However, by using 20 nm SiN or thinner, force sensitivity could be significantly enhanced thanks to the scaling of $\sqrt{S_{FF}}\propto h^{-3/2}$. Reducing the thickness of device A to 20 nm, would optimally improve its sensitivity to 1 $\mathrm{aN}/\sqrt{\mathrm{Hz}}$ force sensitivity at 4 K. Furthermore, employing SiN resonators below 100 mK, should significantly enhance force sensitivities as $E_{2}$ of silicon nitride decreases even more under these conditions yuan2015silicon . At these temperatures, structures with low bending loss may become radiative loss limited. In this regime, it may be possible to systematically study the mechanisms of radiative loss in silicon nitride and its transition to internal loss, as has been performed for other mechanical systems maccabe2019phononic . Another consideration for future studies will be the effect of tether size on thermalization of different PnC designs to the cryogenic environment.

Acknowledgements.

We acknowledge funding from AFOSR PECASE, the NSF under grant number PHYS 1734006, ARO-LPS Cross-Quantum Systems Science & Technology program (grant W911NF-18-1-0103), and a Cottrell Scholar award. We thank Yeghishe Tsaturyan for helpful discussions.

V References

VI Appendix A

We present here the method for converting a 2D unit cell geometry to an effective 1D geometry. This approach involves first considering a wave direction $\hat{k}$ along which we will convert the 2D geometry. We wish to encode the geometry transverse to the wave direction into the resulting 1D geometry. It is known that 1D strings with variable width will have lower wave velocities in the wider regions. We imagine that with more transverse?, the plane wave in the 2D geometry will also have lower wave velocity. Our conversion will then be as follows:

[TABLE]

where $x_{\mathrm{l}}$ and $x_{\mathrm{t}}$ refer to the coordinates along the longitudinal and transverse wave directions respectively, and $G$ refers to the domain of a single unit cell.

Upon performing the conversion, $g_{\mathrm{1D}}(x)$ gives a continuous variable width of the beam from which a variable wave velocity can be derived. To derive the bandstructure, the 1D wave equation with Floquet boundary conditions can be numerically solved for an arbitrary $g_{\mathrm{1D}}(x)$. However, we find it fruitful to further simplify the geometry where possible. Namely, we wish to convert wave velocity to match that of the Kronig-Penney model which can be solved analytically. This involves defining regions of high mass and low velocity (pads) and low mass and high velocity (tethers) from an arbritrary $g_{\mathrm{1D}}(x)$.

For unit cells like those in device A, the selection of the pad and tether regions follows naturally from the 2D geometry. Furthermore, the integration carried out in Eq. 12 can be carried out using geometric properties. This procedure will perform the transformation $(w_{p},w_{t},l_{p},l_{t})\rightarrow(\tilde{w}_{p},\tilde{w}_{t},\tilde{l}_{p},\tilde{l}_{t})$ of the geometric parameters.

Following this procedure, we derive the following:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

In this work we also study unit cells that do not fit the pad-tether model at the outset tsaturyan2017ultracoherent . Fig. 10 shows an example collapse in relation to the 2D unit cell. $g_{\mathrm{1D}}(x)$ exhibits pronounced dips over short regions which we will call the tether regions; all other regions will be considered pads. We can then extract the pad width and tether width by taking the mean width over both regions. It is of note that the definition of the pad and tether regions is an arbitrary choice. However, we see that over a reasonable and wide range of tether definitions, $V$ stays roughly constant between 1.4 and 1.5 (Fig. 11).

VII Appendix B

In this section we will present a detailed derivation of bandstructure in the 1D model. We assume that the displacement field has the form

[TABLE]

In this model we view the PnC structure as a modulation in the group velocity of the this acoustic wave. Interfaces between pads and tethers will then have reflection and transmission coefficients given as the familiar Fresnel coefficients

[TABLE]

[TABLE]

[TABLE]

[TABLE]

To continue with our analysis, we can derive the motion across as a unit cell as a transfer matrix $M$, which describes the transformation of the amplitude coefficients $A_{k}$ and $B_{k}$ across each unit cell of the PnC. The total matrix $M$ will then just be a product of the matrices for each subsection

[TABLE]

Here $M_{hp}$ captures the accumulated phase of the plane waves across the half pads at either end of the unit cell. $M_{t}$ does the same for the tether section, while $M_{tp}$ and $M_{pt}$ account for transmission and reflection at each interface.

[TABLE]

[TABLE]

[TABLE]

[TABLE]

If we assume that our structure is infinite, then the Bloch condition implies

[TABLE]

where

[TABLE]

Therefore we have the following eigenvalue problem

[TABLE]

We note that $M$ is a product of unitary matrices, which places a contraint on its components:

[TABLE]

Solving the eigenvalue problem gives us the implicit bandgap equation

[TABLE]

[TABLE]

or written more explicitly as:

[TABLE]

VIII Appendix C

Devices were fabricated on a $375\text{\,}\mathrm{\SIUnitSymbolMicro m}$ thick 3 inch diameter silicon wafer with 100 nm of grown stoichiometric LPCVD silicon nitride on either side. Designs were patterned using a direct write photolithography system after spinning $1\text{\,}\mathrm{\SIUnitSymbolMicro m}$ of SPR-660 photoresist onto either side. The top of the wafer was patterned with the PnC designs, while the back was patterned with rectangular windows aligned to each PnC. The patterning was done on both sides with $300\text{\,}\mathrm{mJ}\text{\,}{\mathrm{cm}}^{-2}$ of $405\text{\,}\mathrm{nm}$ light. During these steps, the wafer was affixed to a sapphire carrier wafer with Crystalbond 509 in order to protect the bottom side from unwanted processing. Patterning of the silicon nitirde was completed via a $\mathrm{CF_{4}}$ reactive ion etch. The wafer was then cleaned with $\mathrm{O_{2}}$ plasma followed by ultrasound cleaning in an acetone bath. Additional cleaning was performed with isopropyl alcohol and water. To suspend the PnC structures, the window side of the wafer was etched using a 80 C KOH bath. The PnC side was protected via a PEEK wafer holder. Following wet etching, the wafer was cleaned in a Nanostrip bath, acetone and isopropyl alcohol.

IX Appendix D

FEM simulations in this work were performed using COMSOL Multiphysics. Simulations were performed in 2 steps; the stationary stress redistribution was calculated, followed by a eigenfrequency analysis to determine the mode spectrum of each PnC. The simulations were performed with increasingly fine meshes until both the mode frequencies and loss factors were found to converge to a change less than 5 percent. These simulations assumed that the density of silicon nitride was $3100\text{\,}\mathrm{kg}\text{\,}{\mathrm{m}}^{-3}$, the film stress was $1.13\text{\,}\mathrm{GPa}$, the poisson ratio was .27, and the Youngs’ modulus was $250\text{\,}\mathrm{GPa}$.