## Abstract

We present the optical and mechanical design of a mechanically compliant quasi-two-dimensional photonic crystal cavity formed from thin-film silicon in which a pair of linear nanoscale slots are used to create two coupled high-*Q* optical resonances. The optical cavity supermodes, whose frequencies are designed to lie in the 1500 nm wavelength band, are shown to interact strongly with mechanical resonances of the structure whose frequencies range from a few MHz to a few GHz. Depending upon the symmetry of the mechanical modes and the symmetry of the slot sizes, we show that the optomechanical coupling between the optical supermodes can be either linear or quadratic in the mechanical displacement amplitude. Tuning of the nanoscale slot size is also shown to adjust the magnitude and sign of the cavity supermode splitting 2*J*, enabling near-resonant motional scattering between the two optical supermodes and greatly enhancing the *x*^{2}-coupling strength. Specifically, for the fundamental flexural mode of the central nanobeam of the structure at 10 MHz the per-phonon linear cross-mode coupling rate is calculated to be ${\tilde{g}}_{+-}/2\pi =1\phantom{\rule{0.2em}{0ex}}\text{MHz}$, corresponding to a per-phonon *x*^{2}-coupling rate of ${\tilde{g}}^{\prime}/2\pi =1\phantom{\rule{0.2em}{0ex}}\text{kHz}$ for a mode splitting 2*J*/2*π* = 1 GHz which is greater than the radiation-limited supermode linewidths.

© 2016 Optical Society of America

## 1. Introduction

Multimode optomechanical systems consisting of three or more optical and mechanical modes have recently received growing interest within the field of cavity optomechanics [1]. One of the main original motivations for multimoded systems was the prospect of realizing a quantum non-demolition (QND) measurements of mechanical energy or phonon number using position-squared optomechanical coupling arising from the novel dispersion of a multimode optical system [2–4]. While QND measurements turn out to be very challenging in such a set-up [4, 5], multimode optomechanical systems are promising for studying a variety of other interesting phenomena, such as three-mode parametric instability [6–8], multimode optomechanically induced transparency [9,10], synchronization of mechanical oscillators [11–14] and the generation of sub-Poissonian statistics [15–18]. Multimode optomechanical systems may also find application in displacement sensing [19, 20] and optical information processing [21], where they can be used as optical filters [22], switches or delay lines [23–26]. Given the strong interest and opportunity in multimode cavity optomechanics, a number of new experimental platforms have been developed. These include membrane-in-the-middle setups [3, 27–29], nanofabricated chip-scale resonators [10, 29–31], hybrid multimode microwave circuits [32], and gravitational wave detectors [33].

In a recent experimental work, we realized a silicon photonic crystal optomechanical cavity capable of very large position-squared optomechanical coupling [34]. As shown schematically in Fig. 1 the structure consists of a double-waveguide photonic crystal cavity in which two individual waveguide cavity modes are coupled by photon tunneling through a mechanically compliant element. The double-waveguide cavity localizes two optical resonances at telecommunication wavelengths *λ* =1550 nm with high quality factors *Q* > 10^{5}. The optical modes couple efficiently to mechanical modes with frequencies *ω*_{m}/2*π* ranging from 6 MHz up to 1 GHz. A tuning scheme based on electrostatic actuation [22] is also employed to tune the optical tunneling rate between the two slotted waveguide modes as well as the supermode frequencies. This provides a direct dynamical control over the linear and quadratic optomechanical coupling strengths. Here, we present details of this double-slotted photonic crystal cavity design and explore the range of possible optomechanical interactions in such a device. In particular, we focus on the presentation of a unifying description and calculation method for the optomechanical coupling coefficients, the scaling of the coupling to fundamental and higher-order modes of the structure, and the dependence of the optomechanical properties of the photonic crystal cavity on slot width, separation, and asymmetry.

The paper is organized as follows. We start in Sec. 2 with a brief introduction of the theory of multimode optomechanics. We consider a generic Hamiltonian describing a system of two optical modes and a series of mechanical modes. Coupling between the optical modes occurs via mechanically independent photon tunneling, and via absorption or emission of mechanical phonons. By considering specific symmetries of the optical and mechanical modes we show how different orders of coupling in mechancial amplitude can occur. In Sec. 3 we detail the design of the specific double-slotted photonic crystal cavity of interest to our work. We show that by controlling the widths of the waveguide slots of the structure it is possible to adjust the photon tunneling rate of the optical modes, and to even completely suppress this tunneling. In Sec. 4 we analyze the mechanical resonances of the photonic crystal cavity and identify several mechanical modes with frequencies ranging from MHz to GHz that have the appropriate symmetry for realizing large optomechanical coupling. Finally, in Sec. 5 we employ a perturbation theory to calculate the strength of the linear self-mode, linear cross-mode, and position-squared coupling a number of the different mechanical resonances of the photonic crystal structure.

## 2. Multimode Hamiltonian

We consider a multimode cavity optomechanical system consisting of two spatially separated optical modes *a*_{1,2} of frequencies *ω*_{1,2} and independent mechanical modes *b _{k}* of frequencies

*ω*. The individual optical cavity modes are coupled to each other through photon tunneling at a rate

_{m,k}*J*. Conceptually, our structure can be viewed as an on-chip generalization of the membrane-in-middle setups [3, 4] in a chip-scale architecture, as it is schematically represented in Fig. 1. The corresponding Hamiltonian can be written as,

*g*are the coupling strengths between the optical modes

_{ij,k}*a*,

_{i}*a*and the mechanical mode

_{j}*b*. In the particular case of symmetric optical cavities,

_{k}*a*

_{1}and

*a*

_{2}are degenerate with

*ω*

_{1}=

*ω*

_{2}=

*ω*

_{0}. Since the individual optical modes are spatially separated, we can to a good approximation neglect the terms proportional to

*g*with

_{ij,k}*i*≠

*j*and therefore write

*g*simply as

_{ij,k}*g*. By introducing the supermode basis ${a}_{\pm}=\left({a}_{1}\pm {a}_{2}\right)/\sqrt{2}$, we can diagonalize ${\widehat{\mathscr{H}}}_{\text{opt}}$ which yields for the total Hamiltonian,

_{i,k}*ω*

_{+}(0) −

*ω*

_{−}(0) = 2

*J*. The first term inside the brackets of Eq. (2) describes the

*linear self-mode*optomechanical coupling of the

*a*

_{±}supermodes to the mechanical mode of interest. The last term describes the

*linear cross-mode*optomechanical coupling,

*i.e.*the coupling between the

*a*

_{±}supermodes mediated by the mechanical vibrations. In the

*a*

_{±}basis we have for the

*linear self-mode*optomechanical coupling to mechanical mode

*b*

_{k}, and the

*linear cross-mode*optomechanical coupling,

