Two-level system as topological actuator for nanomechanical modes
C. Dutreix R. Avriller B. Lounis F. Pistolesi
cond-mat.mes-hall quant-ph




Université de Bordeaux and CNRS, LOMA, UMR 5798, F-33400 Talence, France


Université de Bordeaux and CNRS, LOMA, UMR 5798, F-33400 Talence, France


Institut d'Optique & CNRS, LP2N UMR 5298, F-33400 Talence, France
Université de Bordeaux, LP2N, F-33400 Talence, France


Université de Bordeaux and CNRS, LOMA, UMR 5798, F-33400 Talence, France


We investigate theoretically the dynamics of two quasidegenerate mechanical modes coupled through an open quantum two-level system.

A mean-field approach shows that by engineering the retarded response of the two-level system with a coherent drive, the non-Hermitian mechanical spectrum exhibits an exceptional degeneracy point where the two modes coalesce.

We show that this degeneracy can be exploited to manipulate the vectorial polarization of the mechanical oscillations.

We find that adiabatically varying the detuning and the intensity of the drive induces a rotation of the mechanical polarization, which enables the topological and chiral actuation of one mode from the other.

This topological manifestation of the degeneracy is further supported by quantum-jump Monte Carlo simulations to account for the strong quantum fluctuations due to the spontaneous emission of the two-level system. 

Our presentation focuses on a promising realization based on flexural modes of a carbon-nanotube cantilever 
coupled to a single-molecule electric dipole irradiated by a laser. 



Two-level system as topological actuator for nanomechanical modes
    F. Pistolesi
    December 30, 2023
=================================================================










§ INTRODUCTION


The manipulation and detection of nanometer oscillators are important challenges in nanomechanics <cit.>, and recent progress has led to unprecedented high-resolution sensors
<cit.>.

Most nano-oscillators exhibit a multimode dynamics
<cit.>.

In particular, the flexural dynamics of suspended nanowires involves nearly degenerate orthogonal modes, which enables the detection of anisotropic and nonconservative force fields <cit.>. 

Such advances in vectorial force microscopy rely on vectorial oscillations whose control is crucial to scan surfaces <cit.>.

Strategies to accurately monitor ultralight cantilevers, such as carbon nanotubes that allow zN/√(Hz) force sensitivity <cit.>, are then highly desirable to develop more sensitive vectorial probes.

One interesting possibility that we investigate here is to exploit the nonconservative force induced by the detection system.

Indeed, open systems can exhibit intriguing degeneracy points in the analytic continuation of their spectra, associated with the coalescence of the eigenstates and known as exceptional points (EPs) <cit.>.

EPs have recently allowed efficient topological energy transfers between two harmonic modes of a membrane placed in the middle of an optical cavity <cit.>.



Two-level systems (TLSs) constitute minimal quantum systems to detect and manipulate mechanical motions <cit.>. 

For example, their coupling to flexural modes allowes the localization of emitters randomly 
distributed in micropillars <cit.>. 

It was also shown that single-molecule TLSs are sensitive local probes to measure, through the Stark effect, the small 
displacements of a charged nanotube <cit.>.

Using a TLS to detect and actuate nanomechanical oscillators brings two fundamental differences 
with respect to the use of optical or electromagnetic cavities.

(i) A strongly pumped optical cavity has a linear behavior, whereas the quantum nature of a TLS is intrinsically nonlinear.

(ii) In highly populated optical cavities the Poissonian fluctuations are negligible, whereas the TLS experiences strong quantum fluctuations due to spontaneous emission.

Whether it is possible to observe EPs in the electromechanical spectrum of a cantilever coupled to a TLS naturally appears as a fundamental question.



In this paper, we show that adiabatically by varying along a closed path the frequency and the intensity of a coherent field driving a TLS coupled to two quasi-degenerate mechanical modes it is possible to induce a change in the state of the mechanical oscillator that depends on the topology of the path.

We show that this behavior is due to the presence of an EP in the mean-field description of the electromechanical spectrum.

Using quantum-jump Monte Carlo simulations, we prove that this property holds in the presence of strong TLS fluctuations.

The  presence of the EP allows one to generate elliptic mechanical eigenmodes, whose axis angles can be controlled. 

Finally, we propose a detection scheme to probe the topological switch between the two quasidegenerate flexural modes in single-molecule spectroscopy. 







§ MECHANICAL MODES COUPLED TO A DRIVEN TLS






 §.§ The system Hamiltonian


We consider the generic Hamiltonian of a TLS driven by a coherent field, linearly coupled to two nearly degenerate mechanical modes,


    ℋ =
        -ω_ TLS/2σ̂_z + Ω_ Lcos(ω_ Lt)σ̂_x + 
    ω_1 b^†_1b_1 + ω_2 b^†_2b_2
        - ∑_i=1,2  
    g_i (b_i^†+b_i) σ̂_z .


The operators σ_x and σ_z are Pauli matrices and describe the
TLS of energy splitting ω_ TLS (ħ=1).

The TLS is driven by a coherent field of frequency ω_ L and intensity (Rabi frequency) Ω_ L.

It couples with strength g_i, the mechanical modes
of frequencies ω_i, and destruction (creation) operators b_i (b^†_i)<cit.>.

The Hamiltonian ℋ describes the unitary evolution.

We further consider dissipation processes.

The driven TLS can spontaneously emit photons toward the electromagnetic environment with decay rate Γ.

The mechanical modes are coupled to thermal baths and have damping rates γ_i <cit.>.

Both baths are assumed to have the same temperature T_0.






 §.§ Electronanomechanical system


The generic model introduced above can describe various physical systems <cit.>.

We will focus our presentation on the system proposed in Ref.puller2013single  and shown schematically in 
Fig. <ref>.







In this case, the TLS is given by the electronic doublet in organic molecules embedded in a solid-state matrix <cit.>. 

As discussed in Ref. puller2013single, when a suspended carbon nanotube kept at a fixed difference of potential from
the substrate oscillates, it modulates the electric field on the electronic doublet, and by Stark effect it modulates the two-level system energy splitting. 

The coupling is proportional to the permanent electric dipole moment of the molecule, which can reach values up to two debyes <cit.>.

Performing molecular spectroscopy, then, allows one to detect the displacement of the oscillator.

In this paper we consider the presence of two  flexural modes that are quasidegenerate for symmetry reasons. 



Concerning the typical parameters one has that the carbon nanotube mass is m≃10^-20 kg, 
 the fundamental frequencies satisfy ω_i/(2π)≃1–10 MHz with 
 quality factors Q_i=γ_i/ω_i≃10^3–10^5 in the underdamped regime <cit.>.

In usual single-molecule experiments performed at liquid-helium temperatures, the TLS exhibits a lifetime limited dephasing rate Γ/2, with Γ/(2π)≃8–10MHz.

Realistic values of the coupling strengths g_i can be as large as g_c≃ 1GHz, which corresponds to a discharged electric field of about 10 mV/nm between the nanotube tip and the molecule substrate <cit.>.






§ MEAN-FIELD MECHANICAL DYNAMICS






 §.§ Langevin equation of motion


Now we aim to describe the non-Hermitian dynamics of the mechanical modes
that occurs when the TLS and the environment are traced out. 

