Charge-vibration interaction effects in normal-superconductor quantum dots

TL;DR

This paper investigates charge-vibration interactions in a normal-superconductor quantum dot system, revealing how inelastic tunneling processes influence current, vibrational states, and enable vibrational spectroscopy at high temperatures.

Contribution

It provides a detailed analysis of inelastic tunneling effects and vibrational state control in quantum dots coupled to superconductors using Keldysh Green's functions.

Findings

Inelastic tunneling peaks are sharp and temperature-independent, enabling vibrational spectroscopy.

Inelastic processes induce diverse vibrational states, including cooling, heating, and instabilities.

Ground-state cooling can occur for multiple vibrational modes simultaneously.

Abstract

We study the quantum transport and the nonequilibrium vibrational states of a quantum dot embedded between a normal and a superconducting lead with the charge on the quantum dot linearly coupled to a harmonic oscillator of frequency $\omega$. To the leading order in the charge-vibration interaction, we calculate the current and the nonequilibrium phonon occupation by the Keldsyh Green's function technique. We analyze the inelastic, vibration-assisted tunneling processes in the regime $\omega <\Delta$, with the superconducting energy gap $\Delta$, and for sharp resonant transmission through the dot. When the energy $\varepsilon_0$ of the dot's level is close to the Fermi energy $\mu$, i.e. $|\varepsilon_0-\mu|\ll \Delta$, inelastic Andreev reflections dominate up to voltage $eV\gtrsim\Delta$. The inelastic quasiparticle tunneling becomes the leading process when the dot's level is close to the gap $|\varepsilon_0-\mu|\sim \Delta \pm \omega$. In both cases, the inelastic tunneling processes appear as sharp and prominent peaks - not broadened by temperature - in the $I$-$V$ characteristic and pave the way for inelastic spectroscopy of vibrational modes even at temperatures $T \gg \omega$. We also found that inelastic Andreev reflections as well as quasiparticle tunneling induce a strong nonequilibrium state of the oscillator. In different ranges on the dot's level, we found that the current produces: (i) ground-state cooling of the oscillator with phonon occupation $n\ll 1$, (ii) accumulation of energy in the oscillator with $n\gg 1$ and (iii) a mechanical instability which is a precursor of self-sustained oscillations. We show that ground-state cooling is achieved simultaneously for several modes of different frequencies. Finally, we discuss how the nonequilibrium vibrational state can be detected by the asymmetric behavior of the inelastic current peaks respect to the gate voltage.

Equations

$$\hat{H}={\hat{H}}_N+{\hat{H}}_S+{\hat{H}}_T+{\epsilon }_0{\hat{n}}_d+\lambda {\hat{n}}_d({\hat{b}}^{†}+{\hat{b}}^{})+\omega {\hat{b}}^{†}{\hat{b}}^{},$$

$$\check{G}(t,t^{\prime })=\left(\begin{array}{cc}{\hat{G}}^R(t,t^{\prime })& {\hat{G}}^K(t,t^{\prime })\\ 0& {\hat{G}}^A(t,t^{\prime })\end{array}\right),$$

$$\hat{F}(\epsilon )=\left(\begin{array}{cc}{f}_{+}(\epsilon )& 0\\ 0& f_{-}(\epsilon )\end{array}\right),$$

$${\hat{\Sigma }}_S^R(\epsilon )=-\frac{{\Gamma }_S}{\Omega (\epsilon )}\left(\begin{array}{cc}\epsilon +i\eta & -\Delta \\ -\Delta & \epsilon +i\eta \end{array}\right)$$

$$\hat{G}^{R}(\varepsilon)=\left(\begin{array}[]{cc}{G}_{+}(\varepsilon)&{F}(\varepsilon)\\ {F}(\varepsilon)&{G}_{-}(\varepsilon)\end{array}\right)\,,$$

$$G_{\pm }(\epsilon )$$

$$F(\epsilon )$$

$$D(\epsilon )=\frac{{\Delta }^2{\Gamma }_S^2}{{\Omega }^2(\epsilon )}-\left[\epsilon +{\epsilon }_0+i{\Gamma }_N+\frac{\epsilon {\Gamma }_S}{\Omega (\epsilon )}\right]\left[\epsilon -{\epsilon }_0+i{\Gamma }_N+\frac{\epsilon {\Gamma }_S}{\Omega (\epsilon )}\right]$$

$${\rho }_d(\epsilon )=-\text{Im}\text{Tr}[{\hat{G}}^R(\epsilon )],$$

$$\check{D}(t,t^{\prime })=\left(\begin{array}{cc}{D}^R(t,t^{\prime })& D^K(t,t^{\prime })\\ 0& D^A(t,t^{\prime })\end{array}\right)$$

$${\Pi }^{R,A}(\epsilon )$$

$$+$$

$${\Pi }^K(\epsilon )$$

$$+$$

$$\mathrm{Im}{\Sigma }_0^R(\epsilon )$$

$${\Sigma }_0^K(\epsilon )$$

$$D^{R,A}(\epsilon )$$

$$\simeq \frac{1}{\epsilon -\omega \pm i{\gamma }_{\mathrm{tot}}}-\frac{1}{\epsilon +\omega \pm i{\gamma }_{\mathrm{tot}}},$$

$$D^K(\epsilon )$$

$$\simeq i\pi \frac{{\Pi }^K(\epsilon )+{\Sigma }_0^K(\epsilon )}{{\gamma }_{\mathrm{tot}}(\epsilon )}[\delta (\epsilon -\omega )+\delta (\epsilon +\omega )].$$

$${\gamma }_{\mathrm{tot}}(\omega )=-\text{Im}[{\Pi }^R(\omega )+{\Sigma }_0^R(\omega )]=\gamma +{\gamma }_0.$$

$$\gamma ={\gamma }_{\mathrm{NN}}+{\gamma }_{\mathrm{NS}}+{\gamma }_{\mathrm{SS}},$$

$${\gamma }_{\mathrm{NN}}={\gamma }_{\mathrm{AR}}+{\gamma }_{\mathrm{NR}}.$$

$${\gamma }_{\mathrm{AR}}$$

$${\gamma }_{\mathrm{NR}}$$

$${\gamma }_{\nu {\nu }^{\prime }}^s(\omega )={\lambda }^2\frac{{\Gamma }_N^2}{2}\int \frac{d\epsilon }{2\pi }{\left|T_{\nu {\nu }^{\prime }}^s(\epsilon )\right|}^2{f}_{\nu }(\epsilon )[1-f_{{\nu }^{\prime }}(\epsilon +s\omega )],$$

$$T_{\nu \bar{\nu }}^s(\epsilon )$$

$$T_{\nu \nu }^s(\epsilon )$$

$${\gamma }_{\mathrm{NS}}$$

$${\gamma }_{\mathrm{NS}}^{s\nu }=16{\lambda }^2{\Gamma }_N{\Gamma }_S\int \frac{d\epsilon }{2\pi }{\mathcal{T}}_{\nu }^s(\epsilon ){\rho }_S(\epsilon +s\omega )f_{\nu }(\epsilon )[1-f_0(\epsilon +s\omega )].$$

$${\mathcal{T}}_{\nu }^s(\epsilon )=|T_{\nu \nu }^s(\epsilon ){|}^2+|T_{\nu \bar{\nu }}^s(\epsilon ){|}^2+\frac{2\Delta }{\epsilon +s\omega }\text{Re}[T_{\nu \nu }^s(\epsilon ){T_{\nu \bar{\nu }}^s(\epsilon )}^{\ast }].$$