Following the approach in [4], we further diagonalize the full Hamiltonian assuming a quasistatic approximation for the mechanical motion. The resulting eigenfrequencies of the *a*_{±}({*x _{k}*}) supermodes are,

*b*and Taylor expanding the optical supermode frequencies as a function of small displacement

_{k}*x*around equilibrium position

_{k}*x*yields,

_{k, eq}Here *g*_{±,}* _{k}*({

*x*

_{k}_{, eq}}) is the linear (self-mode) coupling coefficient and

*g*′

_{±,}

*({*

_{k}*x*

_{k,}_{eq}}) is the quadratic coupling coefficient of the

*a*

_{±}({

*x*}) supermode to the

_{k}*k*

^{th}mechanical mode

*b*. In what follows we will be primarily interested in the linear and quadratic coupling coefficients around the symmetric equilibrium position {

_{k}*x*

_{k,}_{eq}} = {0}, where the linear coupling is trivially

*g*

_{±,}

*({0}) =*

_{k}*g*

_{±,}

*, = (*

_{k}*g*

_{1,}

_{k}+ g_{2,}

*)/2 and the quadratic coupling can be related to the linear cross-mode coupling of the*

_{k}*a*

_{±}supermodes,

We return to the more general result for *x _{k}*

_{, eq}≠ {0} in Sec. 5.4.

From the above expressions, mechanical modes such that *g*_{1,}* _{k}* = −

*g*

_{2,}

*will have vanishing*

_{k}*linear self-mode*couplings to the

*a*

_{±}supermodes, while mechanical modes such that

*g*

_{1,}

*=*

_{k}*g*

_{2,}

*will have vanishing*

_{k}*linear cross-mode*coupling strength. All intermediate cases such as |

*g*

_{1,}

*|≠ |*

_{k}*g*

_{2,}

*| can of course occur in general. In the following, we introduce a photonic crystal optomechanical resonator supporting multiple mechanical modes of different symmetries and with optimized overlap with the optical modes. We show that by engineering the bandstructure and defect of the photonic crystal the splitting 2*

_{k}*J*between the optical supermodes can be tuned to arbitrarily small values, which greatly enhances the

*x*

^{2}-coupling strength.

## 3. Multimode photonic crystal optomechanical cavity

#### 3.1. Bandstructure properties

In this paper, we design our structure assuming a silicon thin film device layer of thickness *t* = 220 nm, Young’s modulus *Y* = 169 GPa, mass density *ρ* = 2329 kg/m^{3} and refractive index *n* = 3.42. Our initial geometry is a quasi-two-dimensional (quasi-2D) periodic photonic crystal membrane patterned with a triangular lattice of holes (lattice constant *a* = 480 nm, circular hole radius *r* = 0.3*a*) This structure has an in-plane photonic bandgap (a pseudo-bandgap) for guided slab modes of predominantly TE polarization (electric field polarized in the plane of the slab) around a free-space wavelength of *λ* ≈ 1550 nm. Referring to Fig. 1(b), addition of the slots breaks the translational periodicity in the transverse direction (*y*-direction), leaving a periodic structure of lattice constant *a* in the longitudinal (*x*) direction. Introducing two air slots of width *s*_{1} = *s*_{2} = *s* ≈ 100 nm into the photonic crystal membrane results in a pair of optical waveguides with guided modes localized to each of the individual air slots due to the pseudo-bandgap. As shown in Fig. 1(b) this splits the triangular lattice into two outer slabs and a central beam.

The optical bandstructure of the guided modes of the double-slotted waveguide structure with a central beam consisting of a single row of air holes is plotted in Fig. 2. The optical bandstructures here and in what follows are computed using the MIT Photonic Bands package [35]. To simplify the bandstructure we only plot the optical bands with TE-like polarization (more accurately we plot those modes with even vector symmetry about the mid-plane of the photonic crystal thin-film slab). In addition to the symmetry about the mid-plane of the silicon thin-film slab – corresponding to the *σ _{z}* mirror operator – the double-slotted waveguide structure also has a symmetry plane about the

*y*= 0 plane as indicated by the green line in Fig. 1(b). The mirror operator corresponding to this symmetry we label

*σ*, and the modes of the double-slotted photonic crystal waveguide can be categorized by their even and odd parity under

_{y}*σ*. Here we use the labeling convention that the waveguide supermode with even

_{y}*E*electric field profile is called the “even” mode, while the waveguide supermode with odd

_{y}*E*electric field profile is termed the “odd” mode (note that classifying the modes by their vector symmetry would swap the mode labels). The even and odd waveguide bands are shown as solid black and red lines in Fig. 2. The back and red dashed curves correspond to the unguided (in the transverse

_{y}*y*-direction) modes of the triangular photonic crystal slabs surrounding each air slot. The shaded grey region corresponds to the region above the light cone of the air cladding surrounding the silicon slab, in which a continuum of radiation modes freely propagate out of the plane of the slab (

*z*direction). The insets to the right of the bandstructure plot in Fig. 2 show the

*E*field profiles of the even and odd waveguide modes at the edge of the first Brillouin zone (

_{y}*X*-point).

In the hopes of achieving large *x*^{2}-coupling strengths, our focus will be to design photonic crystal cavity modes with minimal optical splitting 2*J* (see Eq. (9)). An obvious way to minimize the optical coupling between the slot waveguide modes would be to increase the separation between the slots. This decreases the photon tunneling rate *J*, and hence decreases the frequency splitting between the waveguide modes. For instance, we show in Fig. 3 that the frequency splitting between the odd and even waveguide modes at the *X*-point can be decreased from ∼ 3 THz in the case of waveguide slots separated by a single row of holes, all the way down to ∼ 68 GHz for waveguide slots separated by five rows of holes. In the latter case the bands become nearly degenerate over a significant fraction of the Brillouin zone, whereas in the more strongly-coupled case of a single row of holes the even and odd waveguide bands have different slopes and even cross near the *X*-point. The mode profiles of the even and odd waveguide supermodes are plotted in Fig. 2 for a separation between slots of a single row of holes, with the odd waveguide modes having a node in the center of the central beam, and thus, more of their energy in the air slots. As described below, for a cavity based upon the coupled waveguides, this feature allows one to control the relative mode frequencies of the even and odd cavity supermodes by changing the size of the air slot gaps in the structure.

In addition to slot width and slot separation, the fact that the even waveguide modes have more of their energy in the central beam may also be exploited to tailor the relative frequencies of the waveguide supermodes. We consider a middle slab consisting of a single row of holes (nanobeam) and analyze the impact of the ellipticity and lattice constant on the frequency of the waveguide modes. We parametrize the ellipticity of the holes by the aspect ratio η = *r _{y}*/

*r*, where

_{x}*r*(

_{y}*r*) are the semi-axis of the ellipse in

