Theory and modeling of a quadratically coupled optomechanical sensor for atomic force microscopy

(1)

Theory and modeling of a quadratically coupled

optomechanical sensor for atomic force microscopy

Ariadna Soro ´ Alvarez

Supervised by:

Ermes Scarano and Prof. David B Haviland

July 2020

SK202X Degree project in Applied Physics MSc Engineering Physics - Nanophysics Track

(2)

(3)

Abstract

This thesis is part of the Quantum-Limited Atomic Force Microscopy (QAFM) project [1], which aims at designing and fabricating an optomechanical device capable of detection at the standard quantum limit. In a first step towards this goal, the present work offers a deeper understanding of the behavior of such device through mathematical modeling and simulations.

In particular, we focus in the case in which the coupling between the optical and the mechanical modes of the system is quadratic. We study the feasibility of such sensor for quantum- limited detection by analyzing the optomechanical coupling in the sideband-resolved regime, the effects of parametric amplification and sideband cooling and the influence of different noise sources. We also discuss the interaction between the device and the environment in order to illustrate the complexity of the situation and provide an outlook for future work.

Acknowledgements

First, I would like to thank the supervisors of my thesis, Prof. David Haviland and Er- mes Scarano, for their incommensurable guidance and dedication. Without them, this work would have never been accomplished. I am particularly grateful to David for introducing me to the field of optomechanics and for offering me a place in his group, where he always made sure I felt welcome.

I would like to offer special thanks to Dr. Riccardo Borgani for shedding some light over the noise-related matters and to Gabriele Baglioni for his contributions regarding parametric amplification. Gabriele also deserves my deep gratitude for his advice and criticism, which constituted a great help in the elaboration of the thesis.

Finally, I would like to thank KTH and the QAFM group for branding me the opportunity to work on a master’s thesis of my choice and for providing me the resources to do it.

(4)

(5)

1. Introduction

1.1 Optomechanics

Optomechanics is a branch of physics that studies the interaction between light and mechanical objects at low energy scales. This interaction arises from the momentum carried by light, which generates radiation-pressure forces. Radiation-pressure force was first pos- tulated by Kepler in the 17th century [2] and the first experimental demonstrations were performed in 1901 [3, 4]. However, it was not until the 1970s when optomechanics started taking shape in the hands of Braginsky, in his effort to detect gravitational waves.

Braginsky studied the dynamical influence of radiation pressure on a harmonically suspended end mirror of a cavity. His analysis revealed that the retarded nature of the force, due to the finite cavity lifetime, can cause either damping or antidamping of mechanical motion, and he demonstrated these two effects using a microwave cavity [5, 6]. The damping of the mechanical motion occurs for red-detuned forces and allows oscillators to be cooled down to their motional ground state. This phenomenon not only applies to cavities, but also to smaller systems, such as atoms or ions, and it can be achieved by numerous techniques known as laser cooling. On the other hand, blue-detuned forces lead to amplification of the mechanical motion. Braginsky called the process parametric oscillatory instability and it is otherwise known as radiation-pressure-induced self-oscillations [7].

Braginsky also addressed the fundamental consequences of the quantum fluctuations of radiation pressure and demonstrated that they impose a limit on how accurately the position of a free mass can be measured [8, 9]. In the 1990s, several other aspects of quantum cavity- optomechanical systems were explored theoretically, such as squeezing of light [10, 11] and quantum nondemolition (QND) detection of light intensity [12, 13]. On the experimental side, optical feedback cooling based on the radiation-pressure force was first demonstrated in 1999 for the vibrational modes of a macroscopic end mirror [14].

Since the mid 2000s, cavity optomechanics has advanced rapidly and optomechanical coupling has been reported in numerous novel systems. These include microtoroid resonators [15–18], suspended micromirrors [19, 20], membranes [21] and nanorods [22] inside Fabry- P´erot resonators, whispering gallery microdisks [23, 24] and microspheres [25–27], photonic crystals [28, 29], and evanescently coupled nanobeams [30].

There are several different motivations that drive the rapidly growing interest into cavity optomechanics. On the one hand, there is the highly sensitive optical detection of small forces, displacements, masses, and accelerations. On the other hand, cavity quantum optomechanics promises to manipulate and detect mechanical motion in the quantum regime

(8)

using light, creating nonclassical states of light and mechanical motion. These tools form the basis for applications in quantum information processing, in building hybrid quantum devices, and in testing fundamentals of quantum mechanics [7].

1.2 Atomic Force Microscopy

In 1986, in the quest for measuring ultrasmall forces on particles as small as single atoms, Binnig, Quate and Gerber presented a device called atomic force microscope (AFM) [31]

that was developed after the scanning tunneling microscope (STM) . It was designed to investigate both conductors and insulators on an atomic scale and was able to detect inter- atomic and electromagnetic forces as small as 10⁻¹⁸ N.

The AFM consists of a microcantilever clamped at one end and free at the other, as showed in Figure 1.1. At the free end, a sharp tip is attached whose position is usually measured with an optical lever system. The cantilever behaves as a spring, sensitive to the forces acting between the tip and the sample surface, and changes of the deflection of the cantilever are measured to form an image of the sample. A critical component needed in AFM is a feedback loop controlling the tip-surface separation, to avoid the tip drifting away and losing contact with the surface, or worse, crashing into the surface. This system allows the cantilever to be positioned with sub-nanometer resolution using a piezo-electric scanner. An additional piezo shaker is typically used to excite the cantilever oscillations.

Modifications of the scheme in Figure 1.1 have been explored and implemented throughout the years, but they operate on the same principles. The most well-known variation is the qPlus [32], made out of a quartz tuning fork which uses the piezoelectricity of quartz to read the velocity of the tip as a small electrical current.

Although the AFM was invented more than 30 years ago and it is widely used, there is still a lack of understanding of some aspects of the instrument, mainly due to the nonlinear nature of the tip-surface interaction. Analysis and modeling of AFM is often based on lineariza- tion, which can notably reduce the image contrast or create image artifacts [33]. Another limitation of the AFM is the scanning speed –in order to get images with good resolution, measurements can take from several minutes to several days.