$${\gamma }_{\mathrm{SS}}=\sum _{\nu s}s{\gamma }_{\text{SS}}^{s\nu }.$$

$${\gamma }_{\text{SS}}^{s\nu }=8{\lambda }^2{\Gamma }_S^2\int \frac{d\epsilon }{2\pi }{\mathcal{R}}_{\nu }^s(\epsilon ){\rho }_S(\epsilon ){\rho }_S(\epsilon +s\omega )f_0(\epsilon )[1-f_0(\epsilon +s\omega )].$$

$${\mathcal{R}}_{\nu }^s(\epsilon )$$

$$+\frac{{\Delta }^2}{\epsilon (\epsilon +s\omega )}(T_{\nu \nu }(\epsilon )T_{\bar{\nu }\bar{\nu }}^{\ast }(\epsilon )+T_{\nu \bar{\nu }}(\epsilon )T_{\bar{\nu }\nu }^{\ast }(\epsilon ))].$$

Full text

I Introduction

Electronic transport through nanoscale devices is characterized by a variety of physical phenomena [Weiss, 2009,Cuevas and Scheer, 2010]. One of them is the interplay between quantum transport at the level of single electron and the mechanical motion of localized vibrations in various nanoscale devices. These include single-molecule junctions [Galperin et al., 2007; Aradhya and Venkataraman, 2013; Park et al., 2000; Yu et al., 2004; Qiu et al., 2004; Wu et al., 2004; Pasupathy et al., 2005; Thijssen et al., 2006; Böhler et al., 2007; de Leon et al., 2008; Franke and Pascual, 2012], suspended carbon nanotube quantum dots (CNT-QDs)[Nygård and Cobden, 2001; LeRoy et al., 2004; Sapmaz et al., 2005; LeRoy et al., 2005; Sapmaz et al., 2006; Leturcq et al., 2009; Hüttel et al., 2009; Steele et al., 2009; Lassagne et al., 2009; Utko et al., 2010; Island et al., 2011; Ganzhorn and Wernsdorfer, 2012; Schmid et al., 2012; Meerwaldt et al., 2012; Island et al., 2012; Benyamini et al., 2014; Weber et al., 2015; Deng et al., 2016] and several types of nanoelectromechanical systems (NEMS) [Blencowe, 2005; Greenberg et al., 2012; Poot and van der Zant, 2012] as single-electron transistors [Knobel and Cleland, 2003; Armour et al., 2004; Blanter et al., 2004; Chtchelkatchev et al., 2004; Doiron et al., 2006], superconducting single-electron transistors [LaHaye, 2004; Naik et al., 2006; Clerk and Bennett, 2005; Blencowe et al., 2005; Koerting et al., 2009], single-electron shuttles [Erbe et al., 2001; Gorelik et al., 2001; Novotný et al., 2003; Fedorets, 2003; Gorelik et al., 2006; Smirnov et al., 2004; Fedorets et al., 2004; Pistolesi, 2004; Isacsson and Nord, 2007], shuttling nanopilars [Kim et al., 2013,Prada et al., 2014] and quantum dots in suspended nanostructures [Weig et al., 2004,Okazaki et al., 2016]. In such nanodevices the motion of the resonator can be not only detected but also manipulated via electron transport. At the same time, the mechanical motion strongly influences the transport itself.

For instance, in quantum dots, the interaction between the electrons and the harmonic vibrational modes gives rise to inelastic tunneling processes appearing when consecutive vibrational levels enter the bias window. Such an inelastic spectroscopy has been realized in experiments in single molecules [Yu et al., 2004; Wu et al., 2004; Thijssen et al., 2006; Böhler et al., 2007; de Leon et al., 2008] and in suspended CNT-QDs [Sapmaz et al., 2005; LeRoy et al., 2005; Sapmaz et al., 2006; Leturcq et al., 2009; Hüttel et al., 2009] at least for the high-frequency vibrational modes $\omega\gg T$ ($\hbar=k_{B}=1$), viz. longitudinal or radial modes, but not for the low-frequency flexural modes. As increasing the bias voltage, tunneling electrons can affect strongly the harmonic vibration leading to current-induced nonequilibrium vibrational states, in some cases up to the threshold of molecular dissociation [Franke and Pascual, 2012].

More generally, charge-vibration interaction leads to a plethora of novel and unexpected behaviors of the nanoresonators. Electromechanical backaction effects - as oscillator’s frequency shift and electromechanical damping - have been reported in experiments for the flexural modes in suspended CNT-QDs [Steele et al., 2009,Lassagne et al., 2009,Meerwaldt et al., 2012,Benyamini et al., 2014], in quantum dots coupled to a piezoelectric nanoresonator [Okazaki et al., 2016] and in a dc-SQUID integrating a flexural resonator [Poot et al., 2010]. Further increasing the coupling strength, current suppression is expected due to the Franck-Condon blockade mechanism which was experimentally observed in CNT-QDs [Leturcq et al., 2009]. Mechanical bistabilities and blocked-current states are also theoretically predicted beyond a critical threshold of the charge-vibration coupling strength for low frequency (classical) nanoresonators [Galperin et al., 2005; Mozyrsky et al., 2006; Pistolesi and Labarthe, 2007; Pistolesi et al., 2008; Micchi et al., 2015, 2016]. Other theoretical studies reported negative electromechanical damping [Usmani et al., 2007], self-sustained mechanical oscillations [Sonne et al., 2008; Radić et al., 2011; Atalaya and Gorelik, 2012; Vikstr om et al., 2016] and charge-induced non-linear effects [Meerwaldt et al., 2012,Nocera et al., 2012,Nocera et al., 2013].

The majority of experimentally realized single-molecule junctions - as well as suspended CNT-QDs - can be properly modeled as single quantum dots [Galperin et al., 2007]. Hence, the single-impurity Holstein model [Glazman and Shekhter, 1988; Wingreen et al., 1989; König et al., 1996; Boese and Schoeller, 2001; Lundin and McKenzie, 2002; Kuo and Chang, 2002; McCarthy et al., 2003; Braig and Flensberg, 2003; Flensberg, 2003; Galperin et al., 2004; Mitra et al., 2004; Wegewijs and Nowack, 2005; Kaat and Flensberg, 2005] has become the paradigmatic model to discuss the effects of charge-vibration interaction in such systems. In its simple formulation, one assumes a linear coupling between the electron occupation on the dot and the oscillation’s amplitude of one or more harmonic oscillators representing the local vibrations. This model has been theoretically investigated widely in literature with different variations, extensions, in different regimes and with various theoretical approaches [Koch and von Oppen, 2005a, b; Zazunov et al., 2006a; Koch et al., 2006a, b; Hwang et al., 2007; Lüffe et al., 2008; Tahir and MacKinnon, 2008; Entin-Wohlman et al., 2009; Haupt et al., 2009, 2010; Schultz, 2010; Cavaliere et al., 2010; Yar et al., 2011; Piovano et al., 2011; Novotný et al., 2011; Fang et al., 2011; Li et al., 2012; Skorobagatko, 2012; Rastelli et al., 2012; Santamore et al., 2013; Härtle et al., 2013; Knap et al., 2013; Walter et al., 2013; Agarwalla et al., 2015; Schinabeck et al., 2016; Sowa et al., 2017], in particular using diagrammatic techniques [Ryndyk et al., 2006; Zazunov and Martin, 2007; Egger and Gogolin, 2008; Dash et al., 2010; Rastelli et al., 2010; Dash et al., 2011; Seoane Souto et al., 2014; Chen et al., 2016] and numerical and non-perturbative methods [Mühlbacher and Rabani, 2008; Hützen et al., 2012; Eidelstein et al., 2013; Jovchev and Anders, 2013; Seoane Souto et al., 2015].