_{x}*y*(

*x*) direction, so that η = 1 corresponds to circular holes. Further, by setting ${r}_{y}={r}_{0}/\sqrt{\eta}$ and ${r}_{x}={r}_{0}/\sqrt{\eta}$ where

*r*is circular radius of the unperturbed cells, we keep the air filling fraction invariant between the elliptical holes and the circular holes. Fig. 4(a) shows the frequency shift of the bands at the

*X*-point due to the variation in the aspect ratio of the holes in nanobeam. Only the even band is influenced while the odd band does not shift. For completeness we show in Fig. 4(b) the scaling of the even and odd waveguide modes at the

*X*-point versus a scaling of the in-plane lattice constant (i.e., slab thickness constant). As one might expect, the splitting between the waveguide supermodes of different parity are relatively unaffected by the lattice scaling.

#### 3.2. Optical cavity

To form an optical cavity from the coupled waveguide system described above we need to find a way of closing the ends of the waveguides. To this end, we introduce a “defect” in the waveguide structure by modifying the geometrical parameters of the waveguide along its propagation (*x*) axis. This defect region of the waveguide is then embedded between two “mirror” sections of the waveguide as depicted in Fig. 1. Localized cavity modes result for a waveguide modification that pushes the *X*-point waveguide supermode frequencies in the defect unit cells inside the pseudo-bandgap of the unperturbed mirror unit cells [36,37].

Here the cavity is designed by combining the two defects shown in Fig. 4. The lattice constant is decreased from *a* = 480 nm to *a _{D}* = 468 nm quadratically over

*N*= 7 cells at the center of the structure. Simultaneously with the scaling of the lattice constant, the aspect ratio of the air holes in the central nanobeam are increased from

_{D}*η*= 0.55 to circular holes with

*η*= 1 while keeping a constant air filling fraction. Fig. 5(a) shows the modification of the waveguide mode frequencies from the unperturbed mirror cells to the central defect cell. The change in the lattice constant pushes the mode frequencies inside the bandgap and the change of hole ellipticity reduces the frequency splitting between the even and odd supermodes. We simulated the full cavity structure using the COMSOL Multiphysics [38] finite-element-method (FEM) solver. The length of the entire cavity consists of a total of

*N*= 42 waveguide unit cells along the

_{x}*x*-axis, with

*N*= 7 central defect unit cells. Each of the outer photonic crystal slabs have

_{D}*N*= 9 rows of air holes in the transverse

_{y}*y*-direction. Fig. 5(b) and Fig. 5(c) show the

*E*field profiles for the odd and the even cavity supermodes, respectively. With the given parameters, the odd mode has a simulated free-space wavelength of

_{y}*λ*= 1539 nm and a radiation-limited optical

*Q*-factor of

*Q*= 4.2 × 10

^{6}. The even mode has a free-space wavelength of

*λ*= 1539.8 nm and a substantially lower radiation-limited

*Q*-factor of

*Q*= 4 × 10

^{5}. This asymmetry in

*Q*-factor results from both the larger coupling to even waveguide modes away from the

*X*-point in the mirror section, a result of the non-monotonic dispersion of the even band, and the reduced vertical radiation loss for a mode of odd in-plane symmetry [39].

#### 3.3. Dependence of the cavity supermode frequency splitting on the slot size

Since the optical modes of interest are mostly localized in the air slots, their frequencies are strongly impacted by the width of the air slots *s*. Increasing *s* causes the effective refractive index of the optical waveguide modes to decreases, yielding a blue shift of the resulting optical cavity frequencies. FEM simulations of the waveguide supermodes of the double-slotted photonic crystal structure for slot widths ranging from *s* = 90 nm to *s* = 100 nm are plotted in Fig. 6(a). The frequencies increase approximately linearly with *s* with slightly different slopes, resulting in a crossing of the cavity supermodes at *s* ∼ 95 nm. This crossing is made more apparent in Fig. 6(b), where the relative splitting between the even (*ω*_{+}) and odd (*ω*_{−}) supermode frequencies, Δ*ω* = *ω*_{+} − *ω*_{−} = 2*J*, is plotted.

As we see, by solely modifying the slots width, it is possible to change the splitting between the cavity supermodes (Δ*ω* = *ω*_{+} − *ω*_{−}) from positive to negative. This is again explained by looking at the even and odd field profiles shown in Fig. 5(b) and 5(c). Since the odd mode has a node in the middle of the central nanobeam, and thus a larger energy density in the air slots, its effective refractive index is more sensitive to a change in the air slot size than the even cavity supermode. Therefore, for an equal change of the slot widths the frequency of the odd cavity supermode shifts more than the frequency of the even cavity supermode. It is also worth noting that the crossing of the *cavity* supermodes is only possible because the odd and even *waveguide* bands cross (see Fig. 5(a)), with the even waveguide band having a *higher* frequency than the odd waveguide band slightly away from the *X*-point.

Figure 7(a–c) shows anti-crossing curves obtained where one slot width (*s*_{1}) is fixed and the other is swept from *s*_{2} = *s*_{1} − 5 nm to *s*_{2} = *s*_{1} + 5 nm. In this case, the cavity mode eigenfrequencies anti-cross with a minimal splitting 2*J* achieved for the symmetric *s*_{2} = *s*_{1} situation. The odd and even supermode branches of the anti-crossing curves were identified in Fig. 7 by looking at the parity of electric field profile of the simulated eigenmodes at the center of the anti-crossing (here we label the branches by the parity of the supermodes at the center of the anti-crossing curve). Notice that the lower frequency branch for slot width *s*_{1} = 90 nm (Fig. 7(a)) is the odd cavity supermode, whereas the higher frequency branch for slot width *s*_{1} = 100 nm (Fig. 7(c)) has odd parity. For the in-between case of *s*_{1} = 95 nm, the splitting at the center of the anti-crossing curve is approximately zero (2*J*/2*π* = 17 GHz), as it must be if the parity of the upper and lower frequency supermode branches swap. The possibility to control the supermode frequency splitting in this way provides a unique way to tune the strength of the *x*^{2}-coupling according to Eq. (9).

In addition to tuning the cavity supermode splitting, adjusting the air slot sizes may also be used to tune the absolute and relative optical *Q*-factor of the two cavity supermodes. Figure 7(d) plots the simulated radiation-limited *Q*-factor for the two cavity supermodes as the slot widths are adjusted and the cavity modes sweep through the anti-crossing point. At the center of the anti-crossing curve the two cavity supermodes are to a good approximation even and odd parity modes around the *y* = 0 plane, which as we pointed out earlier results in a significant difference in their *Q*-factor. Motion of the outer photonic crystal slabs or the central nanobeam will also change the relative air slot sizes, thus changing the radiation damping of the optical cavity supermodes. As noted in recent theoretical work [40], this *dissipative* optomechanical coupling can lead to interference of the quantum noise entering the cavity mode system. This in turn can be used to cool the coupled mechanical resonator to the quantum ground-state even in the unresolved sideband regime (bad-cavity limit). More relevantly, in the case of a predominantly *x*^{2} coupled system, this effect can be used to reduce the parasitic linear back-action and enable continuous *x*^{2} measurements of the mechanical motion.