In usual cavity optomechanics, the coupling between the mechanical oscillator 
and the cavity field can be linearized to solve the problem exactly. 

This is not possible with a TLS due to its intrinsic nonlinear nature. 

On the other side, in our case, we can exploit the timescale separation between 
the mechanical and the TLS dynamics, γ_i ≪Γ.

We then follows  Refs. clerk2004quantum,clerk2010introduction 
and derive a Langevin equation for the mechanical degree of freedom 
tracing out the TLS quantum degree of freedom.

This approach captures the Gaussian contributions and, in this sense, it is not limited to weak 
coupling. 

We obtain that the  expectation value of the displacement 
x_i= x^ zpf_i⟨ b_i+b^†_i⟩, where x^ zpf_i

denotes the zero-point fluctuations, satisfies the Langevin equation,


    ẍ_i(t) + γ_iẋ_i(t) + ω_i^2 x_i(t) 
    = 
    g_i/m x^ zpf_i⟨σ̂_z⟩_0 
    + δ F_i(t)/m
          
    + 1/m∑_j=1,2g_i g_j/x^ zpf_ix^ zpf_j∫dt' S_ R(t-t') x_j(t')  .


With ⟨…⟩_0, we indicate quantum averages evaluated in the absence of mechanical 
coupling (g_i=0).

The average force associated with ⟨σ̂_z⟩_0 only shifts the equilibrium position of the mechanical oscillator and we disregard it from now on.

The forces δ F_i denote the Brownian thermal fluctuations, as well as the nonequilibrium stochastic fluctuations due to the spontaneous emission.

The last term describes the TLS-mediated retarded coupling between the mechanical modes.

It involves the retarded response function of the TLS in the presence of the laser field and the electromagnetic environment:

S_ R(t) = 
-iθ(t)⟨ [ δσ̂_z(t), δσ̂_z(0) ] ⟩_0, 

where δσ̂_z(t)=σ̂_z(t)-⟨σ̂_z(t)⟩_0 characterizes the fluctuations of the population difference.







 §.§ Non-Hermitian mean-field dynamics


From ClassicalEQM,  we derive an effective non-Hermitian Hamiltonian that describes the oscillator's mode dynamics. 

We begin by neglecting the fluctuation forces δ F_i.

We then linearize the equation of motion (<ref>) for frequencies close to the 
two mechanical resonances ±ω_i (see Appendix <ref>).

We introduce the positive-frequency complex amplitude  X(t) from which one can obtain
the physical oscillator displacements, x=(x_1, x_2)=2  Re[ X].

The X quantities formally obey a Schrödinger equation

i Ẋ=H  X, 

where the effective Hamiltonian is non-Hermitian,

    H = 
    (
    [  ω_1-iγ_1/2-g_1^2 S_ R          -g_1 g_2 S_ R;          -g_1 g_2 S_ R ω_2-i γ_2/2-g_2^2 S_ R ]) .


Here, we assume that the mechanical frequency splitting is much smaller than the TLS linewidth, ω_1-ω_2≪Γ.

Thus, the coupling between the nearly degenerate mechanical modes depends on the retarded response 
at the mean mechanical frequency ω_0=(ω_1+ω_2)/2, and so S_ R≡ S_ R(ω_0).


To determine the effective coupling between the mechanical modes, 
we need to evaluate S_ R.

This can be done by using a standard Born-Markov approximation for the 
evolution of the TLS reduced density matrix ρ obtained by tracing out the electromagnetic environment.

Defining ρ=(ρ_11,ρ_12,ρ_21,ρ_22), one obtains that it satisfies
the Liouville–von Neumann equation ρ̇=ℒρ associated with the superoperator


    ℒ = 
    ( [        0  iΩ_ L/2 -iΩ_ L/2        Γ;  iΩ_ L/2  -iδ-Γ/2        0 -iΩ_ L/2; -iΩ_ L/2        0   iδ-Γ/2  iΩ_ L/2;        0 -iΩ_ L/2  iΩ_ L/2       -Γ;          ])


in the rotating wave approximation, where δ=ω_ L-ω_ TLS defines the laser detuning.

This corresponds to the optical Bloch equations that are known to provide a realistic description 
of electric dipoles in single-molecule experiments <cit.>.

In the stationary regime ℒρ_0=0, the quantum regression theorem leads to the autocorrelations


    S(t)
           = ⟨σ_z(t) σ_z(0) ⟩
        - ⟨σ_z ⟩^2 
       =⟨ w_0 | ℳ_ze^ℒtℳ_z | ρ_0⟩ - ⟨ w_0 | ℳ_z |ρ_0⟩^2 ,


where w_0=(1,0,0,1) denotes the kernel left-hand eigenvector of ℒ, and ℳ_z=(1,1,-1,-1) 
(see, also, Ref. pistolesi2018bistability).

The power spectral density S(ω)=∫ dt  e^iω t  S(t) characterizes the absorption (ω>0) and emission
 (ω<0) of the TLS.

The frequency asymmetry of the quantum noise relates to the imaginary part of the response function through 2 S_ R(ω)=S(-ω)-S(ω).

We then obtain


    S_ R(ω) = 8δΩ_ L^2/Γ^2+(2δ)^2+2Ω_ L^2iω-Γ/P(-iω) ,


where the polynomial P(z)=∑_n=0^3a_n z^n has real coefficients

a_0 = Γ(Γ^2 + (2δ)^2 + 2Ω_ L^2)/4, 
a_1 = (5Γ^2 + (2δ)^2 + 4Ω_ L^2)/4,
a_2 = 2Γ, and
a_3 = 1.


Figure <ref> shows the coupling strength S_ R≡ S_ R(ω_0) as a function of the laser detuning and the 
Rabi frequency.





The coupling vanishes when the coherent drive is in resonance (δ=0) or strongly detuned (|δ| ≫Γ) with respect to the splitting of the TLS.

It also vanishes when the coherent drive is turned off (Ω_ L=0) or when it is sufficiently strong to 
saturate the populations of the TLS (Ω_ L≫Γ). 

In these situations, the two eigenstates are oscillations along two orthogonal directions.

Otherwise, the TLS mediates an effective coupling that modifies the non-Hermitian dynamics of the mechanical modes 
leading to complex eigenvalues λ_± and eigenstates X_± of H.

In the next section, we show how degeneracies in the complex spectrum of H are associated with singular properties of the eigenstates 
that affect the polarization of the mechanical modes in real space.








§ EXCEPTIONAL DEGENERACY POINTS






 §.§ Electromechanical spectrum


For nonvanishing coupling (g_1 g_2 S_ R≠ 0), the electromechanical spectrum may exhibit exceptional 
degeneracy points in the parameter space (δ,Ω_ L), which describes the laser driving.

They occur when the retarded TLS response satisfies


    S_ R=S_ EP±≡ - Δω/G_± ,


where Δω=ω_1 - ω_2 - i(γ_1-γ_2)/2 characterizes the splitting of the mechanical modes and G_± = -(g_1^2-g_2^2± i 2 g_1g_2) relates to the coupling strengths between the mechanical oscillators and the TLS. 

We focus in particular on the degeneracy point S_ EP≡ S_ EP-, which lies in the positive imaginary plane.

It is then convenient to chose this EP as a new origin of the complex plane, such that S_ R=S_ EP+z.

The mechanical dynamics is described equivalently in terms of an effective Hamiltonian H' similar to H 
(see, also, Appendix <ref>).

We find


    H∼ H' = 1/2(
        [    h_0      0; G_+Δ S    h_0 ])
        +
        z/2(
        [   0 G_-; G_+   0 ]) ,