Current limitations of the AFM call for a new design that incorporates a radical improvement in the speed of acquisition and information content of the images. Ultimately, these advances would lead to the development of a force sensor working at the fundamental quantum limit of action and reaction. With several orders of magnitude improvement in sensitivity beyond the current state-of-the-art sensors, it would enable the acquisition of multi-dimensional data sets in seconds, as opposed to several days as is the current practice [1].

Consequently, we aim at realizing a quantum-limited force sensor suitable for scanning probe applications, by building on ideas from the superconducting circuits community, augmenting them with designs and techniques from MEMS sensors and actuators, and applying them to low temperature AFM [1]. Due to its nature, we named the device quantum-limited atomic

(9)

force microscope (QAFM). As we will see in the following chapters through theory, modeling and simulations, the key to designing such an optimal sensor lies in achieving strong optomechanical coupling and very low noise motion detection.

Figure 1.1: Schematic setup of an AFM, as described in Section 1.2. Figure adapted with permission from [34].

(10)

2. Theoretical background

In order to design a quantum-limited force sensor, we can start by studying the quantum theory of optomechanical systems, since it provides general descriptions that can help us discover which features we need to optimize in the realization of the QAFM.

2.1 Quantum description

Let us consider an uncoupled optomechanical system, which can be described by the following Hamiltonian:

H = ˆˆ H_o+ ˆH_m = 1

2~ωo+ ~ωoˆa^†ˆa +1

2~ωm+ ~ωmˆb^†ˆb −→ H = ~ωˆ oaˆ^†a + ~ωˆ mˆb^†ˆb, (2.1) where the subscripts o and m stand for optical and mechanical respectively, ˆa^† and ˆa are the photon creation and anihilation operators, and ˆb^† and ˆb are the phonon creation and annihilation operators. For simplicity, we can omit the terms ~ω/2 as they are constant and only offset the eigenenergies of the system.

Note that the ladder operators ˆb, ˆb^† are related to the position and momentum of the oscillator like

ˆ

x = x_ZPF(ˆb + ˆb^†), p = −imωˆ _mx_ZPF(ˆb − ˆb^†), (2.2) where

x_ZPF= r

~ 2mω_m =

r

~ωm

2k = s

~ 2√

mk (2.3)

is the zero-point fluctuation (ZPF) amplitude of the mechanical oscillator, i.e. the spread of the coordinate in the ground state h0|ˆx²|0i = x²_ZPF , and where |0i denotes the mechanical vacuum state. In Equation (2.3), we introduced the spring constant k = mω_m².

Now, when we couple the optical and mechanical systems, the optical resonance frequency is modulated by the mechanical oscillations in the following way:

ω_o(x) = ω_o(0) + x∂ω_o

∂x +x² 2

∂²ω_o

∂x² + . . . (2.4)

By introducing Equation (2.4) into the Hamiltonian in (2.1) and using (2.2), we obtain a new term corresponding to the interaction Hamiltonian, which depends on the order of the

(11)

modulation:

H = ˆˆ Ho+ ˆHm+ ˆHint, where ˆHint=











~∂ω_o

∂x x_ZPFˆa^†ˆa (ˆb^†+ ˆb) linear

~ 2

∂²ωo

∂x² x²_ZPFˆa^†ˆa (ˆb^†+ ˆb)² quadratic

(2.5)

By quadratic modulation, we mean that the modulation has no first derivative component, i.e. ωo(x) ≈ x²2 ∂²ωo∂x². This case has been much less studied than the linear one and it provides very interesting features, as any measurement of ˆx² instead of ˆx will be closely connected to the oscillator’s energy, and thus the phonon number. This then permits QND detection of the phonon number in the mechanical resonator [7]. For this reason, we choose to design our device with quadratic coupling, analyze its behavior and study its feasibility for quantum-limited detection.

2.2 Strong coupling

As previously mentioned, one of the key factors in designing a QAFM lies in achieving strong optomechanical coupling, which looking at Equation (2.5) translates into maximizing ˆHint. For quadratic coupling, the vacuum optomechanical coupling strength is defined as

g0= 1 2

∂²ω_o

∂x² x²_ZPF, (2.6)

which leads to the definition of optomechanical coupling strength:

g = g0

√no, (2.7)

where n_o = ˆa^†ˆa is the mean photon number. Strong coupling is achieved when g > Γ_o, where Γ_o = ω_o/Q_o is the damping of the optical mode and Q_o is the quality factor [7].

This implies that, to be in the strong coupling regime, we need to increase the drive power (n) and maximize the optical frequency shift per displacement ( ∂²ω_o∂x²), as well as the mechanical ZPF amplitude (x²_ZPF). Nevertheless, a high value of ZPF means larger quantum fluctuations, which raises the precision limit we can achieve. Hence, when defining the characteristics of our mechanical oscillator, we want to make sure there is a good trade- off between strong optomechanical coupling and small fluctuations. We discuss further the measurement limitations in the following section.

Let us remark that, in optomechanics, it is sometimes more useful to evaluate the system in terms of the so-called optomechanical cooperativity C rather than the coupling strength:

C = 4g²

Γ_oΓ_m, (2.8)

where Γ_m= ω_m/Q_m is the damping of the mechanical mode and Q_m its quality factor [35].

The cooperativity can be understood as the rate at which one gathers information about the oscillator’s motion relative to optical and mechanical decay rates [36].

(12)

2.3 Standard quantum limit

When carrying out a measurement with an optomechanical device, the maximum precision we can achieve is limited by several sources of noise that are present in the system [7].

Therefore, it is often more useful to analyze the mechanical motion in the frequency space, rather than as a time-evolving variable x(t).

Given one particular realization of the trajectory x(t) obtained during a measurement time T , we define the Fourier transform over a finite time interval T :

φ(ω) = FT{x(t)} =˜ 1

√T Z T

0

φ(t)e^iωtdt. (2.9)

Averaging over independent runs, we obtain the spectral density D

|˜x(ω)|² E

, which in the limit T → ∞ is equal to the power spectral density, S_xx(ω) [37]:

T →∞lim

D|˜x(ω)|²E

= S_xx(ω) = Z +∞

−∞

hx(t)x(0)i e^iωtdt. (2.10) From Equation (2.10), we obtain that the area under the mechanical noise spectrum yields the variance of the mechanical displacementx², which for weak damping (Γ_m ω_m) is set by the equipartition theorem:

Z +∞

0

Sxx(ω)dω

π =x² = kBT

mω² = kBT

k , (2.11)

where T is the effective temperature of the oscillator. Let us remark that the variance can- not be smaller than the ZPF, i.e. x² ≥ x²_ZPF.

As depicted in Figure 2.1, we expect the noise spectrum of the mechanical motion to be a Lorentzian curve centered at ω_m, with a FWHM determined by the mechanical damping Γm and area x². Mathematically, this is

S_xx(ω) =x² Γ_m/2π

(ω − Ω_m)²+ (Γ_m/2)². (2.12) The first source of noise that contributes to the power spectral density is the intrinsic noise, of both quantum and thermal nature. As seen in Equation (2.3), the quantum fluctuations correspond to the spread of the mechanical amplitude in the ground state, and are determined by the characteristics of the harmonic oscillator: mass and resonant frequency.

Therefore, the quantum spectral noise is [37]

S_xx^ZPF(ω) = x²_ZPF Γ_m/2

(|ω| − ωm)²+ (Γm/2)². (2.13) On the other hand, the thermal noise, which can be derived from the fluctuation-dissipation theorem, adds a factor equivalent to the Bose-Einstein occupation number nBE [37]:

S_xx^int(ω) = nBE(~ωm) S_xx^ZPF(ω) = x²_ZPF e^~ω^m^/k^B^T − 1

Γm/2

(|ω| − ω_m)²+ (Γ_m/2)². (2.14)

(13)

Figure 2.1: Power spectral density of the mechanical motion, Sxx(ω). The curve presents a Lorentzian shape (see Equation (2.12)) centered at ω_m, with a FWHM equal to Γ_mand areax².

Note that in the classical limit (k_BT ~ωm), integrating the intrinsic noise yields the equipartition theorem result from Equation (2.11).

Another type of noise called backaction noise arises at sufficiently high drive power due to the radiation-pressure forces exerted on the mechanical oscillator. Evading these fluctuations is not as trivial as decreasing the drive intensity, since then imprecision or shot noise becomes relevant due to the quantized nature of light. In this case, the mean number of photons driving the device is so low that fluctuations of a few photons become very large in comparison. In fact, the uncertainty between backaction and imprecision noise is limited by Heisenberg’s principle [37]:

S_{F F}^ba S_xx^imp ≥ ~²

4. (2.15)

As we can see in Figure 2.2, the uncertainty relation is minimal when the measurement imprecision noise and backaction noise each make equal contributions to the added noise.

In such case, we say we have reached the standard quantum limit (SQL) of detection.

We can also define a stronger quantum limit that refers to the total inferred position noise from the measurement, not only the added noise. This limit is reached when the temperature is low enough so that the thermal noise is reduced to level of ZPF, i.e. when the oscillator is cooled to its ground state. This yields a total noise that is twice the oscillator’s zero point noise [37]:

S_xx^tot ≥ 2S_xx^ZPF= x²_ZPF Γm

(|ω| − ωm)²+ (Γm/2)². (2.16) Half the noise here is from the oscillator itself and half is from the added noise of the detector.

(14)

Figure 2.2: Added noise evaluated at the mechanical resonance (pink), plotted in logarithmic scale as a function of the drive power. When T = 0, the imprecision noise (purple) dominates at low powers, while at higher powers the backaction noise (blue) represents the most important contribution. The standard quantum limit (green) is reached at intermediate powers, when S_xx^add= S^ZPF_xx .

Figure 2.3: Full measured noise spectrum, with contributions from the intrinsic fluctuations of the mechanical oscillator (green), the imprecision in the measurement (purple) –typically flat in frequency–, and the backaction heating of the oscillator (blue).

Figures adapted with permission from [38].

Now, we have done a preliminary quantum description of our device that gives us a clearer idea of the parameters we need to optimize when designing the QAFM. In the next chapter, we present our novel design and a classical mathematical model to describe it, as well as an analogy to the quantum theory that we have derived here.

(15)

3. Design and modeling of the QAFM

Having the theory in mind and accounting for the limitations and difficulties of the nanofab- rication, we came up with the design shown in Figures 3.1, 3.2 and 3.3. The optical part of the device consists of an LC circuit, where the inductor is made of superconducting me- anders and the capacitor is interdigitated. One of the plates of the capacitor is fixed and connected to the inductor, whereas the other is part of a flexible grounded cantilever that can move in and out of plane. In this way, the capacitor also acts as the mechanical part of the device, oscillating around the axis marked as a dashed line in Figure 3.3. The tip that interacts with the surface is placed at the opposite end of the rotation axis.

Note that we opted for a microwave circuit instead of an optical cavity because it is easier to drive at a specific frequency and therefore choose the detuning more accurately. However, the disadvantage of microwave circuits is that, to reach the ground state, they need to be operated at cryogenic temperatures.

Figure 3.1: Preliminary design of the QAFM, where the optical circuit is formed by an inductor and a variable capacitor and the mechanical oscillator is the ground plate of the capacitor. The dark green parts will be fabricated with Si, while the lighter green areas will be made of SiN with a thin film of NbTiN on top.

(16)

Figure 3.2: Side view of the QAFM design shown in Figure 3.1. The motion of the capacitor is coupled to the motion of the tip in a torsional/flexural manner.

Figure 3.3: Top view of the QAFM design shown in Figures 3.1 and 3.2, where L and C(φ) denote the inductor and the capacitor, respectively. This design can be modeled as the circuit shown in Figure 3.4.

Vext

R

L

C(φ)

Figure 3.4: Scheme of the optomechanical circuit. Vext denotes the external voltage supplied to the circuit, R = 50Ω denotes the resistance of the transmission line, L is the inductor and C(φ) is the variable capacitor, which also acts as a mechanical oscillator.

(17)