## 4. Mechanical resonances

Having considered the optical cavity modes of the double-slotted planar photonic crystal, we now analyze the mechanical modes of the structure. In order to support mechanical resonances the photonic crystal slabs are suspended. The optical modes can interact with both flexural and localized acoustic modes of the central nanobeam. In Sec. 4.1, we present the flexural modes of the structure. We show that higher orders flexural modes are found to exist up to 1 GHz with significant optomechanical coupling, making them suitable for operation in the resolved sideband regime, where the optical linewidth *κ* is much smaller than the mechanical frequency *ω*_{m} [1]. In Sec. 4.2, we show that the defect developed to form the optical cavity also gives rise to a localized acoustic resonance of a few GHz frequency.

#### 4.1. Flexural mechanical resonances

The outer slabs and nanobeam behave as three independent mechanical resonators supporting various in-plane and out-of-plane flexural mechanical resonances. Here, we focus on the in-plane flexural modes that are asymmetric with respect to the *y* = 0 plane mirror operator (*σ _{y}* = −1), and symmetric with respect to the

*x*= 0 plane mirror operator (

*σ*= +1), about the center of the structure. These flexural modes are represented in Fig. 8 with exaggerated deformation profiles. Since our main focus will be on the fundamental resonances, we denote the fundamental in-plane flexural modes of the two outer slabs and nanobeam as

_{x}*b*

_{1},

*b*

_{2}and

*b*

_{3}respectively. Their respective frequencies are denoted

*ω*

_{b}_{1},

*ω*

_{b}_{2}and

*ω*

_{b}_{3}.

In our design, the outer slabs are suspended by tethers of length *l _{t}* =2.5

*μ*m and width

*w*

_{t}=150 nm, yielding fundamental in-plane flexural resonance of

*ω*

_{m}/2

*π*≃ MHz. As shown in Fig. 8(a), these modes correspond to a uniform displacement of the whole slabs. The displacement of one outer slab causes a uniform change of the width of the adjacent slot, and hence a change of the optical supermode frequencies. The in-plane slab modes provide degrees of freedom for the electromechanical tuning of the slot widths.

In Fig. 8(b), we plot the displacement profiles of the first three lowest frequency (10.8 MHz, 56 MHz and 130 MHz) nanobeam in-plane flexural modes of symmetry {*σ _{x}* = +1,

*σ*= −1}. The

_{y}*y*-polarized electric field profiles |

*E*|

_{y}^{2}are plotted in the inset for both the odd and even optical supermodes

*a*

_{±}. The finite extent of the optical modes along the

*x*-axis of the photonic crystal structure limits the region of the nanobeam that will contribute to the optomechanical interaction. As a result, the nanobeam displacement amplitude

*x*

_{3}can be approximated by a net effective displacement of the whole nanobeam ${\overline{x}}_{3}\approx {x}_{3}$, causing one slot width to change by an amount $+{\overline{x}}_{3}$ and the other to change by $-{\overline{x}}_{3}$. Because of this asymmetric displacement, the optomechanical couplings of these flexural mode to the individual slot modes

*a*

_{1}and

*a*

_{2}are expected to be equal and of opposite sign. This favors the quadratic and linear cross-mode interaction terms introduced in Eqs. (4) and (9).

Higher order flexural modes of the nanobeam have been identified with frequencies up to 1 GHz and are summarized in Tab. 1. Assuming a moderate optical quality factor of *Q* = 5 × 10^{5}, the resolved-sideband regime condition (*κ* < *ω*_{m}) could be met with a flexural mode of frequency *ω*_{m}/2*π* = 400 MHz.

#### 4.2. Localized phononic crystal resonance

In Sec. 3.2 we described how to localize optical waveguide modes of the double-slotted photonic crystal waveguide propagating by engineering a perturbation to the waveguide unit cell in the propagation direction. In particular, we analyzed the photonic bandstructure of the waveguide unit cell and designed a defect based on a combination of change in the lattice constant and change in the central nanobeam hole aspect ratio. Here we study the phononic bandstructure of the nanobeam unit cell and show that our choice of photonic crystal defect parameters makes the nanobeam compatible with the localization of an GHz-frequency acoustic resonance.

Figure 9(a) and 9(b) show the FEM-simulated acoustic bandstructure of the nanobeam unit cell and the frequency shift of the *breathing* mode band at the Γ-point as the nanobeam transitions from the mirror unit cell geometry to the defect unit cell geometry. The breathing mode band is shown as a solid red curve. The nanobeam unit cell and the corresponding normalized displacement field profile **Q**(**r**) of the breathing mode are depicted in Figs. ‘9(e-g). The localized breathing mode is drawn from the Γ-point of the bandstructure in order to have a significant optomechanical coupling to the optical mode [41].

Figure 9(c) and 9(d) detail separately the shifts of Γ-point frequency of the breathing mode caused by the perturbation in the lattice constant and hole ellipticity of the nanobeam. Increasing the aspect ratio *η* will push the breathing mode frequency down into the band gap of the unperturbed mirror cells, while decreasing the lattice constant slightly increases the frequency of the Γ-point mode. As summarized in Fig. 9(a), the same defect used to localize the photonic crystal resonances satisfies the conditions to localize a phononic crystal resonance. Note that, in contrast to the mechanical flexural modes, the displacement of the localized breathing mode is symmetric with respect to *σ _{y}* so the optomechanical coupling to this mode is expected to be the same for both

*a*

_{1}and

*a*

_{2}. Therefore, the breathing mode is expected to have negligible linear cross-mode and quadratic coupling strengths. Nevertheless, the breathing mode presents the advantage of being well in the resolved sideband regime and could potentially be used as an auxiliary mechanical mode in multimode optomechanically induced transparency schemes as proposed in Refs. [9,10].

## 5. Optomechanical coupling relations

With knowledge of the optical and mechanical properties of the double slotted waveguide cavity, we can now turn to the calculation of the different optomechanical coupling factors. We utilize a perturbation theory of Maxwell’s equations [42] suitable for dealing with both spatial shifts in the dielectric boundaries of the cavity structure as well as stress-induced modifications of the local dielectric constant of the deformed structure. Applying first-order perturbation theory to the numerically computed unperturbed optical field profiles and mechanical field profiles allows us to evaluate both the linear self-mode coupling *g*_{±} and the linear cross-mode coupling *g*_{+−} for the nanobeam flexural modes [5,7,43]. By considering the perturbation theory to second-order yields the strength of the quadratic coupling in our structure. Finally, we consider the modification of the coupling strengths due to deviations of the structure from the symmetric equilibrium position of the central nanobeam and outer slabs.