However, all previous mentioned theoretical studies of the Holstein model focused on the case of a quantum dot in contact with two normal-conducting leads (N-QD-N). Other theoretical works investigated the problem of quantum dots coupled to the local vibrations and between two superconductors (S-QD-S) as a Josephson dot coupled to several bosonic modes [Novotný et al., 2005], a suspended carbon-nanotube acting as nanomechanical resonator [Sonne et al., 2008; Zazunov et al., 2006b, c; Sonne et al., 2010a, b; Sadovskyy et al., 2010; Padurariu et al., 2012] or a molecular Josephson junction model [Wu et al., 2012]. The latter was used to explain the resonance peaks experimentally observed in vibrating Nb nanowires [Marchenkov et al., 2007,V. Kretinin et al., ].

On the contrary, the hybrid case of a normal-dot-superconductor (N-QD-S) coupled to the bosonic modes - as photons [Sun et al., 1999,Baranski and Domanski, 2015] or phonons [Zhang et al., 2012; Wang et al., 2013; Dong et al., 2017] - has been less studied in the literature and is theoretically unexplored for the nonequilibrium regime of the oscillator. In this system, one expects that the role of the superconducting lead is unessential for $\Delta\ll\omega$ with the superconducting gap $\Delta$, as the vibrational frequency $\omega$ sets the energy scale of the inelastic tunneling of electrons. The interesting and relevant regime is thus expected for $\Delta\sim\omega$ or $\Delta>\omega$, which corresponds to an unrealizable condition in single-molecule junctions characterized by high-frequency vibrations. This regime can occur readily and naturally in suspended CNT-QDs, see Fig. 1. Moreover, CNT-QDs can be inserted between different types of nanocontacts, as a normal metal and a superconductor [Pillet et al., 2013; Schindele et al., 2014; Gramich et al., 2015, 2016, ].

In a previous recent work [Stadler et al., 2016] we explored the N-QD-S system with charge-vibration interaction in the deep subgap regime corresponding to a gap much larger than all energy scales (viz. $\Delta\rightarrow\infty$). Our main findings was that the vibrational modes can be cooled to the ground state by inelastic vibration-assisted Andreev reflections with phonon number occupation $n\ll 1$ up to the backaction limit dictated by the shot noise and given by $n_{BA}\sim{\left(\Gamma_{N}/\omega\right)}^{2}$ with the tunneling rate $\Gamma_{N}$ to the normal lead. We also showed that such cooling is exploitable even for several mechanical modes. Ground-state cooling by electron transport in suspended CNT-QDs appears as a promising strategy to achieve the quantum regime of low-frequency flexural modes $\omega\ll T$ and it has been the subject of a notable research activity during the past years [Lundin, 2004; Pistolesi, 2009; Zippilli et al., 2009, 2010; Li et al., 2011; Santandrea et al., 2011; Arrachea et al., 2014; Stadler et al., 2014, 2015, 2016]. We recall that cooling by electron transport has been experimentally demonstrated in nanomechanical beams integrating normal metal-insulator-superconductor tunnel junctions [Koppinen and Maasilta, 2009,Pekola and Hekking, 2007].

The objective of this work is to extend the calculations of Ref. [Stadler et al., 2016] and consider a finite superconducting gap $\Delta$. We aim at demonstrating that the N-QD-S is indeed of particular interest as: (i) it is characterized by a rich and interesting behavior of the inelastic transport ruled by inelastic Andreev reflections as well as by inelastic quasiparticle tunneling [see Fig. 2], (ii) it opens the route to inelastic spectroscopy of the low-frequency modes as the flexural ones which fundamentally can not be resolved using normal leads at temperature $T\gg\omega$, (iii) it allows strong manipulation of the vibrational states via electron transport ranging from ground-state cooling up to a mechanical instability, which is possibly a precursor of self-sustained oscillations. Finally, such a system offers intrinsically the possibility to read-out the equilibrium or the nonequilibrium state of the vibration directly by the observation of the dc-current.

The paper is structured as follows. In Sec. II we introduce the Holstein model of a N-QD-S with charge-vibration interaction and discuss the range of validity of our model and of the perturbative diagrammatic expansion.

In Sec. III we discuss the results for the electromechanical damping $\gamma$ showing that it is related to several microscopic processes, in particular the inelastic Andreev reflection and inelastic quasiparticle tunneling between the normal to the superconducting lead which can drive the vibration out of equilibrium [see Fig. 2].

The phonon number occupation is discussed in the Sec. IV. In Sec. IV.1 we discuss the phonon occupation for a single vibrational mode. Remarkably, in the limit $\omega\ll\Delta$, there is a separation of the ranges of the dot’s energy level where the inelastic Andreev reflections dominates over the inelastic quasiparticles tunneling or vice-versa. Both kind of processes can yield ground-state cooling of the oscillator with phonon occupation $n\ll 1$, an accumulation of energy in the oscillator with $n\gtrsim 1$ and a mechanical instability $n\gg 1$ signaled by the negativity of the electromechanical damping. This region, as well as the region of $n\gg 1$, are beyond our perturbative approach and the validity of harmonicity of the resonator. In the Sec. IV.2 we generalize the discussion taking into account many vibrational modes showing that ground-state cooling can be achieved simultaneously for several N modes of different frequencies $\omega_{n}$ via inelastic quasiparticle tunneling [Stadler et al., 2016].

In the Sec. V we show our results for the current. We discuss mainly the inelastic correction to the current and show how the nonequilibrium state of the vibration can be properly read out by the asymmetric behavior of inelastic current peaks respect to the gate voltage. Section VI contains concluding remarks.

II Model and Approximation

In this section, we discuss a model of a quantum dot suspended between a normal-conducting and a superconducting contact which can be realized e.g. with a suspended CNT-QD [see Fig. 1].

We consider the electrons on the quantum dot occupying a spin-degenerate state with energy $\varepsilon_{0}$ and the annihilation and creation operators $\hat{d}^{\phantom{g}}_{\sigma}$ and $\hat{d}^{\dagger}_{\sigma}$ with spin $\sigma$ on the quantum dot. The occupation number reads $\hat{n}_{d}=\hat{d}_{\uparrow}^{\dagger}\hat{d}_{\uparrow}^{\phantom{g}}+\hat{d}_{\downarrow}^{\dagger}\hat{d}_{\downarrow}^{\phantom{g}}$. A single harmonic mode of the resonator with frequency $\omega$ has the bosonic annihilation and creation operators $\hat{b}^{\phantom{g}}$ and $\hat{b}^{\dagger}$. The electrostatic force between the electrons and a gate voltage induces a coupling between the charge and the deflection of the resonator with coupling strength $\lambda$. The full Hamiltonian is given by

[TABLE]