As we have seen in the previous chapter, the mechanical motion of the oscillator is typically denoted with the coordinate x, describing vertical displacement. However, due to the nature of our design, the tip and the capacitor oscillate in a mixture of torsional and flexural movement (see Figure 3.2) that is difficult to describe analytically. Since the tip is closer to the axis of rotation than the capacitor, and the cantilever is free to flex, a very small motion at the tip corresponds to a much larger motion at the capacitor. This mechanical amplification of motion allows us to achieve a larger change of capacitance for small tip movements, typically limited to the order of 1 nm in AFM. This behavior will be studied in more detail in the future through COMSOL^R simulations and laboratory experiments, but for the purpose of this thesis, we will describe the oscillator’s motion through the angular rotation from equilibrium (φ) instead of the vertical displacement (x) and we will assume that this rotation experiences an effective torsional stiffness κ.

Note that in this layout, the capacitance is maximum at the equilibrium point and it depends on the angle of deviation, φ:

C(φ) = C(0) + C⁰φ + C⁰⁰φ²

2 + . . . (3.1)

where C⁰ = 0 and C⁰⁰< 0. This modulation is equivalent to the one in Equation (2.4) and it constitutes the quadratic optomechanical coupling we were aiming for.

In order to study in detail the behavior of the QAFM, we can start by describing it through a lumped element model. According to our design, the capacitor (C) is connected in series with the inductor (L) which is powered by an external voltage (V_ext) through a transmission line. This last element adds a load to the circuit equivalent to a resistance of R = 50Ω, thus turning our model into a series RLC circuit, as shown in Figure 3.4. In this circuit, the optomechanical coupling is introduced through the capacitance: C = C(φ). The series RLC circuit has an optical resonant frequency¹ of ω_o = p 1/LC(0) , while the capacitor plate oscillates at the mechanical resonant frequency denoted by ωm.

It is worth noting that the resonator is overcoupled, meaning that the electrical losses (Γ_o) are dominated by the external coupling to the transmission line, rather than an intrinsic decay in the circuit, i.e. Γex Γ_in ⇔ Γ_o ∼ Γ_ex. This is a desirable effect, since in this way, the losses correspond to signal transferred –and so measured– to the next stages of amplification via the transmission line. On the mechanical side, the losses are represented by the intrinsic damping of the oscillator (Γm).

Mathematically, the RLC circuit is described by V_ext = V_L+ V_R+ V_C = Ld²Q

dt² + RdQ dt + Q

C(φ), (3.2)

where Q is the charge in the capacitor and L = (ω_o²· C(0))⁻¹. We will drive the circuit with

1Note that we use the word optical for consistency with the cavity optomechanics literature. In fact, we are referring to any electromagnetic mode with a resonant frequency much higher than the mechanical resonant frequency, as in practice ωowill be in the microwave regime (∼ 1 − 10 GHz).

(18)

a sinusoidal wave at a frequency ω_d: V_ext = V_dsin(ω_dt).

As for the mechanical oscillator, it follows the equation of motion of a damped harmonic oscillator driven by an external torque τ :

τext= Id²φ

dt² + ηdφ

dt + κφ ^×κ

−1

=⇒ τ_ext κ = 1

ω²_m d²φ

dt² + 1 ωmQm

dφ

dt + φ, (3.3) where I = κ/ω_m² is the moment of inertia, η = κ/(ωmQm) = κΓm/ω²_mis the angular damping constant and κ is the torsional spring constant. In the absence of inertial actuation, the torque is generated by the electrostatic energy E stored in the capacitor and therefore,

τext= −dE

dφ = − d dφ

1 2

Q² C

= Q² 2C²

dC

dφ −

0 Q C

dQ dφ

Eq. (3.1)

= Q²C⁰⁰φ

2C² . (3.4) Note that we are neglecting the charge variation with the angle. We can make this assump- tion because we consider the electrostatic energy stored in the capacitor as if it was isolated from the rest of the circuit. In such circumstances, the charge remains constant with the separation of the capacitor’s plates.

By putting together Equations (3.2)-(3.4), we obtain two coupled differential equations that describe the optomechanical circuit:

Ld²Q

dt² + RdQ dt + Q

C(φ) = Vdsin(ωdt) (3.5)

1 ω²_m

d²φ

dt² + 1 ωmQm

dφ

dt + φ = Q²C⁰⁰φ

2κC²(φ) (3.6)

Note that these equations are coupled by both the variation of the capacitance with the deviation from the equilibrium point, and the variation of the torque (charge) with time.

Therefore, in order to obtain a strong coupling, we need to have C⁰⁰ as large as possible.

This would improve the precision of the measurements, since it would allow the detection of very small changes in φ.

3.1 Strong coupling regime

As we have seen, in our device it is more convenient to measure the displacement of the mechanical oscillator using coordinate φ instead of the more common x. This means that, in our case, the optomechanical strength from Equations (2.6) and (2.7) is defined by:

g = g₀√ n_o= 1

2

∂²ω_o

∂φ² φ²_ZPF√

n_o= C⁰⁰ω_o 4C(0)

~ωm

2κ

√n_o, (3.7)

where we used that ω_o(φ) =

s 1

LC(φ) ' ω_o(0) + φ² 2

∂²ω_o

∂φ² and φ_ZPF= r

~ωm

2κ . (3.8)

(19)

Consequently, we will be in the strong coupling regime when the drive power is high enough to satisfy the following condition:

g > Γo =⇒ √

no > 8κC(0)

~ωmC⁰⁰Q_o. (3.9)

3.2 Implementation

In order to simulate the behavior of our sensor, we can code Equations (3.5)-(3.6) into a script that computes the time evolution of our variables: φ(t), ˙φ(t), Q(t), ˙Q(t). We decide to code in Python and make use of its numerical integrator scipy.integrate.odeint, but since odeint only integrates first-order differential equations, we first need to convert each of our second-order equations of motion into a system of first-order ODEs:

Ld²Q

dt² + RdQ dt + Q

C(φ) = V_dsin(ω_dt) =⇒









 dQ

dt = q dq dt = 1

L

Vdsin(ωdt) − Q

C(φ) − Rq

(3.10)

1 ω²_m

d²φ

dt² + 1 ω_mQ_m

dφ

dt + φ = Q²C⁰⁰φ 2κC²(φ) =⇒









 dφ

dt = Φ dΦ

dt = ω²_m Q²C⁰⁰φ

2κC²(φ) − φ − Φ ω_mQ_m

.