#### 5.1. Linear self-mode optomechanical coupling

Maxwell’s equations in a source-free, linear dielectric medium, yields the following eigenvalue equation for the electric field,

*E*〉 for the electric field

**E**(

**r**) eigenstate with harmonic time dependence

*e*

^{−}

*. Here,*

^{iωt}*c*is the speed of light in vacuum, and

*ϵ*(

**r**) is a dielectric constant which is a function of the spatial coordinate

**r**(and most generally a tensor). We are interested in the change in the modal frequency due to an infinitesimal perturbation δ

*α*to the dielectric structure. The first-order correction term to the mode frequency is expressed as

*ϵ*

^{(0)}(

**r**) is the unperturbed dielectric constant of the structure, |

*E*

^{(0)}〉 and

*ω*(0) are electric field and frequency of the harmonic optical mode of interest, 〈

*E*

^{(0)}|

*ϵ*

^{(0)}|

*E*

^{(0)}〉 ≡ ∫

*d*

_{V}^{3}

*r*(

*E*

^{(0)}·(

*ϵ*

^{(0)}

*E*

^{(0)})*), and

*δα*depends upon the type of perturbation to the dielectric structure.

The change in dielectric constant due to mechanical displacement arises from two main contributions. The first contribution is due to shifting of the interface boundary between two dielectric media. In this case the dielectric function is a high-contrast step function which translates displacements normal to the boundary into local modifications of the dielectric seen by the electric field. As proposed by Johnson, et al., in [42], the appropriate perturbation theory in this so-called *moving boundary* (MB) problem is given by,

*ϵ*

_{1(2)}is the dielectric constant of medium 1(2) at any point in the boundary surface

*A*between two media of differing dielectric constant, Δ

*ϵ*=

*ϵ*

_{1}−

*ϵ*

_{2}, $\mathrm{\Delta}{\mathit{\u03f5}}^{-1}={\mathit{\u03f5}}_{1}^{-1}-{\mathit{\u03f5}}_{2}^{-1}$, $|{E}_{\Vert}^{(0)}|(|{E}_{\perp}^{(0)}|)$ is the magnitude of the unperturbed electric (displacement) field polarized in the plane (out of the plane) of the boundary surface

*A*between medium 1 and medium 2, and

**n**(

**r**) is the outward unit vector normal pointing from medium 1 into medium 2 on boundary

*A*. Here,

**q**(

**r**) is the normalized displacement field of the mechanical mode of interest with maximum displacement equal to unity, max|

**q**(

**r**)|= 1. We can also define an effective mass of the mechanical mode in terms of

**q**,

*ρ*is the mass density of the dielectric material defining the optomechanical structure. This effective mass is the appropriate motional mass for evaluating the zero-point fluctuation amplitude, ${x}_{\text{zpf}}=\sqrt{\hslash /2{m}_{\text{eff}}{\omega}_{m}}$, of the generalized amplitude coordinate corresponding to the point of maximum amplitude of the mechanical mode.

The second contribution to the linear self-mode coupling is due to the *photoelastic e*ff*ect*, resulting from the change of the dielectric constant due to the local strain induced by the mechanical displacement. The first-order perturbation to the dielectric tensor is given by,

*ϵ*^{(0)}is the unperturbed dielectric tensor,

*ϵ*is the permittivity of free space,

**p**is the fourth rank photoelastic tensor, and

**S**is the symmetric strain tensor. For an isotropic medium this simplifies to, in index notation. In matrix form,

The resulting first-order photoelastic (PE) correction to the optical frequency is

In the structures studied here, which are made by etching patterns into a thin-film of silicon, the only two media are silicon and vacuum. As such, for the PE contribution to the linear self-mode coupling we utilize the photoelastic tensor coefficients for silicon in evaluating the integral in the numerator [44]: *p*_{11} = −0.0101, *p*_{12} = 0.009 and *p*_{44} = −0.051.

We begin by considering the calculation of *g*_{+,}_{b}_{1} and *g*_{−,}_{b}_{1} for our double-slotted photonic crystal device, i.e., the linear optomechanical couplings of the optical supermodes *a*_{+} and *a*_{−} to the fundamental in-plane flexural mode of either outer slab (we choose *b*_{1} in this case). Table 2 displays the numerically computed coefficients using the perturbation theory described above in terms of the unperturbed optical and mechanical fields. *g*_{+,}_{b}_{1} and *g*_{−,}_{b}_{1} can also be approximated by fitting the anti-crossing curves of Fig. 7 using the dispersion relation given in Eq. (5). The approximate dispersion relation fit values for the linear couplings are also shown in Tab. 2, and compare well to the exact perturbation theory values despite the fact that the couplings derived from the dispersion relations using Eqs. (3) and (5) neglect cross-coupling terms between *a*_{1} and *a*_{2} mediated by the mechanics.

The localized breathing mode of the central nanobeam was found by FEM simulations to be at a mechanical frequency of *ω*_{m}/2*π* ≈ 4 GHz, with linear coupling rates of ${\tilde{g}}_{+}/2\pi =249\phantom{\rule{0.2em}{0ex}}\text{kHz}$ and ${\tilde{g}}_{-}/2\pi =163\phantom{\rule{0.2em}{0ex}}\text{kHz}$ to the *a*_{±} supermodes, respectively, where we have used the notation ${\tilde{g}}_{\pm}={g}_{\pm}{x}_{\text{zpf}}$.

#### 5.2. Linear cross-mode optomechanical coupling

By analogy with Eq. (11), the first order perturbation term for the linear cross-mode coupling *g _{ij}*, where

*i*≠

*j*, can be written as [24]

In the case of the double-slotted photonic crystal of this work, we have for the shifting boundaries contribution to the cross-mode coupling between the supermodes *a*_{+} and *a*_{−} at the symmetric ({*x*_{k, eq}}=0) equilibrium position (center of the anti-crossing curve of Fig. 7):

Note that for the flexural mechanical modes of the photonic crystal structure (either slab or central nanobeam modes) we expect this to be the dominant contribution to the optomechanical coupling. Expanding ${\mathbf{E}}_{+}=\left({\mathbf{E}}_{1}+{\mathbf{E}}_{2}\right)/\sqrt{2}$ and ${\mathbf{E}}_{-}=\left({\mathbf{E}}_{1}-{\mathbf{E}}_{2}\right)/\sqrt{2}$ in terms of the slot modes, and neglecting the cross terms such as ${\mathbf{E}}_{1}^{*}\cdot {\mathbf{E}}_{2}$ due to the small spatial overlap between the fields of the modes localized in separated slots, we obtain Eq. (4) again, *g*_{+−,}* _{k}* = (

*g*

_{1}

*−*

_{,k}*g*

_{2,}

*)/2.*