where Δ S≡ S_ EP- -S_ EP + depends on the distance between the two EPs and h_0= ω_1+ω_2 -i (γ_1+γ_2)/2-(g_1^2+g_2^2)  S_ R.

This representation explicitly shows that the effective Hamiltonian H supports a Jordan matrix representation at the EP in z=0.


The electromechanical spectrum presents the two eigenvalues


    λ_± = h_0/2±1/2√(G_+G_-(Δ S+z)z) .


Therefore, the EP in z=0 also corresponds to the branch point of the complex square root √(z).

The resonance frequencies Ω_±= [λ_±] and the damping rates Γ_±=- [λ_±] then support a Riemann-surface representation in the vicinity of the EP [Fig. <ref>(a)].


This can be evidenced by varying the detuning δ and Rabi frequency Ω_ L of the TLS drive.











Figure <ref>(b) shows that strongly detuning the TLS drive from the resonance (δ=0) outlines a loop in parameter space.

The response function goes away from and back to S_ R=0, where the mechanical coupling mediated by the TLS vanishes (Fig. <ref>).

This point corresponds to the bare mechanical frequencies and damping rates, which are Ω_±=ω_1,2 and Γ_±=γ_1,2.

One can check in Fig. <ref>(c) that the mechanical frequencies and damping rates go back to their bare initial values when the loop does not enclose the EP in parameter space.

For EP-enclosing loops, however, the mechanical frequencies Ω_± and damping rates Γ_± do not come back on their initial values, 
but are exchanged.

This eigenvalue switch is an evidence of the Riemann-surface topology in the electromechanical spectrum.



For the case of an optical cavity coupled to two eigenmodes of a membrane, this effect has been observed in Ref. xu2016topological.

Here we showed that one can have a similar behavior for a TLS.

Specifically, for the electromechanical system that we propose in Fig. <ref>, the condition of existence for the EP in the positive imaginary plane reads S_ EP≃ iΔω/(2g_i^2), assuming g_1≃ g_2.

Since S_ EP is nearly pure imaginary in the underdamped regime (γ_i≪ω_i), varying the TLS response S_ R around this degeneracy point requires Re[S_ R] to change signs and  |S_ R|≃|S_ EP|.

These two conditions give ω_i>Γ/2 and ΓΔω≃2g_i^2, for the maximum of S_ R is obtained when δ∼Ω_ L∼Γ (Fig. <ref>).

For the TLS of electric dipoles recently observed <cit.>, experiencing the EP then implies a typical coupling g_i/(2π)≈ 0.3 MHz.

It is much smaller than the critical coupling g_c≃10^3 MHz and could be realized, for instance, by positioning the tip of a carbon nanotube 100 nm away from the molecule with a 100 μV bias.

We emphasize that such a coupling is not strong enough to excite higher flexural modes.

According to Euler-Bernoulli beam theory <cit.>, their frequencies are at least six times larger than the fundamental ones.












 §.§ Eigenmode polarization around the EP


The detailed analysis of the eigenvectors unveils a very interesting dynamics of the oscillator tip.



To investigate them, we choose the biorthogonal left and right eigenstates of the effective Hamiltonian H' as 


    Y_± = 1/√(2)( [ ± Z^1/4 ;    Z^-1/4;           ])
        and    X_± = 1/√(2)( [ ± Z^-1/4;    Z^1/4;          ]) ,


where we use the polar representation z=ρ e^iθ around the EP, so that

    Z(θ)=G_+/G_-(Δ S/ρ e^-iθ + 1) .


The exceptional degeneracy point is associated with a phase singularity for the eigenstates, which can also be evidenced by varying S_ R smoothly around the degeneracy.

For a ρ-radius loop that encloses only S_ EP (ρ<|Δ S|), the eigenstates fulfill the condition of parallel transport Y_±(θ)·∇_θ X_±(θ)=0 and are multivalued.

They change as X_±(θ+2nπ) = (-i)^n X_±(θ) if n is even, and X_±(θ+2nπ) = (-i)^n X_∓(θ) if n is odd.

Therefore, additionally to the eigenvalues, the eigenstates also switch after one loop around the EP.



The multivaluation of the eigenstates further affects the polarization of the mechanical oscillations x_±(t) in real space.

The eigenstates are solutions of the effective Schrödinger equation, that is, X_±(t)=x_±(0) e^-iλ_±t X_±, where x_±(0) are initial amplitudes that we assume to be real.

The mechanical dynamics then consists of damped oscillations with elliptical polarization (see, for more details, Appendix <ref>).

We find


    x_±(t) = R(φ)  x_±(0) e^-Γ_±t( [ α cos(Ω_±t); β sin(Ω_±t);             ]) .


The semiaxes α and β of the ellipse determine the eccentricity of the mechanical oscillations, and R(φ)= τ_0 cosφ - iτ_2 sinφ is a rotation matrix where the Pauli matrices τ_i are written in the orthogonal basis (e_1,e_2) of the uncoupled modes.

The rotation angle φ thus characterizes the orientation of the polarization with respect to one of the modes in the absence of coupling (g_1g_2S_R=0).



Figure <ref>(d) presents the dynamics of the mechanical oscillations when the retarded response of the driven TLS performs a loop around the EP in parameter space.

The semiaxes of the elliptical polarization remain unchanged after the loop.

We find that this property generally holds for any loop since the semiaxes vary as α(θ+2π)=α(θ) and β(θ+2π)=β(θ) (Appendix <ref>).

Nevertheless, the mechanical oscillations do not go back to the initial polarization at the end of the loop, for the polarization undergoes a rotation of φ=π/2.

We can more generally show that the polarization rotates as φ(θ-2nπ)=φ(θ)+nπ/2 (Appendix <ref>).

This fourfold invariance is reminiscent of the fourth-root multivaluation of the eigenstates in Eq. (<ref>).

We emphasize that the π/2 rotation of the polarization also comes with the switch of the resonance frequencies and damping rates.

Thus, an eigenmode initially activated as 
x(t) = x(0)e^-γ_2tcos(ω_2t)  e_2 is transferred into the eigenmode of orthogonal polarization, that is,
x(t) ∝ x(0)e^-γ_1tcos(ω_1t+Δϕ)  e_1, 
where Δϕ is the phase accumulated along the EP-enclosing loop. 

This transfer of energy from one eigenmode to the other only depends on whether S_ R encircles S_ EP, regardless of the precise loop geometry.

The actuation between mechanical modes is therefore topological.






 §.§ Chiral nature of the polarization




The topological actuation from one mode to the other is an intrinsic property of the instantaneous eigenstates X_±.

To observe the effects of their multivaluation around the EP, we further investigate their adiabatic transport.

In Hermitian systems, the adiabatic theorem ensures that one can neglect the nonadiabatic transitions over some typical timescale T ≫ 1/|λ_+-λ_-|.