(3.11) By setting initial conditions φ(0), ˙φ(0), Q(0), ˙Q(0), odeint is ready to integrate these four coupled equations. However, there is still one piece missing in our model. As discussed in Section 2.3, an optomechanical device shows different types of noise that so far, we are not considering in the equations of motion.

Physically, noise is a force (or a torque, in our case) that drives the oscillations together with τext. Therefore, we can rewrite Equations (3.6) and (3.11) like:

1 ω²_m

d²φ dt²+ 1

ω_mQ_m dφ

dt+φ = τext+ τnoise

κ =⇒









 dφ

dt = Φ dΦ

dt = ω²_m Q²C⁰⁰φ

2κC²(φ) −τnoise

κ − φ − Φ

ω_mQ_m

. (3.12) We will properly define τ_noise and explain how we can include it into our script in Section 4.3. For now, let us mention that odeint is a variable step integrator, which means that adding random terms into the equations could create unwanted correlations. Thus, we need to be thorough in the implementation of noise into our program or we will end up observing effects in the results that are not actually occurring.

(20)

4. Simulation results and discussion

In this chapter, we present and discuss some key results after simulating the mathematical model developed in the previous section. As shown in Table 4.1, we begin by introducing the parameters we use in the computation.

We want to work with an optical resonant frequency ωo within the typical circuit QED frequency range, and we choose in particular 3 GHz because it suits best the experimental equipment we want to use. Then, by choosing quality factor for the optical resonator in the order of 1000, we obtain a lower bound for the mechanical resonance frequency we need in order to reach the resolved sideband limit (see Section 4.1). Regarding the capacitance at equilibrium, we aim for a value as small as possible but, at the same time, sufficiently large to overcome parasitic effects. By simulating in COMSOL^R and Sonnet, the current design yields a a capacitance in the order of 1 fF. As for the second derivative of the capacitance, it is estimated from Python simulations, assuming that the distance between the capacitor gap and the rotation axis is in the order of 10-100µm .

Table 4.1: Values used in the simulation of the model described by Equations (3.5)-(3.6).

Parameter Symbol Value Calculated from

Optical frequency ωo/2π 3 GHz Experimental equipment

Resistance R 50 Ω Transmission line

Capacitance in equilibrium C(0) 1 fF COMSOL^R simulations Capacitance second derivative C⁰⁰ −2·(10⁻¹¹−10⁻⁹)

F/rad² Python simulations Inductance L 2.81 · 10⁻⁶ H L = 1/ω_o²C(0) Optical quality factor Q_o 1061.03 Q_o = 1/RpL/C(0) Optical damping coefficient Γo/2π 2.83 · 10⁶ Hz Γo= ωo/Qo

Mechanical frequency ω_m/2π 10 MHz Sideband-resolved limit Torsional spring constant κ 2.2 · 10⁻⁹ Nm/rad COMSOL^R simulations Mechanical quality factor Qm ∼ 10⁵− 10⁶ Cantilever in vacuum Mechanical damping coefficient Γm/2π 10 − 100 Hz Γm= ωm/Qm

Zero-point fluctuations φZPF 1.23 · 10⁻⁹ rad Equation (3.8)

(21)

When it comes to mechanical parameters, an effective torsional spring constant is derived through COMSOL^R simulations, while the mechanical quality factor is set to be reasonable for a cantilever in vacuum (∼ 10⁵− 10⁶), although we might work with lower values to save simulation time. The rest of the parameters are calculated according to the formulas shown in Table 4.1.

4.1 Sideband-resolved regime

As described by Equation (3.5) the oscillations of the mechanical system are modulating the capacitance, thus varying the behavior of current and charge in the RLC circuit. If the coupling strength is big enough so that φ(t) induces significant changes in the capacitance, we expect the current to show additional frequency components related to oscillations in C. In fact, if the quality factor of the optical resonator is large enough (Qo ω_o/ωm), sidebands will arise around the drive frequency, ω_d. In this case, we say we are in the sideband-resolved regime [7].

Particularly, when driving on resonance, peaks can be seen at ω_d± ω_m in linear systems, and at ωd± 2ω_m in devices with quadratic coupling. The origin of this factor 2 is easy to pinpoint if we picture the dependence of the capacitor on the angle φ in our QAFM: C(φ) is maximum every time φ = 0, namely twice per mechanical period, hence the modulation of the current occurs at twice the mechanical frequency.

We commonly refer to the left (−) and right (+) sidebands as red-detuned and blue-detuned respectively. As described in Section 1.1, pumping the red sideband causes extra damping, which leads to cooling, whereas pumping the blue sideband causes antidamping, which amplifies the thermal fluctuations and eventually can create instabilities. The physics behind these phenomena is that, when the optical circuit is detuned, the mechanical position is imprinted on the amplitude of the optical field, which then back-acts through radiation pressure upon the mechanical oscillator. Since the optical circuit also induces a delay in the optical response, this dynamical back-action is retarded, with a component of the optical force being proportional to the velocity of the mechanical oscillator. Depending on the sign of this component, it either damps/cools or amplifies/heats the mechanical motion [35].

This effect can also be explained in the context of Raman scattering. Photons impinging at a frequency red detuned from the optical resonance will, via the optomechanical interaction, preferentially scatter upwards in energy in order to enter the optical resonance (anti-Stokes process). They will do so by absorbing phonons from the mechanical oscillator and thus cooling it down. On the other hand, blue-detuned photons will preferably scatter down- wards in energy (Stokes process), leaving phonons behind in the oscillator and therefore heating it [7]. We will study these effects in detail in the following sections.

In Figure 4.1, we show the frequency of the current that circulates through the RLC circuit ( dQ/dt ) by plotting the absolute value of its Fourier transform and we observe the sidebands as expected.

(22)

Figure 4.1: Simulation of the lumped element model showing the response around the electromagnetic resonant frequency in the resolved sideband regime. The electromagnetic mode is driven at ωd= ωo= 3 GHz (delta-function peak in the center). Due to the coupling to the mechanical mode, the two sidebands occur at ωd± 2ωm, with the factor 2 due to the quadratic coupling.

4.2 Optical spring effect