_{k}Consider now the flexural modes of the central nanobeam. At the symmetric ({*x _{k}*

_{,eq}} = 0) equilibrium position, the nanobeam’s in-plane flexural modes are such that

*g*

_{1}= −

*g*

_{2}. Therefore, at the center of the anti-crossing curve

*g*

_{+−}is maximal and equal to the linear coupling of the nanobeam mode to the

*a*

_{1,2}slot modes. Table 3 shows the numerically computed linear cross-mode coupling rate ${\tilde{g}}_{ij}={g}_{ij}{x}_{\text{zpf}}$ for the fundamental and higher order nanobeam in-plane flexural modes along with their respective frequencies simulated for slot sizes

*s*

_{1}=

*s*

_{2}= 95 nm using Eq. (19). Due to the tight localization of optical modes (see Fig. 8) we find there is still significant coupling to higher order flexural modes of frequencies all the way up to 1 GHz. For the numerical simulations of the

*a*

_{±}optical supermodes of the double-slotted photonic crystal structure we also find that the radiation-limited optical quality factor is theoretically equal to 5 × 10

^{6}and 3 × 10

^{5}for the odd and even modes respectively. Therefore, as noted earlier, we can expect mechanical modes of frequencies

*ω*

_{m}/2

*π*> 300 MHz to be in the resolved-sideband regime.

#### 5.3. Quadratic optomechanical coupling

By extending the perturbation theory to the second order, it is also possible to calculate the *x*^{2}-coupling strength [45–47]. We obtain, for a given optical mode *a _{i}*

In the case of the supermodes *a*_{+} and *a*_{−} of the symmetric double-slotted photonic crystal structure ({*x _{k}*

_{,eq}} = 0), the first term vanishes and the only contribution to the

*x*

^{2}-coupling comes from the second term. For optical splittings such that 2

*J*≪

*ω*

_{0},

*x*

^{2}-coupling strengths. The values of ${\tilde{g}}^{\prime}={g}^{\prime}{\left({x}_{\text{zpf}}\right)}^{2}$ are summarized in Table 3 for the nanobeam in-plane flexural modes up to 884 MHz. Here we assume an optical

*a*

_{±}supermode frequency splitting of 2

*J*/2

*π*= 1 GHz, which is close to the minimum splitting based on the estimated optical quality factors which allows the optical supermodes to be selectively excited and interrogated.

#### 5.4. Coupling coefficients as a function of a static displacement of the nanobeam

In the analysis of the optomechanical coupling coefficients described above in Section 5.1–5.3 we considered a symmetric double-slotted structure with equal slot widths *s*_{1} = *s*_{2} (equilibrium position {*x _{k}*

_{,eq}} = 0), and thus the calculations were done for the optical supermodes

*a*

_{±}. In practice, we could find this symmetric condition by tuning one of the slabs to adjust the relative slots sizes until the optical frequency spectrum was at the center of the anti-crossing curve as shown in Fig. 7. Far away from the center of the anti-crossing, however, the optical supermodes correspond more closely to the individual slot modes

*a*

_{1}and

*a*

_{2}, and we expect different optomechanical coupling strengths. Here we describe how the optomechanical coupling coefficients change upon a large, static displacement of the central nanobeam which takes us far from the symmetric condition near the center of the anti-crossing curve.

From the approximate analytical expression of the supermode dispersion (see Eq. 5), we derive here approximate expressions for ${g}_{\pm ,{b}_{3}}\left({x}_{3,\text{eq}}\right)$, ${g}_{+-,{b}_{3}}\left({x}_{3,\text{eq}}\right)$ and ${{g}^{\prime}}_{\pm ,{b}_{3}}\left({x}_{3,\text{eq}}\right)$ as a function of the static displacement amplitude *x*_{3,eq} of the fundamental nanobeam mode:

As discussed in Sec. 4.1, *x*_{3,eq} can be approximated by a static displacement ${\overline{x}}_{3}$ of the whole nanobeam. Using the perturbative calculation of the optomechanical coupling coefficients from the numerically simulated optical and mechanical fields one can then obtain ${g}_{\pm ,{b}_{3}}\left({\overline{x}}_{3}\right)$, ${g}_{+-,{b}_{3}}\left({\overline{x}}_{3}\right)$ and ${{g}^{\prime}}_{\pm ,{b}_{3}}\left({\overline{x}}_{3}\right)$ by simulating a structure with the nanobeam displaced from its equilibrium position by ${\overline{x}}_{3}$ (this becomes the new “unperturbed” structure in our perturbative calculations). Here we consider a structure with nominal slot widths of *s* = 90 nm at equilibrium, and the displacement of the nanobeam from its equilibrium position is swept from ${\overline{x}}_{3}=-3\phantom{\rule{0.2em}{0ex}}\text{nm}$ to ${\overline{x}}_{3}=3\phantom{\rule{0.2em}{0ex}}\text{nm}$ in steps of 0.5 nm. At each position we calculate the coupling coefficients between the optical supermodes and the fundamental in-plane flexural mode of the nanobeam. The results of these simulations and calculations are plotted in Fig. 10. A fit to the numerically calculated coefficients using Eqs. (22–24) show very good agreement.

These results confirm our previous speculations that far from the symmetric equilibrium position the linear optomechanical coupling is the dominant optomechanical interaction between the optical supermodes and the flexural modes of the central nanobeam, while at symmetric equilibrium position of the beam $\left({\overline{x}}_{3}=0\right)$ the optomechanical interaction is predominantly *x*^{2}-coupling. The linear cross-mode coupling between the optical supermodes is also maximal at ${\overline{x}}_{3}=0$. Experimentally, a static displacement of the whole nanobeam by ${\overline{x}}_{3}$ can be mimicked by displacing both outer slabs in the same direction by an amount ${\overline{x}}_{3}$. Conversely, a change in the equilibrium slot size can be achieved by displacing both outer slabs in opposite directions. In our recent experimental realization of the double-slotted photonic crystal cavity structure we used this tuning degree of freedom by integrating a set of independent capacitors on the outer slabs [34]. Note that the magnitude of the *x*^{2}-coupling in these simulations is limited by the large optical mode splitting (2*J*/2*π* = 117 GHz) of the simulated cavity. With the use of capacitive tuning methods much finer splitting can be realized in practice, however, a larger slot structure was chosen for the simulations to reduce the required mesh size of the structure and to allow for more coarse sweeping of the slot size.

## 6. Summary