with the last two terms describing the charge-vibration coupling and the resonator. The Hamiltonians of the normal and the superconducting leads are given by $\hat{H}_{N}=\sum_{k\sigma}\xi_{N,k}\hat{c}^{\dagger}_{k\sigma}\hat{c}^{\phantom{\dagger}}_{k\sigma}$ and $\hat{H}_{S}=\sum_{k\sigma}(\xi_{S,k}^{\phantom{\dagger}}\hat{a}_{k\sigma}^{\dagger}\hat{a}_{k\sigma}+\Delta\hat{a}_{-k\downarrow}^{\dagger}\hat{a}_{k\uparrow}^{\dagger}+\Delta\hat{a}_{k\uparrow}\hat{a}_{-k\downarrow})$ with the real pairing potential $\Delta$ (superconducting gap) and the energy referring to the chemical potential $\xi_{\alpha,k}=\varepsilon_{k}-\mu_{\alpha}$ with $\alpha=(N,S)$. The annihilation operators with spin $\sigma$ and momentum $k$ are given by $\hat{a}_{k\sigma}$ for the superconducting and $\hat{c}_{k\sigma}$ for the normal-conducting lead. The tunneling Hamiltonian is $\hat{H}_{T}=\sum_{k\sigma}(t_{N}^{\phantom{\dagger}}\hat{c}^{\dagger}_{k\sigma}\hat{d}^{\phantom{\dagger}}_{\sigma}+t_{S}^{\phantom{\dagger}}\hat{a}^{\dagger}_{k\sigma}\hat{d}^{\phantom{\dagger}}_{\sigma}+\mathrm{H.c.})$ with the tunneling energies between the leads and the dot $t_{N}$ and $t_{S}$. In the rest of our analysis we assume the wide-band approximation and write the tunnel couplings between the quantum dot with the normal contact and with the superconducting contact as $\Gamma_{N}=\pi{\left|t_{N}\right|}^{2}\rho_{N}$ and $\Gamma_{S}=\pi{\left|t_{S}\right|}^{2}\rho_{S}$, respectively. The electronic density of states at the Fermi level of the normal lead and of the superconducting lead in the normal (non superconducting) phase are denoted by $\rho_{N}$ and $\rho_{S}$.

Before we conclude this section, we discuss the role of the Coulomb interaction in the dot. In the past, the Anderson model with local Hubbard repulsion $U$ in a N-QD-S system without charge vibration interaction was studied to understand the competition between superconducting pairing, Coulomb repulsion and Kondo correlations [Fazio and Raimondi, 1998; Schwab and Raimondi, 1999; Clerk et al., 2000; Cuevas et al., 2001; Tanaka et al., 2007; Martín-Rodero and Levy Yeyati, 2011]. Analytic methods were restricted to truncated approximated expansions in the regime $U\rightarrow\infty$ and a finite superconducting gap [Fazio and Raimondi, 1998; Schwab and Raimondi, 1999; Clerk et al., 2000]. An interpolation approach was used in Ref. [Cuevas et al., 2001] for a wide range of the ratio of $U/\Gamma$, showing that, for $U/\Gamma<1$, superconducting pairing dominates the behavior of the quantum dot. Tanaka et al. [Tanaka et al., 2007] employed a numerical renormalization group method in the limit $\Delta\rightarrow\infty$ which allows to map the problem to the case of a QD with a local pairing amplitude $\Delta_{D}=\Gamma_{S}$. We recovered the same effective Hamiltonian in our previous work [Stadler et al., 2016], viz. $\hat{H}_{dS}=\varepsilon_{0}\sum_{\sigma}\hat{d}_{\sigma}^{\dagger}\hat{d}_{\sigma}^{\phantom{g}}-\Gamma_{S}(\hat{d}_{\downarrow}^{\phantom{g}}\hat{d}_{\uparrow}^{\phantom{g}}+\mbox{h.c.})$. Reference [Tanaka et al., 2007] - limited in the linear transport regime - points out that the effect of the interaction is to renormalize the coupling $\Gamma_{S}$ and the dot level $\varepsilon_{0}$ in the regime $U\lesssim\Gamma_{S}$.

Tackling the problem of systems with Coulomb interaction and charge-vibration shows a formidable task beyond the scope of this work. Here we assume the tunneling rate with the normal lead $\Gamma_{N}$ as the lowest energy involved in the problem and the condition for the energy scales $U\lesssim\Gamma_{S}<\Delta$ which sets a conservative estimation for the range of validity of our approach. Finally, we remark that experiments with $U\lesssim\Delta$ as in Refs. [Gramich et al., 2015,Gramich et al., 2016,Deacon et al., 2010] focused on devices with large asymmetry $\Gamma_{S}\gg\Gamma_{N}$ basically providing a spectroscopic probe of the Andreev spectra. The results can be qualitatively explained in terms of a non-interacting model, viz. the Kondo effect is suppressed by the strong proximity effect and correlation effects seems to play no role.

II.1 Exact results for the quatum dot for $\lambda=0$

In this subsection we recall the exact results for the electronic Green’s function on the dot without charge-vibration interaction and the behavior of the density of state of the quantum dot [Cuevas et al., 1996,Blonder et al., 1982]. These Green’s function are the building block by which we can express the electromechanical damping, the phonon occupation, the zero-order current and the current corrections in presence of the charge-vibration interaction.

In Keldysh space, we write the electronic Green’s function as

[TABLE]

with the retarded (R), advanced (A) and Keldysh (K) Green’s functions defined as $\hat{G}^{R}(t,t^{\prime})=-i\theta(t-t^{\prime})\langle\{\Psi_{d}(t),\Psi^{\dagger}_{d}(t^{\prime})\}\rangle$, $\hat{G}^{A}(t,t^{\prime})=i\theta(t^{\prime}-t)\langle\{\Psi_{d}(t),\Psi^{\dagger}_{d}(t^{\prime})\}\rangle$, and $\hat{G}^{K}(t,t^{\prime})=-i\langle[\Psi_{d}(t),\Psi^{\dagger}_{d}(t^{\prime})]\rangle$. Here, we introduced the spinor in electron-hole space as $\Psi^{\dagger}_{d}(t)=(d^{\dagger}_{\uparrow}(t)\,d_{\downarrow}(t))$ acting in the Nambu space and denote the commutator and anti-commutator with $[\,,\,]$ and $\{\,,\,\}$, respectively. Notice the notation $\check{}$ for the matrix defined on the Keldysh contour space and $\hat{}$ for the matrix defined on the Nambu space.

The Dyson equation for the Green’s function on the dot is $\check{G}(\varepsilon)=\check{g}(\varepsilon)+\check{g}(\varepsilon)[\check{\Sigma}_{N}(\varepsilon)+\check{\Sigma}_{S}(\varepsilon)]\check{G}(\varepsilon)\,,$ with the unperturbed Green’s function of the dot $\check{g}(\varepsilon)$ and the self-energies $\check{\Sigma}_{\alpha}(\varepsilon)$ of the normal and superconducting lead [$\alpha=(N,S)$].

In the following, we apply the voltage on the normal-conducting lead $\xi_{N,k}=\varepsilon_{k}-\mu-eV$, $\xi_{S,k}=\varepsilon_{k}-\mu$ and let $\mu=0$. Using the wide-band approximation, the retarded and Keldysh self-energies of the normal metal are $\hat{\Sigma}^{R}_{N}=-i\Gamma_{N}$ and $\hat{\Sigma}_{N}^{K}(\varepsilon)=-2i\Gamma_{N}[1-2\hat{F}(\varepsilon)]$ with the matrix function

[TABLE]

in which we define the two functions $f_{\nu}(\varepsilon)=\{1+\mathrm{exp}[(\varepsilon-\nu\,eV)/T]\}^{-1}$ with the index $\nu=\pm=(e,h)$ to account for the sign of the charge for the electrons and holes, respectively. The retarded and Keldysh self-energies of the superconductor are

[TABLE]