One thing we notice when analyzing the results of our simulations is that the sidebands are not static, they shift under certain conditions. What we observe is the optical spring effect, a phenomenon in which the radiation pressure interaction prompts a shift of the mechanical resonance frequency [35]. In particular, our device behaves like a parametric oscillator due to the displacement dependence of the capacitance and therefore it is useful to define an effective frequency of the oscillations, ωm,eff.

From Equations (3.3) and (3.4), we have τ_ext= Id²φ

dt² + ηdφ

dt + κφ = Q²(t)C⁰⁰

2C²(φ) φ =⇒ Id²φ dt² + ηdφ

dt +

κ − Q²(t)C⁰⁰ 2C²(φ)

| {z }

κeff

φ = 0. (4.1)

Now, as we can see, the φ-dependence of the torque allows us to define an effective spring constant, κ_eff. Note that its value is not constant anymore, but time-dependent (κ_eff = κeff(t)). However, we can use its mean value to calculate the effective frequency of the oscillations:

hκ_effi^φ1' κ −Q²(t) C⁰⁰

2C²(0) , (4.2)

where we used that the φ-dependence of C²(φ) is negligible, since 1

C²(φ) = 1

C²(0)− C⁰⁰φ²

C³(0) + · · ·^φ1≈ 1

C²(0). (4.3)

(23)

Figure 4.2: Response of the optical circuit when driven at frequency ωo = 3 GHz and different voltage amplitudes. By increasing the drive voltage, the mechanical frequency shifts according to Equation (4.4), which shows in the Fourier transform of the current as a displacement of the sidebands to ωo± 2ωm,eff.

Hence, this new κ_eff induces a shift in frequency as follows:

ω_m,eff² = hκ_effi I = κ

I − Q²(t) C⁰⁰

2IC²(0) = ω²_m−Q²(t) C⁰⁰

2IC²(0) . (4.4)

In this expression, ω_m,effdescribes the actual frequency at which the capacitor oscillates, and therefore, as shown in Figure 4.2, in the resolved sideband regime we find the sidebands at ω_d± 2ω_m,eff. Let us point out that in linear systems, this shift in mechanical frequency can only be achieved by driving the system detuned from the optical resonant frequency [35].

However, due to the quadratic coupling, our system shows the optical spring effect by driving either on or off resonance, as long as the drive power is sufficient to make the charge oscillations relevant.

4.3 Sideband cooling

As we explained in Section 4.1, it is possible to cool down the mechanical oscillations of an optomechanical device through back-action. This effect shows as a decrease of the amplitude of the oscillations further than the one would expect from intrinsic damping. Therefore, in order to observe it in the simulations, we introduce thermal noise that can be reduced if the cooling is large enough.

(24)

Let us recall from Equation (3.12) that noise is implemented in our equations as a torque:

1 ω²_m

d²φ

dt² + 1 ω_mQ_m

dφ

dt + φ = τ_ext+ τ_noise κ

Eq. (4.2)

=⇒

1 ω²_m

d²φ

dt² + 1 ω_mQ_m

dφ dt +κ_eff

κ φ = τ_noise

κ , (4.5)

where we used the definition of κeff to account for the drive. Now, if we Fourier transform Equation (4.5), we can find an expression for the mechanical susceptibility, χ(ω):

− ω² ω_m²

φ(ω) +˜ iω ω_mQ_m

φ(ω) +˜ κeff

κ

φ(ω) =˜ τ˜noise(ω)

κ =⇒

χ(ω) = φ(ω)˜

˜

τnoise(ω) = 1/κ κ_eff

κ − ω²

ω_m² + i ω ω_mQ_m

= ω_m² κ

1

ω²_m,eff− ω²+ iωΓm

. (4.6)

Here, χ(ω) is the linear response function of the mechanical oscillator to the noise. The nature of this τ_noise is a stochastic torque that represents Gaussian white noise forcing the harmonic oscillator. This means that for a fixed time, different realizations of the noise are drawn from the Gaussian distribution

g(τ ) = 1

√2πσ exp τ² 2σ²

, (4.7)

where the mean of the distribution is zero and the variance is σ². For a single realization of the noise, its time samples are also drawn from (4.7), and σ² =τ_noise² (t).

We refer to white noise as a random signal that has equal intensity at different frequencies, giving it a constant power spectral density (PSD), i.e. S_{τ τ}^white(ω) = c. This is purely a theoretical construction, as its variance is infinite:

σ² = 1 2π

Z _+∞

−∞

S_{τ τ}^white(ω) dω = ∞. (4.8)

To get around this, we consider a Gaussian force whose power spectral density is constant within some bandwidth B and then zero:

S_{τ τ}(ω) =