We presented a general formalism for studying linear and quadratic coupling of mechanical motion with optical fields in a multimoded optical cavity. This formalism is used to model and design a double-slotted photonic crystal cavity structure previously studied by us experimentally in [34]. The device supports two high-*Q* optical resonances at telecommunication wavelengths and mechanical resonances spanning both the unresolved sideband and resolved sideband regimes. We find that depending on the symmetry of the photonic crystal structure, the optical supermodes can be made to interact either linearly or quadratically with mechanical motion. In particular, we show that the splitting between the optical supermodes can be strongly suppressed, which provides a significant enhancement of the quadratic (*x*^{2}) coupling. The linear and quadratic optomechanical coupling coefficients are calculated using perturbation theory and finite-element-method simulations of the (unperturbed) optical modes and mechanical displacement fields. The simulated optical quality factors of the photonic crystal cavity are calculated to be of order *Q* ≥ 5 × 10^{5}, placing mechanical modes of frequencies *ω*_{m}/2*π* ≥ 400 MHz in the resolved sideband regime. With such parameters, it is shown that the designed photonic crystal device can achieve zero-point *x*^{2}-couplings as large as ${\tilde{g}}^{\prime}/2\pi =1\phantom{\rule{0.2em}{0ex}}\text{kHz}$ for the 10.8 MHz fundamental in-plane flexural mode of the structure, and ${\tilde{g}}^{\prime}/2\pi =1\phantom{\rule{0.2em}{0ex}}\text{Hz}$ for a higher order flexural mode of frequency 884 MHz in the resolved sideband regime. Notably, the splitting of the photonic crystal cavity supermodes may be tuned by adjusting the cavity slot widths to a minimal value approaching that of the cavity mode linewidths (~ 1 GHz), greatly enhancing the achievable *x*^{2} coupling. Although QND measurements of the stored mechanical energy will require substantially lower optical cavity losses [4, 5, 34], with this scale of *x*^{2}-coupling it is feasible to consider a number of other interesting experiments with slightly less restrictive parameters. These include the quantum measurement of phonon or photon shot noise [48], or utilizing interference of quantum noise in the bad-cavity limit [40], a continuous position measurement of *x*^{2}. The latter measurement requires additional dissipative optomechanical coupling, which the photonic crystal cavity structure of this work also possesses.

## Funding

This work was supported by the AFOSR Hybrid Nanophotonics MURI (FA9550-12-1-0024), the Institute for Quantum Information and Matter, an NSF Physics Frontiers Center (NSF Grant PHY-1125565) with support of the Gordon and Betty Moore Foundation (GBMF-2644), the Alexander von Humboldt Foundation, and the Max Planck Society.

## References and links

**1. **M. Aspelmeyer, T. J. Kippenberg, and F. Marquardt, “Cavity optomechanics,” Rev. Mod. Phys. **86**, 1391–1452 (2014).

**2. **A. Jayich, J. Sankey, B. Zwickl, C. Yang, J. Thompson, S. Girvin, A. Clerk, F. Marquardt, and J. Harris, “Dispersive optomechanics: a membrane inside a cavity,” New Journal of Physics **10**, 095008 (2008).

**3. **J. Thompson, B. Zwickl, A. Jayich, F. Marquardt, S. Girvin, and J. Harris, “Strong dispersive coupling of a high-finesse cavity to a micromechanical membrane,” Nature **452**, 72–75 (2008).

**4. **M. Ludwig, A. H. Safavi-Naeini, O. Painter, and F. Marquardt, “Enhanced quantum nonlinearities in a two-mode optomechanical system,” Phys. Rev. Lett. **109**, 063601 (2012).

**5. **H. Miao, S. Danilishin, T. Corbitt, and Y. Chen, “Standard quantum limit for probing mechanical energy quantization,” Phys. Rev. Lett. **103**, 100402 (2009).

**6. **V. Braginsky, S. Strigin, and S. P. Vyatchanin, “Parametric oscillatory instability in fabry–perot interferometer,” Phys. Lett. A **287**, 331–338 (2001).

**7. **V. B. Braginsky, S. E. Strigin, and S. P. Vyatchanin, “Analysis of parametric oscillatory instability in power recycled ligo interferometer,” Phys. Lett. A **305**, 111–124 (2002).

**8. **V. Braginsky and S. Vyatchanin, “Low quantum noise tranquilizer for fabry–perot interferometer,” Phys. Lett. A **293**, 228–234 (2002).

**9. **T. Ojanen and K. Børkje, “Ground-state cooling of mechanical motion in the unresolved sideband regime by use of optomechanically induced transparency,” Phys. Rev. A **90**, 013824 (2014).

**10. **L. Fan, K. Y. Fong, M. Poot, and H. X. Tang, “Cascaded optical transparency in multimode-cavity optomechanical systems,” Nat. Commun. **6**, 5850 (2015).

**11. **G. Heinrich, M. Ludwig, J. Qian, B. Kubala, and F. Marquardt, “Collective dynamics in optomechanical arrays,” Phys. Rev. Lett. **107**, 043603 (2011).

**12. **M. Zhang, G. S. Wiederhecker, S. Manipatruni, A. Barnard, P. McEuen, and M. Lipson, “Synchronization of micromechanical oscillators using light,” Phys. Rev. Lett. **109**, 233906 (2012).

**13. **M. Bagheri, M. Poot, L. Fan, F. Marquardt, and H. X. Tang, “Photonic cavity synchronization of nanomechanical oscillators,” Phys. Rev. Lett. **111**, 213902 (2013).

**14. **M. H. Matheny, M. Grau, L. G. Villanueva, R. B. Karabalin, M. C. Cross, and M. L. Roukes, “Phase synchronization of two anharmonic nanomechanical oscillators,” Phys. Rev. Lett. **112**, 014101 (2014).

**15. **X.-W. Xu and Y.-J. Li, “Antibunching photons in a cavity coupled to an optomechanical system,” J. Phys. B At. Mol. Opt. Phys. **46**, 035502 (2013).

**16. **H. Flayac, D. Gerace, and V. Savona, “An all-silicon single-photon source by unconventional photon blockade,” Sci. Rep. **5**, 11223 (2015).

**17. **M.-A. Lemonde, N. Didier, and A. A. Clerk, “Antibunching and unconventional photon blockade with Gaussian squeezed states,” Phys. Rev. A **90**, 063824 (2014).

**18. **N. Lörch and K. Hammerer, “Sub-Poissonian phonon lasing in three-mode optomechanics,” Phys. Rev. A **91**, 061803 (2015).

**19. **T. Kippenberg, H. Rokhsari, T. Carmon, A. Scherer, and K. Vahala, “Analysis of radiation-pressure induced mechanical oscillation of an optical microcavity,” Phys. Rev. Lett. **95**, 033901 (2005).

**20. **A. G. Krause, M. Winger, T. D. Blasius, Q. Lin, and O. Painter, “A high-resolution microchip optomechanical accelerometer,” Nat. Photon. **6**, 768–772 (2012).

**21. **M. Schmidt, M. Ludwig, and F. Marquardt, “Optomechanical circuits for nanomechanical continuous variable quantum state processing,” New Journal of Physics **14**, 125005 (2012).