and $\hat{\Sigma}_{S}^{K}(\varepsilon)=[1-2f_{0}(\varepsilon)][\hat{\Sigma}_{S}^{R}(\varepsilon)-\hat{\Sigma}_{S}^{A}(\varepsilon)]$ with $\Omega(\varepsilon)=\sqrt{\Delta^{2}-(\varepsilon+i\eta)^{2}}$, an infinitesimal small real part $\eta$ and the Fermi function of the superconducting lead $f_{0}(\varepsilon)$ $(\nu=0)$. Solving the Dyson equation, we obtain the retarded Green’s function of the dot

[TABLE]

with

[TABLE]

and

[TABLE]

The Keldysh Green’s function of the dot is obtained from the Dyson equation and is given by $\hat{G}^{K}(\varepsilon)=\hat{G}^{R}(\varepsilon)(\hat{\Sigma}_{N}^{K}(\varepsilon)+\hat{\Sigma}_{S}^{K}(\varepsilon))\hat{G}^{A}(\varepsilon)$.

The density of states $\rho_{d}(\varepsilon)$ on the quantum dot is calculated by

[TABLE]

with the retarded Green’s function given by Eq. (5). An example of the density of states as a function of energy is shown in Fig. 3(a). For $|\varepsilon|\gg\Delta$, $\rho_{d}(\varepsilon)$ shows a similar behavior as for a quantum dot coupled to two normal leads. At low energy, the superconducting correlation penetrate on the quantum dot and lead to a suppression of the density of states and to the formation of Andreev levels pinned inside the superconducting gap and signaled by the two peaks. These Andreev levels $E_{A}$ have an energy-dependence given by the zeroing of $D(\varepsilon)$ in Eq. (8) [see Fig. 3(b)]. For large superconducting gap, the Andreev levels are approximated by $E_{A}\simeq\pm\sqrt{\varepsilon_{0}^{2}+\Gamma_{S}^{2}}$ valid for $|E_{A}|\ll\Delta$.

II.2 Phonon Green’s function with charge-vibration interaction

In this section we discuss the phonon Green’s function for one harmonic oscillator. We use the perturbative expansion to calculate the phonon self-energy $\check{\Pi}(\varepsilon)$ to lowest order in the charge-vibration coupling [Mitra et al., 2004]. In this way, we can obtain the electromechanical damping and, by solving the Dyson equation, the nonequilibrium phonon occupation. Analytic formulas are possible in the limit of small damping $\gamma\ll\omega$ and a high quality factor of the resonator. Such an approach is equivalent to solve the quantum Langevin equation for the harmonic oscillator assuming as input the bare quantum noise correlators of the dot’s charge in the N-QD-S without charge-vibration interaction [Clerk, 2004] and using the Markovian approximation. The latter approximation holds in our case as the typical frequency scale $\delta\omega$ of the frequency depending electromechanical damping and of the noise scales as $\delta\omega\sim\Gamma_{N}$ and is still larger than then linewidth of the oscillator’s susceptibility which scales as $\sim\lambda^{2}\Gamma_{N}/\omega_{0}^{2}$ in the regime $\Gamma_{N}<\omega$.

We introduce the symmetrized bosonic operator $\hat{A}(t)=\hat{b}^{\dagger}(t)+\hat{b}(t)$, the contour ordering operator $T_{c}$ on the Keldysh contour and define the phonon propagator on the Keldysh contour as $D(\tau,\tau^{\prime})=-i\langle T_{c}\hat{A}(\tau)\hat{A}^{\dagger}(\tau^{\prime})\rangle$ with the contour variable $\tau$ on the Keldysh contour. Then, we transform the contour variable $\tau$ to the real time and use the Larkin-Ovchinnikov rotation to write the Dyson equation in the triangular form. The phonon Green’s function in Keldysh space are given by

[TABLE]

with the phonon Green’s functions $D^{R}(t,t^{\prime})=-i\theta(t-t^{\prime})\langle[\hat{A}(t),\hat{A}(t^{\prime})]\rangle$, $D^{A}(t,t^{\prime})=i\theta(t^{\prime}-t)\langle[\hat{A}(t),\hat{A}(t^{\prime})]\rangle$, and $D^{K}(t,t^{\prime})=-i\langle\{\hat{A}(t^{\prime}),\hat{A}(t)\}\rangle$. Finally, we obtain the Dyson equation in Fourier space as $\check{D}(\varepsilon)=\check{D}_{0}(\varepsilon)+\check{D}_{0}(\varepsilon)\check{\Pi}(\varepsilon)\check{D}_{0}(\varepsilon).$ The bare phonon Green’s functions in $\check{D}_{0}(\varepsilon)$ are given by $D^{R,A}_{0}(\varepsilon)=2\omega/[(\varepsilon\pm i\eta)^{2}-\omega^{2}]$ and $D^{K}_{0}(\varepsilon)=-2\pi i\left[\delta(\varepsilon\mathord{-}\omega)\mathord{+}\delta(\varepsilon\mathord{+}\omega)\right]\mathrm{coth}[\omega/(2T)]$. The polarization $\check{\Pi}(\varepsilon)$ corresponds to the phonon self-energy (polarization diagram) to the leading order in the electron-vibration interaction $\sim\lambda^{2}$ between the resonator and the electrons. These self-energies are given by

[TABLE]

and

[TABLE]

in which the trace $\mathrm{Tr}$ and the Pauli matrix $\hat{\tau}_{z}$ acts on the Nambu space.

To take into account the intrinsic damping of the oscillator, we additionally include a self-energy $\check{\Sigma}_{0}(\varepsilon)$ in Keldysh space. From the phenomenological model [Stadler et al., 2015], one obtains the expressions

[TABLE]

in which the coefficient $Q$ corresponds to the quality factor of the resonator.

Using the Dyson equation, we obtain the retarded, advanced and Keldysh equations

[TABLE]

with the total damping rate given by

[TABLE]

which can be separated into an intrinsic damping $\gamma_{0}=-\mathrm{Im}\,\Sigma_{0}^{R}(\omega)=\omega/Q$ and the damping $\gamma=-\mathrm{Im}\,\Pi^{R}(\omega)$ due to the charge-vibration interaction. In Eq. (17) we expanded the retarded and advanced phonon Green’s function close to $\varepsilon\simeq\omega$ taking into account the frequency shift of the resonator $\Delta\omega=\textrm{Re}\,\Pi^{R}(\omega)$.

III Results for the electromechanical damping

We found that the electromechanical damping $\gamma$ in Eq. (19) due to the electron-vibration interaction can be divided into basic reflection and transmission processes of charge transport occurring in N-QD-S systems [Cuevas et al., 1996,Blonder et al., 1982]. After some algebra, we cast the electromechanical damping as

[TABLE]

in which the damping coefficients $\gamma_{\alpha\beta}$ with the labels $\alpha=(\mathrm{N},\mathrm{S})$ and $\beta=(\mathrm{N},\mathrm{S})$ are associated to an inelastic process with an absorption or emission of one phonon with an incoming particle and an outgoing particle involving the leads $\alpha$ and $\beta$. Examples of such processes are given in Fig. 2.

The first term in Eq. (20) corresponds to reflections at the normal lead and can be divided into inelastic Andreev reflection (AR) and normal inelastic reflection (NR),

[TABLE]

An schematic view of inelastic AR and NR with an absorption of a phonon before an AR or NR is shown in Fig. 2(a) in which the excitation spectra of the normal and the superconducting lead around the Fermi momenta $k_{F}$ are shown as gray lines. Below the superconducting gap, these processes are the only possible processes as an incoming electron can only be Andreev reflected or normal reflected.