In open systems, however, this is no longer true for all the eigenstates <cit.>.

For a two-state system, in particular, only the least dissipative state is expected to be transported adiabatically around the EP <cit.>.

Here we study this issue by solving numerically the Schrödinger-like equation iẊ=H  X.



We focus on the EP-enclosing loop associated with the Rabi frequency Ω_ L=0.3Γ in Fig. <ref>(b).

We ramp the detuning linearly over the time scale T=10^5/Γ between δ=0 and δ=-8Γ.

Then, | X_±(t)|^2=| Y_±(t)𝒰(t,0) X(0)|^2 provides the exact dynamics, where 𝒰 denotes the time-evolution operator for the Hamiltonian H.

We assume the initial eigenstates are either X(0)= X_+= X_1 or X(0)= X_-= X_2 associated with the frequencies
ω_1=Γ and ω_2=ω_1-10^-3Γ.

We then compare the exact dynamics to their adiabatic evolutions |Ad_±(t)|^2=|exp[-i∫_0^t λ_±(τ) dτ] Y_±(0)  X(0)|^2.

The results are presented in Fig. <ref> for two orientations of the EP-enclosing loop in parameter space.





For a clockwise loop [Fig. <ref>(a)], the eigenstate X_- is the least dissipative one [Fig. <ref>(b)]


    ∫_0^Tdt [Γ_+(t)-Γ_-(t) ]>0  .


We find that only this state can experience the adiabatic transport [Figs. <ref>(c) and (d)].

When reversing the orientation of the loop, the situation is reversed too.

The eigenstate X_+ becomes the least dissipative one and follows the adiabatic evolution, whereas X_- does not (Figs. <ref>e-h).

Thus, both modes x_± can experience the topological actuation, but for opposite orientations of the loop.

This asymmetry with respect to the orientation of the loop reveals the chiral nature of the π/2 rotation of the eigenmode polarization around the EP.












 §.§ Eigenmode polarization at the EP


We can also investigate the fate of the mechanical oscillations when approaching the EP.

Figure <ref>(e) illustrates this situation.

We find that the two mechanical modes coalesce into a single mode of circular polarization.

The coalescence results from the eigenstates that become collinear in the limit of small ρ in Eq. (<ref>).

The circular polarization is reminiscent of the Jordan matrix representation of H and is then a manifestation of the EP.

Nevertheless, it is a local property in parameter space for which noise may be detrimental, especially because noise is known to be enhanced by the extreme nonorthogonality of the eigenstates close to the EP <cit.>.

Another drawback is that the EP is a degeneracy point where, by definition, the band gap closes and around which nonadiabatic transitions then become unavoidable <cit.>.

In contrast, the π/2 rotation of the polarization that we have introduced above evidences a global (topological) property of the EP for which the adiabatic transport is possible.

It does not require to approach the EP and, therefore, should be a more robust manifestation.







§ QUANTUM NOISE AND BACK-ACTION






 §.§ TLS quantum fluctuations


The mean-field picture discussed above is appealing, for it introduces a topological property of the EP to manipulate the vectorial polarization of mechanical oscillations.

Now, we show that this manifestation of the EP can also survive the fluctuations of the open quantum system that dresses the mechanical modes. 

In optical cavities, the number of circulating photons n_c obeys a Poisson distribution, and so ⟨δ n_c^2⟩=⟨ n_c⟩.

The fluctuations of the radiation pressure force F then satisfy ⟨δ F^2⟩/⟨ F ⟩^2=1/⟨ n_c⟩, which becomes negligible for the usual cavities that are strongly populated.

For the TLS, however, the force scales as F∝n̂_g-n̂_e, where the populations of the ground and excited states satisfy the constraints n̂_g+ n̂_e=1 and n̂_g - n̂_e=σ_z.




















The force fluctuations verify


    ⟨δ F^2⟩/⟨ F ⟩^2=⟨n̂_e⟩(1-⟨n̂_e⟩)/1/4-⟨n̂_e⟩(1-⟨n̂_e⟩) ,


where we used ⟨n̂_e^2⟩=⟨n̂_e⟩.

Since ⟨n̂_e⟩≤1/2, the fluctuations can be arbitrarily large compared to the mean force, and eventually ⟨δ F^2⟩/⟨ F ⟩^2→∞.



To describe the TLS force fluctuations due to spontaneous emission and test the mean-field description based on Eq. (<ref>), we perform 
quantum-jump Monte Carlo simulations
<cit.>.

This implies solving, on short timescales Δ t≪ 1/Γ, the differential equations

    i ȧ_ g   =   
    δ[ x(t)] a_ g/2 + Ω_ La_ e/2 
    
    i ȧ_ e   = (δ[ x(t)]/2-iΓ) a_ e + Ω_ La_ e/2 
    
    mẍ_i    = -mω_i^2 x_i-mγ_iẋ_i + 
    δ F_ th +g_i/x_i^ zpf⟨σ_z⟩ ,
 
where we introduce δ[ x(t)]=δ-∑_i g_i x_i(t)/x_i^ zpf, the Brownian thermal force δ F_ th with variance √(2mγ_i k_B T_i),
and the wavefunction of the TLS, |ψ(t)⟩=a_ g(t)| g⟩+a_ e(t)| e⟩, so that
⟨σ_z(t)⟩=|a_ g(t)|^2-|a_ e(t)|^2.

At each time step, one randomly either allows a transition to the ground state with probability 
ΓΔ t |a_e|^2, or normalizes the wavefunction and proceeds with the time evolution. 












 §.§ Mechanical cooling via the TLS noise


The mean-field description introduced above involves the response function of the bare TLS, that is, in the absence of electromechanical coupling (g_i=0).

This neglects the effects of the detuning shift δ[ x(t)] induced by the mechanical displacement.

In particular, large oscillations may lead to important effects when g_i x_i/x_i^ zpf∼Γ, since the oscillations can effectively change the laser-TLS detuning.


It is then interesting to work at low temperature. 

We investigate here how the system can be cooled down by coupling to the TLS <cit.>.


The power spectrum S(ω) of the TLS noise can be obtained from Eq. (<ref>). We denote its symmetric- and asymmetric-in-frequency parts, S_±(ω)=S(ω)± S(-ω).

The optical damping induced by the TLS noise relates to the asymmetry between emission and absorption as γ_ TLS=g_0^2 S_-(ω_0), where we assume g_1,2= g_0.

Cooling the mechanical oscillator via the TLS then requires the optical damping to be stronger than the intrinsic damping of the thermal baths, γ_0≪γ_ TLS, where γ_1,2 = γ_0.

The natural frequency scale of the noise must relate to the spontaneous emission rate Γ, so one can expect S_-=-2 S_ R∼ 1/Γ.

One can check that this is indeed the case for δ=Ω_L=Γ=ω_0 in Eq. (<ref>).

As the coupling cannot be larger than the critical one of the discharge electric field between the TLS and the carbon nanotube, 
we find that mechanical cooling is possible when


    ω_0γ_0 ≪ g_0^2 < g_c^2  .


For the electromechanical system in Fig. <ref>, one has ω_0γ_0 that is at least three orders of magnitude smaller than the critical coupling g_c, which leads to a large range of possible coupling strength g_0.