**22. **M. Winger, T. D. Blasius, T. P. M. Alegre, A. H. Safavi-Naeini, S. Meenehan, J. Cohen, S. Stobbe, and O. Painter, “A chip-scale integrated cavity-electro-optomechanics platform,” Opt. Express **19**, 24905–24921 (2011).

**23. **K. Stannigel, P. Komar, S. J. M. Habraken, S. D. Bennett, M. D. Lukin, P. Zoller, and P. Rabl, “Optomechanical quantum information processing with photons and phonons,” Phys. Rev. Lett. **109**, 013603 (2012).

**24. **D. Chang, A. H. Safavi-Naeini, M. Hafezi, and O. Painter, “Slowing and stopping light using an optomechanical crystal array,” New Journal of Physics **13**, 023003 (2011).

**25. **A. H. Safavi-Naeini, T. P. Mayer Alegre, J. Chan, M. Eichenfield, M. Winger, Q. Lin, J. T. Hill, D. E. Chang, and O. Painter, “Electromagnetically induced transparency and slow light with optomechanics,” Nature **472**, 69–73 (2011).

**26. **S. Weis, R. Rivière, S. Deléglise, E. Gavartin, O. Arcizet, A. Schliesser, and T. J. Kippenberg, “Optomechanically induced transparency,” Science **330**, 1520–1523 (2010).

**27. **J. C. Sankey, C. Yang, B. M. Zwickl, A. M. Jayich, and J. G. Harris, “Strong and tunable nonlinear optomechanical coupling in a low-loss system,” Nature Physics **6**, 707–712 (2010),

**28. **D. Lee, M. Underwood, D. Mason, A. B. Shkarin, S. W. Hoch, and J. G. E. Harris, “Multimode optomechanical dynamics in a cavity with avoided crossings,” Nat. Commun. **6**, 6232 (2015).

**29. **X. Chen, C. Zhao, S. Danilishin, L. Ju, D. Blair, H. Wang, S. P. Vyatchanin, C. Molinelli, A. Kuhn, S. Gras, T. Briant, P.-F. Cohadon, A. Heidmann, I. Roch-Jeune, R. Flaminio, C. Michel, and L. Pinard, “Observation of three-mode parametric instability,” Phys. Rev. A **91**, 033832 (2015).

**30. **I. S. Grudinin, H. Lee, O. Painter, and K. J. Vahala, “Phonon laser action in a tunable two-level system,” Phys. Rev. Lett. **104**, 083901 (2010).

**31. **C. Doolin, B. D. Hauer, P. H. Kim, a. J. R. MacDonald, H. Ramp, and J. P. Davis, “Nonlinear optomechanics in the stationary regime,” Phys. Rev. A **89**, 1–6 (2014).

**32. **F. Massel, S. U. Cho, J.-M. Pirkkalainen, P. J. Hakonen, T. T. Heikkilä, and M. A. Sillanpää, “Multimode circuit optomechanics near the quantum limit,” Nat. Commun. **3**, 987 (2012).

**33. **M. Evans, S. Gras, P. Fritschel, J. Miller, L. Barsotti, D. Martynov, A. Brooks, D. Coyne, R. Abbott, R. X. Adhikari, K. Arai, R. Bork, B. Kells, J. Rollins, N. Smith-Lefebvre, G. Vajente, H. Yamamoto, C. Adams, S. Aston, J. Betzweiser, V. Frolov, A. Mullavey, A. Pele, J. Romie, M. Thomas, K. Thorne, S. Dwyer, K. Izumi, K. Kawabe, D. Sigg, R. Derosa, A. Effler, K. Kokeyama, S. Ballmer, T. J. Massinger, A. Staley, M. Heinze, C. Mueller, H. Grote, R. Ward, E. King, D. Blair, L. Ju, and C. Zhao, “Observation of Parametric Instability in Advanced LIGO,” Phys. Rev. Lett. **114**, 161102 (2015).

**34. **T. K. Paraïso, M. Kalaee, L. Zang, H. Pfeifer, F. Marquardt, and O. Painter, “Position-Squared Coupling in a Tunable Photonic Crystal Optomechanical Cavity,” Phys. Rev. X **5**, 041024 (2015).

**35. **S. G. Johnson and J. D. Joannopoulos, “Block-iterative frequency-domain methods for maxwell’s equations in a planewave basis,” Opt. Express **8**, 173–190 (2001).

**36. **M. Eichenfield, R. Camacho, J. Chan, K. J. Vahala, and O. Painter, “A picogram-and nanometre-scale photonic-crystal optomechanical cavity,” Nature **459**, 550–555 (2009).

**37. **J. Chan, M. Eichenfield, R. Camacho, and O. Painter, “Optical and mechanical design of a "zipper" photonic crystal optomechanical cavity,” Opt. Express **17**, 3802 (2009).

**38. **COMSOL Multiphysics, version 3.5a.

**39. **K. Srinivasan and O. Painter, “Momentum space design of high-*q* photonic crystal optical cavities,” Opt. Express **10**, 670–684 (2002).

**40. **Y. Yanay, J. C. Sankey, and A. A. Clerk, “Quantum backaction and noise interference in asymmetric two-cavity optomechanical systems,” Phys. Rev. A **93**, 063809 (2016).

**41. **M. Eichenfield, J. Chan, R. Camacho, K. Vahala, and O. Painter, “Optomechanical crystals,” Nature **462**, 78–82 (2009).

**42. **S. G. Johnson, M. Ibanescu, M. Skorobogatiy, O. Weisberg, J. Joannopoulos, and Y. Fink, “Perturbation theory for maxwell’s equations with shifting material boundaries,” Phys. Rev. E **65**, 066611 (2002).

**43. **C. Zhao, L. Ju, H. Miao, S. Gras, Y. Fan, and D. G. Blair, “Three-mode optoacoustic parametric amplifier: a tool for macroscopic quantum experiments,” Phys. Rev. Lett. **102**, 243902 (2009).

**44. **D. Biegelsen, “Frequency dependence of the photoelastic coefficients of silicon,” Phys. Rev. B **12**, 2427 (1975).

**45. **H. Kaviani, C. Healey, M. Wu, R. Ghobadi, A. Hryciw, and P. E. Barclay, “Nonlinear optomechanical paddle nanocavities,” Optica **2**, 271–274 (2015).

**46. **A. W. Rodriguez, A. P. McCauley, P.-C. Hui, D. Woolf, E. Iwase, F. Capasso, M. Loncar, and S. G. Johnson, “Bonding, antibonding and tunable optical forces in asymmetric membranes,” Opt. Express **19**, 2225–2241 (2011).

**47. **A. H. Safavi-Naeini, “Quantum Optomechanics with Silicon Nanostructures,” Ph.D. thesis, California Institute of Technology (2013).

**48. **A. A. Clerk, F. Marquardt, and J. G. E. Harris, “Quantum measurement of phonon shot noise,” Phys. Rev. Lett. **104**, 213603 (2010).