The damping rate $\gamma_{\mathrm{AR}}$ and $\gamma_{\mathrm{NR}}$ can be written as the sum of elemental processes

[TABLE]

with the index $\nu=\pm=(e,h)$ for electrons and holes, respectively, and the notation $\bar{\nu}=-\nu$. The index $s=\pm$ refers to an inelastic process with an absorption ($s=+$) and the emission ($s=-$) of one phonon. As an example, $\gamma_{eh}^{s}$ corresponds to an inelastic AR with an absorption or emission of a phonon for an incoming electron and an outgoing hole at the normal lead. The damping rates $\gamma_{\nu\nu^{\prime}}^{s}$ in Eqs. (22) and (23) can be written as

[TABLE]

with the Fermi functions $f_{\nu}(\varepsilon)=1/[1+\mathrm{exp}((\varepsilon-\nu\,eV))/T]$ for electron and holes, respectively, and the transmission functions

[TABLE]

As discussed in Ref. [Stadler et al., 2016] the transmission functions consist of a coherent sum of two amplitudes that are associated to the two possible paths in which the phonon is emitted or absorbed before or after a reflection. To the lowest order in the tunneling rates, the electromechanical dampings due to AR and NR are proportional to $\sim\Gamma^{2}_{S}\Gamma_{N}^{2}$ and $\sim\Gamma^{2}_{N}$, respectively. An example of the behavior of $\gamma_{\mathrm{AR}}$ and $\gamma_{\mathrm{NR}}$ is shown in Fig. 4. We analyzed in detail these processes in Ref. [Stadler et al., 2016].

We now turn to discuss the term $\gamma_{\mathrm{NS}}$ in which a particle is exchanged between the normal contact (electron or hole) and the superconducting contact (quasiparticle). We can break up these terms in absorption or emission processes with a initial electron or hole, namely

[TABLE]

with

[TABLE]

In Eq. (28) we have introduced the superconducting density of states defined as ${\rho}_{S}(\varepsilon)=\theta(\varepsilon^{2}-\Delta^{2})\varepsilon/\sqrt{\varepsilon^{2}-\Delta^{2}}$ normalized to the one of the normal phase, and the transmission function

[TABLE]

The processes described by Eq.(28) admit two equivalent descriptions: for instance the term for $\nu=e$ can be seen as an impinging electron from the normal lead $N$ to the superconductor $S$, leaving an hole in the lead $N$, or it can be seen as an impinging quasiparticle from the superconductor $S$ transmitted as an hole in the normal lead $N$.

Comparing the results Eqs. (28) and (29) with the ones reported in literature [Blonder et al., 1982,Cuevas et al., 1996], it turns out that the term in Eq. (29) proportional to $|T_{\nu\nu}(\varepsilon)|^{2}$ corresponds to the inelastic direct tunneling of quasiparticles (DT) and the term with $|T_{\nu\bar{\nu}}(\varepsilon)|^{2}$ corresponds to branch crossing (BC), according to the terminology introduced in the BTK model [Blonder et al., 1982]. When an electron is transmitted via inelastic DT, it occupies an electron-like branch on the superconductor. Branch crossing refers to the transmission of the electron on the hole-like branch of the superconductor. Examples of processes corresponding to inelastic DT and BC for an incoming electron from the normal metal are sketched in Fig. 2(b), in which a phonon is absorbed when the electron is transmitted to the superconductor. To the leading order in the tunneling rates, DT and BC are proportional to $\Gamma_{S}\Gamma_{N}$ and $\Gamma_{S}^{2}\Gamma_{N}^{2}$, respectively. The last term in Eq. (29) is proportional to the product of two transmission functions and can be interpreted as tunneling via an intermediate state [Cuevas et al., 1996] and scales as $\Gamma_{N}\Gamma_{S}^{2}$.

An example of the behavior of $\gamma_{\mathrm{NS}}^{s}=\sum_{\nu}\gamma_{\mathrm{NS}}^{s\nu}$ is shown in Fig. 5. By comparison of Fig. 5 with Fig. 4 we notice that the damping associated to the inelastic tunneling of quasiparticles is negligible when the dot’s level is well inside the gap $|\varepsilon_{0}|\ll\Delta$ and becomes relevant when $|\varepsilon_{0}|\sim\Delta$ whereas the damping associated to inelastic Andreev reflections dominates at low energy $|\varepsilon_{0}|\ll\Delta$. We will explain this behavior in the next section.

Finally, for completeness, we report the expression for the reflections involving the quasiparticles of the superconductor. This term is given by

[TABLE]

with each term referring to individual process of emission or absorption for electron/hole. They read

[TABLE]

and the reflection function

[TABLE]

The damping contribution associated to the quasiparticle reflections in the superconductor $\gamma_{\mathrm{SS}}$ is independent of the applied voltage. Moreover, for low temperatures compared to the superconducting gap $T\ll\Delta$, the population of quasiparticles above the gap is exponentially small such that we can disregard $\gamma_{\mathrm{SS}}$ for the rest of our discussion.

IV Results for the phonon occupation

Applying a bias voltage, the electron current drives the oscillator to a nonequilibrium state with phonon occupation

[TABLE]

In the limit $\gamma_{\textrm{tot}}\ll(\omega,\Gamma_{N},\Gamma_{S},T,eV)$ and separating the contribution of the intrinsic damping and the damping due to the charge-vibration interaction, the phonon occupation can be written as

[TABLE]

with the Bose function $n_{B}(\omega)=[\mathrm{exp}(\omega/T)-1]^{-1}$. The nonequilibrium occupation due to inelastic ARs reads [Stadler et al., 2016]

[TABLE]

whereas the nonequilibrium occupation $n_{\textrm{NS}}$ is given by

[TABLE]

with the Bose function shifted by $eV$ due to the tunneling of quasiparticles. Notice that the normal reflection (NR) and the quasiparticle reflections (SS) can drive the oscillator only to the thermal equilibrium.

For low temperature compared to the gap $T\ll\Delta$ and for the case $\gamma_{0}\ll\gamma_{\mathrm{NR}}$, we can approximate $\gamma_{\mathrm{NR}}+\gamma_{\textrm{SS}}+\gamma_{0}\simeq\gamma_{\mathrm{NR}}$. Moreover, in the high voltage limit $eV\gg T$ and for positive bias $eV>0$, we have essentially that electrons from the normal lead are moving to the right superconducting lead (or equivalently the holes move from the right to the left). In this case, we can approximate the expression for $\bar{n}$ as

[TABLE]

The Eq. (37) generalizes the main result in Ref. [Stadler et al., 2016] as it takes into account the effects of the quasiparticles on the nonequilibrium phonon occupation of the mechanical oscillator.

IV.1 Single mode phonon occupation

Figure 6 shows an example of the phonon occupation $\bar{n}$ as a function of $\varepsilon_{0}$ and source-drain voltage $eV$ for a single mode with frequency $\omega\ll\Delta$ and a temperature $T=10\,\omega$. As anticipated in the introduction, we have: (i) regions of ground state cooling $\bar{n}\ll 1$, (ii) regions of strong energy pumping into the resonator $\bar{n}\gg 1$ and (iii) regions of a mechanical instability in which the total electromechanical damping becomes negative $\gamma+\gamma_{0}<0$ and our perturbative method breaks down. In analogy to previous results in Ref. [Usmani et al., 2007], one expects that, going beyond the perturbative approach for the charge-vibration interaction and taking into account the anharmonicity of the resonator, such regions can host self-sustained oscillations.