(α² if |ω| < 2πB/2

0 otherwise. (4.9)

From here we can derive the definition of noise variance:

σ² = 1 2π

Z +∞

−∞

S_{τ τ}(ω) dω = Z B/2

−B/2

α²df = α²B, (4.10)

where f = ω/2π.

(25)

The noise torque τ_noise induces a noise on the position of the mechanical oscillator whose PSD satisfies

S_φφ(ω) = |χ(ω)|²Sτ τ(ω) = α²

ω²_m κ

1

ω²_m,eff− ω²+ iωΓ_m

2

. (4.11)

Let us remark that this expression is equivalent to Equation (2.12), where we used variables x, k instead of φ, κ, and we approximated χ² by a Lorentzian. In the same way we did in Section 2.3, we can obtain the variance of the position as the area under the curve S_φφ:

φ² = 1 2π

Z +∞

−∞

S_φφ(ω) dω = 1 2π

Z +∞

−∞

|χ(ω)|²S_{τ τ}(ω) = α²ω_mQ_m

2κ² . (4.12) Note that this result does not depend on the bandwidth, as long as B ω_m/π.

We can now relate Equation (4.12) to the effective temperature T_effof the oscillator through the equipartition theorem:

φ² = kBT_eff

κ = α²ωmQm

2κ² =⇒ α² = 2k_BT_eff κ

ω_mQ_m. (4.13) If the oscillator is its ground state, thenφ² = φ²_ZPF and

φ²_ZPF= h0|φ²|0i = ~ωm

2κ

Eq. (4.13)

=⇒ α²_ZPF= 2k_BT_ZPF κ

ω_mQ_m = ~κ

Q_m. (4.14) Using all these relations, we can set an equilibrium temperature in our system and introduce the corresponding thermal fluctuations into our equations of motion to simulate noise. In this way, by tracking the value of the effective temperature (or equivalently, the value of φ²), we can observe cooling or heating even if the drive voltage is zero and the only force powering the optomechanical device is the noise.

4.3.1 Noise implementation

Let us recall from Section 3.2 that the implementation of noise into the simulations is not trivial. When using a variable-step integrator such as odeint, we can create unwanted correlations that generate artifacts in the results, so we have to find a way to get around it.

We begin by generating an array of random numbers from a Gaussian distribution, with mean zero and σ² = α²B. For a fixed-step integrator, setting this array equal to τnoise

would be enough to compute the equations of motion, but in our case, we need to inter- polate the array into a continuous function that can be evaluated at any time (see Figure 4.3). This is the function that will act as τnoise in Equation (3.12). In this way, odeint can evaluate the function even if it chooses smaller time steps than the original array.

(26)

Figure 4.3: Array of random noise points (blue) and their cubic interpolation (green), extracted from a Gaussian distribution of mean 0 and σ² = α²B, where α is calculated as in Equation (4.13) and B = 1.6 · 10⁸ Hz.

Note that the interpolation will create correlations in the noise such that it will not be white for higher frequencies than the inverse time step. Therefore, for this method to work, we have to generate the random array with a time step sufficiently small so that this threshold frequency is well above ωm.

4.3.2 Simulation results

Once we know how to add noise into our system, we are ready to simulate it, but it might be useful to make a few remarks first.

By increasing the drive voltage, we increase the amplitude of the charge oscillations, and therefore in general, it is satisfied that

κeff= κ − Q²C⁰⁰

2C² ≥ κ, (4.15)

since Q², C²≥ 0 and C⁰⁰< 0, which will induce a frequency shift towards higher frequencies:

ω_eff= ω_mr κ_eff

κ ≥ ω_m. (4.16)

This means that by increasing the drive voltage, κ_eff and ω_eff will be larger, which will, in turn, decrease the mechanical oscillations, and consequently reduce S_φφ. However, this phenomenon in itself does not necessarily imply cooling, as the effective temperature of the system might remain constant. In Figure 4.4, we show a decrease in the noise level caused by an increase in the effective stiffness of the oscillator. This increase in κ_eff reduces the zero-point fluctuations, which lower the noise level, but does not affect the thermal fluctuations.

For linear systems, when driving at the red sideband and increasing the drive power, the cooling of the mechanical oscillator shows as a shift in resonant frequency towards higher

(27)

Curve φ²

Teff

Blue (5.4 ± 0.9) · 10⁻¹⁷ rad² (9 ± 1) mK Purple (1.6 ± 0.3) · 10⁻¹⁷ rad² (10 ± 1) mK

Green (5 ± 1) · 10⁻¹⁸rad² (9 ± 1) mK

Figure 4.4: Reduction of the noise level due to effective stiffening in a system initially set at T = 10 mK. The plot shows the power spectral density calculated as S_φφ(ω) = |F T (φ)|² for three different drive voltages. The effective temperature is obtained from the variance of the oscillations according to the equipartition theorem (Teff=φ² κeff/kB).

values and a broadening and flattening of the power spectral density curve (noise squashing).

In the strong coupling regime, we can even observe normal mode splitting [39]. However, in our quadratically coupled system we have not yet managed to observe cooling of the mechanical oscillator. It is left for future work to derive a mathematical expression for Γm(ω, g) and show if it is indeed possible to cool down the oscillator by driving at the proper frequency with sufficient power.

4.4 Parametric amplification

In a parametric oscillator, there is a regime known as parametric resonance or parametric amplification in which the solution of the differential equation becomes unstable and grows exponentially. This behavior occurs at certain frequencies and can sometimes be unwanted, but it has interesting applications in coupled systems. For instance, in our device, parametric amplification allows driving the mechanical oscillator indirectly, only by providing a voltage to its coupled RLC circuit. The instability can also be used to amplify the response of a weak oscillating signal.

(28)

Let us study in detail the effect of the parametric resonance in our sensor, assuming first that there is no damping. From Equation (4.1), it follows the expression of κ_eff(t):

κ_eff(t) = κ − C⁰⁰Q¯²

2C²(0)sin²(ω_dt) = κ − C⁰⁰Q¯²

4C²(0)(1 − cos(2ω_dt)), (4.17) where ¯Q is the amplitude of oscillation in the charge, and ωd is the drive frequency. The- ory predicts that, in linear systems, parametric resonances occur near drive frequencies of ω_d= 2ω_m/n, where n ∈ Z [40–42]. However, in nonlinear systems, achieving parametric excitation is not always as trivial. For example, in our case, the nonlinearity is introduced in κ_eff(t) through the charge squared (Q²(t)), which adds a constant term that shifts us from resonance:

κ_eff= κ − C⁰⁰Q¯² 4C²(0)

| {z }

ξ = constant

+ C⁰⁰Q¯² 4C²(0)

| {z }

ξ

cos(2ω_dt), (4.18)

where ξ < 0 since C⁰⁰< 0. Now, by making a the change of variable τ = 2ω_dt, we can rewrite Equation (4.1) in the form of the well-known Mathieu’s equation ( ¨φ+(δ −ε cos τ )φ = 0) [43], used for describing parametric oscillators:

Id²φ

dt² + 0 ηdφ

dt + κ_eff(t)φ = 0 ^Eq.(4.18)=⇒

=⇒ Id²φ

dt² + (κ − ξ + ξ cos 2ω_dt)φ = 0 ^×I

−1

=⇒

=⇒ d²φ dt² + ω²_m

1 −ξ

κ+ ξ

κcos 2ω_dt

φ = 0 ^{τ = 2ω}=⇒^d^t

=⇒ d²φ dτ² +ω²_m

ω²_d

1 4 − ξ

4κ+ ξ 4κcos τ

φ = 0 =⇒

=⇒ φ + (δ − ε cos τ )φ = 0,¨ (4.19)

where, in our case,

ε = − C⁰⁰Q¯²ω_m²

16κC²(0)ω_d², δ = ω²_m

4ω_d² + C⁰⁰Q¯²ω_m²

16κC²(0)ω_d² = ω_m²

4ω²_d − ε. (4.20) Once we have established the equivalence between our equation of motion and Mathieu’s equation, it becomes much easier to find the parametric resonance regime, as Mathieu’s stability problem has a known analytical solution. In particular, we can draw a stability diagram in the δ-ε plane (see Figure. 4.5a) thanks to the expressions for the transition curves that can be calculated by harmonic balance [42].

In the diagram, we can clearly distinguish the unstable regions where the solution grows exponentially. It is very important to remark that in Figure 4.5a the instability areas reach values as low as ε = 0 because we assume that there is no damping. However, when the oscillations are limited by a finite quality factor, parametric resonance only occurs after a certain threshold has been surpassed, namely when the amplification is big enough to coun- teract the damping. In the diagram, this threshold is represented by a rounding of the sharp

(29)

(a) (b)

Figure 4.5: (a) Stability diagram in the δ-ε plane for the solution of Mathieu’s equation. No damping considered. (b) Effect of the damping η on the stability lines.

Figure 4.6: Stability diagram for different values of ωd and Vd. By setting ωd, we determine where the line δ = ω_m²/4ω²_d− ε stems from in the horizontal axis. In the plot, we differentiate two characteristic drive frequencies: ωo (pink dotted line) and ωm (green dotted line). On the other hand, by tuning Vdwe move up and down the line, which allows us to choose whether the system is stable or unstable. Depicted here are voltages Vd1 Vd2 and Vd3 Vd4.

minima of the instability regions, as depicted in Figure 4.5b.

As previously stated, our system does not behave like a linear parametric oscillator. In a simple parametric oscillator, the system moves in the δ-ε space along a vertical line: when driving at the parametric resonant frequencies, i.e. driving at the δ-points where we find the minima of the instability regions, the amplification is reached for all ε > 0. On the other hand, as shown in Figure 4.6, our system is bound to the line δ = ω_m²4ω²_d− ε, which means we can reach different regions of instability by tuning the parameters included in the definition of ε. In particular, the drive frequency (ω_d) determines where the line stems from in the δ-axis, while the amplitude of the external voltage (Vd) sets the vertical position in said line.

(30)

Let us show some results from the simulation illustrating the three different regimes: stability, boundary and instability. Firstly, in Figure 4.7, we observe how the system behaves in a stable region when approaching the boundary. As we can see, the angle oscillates randomly due to the thermal fluctuations in the order of 10⁻¹⁰ rad. Then, by increasing the drive voltage, we situate the system at the boundary. In Figure 4.8, we see how the oscillations increase exponentially until they are counteracted by the mechanical damping. Finally, by increasing the drive a bit more, we are fully in the unstable region, where the solution grows exponentially until it breaks down (see Figure 4.9). This occurs because the larger the oscillations are, the smaller the capacitance becomes, so when C(φ) → 0, the optomechanical coupling is lost and the system comes back to its equilibrium position.

(a)

(b)

Figure 4.7: (a) Stability diagram and (b) time evolution of the oscillations for a system in the stable regime, driven by thermal noise and Vd = 35V .

(a)

(b)

Figure 4.8: (a) Stability diagram and (b) time evolution of the oscillations for a system at the instability boundary, driven by thermal noise and Vd= 38V .

(31)

(a)

(b)

Figure 4.9: (a) Stability diagram and (b) time evolution of the oscillations for a system in the unstable regime, driven by thermal noise and Vd= 40V .

The drive voltages in Figures 4.7 - 4.9 may appear unreasonably high, but there are several ways around this issue. Firstly, these simulations are using the lowest value within the predicted range for C⁰⁰ (see Table 4.1). By using a higher value of C⁰⁰, we would be able to reach the instability region with much lower voltages. Secondly, we could drive the circuit at a lower frequency, closer to ω_mthan ω_o, which would make the line δ = ω_m²4ω_d²− ε stem closer to the instability zone, thus requiring less power for our system to become unstable.

Note that this last scenario is experimentally feasible, as a lower drive frequency will not be blocked by the series inductor.

(32)

5. Tip-surface interaction

Up until this point, we have studied how the sensor behaves in vacuum when it is isolated, however the main purpose of the device is not to operate alone, but to analyze a sample.

When we place a sample close to the mechanical oscillator, it interacts with the tip by changing the amplitude and the frequency of its motion. Ultimately, we want to be able to read these variations and decipher the information they contain about the material, but in order to achieve such deep understanding, we will begin with the inverse process. In this chapter, we will apply a known attractive/repulsive force between the tip and the surface and we will study its effect on the oscillations of angle, charge and current.

5.1 Tip-surface model

The total tip-surface force (F_ts) can be viewed as an additional force driving the system, and therefore can be modeled in our equation of motion (3.6) like

1 ω_m²

d²φ

dt² + 1 ω_mQ_m

dφ

dt + φ = 1

κ(τ_ext+ τ_noise+ τ_ts), (5.1) where τ_ts denotes the tip-surface torque. This force has conservatives and dissipative components and it depends on the position, the velocity and the past trajectory of the tip [44].

However, for simplicity, we will assume a conservative force only dependent on the instant position, as most models do.

The forces relevant to AFM are usually of electromagnetic origin, but different intermolecu- lar, surface and macroscopic effects give rise to interactions with distinctive distance depen- dencies. For this reason, added to the fact that the sample’s surface is not necessarily flat, it is important to properly define the position coordinates we are going to use. Due to the nature of the experimental setup, we typically define the coordinates in the laboratory frame of reference, to which the sample is rigidly attached, and we denote the equilibrium position of the oscillator by h(t). This height is continuously controlled with the piezo scanner to prevent the tip from crashing into the surface. By maintaining the notation φ(t) for the deviation from equilibrium, we can name the deflection of the tip d(t) and its instantaneous position with respect to the reference point z(t). In this way, it is satisfied that z = h + d.

Finally, if the sample’s surface is oriented in the uv-plane, then its height is denoted by z0(u, v). We display a scheme of the coordinate system in Figure 5.1.

For simplicity in the mathematical model, we will assume z₀ = 0 for all u and v, which is equivalent to bringing the reference point to the sample’s surface.

Theory and modeling of a quadratically coupled optomechanical sensor for atomic force microscopy