Besides, the TLS temperature verifies 2k_ B T_ TLS/(ħω_0) = S_+/S_-.

This leads to the effective temperature T_ eff for the mechanical modes dressed by the TLS,


    2k_ B T_ eff/ħω_0
        = γ_02k_ B T_0/ħω_0+γ_ TLS2k_ B T_ TLS/ħω_0/γ_0+γ_ TLS .


Figure <ref> (a) represents the map of the effective temperature as a function of TLS drive parameters δ and Ω_ L.

It shows that the mechanical system can be  cooled below ħΓ (≃0.5 mK).



We further verify this possibility of mechanical cooling by means of quantum-jump Monte Carlo simulations.

Figure <ref>(b) presents the evolution of the mechanical energies averaged over 1000 dynamics for δ=-Γ and Ω_ L=Γ, which corresponds to the region of maximum cooling predicted in Eq. (<ref>) and Fig. <ref>(a).

We study the mechanical energy through the dimensionless parameter ℰ_i defined as

    ℰ_i=4ω_i/Γ1/ħΓ(p_i^2/2m+1/2m ω_i^2 x_i^2)  .


Each mechanical mode is initially assumed to be in thermal equilibrium with a bath of typical dilution-fridge temperature T_0=10 mK.

Thus, we begin with randomly generating the position and momentum of each mode according to a Boltzmann distribution.

The equipartition theorem implies that the initial mean mechanical energy in the figures corresponds to ℰ_i=4ω_n/Γk_ BT_n/ħΓ=80.

We then simulate the stochastic dynamics based on Eq. (<ref>) with the usual parameters 
ω_1=1Γ, ω_2=ω_1-10^-3Γ, γ_i=10^-5Γ, and g_i=0.1Γ.

Figure <ref>(c) shows from the mean dynamics that the system can be cooled down ħΓ (∼0.5 mK), to about 50 μK, in agreement the mechanical cooling predicted in Fig. <ref>(a). 






 §.§ Mean rotation of the mechanical polarization


We come to the main question we want to address with the simulations: does the EP remain stable in the presence of 
fluctuations and backaction?

To test the mean-field manifestation of the EP, we study the evolution of the mechanical mode 1 or 2 with initial amplitude 10x^ zpf and ramp adiabatically δ from 0 to -8Γ over the timescale T=10^4/Γ for various values of Ω_ L.

Such an initial actuation can be achieved by forcing the charged nanotube tip with a transverse electric field oscillating at a frequency tuned on the selected mode frequency.

The mean dynamics of the oscillator, averaged over 10^2 simulated trajectories, consists of quasiperiodic elliptical oscillations in real space.

We thus identify the semiaxes in each quasiperiod and determine their rotation angle φ at a given value of δ and Ω_ L.



Figure <ref>(a) presents the prediction of the mean-field description, clearly showing the presence of the EP at δ≃-0.85Γ.

In comparison, Fig. <ref>(b) shows the rotation angle for the same parameters, but obtained from the Monte Carlo
simulations. 

Though fluctuations seem to blur the EP, where the noise is enhanced by the coalescence of the nonorthogonal eigenstates <cit.>, the axis of the oscillations performs on average a rotation in agreement with the mean-field description.

Starting with an oscillation in the horizontal direction (φ=0), one ends with a perpendicular oscillation (φ=π/2), showing that the adiabatic picture remains valid.

If the rotation implies a transfer of energy from mode 1 to mode 2, continuously following the evolution of the oscillation axis proves that this transfer is purely due to the adiabatic evolution predicted by the mean-field approach, and not to dissipation or stochastic effects.

Figure <ref>(c) shows the same evolution from mode 1.

No energy transfer is realized, as expected from the mean-field adiabatic transport.

We would like to stress that the numerical calculation fully takes into account the nonlinearity and non-Gaussian behavior of the TLS, 
thus confirming the, at least approximate, validity of the mean-field approach.

In the next section, we finally discuss a possible experimental implementation for the direct detection of the energy transfer. 




§ DETECTION VIA FREQUENCY MODULATIONS





After the adiabatic loop in parameter space, the mechanical energy of one mode is transferred to 
the other one.

This implies that the mechanical oscillation frequency has changed, but the variation is extremely small.  

Here we show that detecting such a small change can be performed by modulation of the frequency of the coherent field that drives the TLS.


At the semiclassical level, the Hermitian Hamiltonian of the driven TLS  can be written as


    H_ TLS(t) = 
    ( [       -ω_ TLS/2 -g_i x_i(t) Ω_ L( e^iF(t) + e^-iF(t))/2; Ω_ L( e^iF(t) + e^-iF(t))/2       ω_ TLS/2 + g_i x_i(t);                             ])  ,


where F(t)= ∫_0^tdτ  [ω_ L + f(τ)] and f(t) = a_λcos(ω_λt+ϕ_λ) describes the frequency modulation around frequency ω_ L. 

Let ψ(t) be a wave function satisfying the time-dependent Shrödinger equation i ψ̇(t) = H_ TLS(t)  ψ(t).

We can then introduce the gauge transformation ϕ(t) = U^†(t) ψ(t) based on the unitary operator


    U = 
    ( [  e^iF(t)/2          0;          0 e^-iF(t)/2;            ]) .


This leads to i ϕ̇(t) = [ U^†(t)H_TLS(t)U(t) - i U^†(t)U̇(t) ] ϕ(t), where the effective Hamiltonian is


    H̃_ TLS(t) = 
    ( [   δ(t)/2 -g_i x_i(t)               Ω_ L/2;               Ω_ L/2 -δ(t)/2 + g_i x_i(t);                      ]) ,


and δ(t)=δ-(β_λω_λ)cos(ω_λ t)-g_i x_i(t).

Thus, the frequency modulation induces an additional detuning shift, which adds up to the one of the mechanical oscillator in the modulated rotating frame.

The frequency modulation and the mechanical oscillations may then lead to interference when ω_λ=ω_i, since only the mode i is exited at the end of the loop.

We can show, in particular, that the interference affects the excited-state population of the driven TLS, and so its
luminescence excitation spectrum (Appendix <ref>).

This is illustrated in Fig. <ref>(d).

The Lorentzian background of width Γ and centered on ω_λ=0 already exists in the absence of mechanical mode (g_i=0) and, thus, is not due to any interference.

The interference, however, appears through a narrow Lorentzian peak of width γ_i centered on the mechanical frequency ω_i in the figure.

The width of the peak verifies γ_i≪ω_1-ω_2, so that the two quasidegenerate modes could be resolved clearly in the experiments.






§ CONCLUSION


Manipulating mechanical systems at the nanometer scale is an important challenge of present research.

The possibility of using EP in the excitation spectrum to transfer energy from one mechanical mode to the other had been 
proposed and observed in the past for mechanical modes coupled to optical cavities. 

In this paper, we have shown that EP can be equally generated by coupling mechanical modes to TLSs. 

Specifically, we considered a concrete example of single molecules coupled to flexural modes of carbon nanotubes, for which we performed detailed simulations.

We have shown quite generally that the topological and chiral energy transfer is possible.