The first striking observation in Fig. 6 is the weakly dependence of the phonon occupation $\bar{n}$ on the bias-voltage whereas $\bar{n}$ is strongly dependent on $\varepsilon_{0}$. In other words, the phonon occupation is determined by the bare dot’s level $\varepsilon_{0}$, or the gate voltage, whereas no additional features appear at $eV\sim\Delta$. This behavior can be explained as we assume a sharp resonant regime with a small coupling to the normal lead $\Gamma_{N}=0.1\omega\ll\Delta$. As discussed in Sec. II.1, the Andreev levels $E_{A}$ are formed within the superconducting gap whose position depends on $\varepsilon_{0}$ (and the coupling to the superconductor $\Gamma_{S}$) and whose small broadening is set by $\Gamma_{N}$. This translates into sharp peaks entering in the transmission functions in Eqs. (25), (26), and (29), and yield a saturation for the damping rates $\gamma_{\mathrm{AR}},\gamma_{\mathrm{NR}}$ and $\gamma_{\textrm{NS}}$ as increasing the bias voltage. Thus we divide the following discussion about the phonon occupation of Fig. 6 into the deep subgap regime with $|\varepsilon_{0}|\ll\Delta$, and a regime close to the gap $|\varepsilon_{0}|\sim\Delta$. This scale separation is valid for $\omega\ll\Delta$.

In the first region at low energy, $|\varepsilon_{0}|\ll\Delta$, inelastic AR dominates [see Fig. 4]. This subgap regime was discussed in Ref. Stadler et al., 2016 in which we showed that resonances in both amplitudes of the transmission function in Eq. (25) as well as interference between the two amplitudes result in: (i) an enhancement of the absorption rates compared to the emission rates for $\varepsilon_{0}<0$ leading to ground-state cooling, (ii) a large phonon occupation eventually approaching a mechanical instability for $\varepsilon_{0}>0$. Here, we extend the previous result showing that $\bar{n}$ is controlled by the gate voltage even for bias voltages larger than the superconducting gap $eV>\Delta$. If additionally $|\varepsilon_{0}|\ll\Delta$, we can neglect the contribution of quasiparticles since the energy levels of the Andreev bound states lie deep inside the gap and quasiparticles inelastic transport occurs out of resonance.

By contrast, the inelastic tunneling of quasiparticles strongly controls the phonon occupation when the gate voltage approaches the superconducting gap $|\varepsilon_{0}|\simeq\Delta$ [see Fig. 6]. Indeed, in this region, inelastic quasiparticle tunneling are the leading processes as shown in Fig. 5. These processes, similarly to the inelastic ARs, yield a similar behavior of the nonequilibirium occupation with alternating regions of ground-state cooling, energy “heating” and instability. The line separating the ground state cooling $\bar{n}\ll 1$ from the region $\bar{n}\gg 1$, is a narrow range in which the resonator is close to the thermal equilibrium. These results can be explained by considering the subtle interplay between the density of states $\rho_{d}(\varepsilon)$ of the quantum dot [see Fig. 3(a)], that acts as energy filter via the transmission coefficient $\mathcal{T}^{s}_{\nu}(\varepsilon)$ in Eq. (26), and the behavior of the superconducting density of states with an energy gap and a divergence at $\varepsilon=\Delta$. The basic idea is schematically reported in Fig. 7. Hereafter we discuss in detail the positive part $\varepsilon_{0}>0$ at $\varepsilon_{0}\simeq\Delta$. Similar arguments apply for the negative part $\varepsilon_{0}<0$ at $\varepsilon_{0}\simeq-\Delta$.

We start the discussion with the case shown in Fig. 7(a) when the Andreev states are well formed inside the dot and the energy level $\varepsilon_{0}<\Delta$. Tunneling from $N$ to $S$ with phonon absorption or phonon emission are possible for electrons with the appropriate energy threshold, namely $\varepsilon\geq\Delta-\omega$ or $\varepsilon\geq\Delta+\omega$, respectively. However, the density of states $\rho_{d}(\varepsilon)$ (and hence the transmission function $\mathcal{T}^{s}_{\nu}(\varepsilon)$) is sharply peaked at the Andreev levels such that one can assume that electron are essentially injected at energy $\varepsilon=|E_{A}|$. In this case shown in Fig. 7(a), the emission of a phonon is simply blocked due the superconducting gap.

Increasing $\varepsilon_{0}$, we approach the region $\varepsilon_{0}>\Delta$. Here, beyond Andreev levels inside the gap, a broadened peak appears in $\rho_{d}(\varepsilon)$

• and in the transmission $\mathcal{T}^{s}_{\nu}(\varepsilon)$ - above the gap, reminiscent of the bare dot’s level [see Fig. 3(a) and Fig. 7(b)]. In this case, beyond the fact that phonon emission is still blocked due the superconducting gap, the absorption rate is in addition enhanced due to the increase in $\rho_{d}(\varepsilon)$ and $\mathcal{T}^{s}_{\nu}(\varepsilon)$ at $\varepsilon>\Delta$.

Finally, further increasing $\varepsilon_{0}$, the peak of the functions ${\rho}_{S}(\varepsilon)$ or $\mathcal{T}^{s}_{\nu}(\varepsilon)$ moves far away from the gap edge such that the condition $\varepsilon_{0}=\Delta+\omega$ is reached. In this case, electrons from the normal lead N are again filtered by the dot at energy $\varepsilon=\varepsilon_{0}$ and phonon emission are greatly enhanced in this case due to the divergence of ${\rho}_{S}(\varepsilon)$ at the gap edge, viz. $\gamma^{-}_{\textrm{NS}}>\gamma^{+}_{\textrm{NS}}$, leading to the instability.

IV.2 Multimode phonon ground-state cooling

In this section, we illustrate that several nondegenerate mechanical modes can be cooled to the ground state in the N-QS-S system.

We generalized the Hamiltonian Eq. (1) considering several vibrational eigenmodes of the resonator. The frequencies of these vibrational modes are denoted by $\omega_{k}$ and the bosonic annihilation and creation operators for each mode by $\hat{b}_{k}^{\phantom{g}}$ and $\hat{b}_{k}^{\dagger}$. The Hamiltonian reads:

[TABLE]

To illustrate the multimode cooling, we assume a linear low-frequency spectrum by $\omega_{k}=k\omega$ with $k=1,2,\dots,N$ although this assumption is not essential. The calculation of the multimode phonon occupation can be limited by considering a finite number of modes, as we have a natural cut-off given by the temperature. For the case $T=10\omega$ we choose $N=7$ modes. High-frequency modes with $\omega\gtrsim T$ are close to the ground-state. The nonequilibrium value $\bar{n}_{k}$ for each mode is calculated by Eq. (34) from which we can plot the total mechanical energy defined by $E_{\mathrm{tot}}=\sum_{k=1}^{N}\omega_{k}\bar{n}_{k}$. Such a approach holds since, in the regime of weak charge-vibration coupling, we can neglect the intercorrelation effect between the different modes.

In Fig. 8 we show the multimode phonon occupation at bias voltage larger than the gap $eV>\Delta$ and as a function of $\varepsilon_{0}$. Similarly to the discussion in Sec. IV.1 of the phonon occupation for a single mode, we found that the total energy $E_{\mathrm{tot}}$ has a weak bias-voltage dependence and its behavior is controlled by the gate voltage $\varepsilon_{0}$.

As can be seen from Fig. 8, the thermal equilibrium value of the total mechanical energy $E_{\mathrm{tot}}\approx 57.2\omega$ is approached if the gate voltage is much larger than the gap, $|\varepsilon_{0}|\gg\Delta$. Figure 8 has a similar behavior as the single mode results $\bar{n}$ with features in the deep subgap region $|\varepsilon_{0}|\ll\Delta$ and close to the gap edges $|\varepsilon_{0}|\sim\Delta$. Around these ranges, the total energy shows dips below the thermal equilibrium value followed by sharp increases of the multimode phonon occupation.

When the gate voltage is close to $\varepsilon_{0}\simeq-\Delta$ and for $eV=3\Delta$, the total energy has a small dip (weak cooling) followed by a increase due to inelastic tunneling of quasiparticles. In this case, since the temperature is small compared to the superconductor gap $T\ll\Delta$, the states in the superconductor are almost filled and the phonon occupations of the several modes are only slightly reduced.

As discussed in a previous work in Ref. [Stadler et al., 2016], for negative $\varepsilon_{0}<0$ and $|\varepsilon_{0}|\ll\Delta$ , all the vibrational modes are strongly cooled by inelastic-vibration assisted Andreev reflections. Due to an interference mechanism, such a cooling is possible without matching a resonant condition between levels of the quantum dot and the vibrational frequency. Here we demonstrate that this mechanism is working even at bias voltage larger than the gap.

Interestingly, the phonon occupation of several modes is strongly reduced even closely below the superconducting gap, for $\varepsilon_{0}\simeq\Delta$, where the dominant processes are the inelastic tunneling of quasiparticles, similarly as discuss in the previous section.

V Results for the current

In this section, we calculate the current to the leading order in the charge-vibration coupling for a single mode with frequency $\omega$.

Using the spinor notation in Nambu space in the normal lead $\Psi_{k}^{\dagger}=(c_{k\uparrow}^{\dagger}\,c_{-k\downarrow})$ and the number operator $\hat{N}=\frac{1}{2}\sum_{k}\Psi_{k}^{\dagger}\hat{\tau}_{z}\Psi_{k}$, with the Pauli matrix $\hat{\tau}_{z}$, the current through the normal lead $I=-e\langle\mathrm{d}\hat{N}_{n}/\mathrm{d}t\rangle$ can be expressed as

[TABLE]

in which $\hat{\mathcal{G}}^{K}(\varepsilon)$ and $\hat{\mathcal{G}}^{R}(\varepsilon)$ refer to the dot Green’s functions and they are defined similarly as in Sec. II.1. They correspond to the exact Green’s function including the coupling to the leads and the charge-vibration interaction. To calculate the current in Eq. (39) we perform a perturbation expansion of these Green’s functions to the leading order in the charge-vibration coupling $\lambda$ [Egger and Gogolin, 2008,Stadler et al., 2015]. Hence the current can be decomposed into a elastic current, an elastic correction and an inelastic term as

[TABLE]

with

[TABLE]

and the self-energies associated to the charge-vibration interaction $(\sim\lambda^{2})$

[TABLE]

The elastic current $I_{0}$ shows a peak at $\varepsilon_{0}=0$ as a function of the gate voltage which is associated to the Andreev reflections. In the deep subgap regime, this current can be simplified to

[TABLE]

In addition, $I_{0}$ increases close to the gap edges at $\varepsilon_{0}\simeq|\Delta|$ of the gate voltage, where quasiparticle tunneling occurs and contributes to transport.

The term $I_{\mathrm{ec}}^{(\lambda^{2})}$ represents the corrections to the elastic part and shows no additional features with respect to $I_{0}$. On the other hand, the inelastic term of the current $I_{\mathrm{in}}^{(\lambda^{2})}$ has different features. This term can be approximated in the high-voltage limit as

[TABLE]

for $eV\gg T$. The inelastic current is directly related to the average phonon occupation. In the following, we discuss two different cases: First we consider the vibrational states at thermal equilibrium $\bar{n}\simeq n_{B}(\omega)$ and second the case of vibrational states out of equilibrium such that $\bar{n}$ is completely controlled by the charge passing through the dot.

As example of this comparison, we show in Fig. 9 the current as a function of $\varepsilon_{0}$ at large bias voltage $eV=10\Delta$ for a resonator at equilibrium (a) and in a nonequilibrium state (b).

In both cases, we observe that the inelastic current shows peaks which can be divided into contribution due to inelastic Andreev reflection and quasiparticle tunneling. Close to the zero-bias peak $\varepsilon_{0}=0$, inelastic-vibration assisted Andreev reflection at $\varepsilon_{0}=\pm\omega/2$ give rise to transport: The peak at $\varepsilon_{0}=-\omega/2$ is attributed to the absorption of a phonon whereas the peak at $\varepsilon_{0}=\omega/2$ to the emission of one phonon [Stadler et al., 2016].

Additionally to these peaks, we also found peaks close to the gap edges where inelastic quasiparticles tunneling give a contribution to transport at $E_{A}=\pm|\Delta\pm\omega|$, using the same argument discussed at the end of the Sec. IV.1 [see Fig. 7]. The peaks at $E_{A}=|\Delta-\omega|$ corresponds to the absorption of one phonon and $E_{A}=-|\Delta-\omega|$ to the emission of one phonon, namely the role of the negative and positive peaks are inverted compared to peaks associated to the ARs.

The difference behaviors of the peaks in Fig. 9(a) and Fig. 9(b) can be explained by the phonon occupation Fig.6 with the help of the expression of the inelastic current in Eq. (47).

When the vibration is at thermal equilibrium as in Fig. 9(a), both kind of the peaks due to the inelastic processes involving ARs and quasi-particles are almost of the same height. The peak due emission are slightly larger than the peaks due to emission of one phonon due factor $n_{B}(\omega)+1$ in front of emission processes compared to the factor $n_{B}(\omega)$ in front of the absorption processes.

By contrast, in Fig. 9(b) the resonator is in a nonequilibrium state and the paired sidepeaks are of different height: the peaks which are associated to the phonon absorption are suppressed compared to the peaks associated to the phonon emission. This result implies that the heights of the inelastic peaks lying in two different range of the gate voltage are correlated and they can be used to check if the resonator is in a nonequilibrium state. Using the same argument discussed in Ref. [Stadler et al., 2016], the ratio between the brightness of paired peaks - defined as the underlying area for each peak - can be used to extract information about the phonon occupation of the non equilibrated resonator.

VI Conclusion

We studied the phonon occupation and the quantum transport of a quantum dot embedded between a normal-conducting and a superconducting lead showing that the vibrational states of the resonator are controlled by applying a bias-voltage. The enhancement or suppression of the phonon occupation depends strongly on the resonant level of the quantum dot. When the energy of the resonant level is well inside the superconducting gap, inelastic vibration-assisted Andreev reflections drive the resonator to a nonequilibrium state. At gate-voltages close to the superconducting gap edge, inelastic tunneling of quasiparticles control the state of the resonators. As a result we obtained that inelastic vibration-assisted Andreev reflections and inelastic tunneling of quasiparticles cool the vibrational state of the resonator to the ground-state even for several vibrational modes. The current shows characteristic signatures of the inelastic tunneling processes and can be exploited to detect the resonator in a nonequilibrium state. Similar inelastic processes were observed in recent experiments [Schindele et al., 2014; Gramich et al., 2015, 2016, ] suggesting that our proposal in within the reach of current research.

Acknowledgements.

We acknowledge J. C. Cuevas and E. Scheer for interesting discussions and useful comments. This research was supported by the Zukunftskolleg of the University of Konstanz and by the DFG through the collaborative research center SFB 767.