Remarkably, the prediction of the analytical mean-field theory is confirmed by the quantum-jump Monte Carlo approach.

This guarantees that even if a TLS is quite different from an optical cavity, since it is a strongly nonlinear quantum 
system and has strong quantum fluctuations, it can be used to manipulate a mechanical oscillator exploiting the EP. 


From a conceptual point of view, the flexural mode eigenstates allow one to understand, 
in a transparent way, the mechanism of formation of the EP. 

We find that the evolution from the standard orthogonal eigenvectors to the coincident eigenvectors at the EP is 
performed by evolving the eigenvectors into elliptical oscillations, which eventually become circular at the EP.

The energy transfer is then simply obtained by a rotation of the axis of the elliptic oscillation of π/2 when a loop is performed around the EP.

These findings can allow a manipulation of the tip of the nanotube without the addition of any external electric fields. 

This possibility of quantum manipulation could, for instance, find applications in vectorial force microscopy, where monitoring the eigenmode polarization is crucial to scan a surface.




We gratefully acknowledge Conseil Régional de Nouvelle-Aquitaine for financial support (Grant No. 2016-1R60306-00007470) and Idex Bordeaux (Maesim Risky project 2019  of the LAPHIA Program). B.L. acknowledges the Institut Universitaire de France.










§ DERIVATION OF THE NON-HERMITIAN HAMILTONIAN



The expectation value of the displacement operator obeys the Langevin equation (<ref>). 

We begin with neglecting the fluctuating forces δ F_i and ⟨σ̂_z⟩_0, which only shifts the equilibrium position of the oscillator.

In Fourier space, this leads to


    m(-ω^2 - iγ_iω +ω_i^2) x_i(ω)    = ∑_j=1,2g_i/x_i^ zpfg_j/x_j^ zpf S_ R(ω) x_j(ω)  .


Since we consider the underdamped regime (γ_i≪ω_i), the bare mechanical poles verify

ω_i±≃±ω_i - i γ_i/2.

We can linearize Eq. (<ref>) in the vicinity of the two mechanical frequencies ω_i±.

The displacement can then be written as x_i ≃ x_i++x_i-, where x_i+ (x_i-) describes the mechanical oscillations of positive (negative) frequency ω_i+ (ω_i-).

The positive- and negative-frequency oscillations satisfy:


    {[   +i2mω_i(-iω) x_i+(ω) = 2mω_i(ω_i -i γ_i/2) x_i+(ω) - ∑_j=1,2g_i/x_i^ zg_j/x_j^ z S_ R x_j+(ω); -i2mω_i(-iω) x_i-(ω) = 2mω_i(ω_i +i γ_i/2) x_i-(ω) - ∑_j=1,2g_i/x_i^ zg_j/x_j^ z S_ R^* x_j-(ω);                                                                                                 ]. ,


where S_ R≡ S_ R(ω_0) and ω_0=(ω_1+ω_2)/2 is the mean mechanical frequency.

Since x_i≃ x_i++x_i-=2 Re[x_i+], we only focus on the positive-frequency solution x_i+.

Fourier transforming back to the time domain finally results in


    -iẋ_i+(t)    = (ω_i -i γ_i/2) x_i+(t) + ∑_j=1,2 g_i g_j S_ R x_j+(t)  ,


where we have used 2mω_i≃ x_i^ zpfx_j^ zpf.



Therefore, the positive-frequency oscillations of the mechanical modes are described by the vector X=(x_1,+,x_2,+).

It obeys the Shrödinger-like equation i Ẋ = H  X with the non-Hermitian effective Hamiltonian,


    H    =
    ( [ ω_1-i γ_1/2 - g_1^2 S_ R           - g_1 g_2 S_ R;           - g_1 g_2 S_ R ω_2-i γ_2/2 - g_2^2 S_ R;                          ]) .


This relates to the mechanical displacement vector x=( x_1, x_2) in real space as x=2 [  X ],










§ EXCEPTIONAL POINT PROPERTIES






 §.§ Branch point in the spectrum


We are now interested in the consequences of degeneracies in the non-Hermitian spectrum of H on the dynamics of the mechanical modes.

To make them more explicit, we perform two subsequent π/2 rotations of axes y and z with respect to the Bloch sphere of H eigenstates.

Thus, we introduce the unitary operator


    U =   
    e^-iπ/4τ_y e^-iπ/4τ_z
    =
    1/2( [  1-i -1-i;  1-i  1+i;      ])


and perform the transformation H' = U H U^†.

This results in H' =∑_i h_iτ_i/2 with


    {[                                          h_0 = ω_1+ω_2 -i (γ_1+γ_2)/2-(g_1^2+g_2^2)  S_ R;                                          h_1 =ω_1-ω_2 -i (γ_1-γ_2)/2 -(g_1^2-g_2^2)  S_ R;                                          h_2                             =-2 g_1g_2  S_ R;                                          h_3                                           =0 ].
     .


Degeneracies in the electromechanical spectrum then occur when S_ R=S_ EP±, where


    S_ EP±≡ - Δω/G_± .


Here, Δω=ω_1 - ω_2 - i(γ_1-γ_2)/2 characterizes the splitting of the mechanical modes and G_± = -(g_1^2-g_2^2± i 2 g_1g_2) relates to the coupling strengths between the mechanical oscillator and the TLS. 

We focus, in particular, on the degeneracy point S_ EP≡ S_ EP-, which lies in the positive imaginary plane.

It is then convenient to set it as the new origin of the complex plane by introducing the variable z=S_ R-S_ EP.

Thus, we find that the effective Hamiltonian reads


    H' = 1/2(
        [    h_0      0; G_+Δ S    h_0 ])
        +
        z/2(
        [   0 G_-; G_+   0 ]) ,


where Δ S≡ S_ EP- -S_ EP + depends on the distance between the two EPs.

This explicitly shows that the effective Hamiltonian H supports a Jordan matrix representation at the EP in z=0.

The electromechanical spectrum relies on the two eigenvalues


    λ_± = h_0/2±1/2√(G_+G_-(Δ S+z)z) .


Therefore, the EP in z=0 also corresponds to the branch point of the complex square root √(z), and hence the Riemann surface of the electromechanical spectrum in Fig. <ref>(a).

The same features hold when focusing on the degeneracy point S_ EP+ in the negative imaginary plane.





 §.§ Eigenstate multivaluation


The squareroot behavior of the spectrum near the EP also leads to singular properties for the eigenstates.

To investigate them, we choose the biorthogonal left and right instantaneous eigenstates of the effective Hamiltonian H' as 


    Y_±(θ) = 1/√(2)( [ ± Z^1/4 ;    Z^-1/4;           ])
        and    X_±(θ) = 1/√(2)( [ ± Z^-1/4;    Z^1/4;          ])  .


Using the polar representation z=ρ e^iθ around the EP, we find


    Z(θ)=G_+/G_-(Δ S/ρ e^-iθ + 1)  .


Let us assume that S_ R can be varied smoothly along a ρ-radius loop that encloses only S_ EP (ρ<|Δ S|).

Then, we find

    Z^1/4(θ=2π)=-iZ^1/4(θ=0)  ,

so that the eigenstates, which meet the condition of parallel transport Y_±·∇_θ X_±=0, change as follows:


    {[       X_±(2π) = -i   X_∓(0);       X_±(4π) = -1   X_±(0);       X_±(6π) = +i   X_∓(0);       X_±(8π) = +1   X_±(0) ] .  .


Thus, we expect the two eigenstates to swap after one loop, to pick up a geometrical Berry phase π after two loops, and then come the eigenstates come back onto the initial states without any geometrical phase only after four loops.

The eigenstate swap comes from the multivaluation of the eigenstates around the EP.

The fourfold multivaluation results from the complex fourth root in the parallel-transported eigenstates in Eq. (<ref>).



Note that if the loop encircles the two exceptional points (ρ>|Δ S|), one can choose single-valued eigenstates along the loop, so that X_±(2π n) =  X_±(0) for any integer n.

The eigenstates neither pick a geometrical phase nor swap in this case.






 §.§ Eigenmode polarization


The instantaneous eigenstates obey the Shrödinger-like equation iẊ(t) _±=H' X_±(t), and so X_±(t) = x_±(0)   e^-iλ_±t X_±. 

In real space, the mechanical eigenmodes evolve in time as


    x_±(t)    = 2  Re[U^† X_±(t) ] 
       = x_±(0)  Re( [   e^-iλ_± te^iπ/4(Z^1/4± Z^-1/4); -ie^-iλ_± te^iπ/4(Z^1/4∓ Z^-1/4);                                  ]) .


We can ignore the overall phase shift of π/4 and rewrite it as


    x_±(t) 
           = x_±(0)   e^-Γ_± t( [ a_±cos(Ω_± t+ϕ_±); a_∓sin(Ω_± t+ϕ_∓);                   ]) .


The amplitudes and the phases of the oscillations depend on the polar angle θ around the EP through

    a_±   = |Z^1/4± Z^-1/4| 
    ϕ_±   =  Arg[Z^1/4± Z^-1/4]  ,

where Z has been introduced in Eq. (<ref>).

In particular, if S_ R varies on a clockwise-oriented loop around S_ EP, the multivaluation in Eq. (<ref>) requires the amplitudes and the phases of the oscillations to vary as

    a_±(-2π)   =a_∓(0) 
    ϕ_±(-2π)    = ϕ_∓(0) +π/2 .



In addition, Eq. (<ref>) is the parametric equation of a rotated (damped) ellipse.

To make it explicit, we introduce a counterclockwise rotation matrix of angle φ, namely


    R(φ) =
        (
        [  cosφ -sinφ;  sinφ  cosϕ ]) ,

so that


    ( [ a_±cos(Ω_± t+ϕ_±); a_∓sin(Ω_± t+ϕ_∓);                   ])
        =
        R(φ)
        (
        [ αcos(ω_± t); βsin(ω_± t) ]) .


This straightforwardly leads to the following equalities:


    a_+ cos(ϕ_+)    = αcos(φ) 
    
        a_+ sin(ϕ_+)    = βsin(φ) 
    
        a_- cos(ϕ_-)    = βcos(φ) 
    
        a_- sin(ϕ_-)    = αsin(φ)  .


Along with Eq. (<ref>), this implies that if S_ R varies on a clockwise-oriented loop around S_ EP from θ=0 to θ=-2π,


    α(-2π)cos(φ(-2π))
    
    
    
    
           = α(0)cos(φ(0)+π/2) 
    α(-2π)sin(φ(-2π))
           = α(0)sin(φ(0)+π/2)  ,

and we find similar relations for β.

This demonstrates that the semiaxes α and β have undergone a π/2 rotation at the end of the loop, albeit their lengths remain unchanged.


Therefore, the mechanical oscillations of the two eigenmodes have elliptical polarizations in real space,


    x_±(t)
           = x_±(0)   e^-Γ_± t
        R(φ)
        (
        [ αcos(Ω_± t); βsin(Ω_± t) ]) .


The φ rotation of the polarization semiaxes α and β is controlled by the polar angle θ around the EP in parameter space.

This behavior is illustrated in Fig. <ref>(d).

If the distance ρ to the EP vanishes, then a_+≃ a_- and ϕ_+≃ϕ_- according to Eq. (<ref>).

The ellipse excentricity vanishes (α≃β) and the two mechanical eigenmodes coalesce into a single circularly polarized mode.

This is shown in Fig. <ref>(e).







§ FREQUENCY MODULATION OF THE TLS DRIVE




We start from the effective Hamiltonian of the driven TLS,


    H̃_ TLS(t) = 
    ( [   δ(t)/2 -g x(t)           Ω_ L/2;           Ω_ L/2 -δ(t)/2 + g x(t);                  ])  .


The luminescence excitation spectrum scales linearly with the excited-state population of the TLS.

The equations of motion for the reduced density matrix components σ_12=σ_21^* and σ_22=1-σ_11 are obtained from H̃_ TLS within the Born-Markov approximation.

They correspond to the usual optical Bloch equations,


    {[                                                 σ̇_12 = -i[δ(t) - 2g x(t) -i Γ/2] σ_12 - i Ω_ L (2σ_22-1)/2;                                                 σ̇_22                                    = Ω_ Lσ_12 - Γσ_22 ].  .


They can be rewritten as


    {[                          dX/dτ =-i[D(τ)-i/2] X -iϵ[ Y - 1/2 ];                          dY/dτ                      =-Y + ϵ X ].  .


where we have introduced X=σ_12=σ_21^*, Y=σ_22=1-σ_11, τ = Γ t, and ϵ=Ω_ L/Γ.

We then look for perturbative solutions of the type X = ∑_nϵ^n X_n and  Y = ∑_nϵ^n Y_n in the limit ϵ≪ 1.





 §.§ Order n=0


The Bloch equations lead to


    {[           dX_0/dτ =-i[D(τ)-i/2] X_0;           dY_0/dτ             =-Y_0 ].  ,


and the solutions generically read


    {[                                        X_0 = X_0(τ_0)  e^-(τ-τ_0)/2-i∫_τ_0^τdτ' D(τ');                                        Y_0                     = Y_0(τ_0)  e^-(τ-τ_0) ].  ,


In the limit τ_0→ -∞, they reduce to


    {[ lim_τ_0→ -∞ X_0             = 0; lim_τ_0→ -∞ Y_0             = 0 ].  .






 §.§ Order n=1


The Bloch equations lead to


    {[                            dX_1/dτ =-i[D(τ)-i/2] X_1 -i [ Y_0 - 1/2 ];                            dY_1/dτ                       =-Y_1 +  X_0 ].  ,


and the solutions generically read


    {[                                                                                                    X_1 = X_1(0)  e^-i(τ-τ_0)/2-i∫_τ_0^τdτ' D(τ') - i ∫_τ_0^τdτ' [ Y_0(τ') -1/2 ]  e^(τ'-τ)/2+i∫_τ^τ'dτ” D(τ”);                                                                                                    Y_1                                                 = Y_1(0)  e^-(τ-τ_0) + ∫_τ_0^τ dτ'   X_0(τ')  e^(τ'-τ) ].  .


In the limit τ_0→ -∞, they reduce to


    {[                              lim_τ_0→ -∞ X_1 = i/2 ∫_-∞^τdτ'  e^(τ'-τ)/2+i∫_τ^τ'dτ” D(τ”);                              lim_τ_0→ -∞ Y_1                                          = 0 ].  .






 §.§ Order n=2


The Bloch equations lead to


    {[                        dX_2/dτ = - i [D(τ)-i/2 ] X_2 - i  Y_1;                        dY_2/dτ                = - Y_ 2 +  X_1 ].  ,


and the solutions generically read


    {[                                                                                         X_2 = X_2(0)  e^-(τ-τ_0)2-i∫_τ_0^τdτ' D(τ') - i ∫_τ_0^τdτ' Y_1(τ')  e^(τ'-τ)/2+i∫_τ^τ'dτ” D(τ”);                                                                                         Y_2                                      = Y_2(0)  e^-(τ-τ_0) + ∫_τ_0^τ dτ'   X_1(τ')  e^(τ'-τ) ].  .


In the limit τ_0→ -∞, they reduce to


    {[                                      lim_τ_0→ -∞ X_2 = - i ∫_-∞^τdτ' Y_1(τ')  e^(τ'-τ)/2+i∫_τ^τ'dτ” D(τ”);                                      lim_τ_0→ -∞ Y_2                     = ∫_-∞^τ dτ'   X_1(τ')  e^(τ'-τ) ].  .




At the end of the loop around the EP, either the flexural modes have gotten exchanged or they have not.

In both cases, the flexural dynamics is of the type x(t)=Ae^-γ_itcos(ω_it+ϕ_i).

From now on, we consider ϕ_i=0, which fixes an arbitrary origin of time.

The luminescence excitation spectrum scales linearly with


    σ_22(t)    ≃ϵ^2 Y_2(t) 
       ≃Ω_ L^2/2 e^-Γ t∫_-∞^tdt_1  e^(Γ/2-iδ)t_1-iβ_λsin(ω_λt_1+ϕ_λ)∫_-∞^t_1dt_2  e^(Γ/2+iδ)t_2+iβ_λsin(ω_λt_2+ϕ_λ)
    e^-i2g∫_t_1^t_2dt_3x(t_3) ,


where β_λ=a_λ/ω_λ, and we introduce the dimensionless parameter β=2gA/ω_i.

In the limit β≪1, the luminescence excitation spectrum can be approximated by


    σ_22(t)
       ≃Ω_ L^2/4∑_m,nJ_m(β_λ)J_n(β_λ) [ e^-i(m-n)ϕ_λ (A_0+A_1-A_2-A_3+A_4)]


where the sum involves positive and negative values of the integers m and n.

In addition, J_k denotes the kth-order Bessel function of the first kind, and


    A_0   =e^-i(m-n)ω_λt/(δ+nω_λ-iΓ/2)((m-n)ω_λ+iΓ)
    
    A_1   = β(1-iγ_i/ω_i) e^[-i(m-n)ω_λ+iω_i-γ_i]t/(δ+nω_λ-iΓ/2)((m-n)ω_λ-ω_i+iΓ)
    
    A_2   = β(1-iγ_i/ω_i) e^[-i(m-n)ω_λ+iω_i-γ_i]t/(δ+nω_λ+ω_i-iΓ/2)((m-n)ω_λ-ω_i+iΓ)
    
    A_3   = β(1+iγ_i/ω_i) e^[-i(m-n)ω_λ-iω_i-γ_i]t/(δ+nω_λ-iΓ/2)((m-n)ω_λ+ω_i+iΓ)
    
    A_4   = β(1+iγ_i/ω_i) e^[-i(m-n)ω_λ-iω_i-γ_i]t/(δ+nω_λ-ω_i-iΓ/2)((m-n)ω_λ+ω_i+iΓ) .


The term A_0 does not depend on β and describes the luminescence excitation spectrum in the absence of electromechanical coupling.

We focus on Δσ_22=σ_22⌋_g≠0-σ_22⌋_g=0.

We can then accumulate the luminescence over time, which results in


    ∫_0^∞dt Δσ_22(t)
       ≃βΩ_ L^2/4∑_m,nJ_m(β_λ)J_n(β_λ) [ ie^-i(m-n)ϕ_λ (B_1-B_2-B_3+B_4)]


where


    B_1   = (1-iγ_i/ω_i)
    δ+nω_λ+iΓ/2/(δ+nω_λ)^2+(Γ/2)^2(m-n)ω_λ-ω_i-iΓ/((m-n)ω_λ-ω_i)^2+Γ^2(m-n)ω_λ-ω_i+iγ_i/((m-n)ω_λ-ω_i)^2+γ_i^2
    
    B_2   = (1-iγ_i/ω_i)
    δ+nω_λ+ω_i+iΓ/2/(δ+nω_λ+ω_i)^2+(Γ/2)^2(m-n)ω_λ-ω_i-iΓ/((m-n)ω_λ-ω_i)^2+Γ^2(m-n)ω_λ-ω_i+iγ_i/((m-n)ω_λ-ω_i)^2+γ_i^2
    
    B_3   = (1+iγ_i/ω_i)
    δ+nω_λ+iΓ/2/(δ+nω_λ)^2+(Γ/2)^2(m-n)ω_λ+ω_i-iΓ/((m-n)ω_λ+ω_i)^2+Γ^2(m-n)ω_λ+ω_i+iγ_i/((m-n)ω_λ+ω_i)^2+γ_i^2
    
    B_4   = (1+iγ_i/ω_i)
    δ+nω_λ-ω_i+iΓ/2/(δ+nω_λ-ω_i)^2+(Γ/2)^2(m-n)ω_λ+ω_i-iΓ/((m-n)ω_λ+ω_i)^2+Γ^2(m-n)ω_λ+ω_i+iγ_i/((m-n)ω_λ+ω_i)^2+γ_i^2 .


Every term B_k consists of a product of three Lorentzian functions.

The last one, which has the smallest width γ_i(≪Γ), results from the interferences due to the frequency modulation.

It describes narrow luminescence peaks of width γ_i every time the condition (m-n)ω_λ=ω_i is fulfilled.

For the typical values of the parameters that we are considering here, the nearly-degenerate flexural frequencies verify |ω_1-ω_2|≪γ_0, so the two modes should be well resolved through the interference peaks in experiments.






We can further estimate the characteristic amplitude of the interference peaks.

We consider that the laser frequency is modulated around the resonance (δ=0), where the two flexural modes are not coupled.

The interference condition reads (m-n)ω_λ=ω_i and requires m≠ n.

The main contribution in Eq. (<ref>) involves the zeroth-order Bessel function, and m and n have to be as small as possible.

Thus, the main contributions when ω_λ=ω_i arise from B_1 for (m,n)=(1,0), B_2 for (m,n)=(0,-1), B_3 for (m,n)=(-1,0), and B_4 for (m,n)=(0,1).

This leads to a peak of amplitude


    ∫_0^∞dt Δσ_22(t)
       ≃ - J_0(β_λ)J_1(β_λ)(cosϕ_λ-γ_i/ω_isinϕ_λ)β/γ_i(Ω_ L/Γ)^2 .


This explains the two interference peaks centered on the mechanical frequencies ω_1 and ω_2 in Fig. <ref>(b).





apsrev4-1


