How long does it take to form the Andreev quasiparticles ?

TL;DR

This paper investigates the transient dynamics of Andreev quasiparticle formation in a quantum dot system coupled to metallic and superconducting leads, providing analytic expressions for key time scales and effects of electron interactions.

Contribution

It offers the first analytic derivation of time-dependent charge, pairing amplitude, and currents in a quantum dot with proximity-induced superconductivity, including interaction effects.

Findings

Derived explicit formulas for charge occupancy and pairing amplitude over time.

Identified the characteristic time scale for quasiparticle formation.

Showed the interplay between Coulomb interactions and electron pairing effects.

Abstract

We study transient effects in a setup, where the quantum dot (QD) is abruptly sandwiched between the metallic and superconducting leads. Focusing on the proximity-induced electron pairing, manifested by the in-gap bound states, we determine characteristic time-scale needed for these quasiparticles to develop. In particular, we derive analytic expressions for (i) charge occupancy of the QD, (ii) amplitude of the induced electron pairing, and (iii) the transient currents under equilibrium and nonequilibrium conditions. We also investigate the correlation effects within the Hartree-Fock-Bogolubov approximation, revealing a competition between the Coulomb interactions and electron pairing.

Equations

$$\hat{H}=\sum _{\sigma }{\epsilon }_{\sigma }{\hat{d}}_{\sigma }^{†}{\hat{d}}_{\sigma }+U{\hat{n}}_{\uparrow }{\hat{n}}_{\downarrow }+\sum _{\beta }\left({\hat{H}}_{\beta }+{\hat{V}}_{\beta -QD}\right)$$

$$\displaystyle V_{{\bf k}({\bf q})}(t)=\left\{\begin{array}[]{cc}0&\hskip 14.22636pt\mbox{\rm for }t\leq 0\\ V_{{\bf k}({\bf q})}&\hskip 14.22636pt\mbox{\rm for }t>0,\end{array}\right.$$

$${\hat{d}}_{\sigma }(s)={\int }_0^{\infty }{e}^{-st}{\hat{d}}_{\sigma }(t)dt\equiv \mathcal{L}\left\{{\hat{d}}_{\sigma }(t)\right\}(s).$$

$$n_{\sigma }(t)=⟨{\mathcal{L}}^{-1}\left\{{\hat{d}}_{\sigma }^{†}(s)\right\}(t){\mathcal{L}}^{-1}\left\{{\hat{d}}_{\sigma }(s)\right\}(t)⟩,$$

$$n_{\uparrow }(t)$$

$$n_{\sigma }(t)=n_{\sigma }(0)+\left[1-n_{\sigma }(0)-n_{-\sigma }(0)\right]\text{mbox}{\sin }^2\left(\frac{\sqrt{\delta }}{2}t\right)\frac{{\Gamma }_S^2}{\delta }.$$

$$n_{\sigma }(t)={\cos }^2\left(\frac{{\Gamma }_S}{2}t\right)n_{\sigma }(0)+{\sin }^2\left(\frac{{\Gamma }_S}{2}t\right)\left[1-n_{-\sigma }(0)\right],$$

$$P(t)=\frac{{\Gamma }_S^2}{(E_1-E_2{)}^2+{\Gamma }_S^2}\text{mbox}{\sin }^2\left(\frac{t}{2}\sqrt{(E_1-E_2{)}^2+{\Gamma }_S^2}\right),$$

$$n_{\sigma }(t)$$

$$\chi (t)$$

$$j_{N\sigma }(t)=2\text{mbox}\mathrm{Im}\sum _{\mathbf{k}}{V}_{\mathbf{k}}⟨{\hat{d}}_{\sigma }^{†}(t){\hat{c}}_{\mathbf{k}\sigma }(t)⟩.$$

$${\hat{c}}_{\mathbf{k}\sigma }(t)={\hat{c}}_{\mathbf{k}\sigma }(0)e^{-i{\epsilon }_{\mathbf{k}\sigma }t}-i{\int }_0^{t}d{t}^{\prime }{V}_{\mathbf{k}}{e}^{-i{\epsilon }_{\mathbf{k}\sigma }(t-t^{\prime })}{\hat{d}}_{\sigma }(t^{\prime }),$$

$$j_{N\sigma }(t)$$

$$j_{S\sigma }(t)=\frac{{\Gamma }_S^2}{2\sqrt{\delta }}\text{mbox}\sin \left(\sqrt{\delta }t\right)e^{-{\Gamma }_{N}t}\left[1-\sum _{{\sigma }^{\prime }}{n}_{{\sigma }^{\prime }}(0)\right].$$

$$j_{S\sigma }(t)+j_{N\sigma }(t)=\frac{d}{dt}{n}_{\sigma }(t).$$

$$G_{\sigma }(\mu ,t)$$

$$G(\mu ,t)$$

$$G(\mu ,\infty )=\frac{{\Gamma }_N^2{\Gamma }_S^2}{4\left[{\left(\frac{{\Gamma }_N}{2}\right)}^2+{\mu }_{-}^2\right]\left[{\left(\frac{{\Gamma }_N}{2}\right)}^2+{\mu }_{+}^2\right]}.$$

$${\hat{d}}_{\uparrow }^{†}{\hat{d}}_{\uparrow }{\hat{d}}_{\downarrow }^{†}{\hat{d}}_{\downarrow }$$

$$(s+i{\epsilon }_{\uparrow }){\hat{d}}_{\uparrow }(s)$$

$$(s+i{\epsilon }_{\mathbf{k}}){\hat{c}}_{\mathbf{k}\uparrow }(s)$$

$$(s+i{\epsilon }_{\mathbf{q}}){\hat{c}}_{\mathbf{q}\uparrow }(s)$$

$$(s-i{\epsilon }_{\mathbf{q}}){\hat{c}}_{-\mathbf{q}\downarrow }^{†}(s)$$

$$(s-i{\epsilon }_{\downarrow }){\hat{d}}_{\downarrow }^{†}(s)$$

$$(s-i{\epsilon }_{\mathbf{k}}){\hat{c}}_{\mathbf{k}\downarrow }^{†}(s)$$

$$(s-i{\epsilon }_{\mathbf{q}}){\hat{c}}_{\mathbf{q}\downarrow }^{†}(s)$$

$$(s+i{\epsilon }_{\mathbf{q}}){\hat{c}}_{-\mathbf{q}\uparrow }(s)$$

$${\hat{d}}_{\uparrow }(s)M_{\uparrow }^{(+)}(s)$$

$${\hat{d}}_{\downarrow }^{†}(s)M_{\downarrow }^{(-)}(s)$$

$$M_{\sigma }^{(\pm )}(s)$$

$$K(s)$$

$$\hat{A}(s)$$

$$\hat{B}(s)$$

Full text

I Introduction

Quantum impurity attached to a superconducting bulk material absorbs the Cooper pairs, developing the quasiparticle states in its subgap spectrum $|\omega|\leq\Delta$, where $\Delta$ is the energy gap of superconducting reservoir Balatsky et al. (2006); Martín-Rodero and Levy Yeyati (2011). These bound Andreev (or Yu-Shiba-Rusinov) states have been observed in numerous STM studies, using impurities deposited on superconducting substrates Yazdani et al. (1997) and in tunneling experiments via quantum dots arranged in the Josephson De Franceschi et al. (2010), Andreev Deacon et al. (2010a) and more complex (multi-terminal) configurations Schindele et al. (2014); Gramich et al. (2017). Since measurements can be nowadays done with state-of-art precision probing the time-resolved properties, we address this issue here and determine some characteristic temporal scales of the in-gap quasiparticles.

Any abrupt change of the model parameters (quantum quench) is usually followed by a time-dependent thermalization of the many-body system , where continuum states play a prominent role Souto et al. (2017). Dynamics of these processes has been recently explored in the solid state and nanoscopic physics Eisert et al. (2015). From a practical point of view, especially useful could be nanoscopic heterostructures with the correlated quantum dot (QD) embedded between external (metallic, ferromagnetic or superconducting) leads which enable measurements of the transport properties under tunable nonequilibrium conditions Bidzhiev and Misguich (2017).

Transport phenomena through QD coupled between the normal or superconducting leads have so far explored mainly in the static cases. Since novel experimental methods allow to study the QDs subjected to voltage pulses or abrupt changes of the system parameters, it would be very desirable to calculate the time-dependent currents and their conductances. In particular, one can ask the question: how fast does the QD respond to an instantaneous perturbation. For this purpose analytical estimation of the transient oscillations and long-time (asymptotic) behaviour of the measurable quantities would be very useful. Some early theoretical works have investigated time-dependent transport via QD between the normal and superconducting leads Sun et al. (2000); Zhao and Wang (2001); Wei and Wang (2002); Cao et al. (2015), however, analytic results are hardly available. As regards the QD coupled to both normal leads, the transient current and charge occupancy have been determined for abrupt voltage pulses or after an instantaneous switching of constituent parts of the system Stefanucci and Almbladh (2004); Maciejko et al. (2006); F. M. Souza et al. (2007); Schmidt et al. (2008); Komnik (2008); Werner et al. (2009); Jin et al. (2010); Segal et al. (2011); Joho et al. (2012); Kulkarni et al. (2013); Albrecht et al. (2013a); Tuovinen et al. (2014); Odashima and Lewenkopf (2017).

Time-resolved techniques could provide an insight into the many-body effects. For instance, the pump-and-probe experiments Orenstein (2012) and the time-resolved ARPES Smallwood et al. (2016) have determined life-time of the Bogoliubov quasiparticles in the high temperature superconductors. Transient effects have been investigated in nanoscopic systems, considering mainly the quantum dots hybridized with the conducting (metallic) leads. There has been studied the time-scale needed for the Kondo peak to develop at the Fermi energy Nordlander et al. (1999), dynamical correlations in electronic transport via the quantum dots Michałek and Bułka (2009), or oscillatory behavior in the charge transport through the molecular junctions Seoane Souto et al. (2015).

Dynamical phenomena of the quantum dots attached to superconducting bulk reservoirs have been studied much less intensively. There has been analyzed: photon-assisted Andreev tunneling Sun et al. (1999), response time on a step-like pulse Xing et al. (2007), temporal dependence of the multiple Andreev reflections Stefanucci et al. (2010), time-dependent sequential tunneling Contreras-Pulido et al. (2012), transient effects caused by an oscillating level Albrecht et al. (2013b), time-dependent bias K.J. Pototzky (2014), the waiting time distributions in nonequilibrium transport Rajabi et al. (2013); Michałek et al. (2018), the short-time counting statistics Stegmann and König (2016), metastable configurations of the Andreev bound states in a phase-biased Josephson junction Souto et al. (2016, 2017) etc. None of these studies, however, addressed the time-scale typical for development of the subgap quasiparticle states in a setup, comprising the quantum dot (QD) coupled to the normal lead (N) on one side and to the isotropic ($s$-wave) superconductor (S) on the other side. Our present study reveals, that a continuous electronic spectrum of the metallic lead enables a relaxation of the Andreev states, whereas the superconducting electrode induces the (damped) quantum oscillations with a period sensitive to the energies of the in-gap quasiparticles. In what follows we evaluate the time-scale at which such Andreev quasiparticle start to form and another one, when they are finally established.

The paper is organized as follows. In Sec. II we introduce the microscopic model and discuss the method for the time-dependent phenomena. Sec. III presents analytical results for the uncorrelated quantum dot, such as: (i) charge occupancy, (ii) complex order parameter, and (iii) charge current for the unbiased and biased heterojunction. In Sec. IV we discuss the correlation effects and finally in Sec. V we summarize the main results.

II Microscopic model

For description of the N-QD-S heterostructure (see Fig. 1) we use the single impurity Anderson Hamiltonian

[TABLE]

where $\beta$ refers to the normal ($N$) and superconducting ($S$) electrodes, respectively. As usually $\hat{d}_{\sigma}$ ($\hat{d}^{\dagger}_{\sigma}$) is the annihilation (creation) operator for the quantum dot (QD) electron with spin $\sigma$ and energy $\varepsilon_{\sigma}$. Potential of the Coulomb repulsion between the opposite spin electrons is denoted by $U$. We treat the external metallic lead as free fermion gas $\hat{H}_{N}\!=\!\sum_{{\bf k},\sigma}\varepsilon_{\bf k}\hat{c}_{{\bf k}\sigma}^{\dagger}\hat{c}_{{\bf k}\sigma}$, and describe the isotropic superconductor by the BCS model $\hat{H}_{S}\!=\!\sum_{{\bf q},\sigma}\varepsilon_{{\bf q}}\hat{c}_{{\bf q}\sigma}^{\dagger}\hat{c}_{{\bf q}\sigma}\!-\!\sum_{\bf q}\Delta\left(\hat{c}_{{\bf q}\uparrow}^{\dagger}\hat{c}_{-{\bf q}\downarrow}^{\dagger}+\hat{c}_{-{\bf q}\downarrow}\hat{c}_{{\bf q}\uparrow}\right)$, where $\varepsilon_{{\bf k}({\bf q})}$ is the energy measured from the chemical potential $\mu_{N(S)}$, and $\Delta$ denotes the superconducting energy gap. Hybridization between the QD electrons and the metallic lead is given by $\hat{V}_{N-QD}=\sum_{{\bf k},\sigma}\left(V_{{\bf k}}\;\hat{d}_{\sigma}^{\dagger}\hat{c}_{{\bf k}\sigma}+\mbox{\rm h.c.}\right)$ and $\hat{V}_{S-QD}$ can be expressed by interchanging ${\bf k}\leftrightarrow{\bf q}$.

Since our study refers to the subgap quasiparticle states, we assume the constant couplings $\Gamma_{N(S)}=2\pi\sum_{{\bf k}({\bf q})}|V_{{\bf k}({\bf q})}|^{2}\;\delta(\omega\!-\!\varepsilon_{{\bf k}({\bf q})})$. In the deep subgap regime $|\omega|\ll\Delta$ (so called, superconducting atomic limit) the coupling $\Gamma_{S}/2$ can be regarded as a qualitative measure of the induced pairing potential, whereas $\Gamma_{N}$ controls the inverse life-time of the in-gap quasiparticles. As we shall see, both these couplings play important (though quite different) role in transient phenomena.

We assume that all three constituents of the N-QD-S heterostructure are disconnected from each other until $t\leq 0$. Let us impose the external (N, S) reservoirs to be suddenly coupled to the quantum dot

[TABLE]

inducing the transient effects. Later on, we shall relax this assumption. Our problem resembles the Wiener-Hopf method Janiš (1997) applied earlier in the studies of X-ray absorption and emission of metals Nozierés and De Dominicis (1969).

In what follows, we explore the time-dependence of physical observables $\hat{O}$, based on the Heisenberg equation of motion $i\hbar\frac{d}{dt}\hat{O}=\left[\hat{O},\hat{H}\right]$. In particular, we shall investigate expectation values of the QD occupancy $\langle\hat{d}_{\sigma}^{\dagger}(t)\hat{d}_{\sigma}(t)\rangle$, the induced on-dot pairing $\langle\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\rangle$, and the transient charge currents flowing between the QD and external electrodes (both under equilibrium and nonequilibrium condictions).

Our strategy is based on the following three steps. First, we formulate the differential equations of motion for the annihilation $\hat{d}_{\sigma}(t)$ and creation $\hat{d}_{\sigma}^{\dagger}(t)$ operators of QD and similar ones for the mobile electrons $\hat{c}_{{\bf k}({\bf q})\sigma}(t)$ and $\hat{c}_{{\bf k}({\bf q})\sigma}^{(\dagger)}(t)$, respectively. Next, we solve them using the Laplace transformations, e.g. for $\hat{d}_{\sigma}(t)$ we denote

[TABLE]

For the uncorrelated QD the analytical expressions for $\hat{d}_{\sigma}(s)$ and $\hat{d}_{\sigma}^{\dagger}(s)$ can be obtained (see Appendix A). Finally, using the corresponding inverse Laplace transforms, we compute the time-dependent expectation values of a the QD occupancy, the QD pair amplitude and currents flowing between QD and both leads. For example QD occupancy $n_{\sigma}(t)\equiv\langle\hat{d}_{\sigma}^{\dagger}(t)\hat{d}_{\sigma}(t)\rangle$ is given by

[TABLE]

where ${\cal{L}}^{-1}\left\{\hat{d}_{\sigma}(s)\right\}(t)$ stands for the inverse Laplace transform of $\hat{d}_{\sigma}(s)$.

In our calculations we make use of the wide-band limit approximation ($\Gamma_{\beta}=$const) and set $e=\hbar=k_{B}\equiv 1$, so that all energies, currents and time are expressed in units of $\Gamma_{S}$, $e\Gamma_{S}/\hbar$ and $\hbar/\Gamma_{S}$, respectively. We also treat the chemical potential $\mu_{S}=0$ as a convenient reference energy point and perform the calculations for zero temperature. For experimentally available value $\Gamma_{S}\sim 200\mu eV$ Deacon et al. (2010b); Pillet et al. (2013); Eichler et al. (2007), the typical period of transient oscillations would be $\sim 3.3psec$.

III Uncorrelated QD case

We start by addressing the transient effects of the uncorrelated quantum dot ($U=0$), focusing on the superconducting atomic limit ($\Delta=\infty$) for which analytical expressions can be obtained. More general cosiderations are presented in Appendix A.

III.1 Time-dependent QD charge

Let us inspect the time-dependent occupancy $n_{\sigma}(t)$ driven by an abrupt coupling of the QD to both external leads. This quantity, defined in Eq. (6), can be determined explicitely for arbitratry $\Delta$ (derivation is presented in Appendix A). Here we shall consider the formula (A.2) simplified for the superconducting atomic limit

[TABLE]

where $f_{N}(\omega)$ is the Fermi-Dirac distribution function of the normal lead and parameters $s_{1},s_{2},s_{3},s_{4}$ and $\delta$ are presented in Eq. (35). The occupancy $n_{\downarrow}(t)$ can be obtained from the same expression (III.1) upon replacing the set $(s_{1},s_{2},s_{3},s_{4})$ by $(s_{3},s_{4},s_{1},s_{2})$. Expressions given in the second and third lines of Eq. (III.1) could be presented in the compact analytical form in the case $\varepsilon_{\sigma}=0$ (see Eqs. (A.2-A.2)). Otherwise they are rather lengthy (even though accessible), therefore we skip them.

Another simplification of Eq. (III.1) is possible upon neglecting the normal lead ($\Gamma_{N}=0$). QD occupancy is then characterized by non-vanishing quantum oscillations

[TABLE]

For $\varepsilon_{\sigma}\!=\!0$ Eq. (8) reduces to

[TABLE]

implying the period of transient oscillations $T=2\pi/\Gamma_{S}$, except of the initial conditions $n_{\sigma}(0)=1$ and $n_{-{\sigma}}(0)=0$ when the QD occupancy is preserved.

The formula (9), obtained in the case $\Gamma_{N}=0$, resembles the Rabi oscillations of a typical two-level quantum system. Indeed, the proximitized QD is fully equivalent to such scenario. To prove it, let us consider the effective Hamiltonian $\hat{H}=\sum_{\sigma}\varepsilon_{\sigma}\hat{n}_{\sigma}+\frac{\Gamma_{S}}{2}\left(\hat{d}^{\dagger}_{\uparrow}\hat{d}^{\dagger}_{\downarrow}+\mbox{\rm h.c.}\right)$, assuming that at $t=0$ the QD is empty $n_{\uparrow}(0)=0=n_{\downarrow}(0)$. For arbitrary time $t>0$ we can calculate the probability $P(t)$ of finding the QD in the doubly occupied configuration $n_{\uparrow}(t)=1=n_{\downarrow}(t)$ within the standard treatment of a two-level system Cohen-Tannoudji et al. (1977). This probability is given by

[TABLE]

where $E_{1}=0$ and $E_{2}=\varepsilon_{\uparrow}+\varepsilon_{\downarrow}$ are the energies of empty and doubly occupied configurations, respectively. This result exactly reproduces our expression (9).

For the QD suddenly coupled to both the normal and superconducting leads ($\Gamma_{N,S}\neq 0$) such oscillations become damped (see Fig. 2). This effect comes partly from the exponential factor $\exp(-\Gamma_{N}t)$ appearing in front of the first term in Eq. (III.1) and partly from the second and third contributions. This can be illustrated, by considering the case $\varepsilon_{\sigma}=0$, $\mu_{N}=0$, for which Eq. (III.1) implies

[TABLE]

Under such circumstances, the QD occupancy approaches asymptotically a half-filling, $\lim_{t\rightarrow\infty}n_{\sigma}(t)=\frac{1}{2}$. Fig. 2 displays $n_{\uparrow}(t)$ obtained in absence of external voltage for several values of $\Gamma_{N}$, assuming $\varepsilon_{\sigma}=0$ and $n_{\sigma}(0)=0$ for both spins. The quantum oscillations occur with a period $2\pi/\Gamma_{S}$ and their damping is governed by the envelope function $e^{-\Gamma_{N}t}$ indicating, that a continuous spectrum of the metallic lead is responsible for the relaxation processes. For a weak enough coupling $\Gamma_{N}$ these oscillations could indirectly probe the dynamical transitions between the subgap bound states, as recently emphasized by J. Gramich et al [Gramich et al., 2017].

Fig. 3 shows the QD occupancies obtained for several initial conditions, assuming $\mu_{N}=\mu_{S}=0$ and $\varepsilon_{\sigma}=0$. The case $n_{\downarrow}(0)=0=n_{\uparrow}(0)$ allows quantum oscillations between two eigenstates of the proximitized QD which are damped due to coupling to the normal lead (see Fig. 2). For the initial condition $n_{\sigma}(0)=1$, $n_{-\sigma}(0)=0$, the transient effects are completely different. The first term in Eq. (III.1) for $(n_{\uparrow}(0),n_{\downarrow}(0))=(1,0)$ or $(0,1)$ equals $e^{-\Gamma_{N}t}$ or vanishes and together with the last term they yield $\frac{1}{2}(1-e^{-\Gamma_{N}t})$ - see the upper curve in Fig. 3 or $\frac{1}{2}(1+e^{-\Gamma_{N}t})$ - the lower curve, respectively. This stems from the fact that proximity-induced pairing affects only the empty and doubly occupied configurations and it is inefficient in the case considered here. In consequence the quantum oscillations are absent and the QD occupancy exponentially evolves towards a halfilling. Let us also remark, that for $\Gamma_{S}=0$ Eq. (III.1) simplifies to the standard formula obtained by the non-equilibrium Green’s function method Jauho et al. (1994) (see Eq. 41).

III.2 Development of the proximity effect

Occupancy of the QD only indirectly tells us about emergence of the subgap bound states. To get some insight into the superconducting proximity effect we shall study here the time evolution of the order parameter $\chi(t)=\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$. The general formula is explicitly given by Eq. (A.3). Expressing its first two terms (which depend on the initial QD occupancy) the pair correlation function can be written as

[TABLE]

where $\Phi_{\uparrow}$ is given by Eq. (A.4). In Appendix A we show, that for $\mu_{N}=0$ the real part of $\Phi_{\uparrow}$ vanishes. Let us next analyze Eq. (12) for different initial conditions and values of the QD energy levels. For $n_{\sigma}(0)=0$, $n_{-\sigma}(0)=1$ and $\mu_{N}=0$ the function $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ is real and non-oscillating in time and is equal to $-\frac{\Gamma_{N}\Gamma_{S}}{4\pi}\mbox{\rm Im}\Phi_{\uparrow}$, regardless of $\varepsilon_{\sigma}$. However, for $\mu_{N}\neq 0$ also imaginary part of $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ equals $-\frac{\Gamma_{N}\Gamma_{S}}{4\pi}\mbox{\rm Re}\Phi_{\uparrow}$ and is non-oscillating function. For the initial conditions $(n_{\sigma}(0),n_{-\sigma}(0)=0)=(0,0)$ or $(1,1)$ the picture is completely different. Depending on the value of $\varepsilon_{\uparrow}+\varepsilon_{\downarrow}$ the real part of $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ oscillates for $\varepsilon_{\uparrow}+\varepsilon_{\downarrow}=0$ or is a smooth function of time for $\varepsilon_{\uparrow}+\varepsilon_{\downarrow}\neq 0$. Simultaneously, the imaginary part of the QD on-dot pairing oscillates irrespective of $\varepsilon_{\sigma}$. The oscillatory part of $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ are dumped via $e^{-\Gamma_{N}t}$ factor, emphasizing the crucial role of continuum states of the normal electrode in relaxation processes.

In Fig. 4 we show the imaginary part of the on-dot pairing $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ assuming the initial QD occupancy $n_{\sigma}(0)=0$. Period $T$ of the damped quantum oscillations depends on the excitation energy between the subgap Andreev quasiparticles Gramich et al. (2017) via $T=2\pi/\sqrt{\left(\varepsilon_{\downarrow}+\varepsilon_{\uparrow}\right)^{2}+\Gamma_{S}^{2}}$. For $\mu_{N}=0$ these oscillations are related to the transient current $j_{S\sigma}(t)$ flowing between the proximitized QD and the superconducting lead (see Sec. III.3) in analogy to the Josephson junction comprising two superconducting pieces, differing in phase of the order parameter. On the other hand, the real part (Fig. 5) evolves monotonously to its asymptotic value, except of one particular case $\Gamma_{N}=0$, when the real part of $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ vanishes.

III.3 Transient currents for unbiased system

So far we have discussed the quantities which are important, but unfortunately they are not directly accessible experimentally. Let us now consider the measurable currents $j_{N\sigma}(t)$ and $j_{S\sigma}(t)$, flowing from the QD to the external leads. Formally the transient current is defined by $j_{\beta\sigma}(t)=\left<\frac{d\hat{N}_{\beta}(t)}{dt}\right>$, where $\hat{N}_{\beta}(t)$ counts the total number of electrons in electrode $\beta=N,S$. For instance $j_{N\sigma}(t)$ simplifies to the standard formula Jauho et al. (1994)

[TABLE]

Assuming the energies of itinerant electrons to be static $\varepsilon_{{\bf k}\sigma}(t)=\varepsilon_{{\bf k}\sigma}$ one obtains

[TABLE]

and within the wide-band-limit approximation it yields

[TABLE]

Finally, inserting the time-dependent operator $\hat{d}^{\dagger}_{\sigma}(t)$ [Eq. (31) to Eq. (13)] we obtain

[TABLE]

To compute the transient current of opposite spin electrons, $j_{N-\sigma}(t)$, one should replace the set of auxiliary parameters $(s_{1},s_{2},s_{3},s_{4})$ by the following one $(s_{3},s_{4},s_{1},s_{2})$. In particular, for $\varepsilon_{\sigma}=0$ we get

[TABLE]

where $\omega_{\pm}=\frac{\Gamma_{S}}{2}\pm\omega$. In absence of the superconducting lead this formula is identical with the result obtained by means of the nonequilibrium Green’s function method.

In Fig. 6 we present transient behavior of the current $j_{N\uparrow}(t)$ induced by an abrupt coupling of the QD to both external leads for $\mu_{N}=\mu_{S}=0$ (i.e. without any bias). Similarly to the time-dependent QD occupancy (Fig. 2) we observe the quantum oscillations of the period $2\pi/\Gamma_{S}$ exponentially decaying with the envelope coefficient $e^{-\Gamma_{N}t}$. Large value of the current at $t=0^{+}$ is a consequence of the abrupt switching (4). One may ask whether this instantenous switching could be realistic in experimental situations. To check if any smooth (gradual) coupling process would affect our main conclusions we have computed the transient currents, assuming the sinusoidal switching profile $V_{{\bf k},{\bf q}}(t)=\frac{V_{{\bf k},{\bf q}}}{2}\left(\sin{(\pi\Gamma_{N}t-\pi/2)}+1\right)$ for $0<t\leq 1/\Gamma_{N}$ and keeping constant value $V_{{\bf k},{\bf q}}$ for $t>1/\Gamma_{N}$. We have solved this problem numerically and present some representative results (for $\Gamma_{N}/\Gamma_{S}=0.1$) in the inset in Fig. 6. We noticed, that for $t>1/\Gamma_{N}$ all the time-dependent quantities are not particularly affected. The only difference (in comparison to the abrut coupling) is in the early time region $0<t<1/\Gamma_{N}$. For instance, the transient current $j_{N\uparrow}(t)$ smoothely evolves from zero to its asymptotic behavior with the same period of quantum oscillations.

In similar steps, we have also determined the transient current $j_{S\sigma}(t)=2\;\mbox{\rm Im}\sum_{\bf k}V_{\bf q}\left<\hat{d}^{\dagger}_{\sigma}(t)\hat{c}_{{\bf q}\sigma}(t)\right>$. Effective quasiparticles in superconductors are represented by a coherent superposition of the particle and hole degrees of freedom, so for this reason the time-dependent operator $\hat{c}_{{\bf q}\sigma}(t)$ consists of four contributions (see Eq. 33). Final expression for $j_{S\sigma}(t)$ becomes rather lengthy, therefore we present it in Appendix A.4. However, in absence of external voltage the current (A.4) simplifies to

[TABLE]

When the QD is initially empty/full the transient current $j_{S\sigma}=\pm\frac{{\Gamma}^{2}_{S}}{2\sqrt{\delta}}\mbox{\rm sin}\left(\sqrt{\delta}t\right)e^{-\Gamma_{N}t}$ reveals the damped oscillations. Contrary to this behavior, for the initial occupancies $n_{\sigma}(0)=0$ and $n_{-{\sigma}}(0)=1$ the current (18) vanishes. We assign this feature to inefficiency of the proximity effect whenever the QD is singly occupied, because electron pairing operates only by mixing the empty with the doubly occupied QD configurations. Initial conditions have thus important influence on transient phenomena.

Furthermore, Eq. A.3 for $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$ and Eq. A.4 imply the exact relationship $j_{S\sigma}(t)=-\Gamma_{S}\mbox{\rm Im}\left<\hat{d}_{\downarrow}\hat{d}_{\uparrow}\right>$ which is popular in considerations of charge transport through Josephson junctions Žonda et al. (2015). The transient current $j_{S\sigma}(t)$ can hence be simply inferred from Fig. 4. At this level it is important to remark, that the charge conservation of our system is correctly satisfied, i.e.

[TABLE]

III.4 Transient currents for biased system

We have seen so far, that time-dependent QD occupancy and transient currents provide indirect information about the subgap quasiparticle energies and dynamical transitions between them. In absence of any voltage ($\mu_{N}=\mu_{S}=0$) these transient currents finally vanish, with a rate dependent on the relaxation processes caused by the coupling $\Gamma_{N}$ with a continuum of metallic lead. From the practical point of view, much more convenient way for probing the time-scales characteristic for the Andreev/Shiba quasiparticles could be provided by transient properties of the biased system $\mu_{N}\neq\mu_{S}$. Following the steps discussed in previous Sec. III.3 we shall study here the time-dependent differential conductance $G_{\sigma}(\mu,t)\equiv\frac{d}{d\mu}j_{N\sigma}(t)$ as a function of external voltage $\mu\equiv\mu_{N}$ (throughout this work the superconducting lead is assumed to be grounded $\mu_{S}\!=\!0$). At zero temperature Eq. (16) implies

[TABLE]

where the conductance is expressed in units of $\frac{2e^{2}}{h}$. Expression for $G_{\downarrow}(\mu,t)$ can be obtained by the replacement $(s_{1},s_{2},s_{3},s_{4})\rightarrow(s_{3},s_{4},s_{1},s_{2})$. Using the corresponding inverse Laplace transforms we find (for $\varepsilon_{\sigma}=0$, $G_{\uparrow}=G_{\downarrow}=G$):

[TABLE]

where $F_{1}(\mu,t)$ and $F_{2}(\mu,t)$ are given in Eqs. (A.2,A.2), and $\mu_{+/-}=\mu\pm\Gamma_{S}/2$. In the steady limit, $t\rightarrow\infty$ and for $\varepsilon_{\sigma}=0$, keeping only terms that survive at late times, we obtain the expression identical with the result derived for the same setup within the Büttiker-Landauer approach Domański and Donabidowicz (2008)

[TABLE]

For $\Gamma_{S}\gg\Gamma_{N}$ the local extrema of this expression occur at $\mu=\pm\frac{\Gamma_{S}}{2}$ and they correspond to the energies of subgap bound states. For an arbitrary set of model parameters such information is encoded in Eq. (III.4) which quantitatively specifies development of the in-gap states driven by the sudden switching at $t=0$. In Fig. 7 we present the differential conductance obtained numerically for $\Gamma_{N}/\Gamma_{S}=0.1$ and $0.7$. Let us notice, that differential conductance approaches its steady-limit shape $G_{\uparrow}(\mu,t=\infty)$ characterized by two Lorentzian quasiparticle peaks centered at $\sim\pm\frac{\Gamma_{S}}{2}$. Their broadening $\Gamma_{N}$ is related to the inverse life-time.

More careful examination of $G_{\uparrow}(\mu,t)$ indicates, that development of the subgap quasiparticles proceeds in three steps with typical time-scales $\tau_{1}$, $\tau_{2}$ and $\tau_{3}$, as can be deduced from Fig. 8 and 7. in Fig. 8 we show how the position of the quasiparticle maxima develops in time for different $\Gamma_{N}$. At $t\!=\!\tau_{1}$ there emerge two maxima from the single broad structure where $\tau_{1}$ changes approximately from 5 (for $\Gamma_{N}=0.2$) up to 10 (for $\Gamma_{N}=0.9$) units of time. These maxima move rapidly (essentially during 1-2 units) from $\mu=0$ up to some value of $\mu$ which depends on $\Gamma_{N}$. Next, the position of the quasiparticle peaks evolve continuously to their steady limit position $\mu=\pm\sqrt{\Gamma_{S}^{2}-\Gamma_{N}^{2}}$ with $\tau_{2}$ approximately changing from 15 (for $\Gamma_{N}=0.9$) up to 30 (for $\Gamma_{N}=0.2$) units of time. Finally, the asymptotic quasiparticle feature is achived with the evolve function $1-\mbox{\rm exp}\left(-t/\tau_{f}\right)$ where $\tau_{f}=2/\Gamma_{N}$, see Eqs. (III.4,A.2,A.2) where the terms proportional to ${\rm exp}\left(-\Gamma_{N}t/2\right)$ are responsible for such asymptotic behaviour. We also clearly see, that near the quasiparticle peaks the total differential conductance $\sum_{\sigma}G(\mu,t\rightarrow\infty)$ acquires its optimal value $4e^{2}/h$ known from the previous studies (see e.g. the Ref. Martín-Rodero and Levy Yeyati, 2011).

IV Correlation effects

Local repulsive interactions $U\hat{n}_{\uparrow}\hat{n}_{\downarrow}$ compete with the proximity-induced electron pairing. This issue has been addressed in the steady limit by numerous methods Martín-Rodero and Levy Yeyati (2011). In particular, it has been shown Bauer et al. (2007) that effective pairing (manifested by the in-gap states) is predominantly sensitive to the ratio $U/\Gamma_{S}$ and depends on the energy level $\varepsilon_{\sigma}$. Various experimental realizations of the correlated quantum dot in N-QD-S geometry Deacon et al. (2010b); Lee et al. (2012); Pillet et al. (2013); Žitko et al. (2015) indicated that the Coulomb potential $U$ safely exceeds (at least one order of magnitude) the superconducting energy gap $\Delta$. Under such circumstances the correlation effects show up in the subgap regime $|\omega|<\Delta$ merely by a quantum phase transition (or crossover) from the spinless (BCS-type) state $u\left|0\right>+v\left|\uparrow\downarrow\right>$ to the spinful (singly occupied) configuration $\left|\sigma\right>$. This changeover occurs upon increasing the ratio $U/\Gamma_{S}$ and above some critical value of the Coulomb potential $U_{cr}$ there can be observed the subgap Kondo effect (even in the superconducting atomic limit) Žitko et al. (2015); Domański et al. (2016). We shall briefly analyze some correlation effects, focusing on the transient effects.

IV.1 Competition between pairing and correlations

The aforementioned quantum phase transition can be qualitatively captured already within the lowest order (Hartree-Fock-Bogoliubov) decoupling scheme

[TABLE]

Using this approximation (23) one can incorporate the Hartree-Fock terms into the renormalized energy level $\tilde{\epsilon}_{\sigma}\equiv\epsilon_{\sigma}+Un_{-\sigma}(t)$, whereas the anomalous (pair source and drain) terms rescale the effective pairing potential $\tilde{\Gamma}_{S}/2\equiv\Gamma_{S}/2+U\chi(t)$. This decoupling procedure (23) can give a crossing of the subgap quasiparticle energies at some critical ratio $U/\Gamma_{S}$, dependent also on $\varepsilon_{\sigma}$. In the Josephson junctions such effect would cause a reversal of the d.c. tunneling current, so called, $0-\pi$ transition Žonda et al. (2015); Domański et al. (2017). In our N-QD-S heterostructure its influence is noticeable, but rather less spectacular.

Analytical determination of the dynamical observables (discussed in Sec. III) is unfortunately not feasible in the present case, because the renormalized energy level $\tilde{\varepsilon}_{\sigma}(t)$ and effective pairing potential $\tilde{\Gamma}_{S}(t)$ are time-dependent in nonexplicit way and the method used in the previous section is useful only for consideration of systems with constant QD energy levels and couplings with the leads. Therefore, in what follows, we consider the Coulomb repulsion in the system of the proximitized QD coupled only to the normal lead, applying the Hartree-Fock-Bogoliubov approximation (23). We have computed numerically $n_{\sigma}(t)$, $\langle\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\rangle$ and $j_{N\sigma}(t)$, solving the closed set of differential equations for time-dependent functions $n_{\sigma}(t)$ and $\langle\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\rangle$, respectively (see Appendix B). At intermediate steps we had to compute additionally the expectation values $\langle\hat{d}^{\dagger}_{\sigma}(t)c_{{\bf k}\sigma}(0)\rangle$ and $\langle\hat{d}_{\sigma}(t)\hat{c}_{{\bf k}-\sigma}(0)\rangle$. All these quantities have been determined within the Runge-Kutta numerical algorithm.

Fig. 9 displays influence of the Coulomb potential $U$ on the induced order parameter $\chi(t)$ for the unbiased system. The imaginary part, which is strictly related to the transient current, exhibits the damped quantum oscillations. Their period and amplitude are substantially suppressed by the Coulomb potential. We assign this fact to a competition between the on-dot pairing and local Coulomb repulsion. The real part of $\chi(t)$ is characterized by the same quantum oscillations. The asymptotic value of the complex order parameter $\chi(t\rightarrow\infty)$ with respect to the Coulomb potential $U$ is shown in Fig. 10 for $\varepsilon_{\sigma}=0$, $\Gamma_{N}/\Gamma_{S}=0.2$. Such monotoneously decreasing Re$\chi(\infty)$ confirms a competing relationship between the on-dot pairing and the local repulsion.

In Fig. 11 we show influence of the Coulomb potential $U$ on the QD occupancy $n_{\uparrow}(t)$. Besides the quantum oscillations, similar to the ones observed in the complex order parameter (Fig. 9), we notice a partial reduction of the QD charge upon increasing $U$. Apparently this is caused by the Hartree term $Un_{-\sigma}(t)$, that lifts the renormalized QD level $\tilde{\varepsilon}_{\sigma}(t)$. In Appendix C we briefly discuss the time-dependent subgap Kondo effect.

V Summary

We have investigated transient effects driven by a sudden coupling of the quantum dot to the metallic and superconducting leads. Our study has revealed a gradual buildup of the subgap Andreev quasiparticle states, which is controlled by the coupling $\Gamma_{N}$ to a continuous spectrum of the metallic lead. Depending on the initial quantum dot occupancy, we have also found the damped quantum oscillations of the charge occupancy $n_{\sigma}(t)$, the complex order parameter $\chi(t)$ and the transient currents $j_{N\sigma}(t)$, $j_{S\sigma}(t)$. Period of these oscillations would be sensitive to the Andreev quasiparticle energies which can be indirectly controlled via a coupling $\Gamma_{S}$ to the superconducting reservoir.

Analogous effects (relaxation and quantum oscillations) have been recently reported in Refs [Souto et al., 2017] and [Souto et al., 2016] in studies of the metastable subgap states for the Josephson junction, considering finite value of the superconducting gap. We estimate that in realistic systems, where $\Gamma_{S}\sim 0.2$ meV, the period of quantum oscillations would be a fraction of nanoseconds or in picoseconds regime (hence should be empirically detectable). Buildup of the subgap Andreev quasiparticle states is expected to be formed in N-QD-S junctions on much longer time-scale, corresponding to a microsecond regime. Our estimations seem to reliable, when comparing them with dynamical transitions between the subgap bound states of nanotubes Gramich et al. (2017) and parity switchings observed in the superconducting atomic contacts Janvier et al. (2015).

We also addressed the correlation effects by means of the Hartree-Fock-Bogoliubov approximation, revealing that the repulsive Coulomb potential $U$ suppresses the proximity-induced electron pairing. We have explored some time-dependent signatures of this competition. In particular, we have found that $\Gamma_{N}$ controls the rate at which the stationary limit behavior is achieved, whereas period of the damped quantum oscillations is dependent on the Coulomb potential due to its influence of the Andreev quasiparticle energies Bauer et al. (2007).

Finally, we have tried to evaluate the characteristic time-scale $\tau_{K}$ needed for the subgap Kondo effect to develop. Upon approaching the quantum phase transition from the (spinful) doublet side we predict the strong reduction of this scale $\tau_{K}$, originating from a subtle interplay between the induced on-dot pairing and the Coulomb repulsion Žitko et al. (2015); Domański et al. (2016). We hope that such variety of dynamical effects of the proximitized quantum dots could be verified experimentally.

Acknowledgements.

We acknowledge instructive discussions with V. Janiš and thank A. Baumgartner for useful remarks on observability of the transient effects in multi-terminal heterostructures. We also kindly thank T. Kwapiński for technical assistance. This work is supported by the National Science Centre (Poland) through the grant DEC-2014/13/B/ST3/04451 (TD).

Appendix A

In this Appendix we derive the Laplace transforms $\hat{d}_{\sigma}(s)$ and $\hat{c}_{{\bf q}\sigma}(s)$ which are needed for calculating the statistically averaged physical quantities considered in this work. We present the explict formulas for the QD occupancy, the pair correlation function and the transient currents flowing between QD and external reservoirs for the case $\Delta=\infty$ and $U=0$.

A.1 Laplace transforms

To calculate the expectation values of the quantities studied in this work we need the time-dependent operator $\hat{d}_{\sigma}(t)$. We can find it taking the corresponding inverse Laplace transform: ${\cal{L}}^{-1}\!\left\{\hat{d}_{\sigma}(s)\right\}(t)$. To find $\hat{d}_{\sigma}(s)$ we write, in the first step, the equations of motion for the closed set of eight operators: $\hat{d}_{\sigma}(t)$, $\hat{d}_{-\sigma}^{\dagger}(t)$, $\hat{c}_{{\bf k}\sigma}(t)$, $\hat{c}_{{\bf k}-\sigma}^{\dagger}(t)$, $\hat{c}_{{\bf q}\sigma}(t)$, $\hat{c}_{-{\bf q}-\sigma}^{\dagger}(t)$, $\hat{c}_{{\bf q}-\sigma}^{\dagger}(t)$ and $\hat{c}_{-{\bf q}\sigma}(t)$. Performing the Laplace transforms of these differential equations we have e.g. for $\sigma=\uparrow$ (for arbitrary energy gap $\Delta$ neglecting the correlation effects, $U=0$):

[TABLE]

Here we have assumed $\varepsilon_{{\bf q}\sigma}=\varepsilon_{\bf q}$, $\varepsilon_{\bf q}=\varepsilon_{-\bf q}$ and the subscript ${\bf k}({\bf q})$ corresponds to the normal (superconducting) electrode.

To find $\hat{d}_{\uparrow}(s)$ we calculate $\hat{c}_{-{\bf q}\downarrow}^{\dagger}(s)$ from Eq. (24) and insert it to $\hat{c}_{{\bf q}\uparrow}(s)$ obtained from Eq. (24). In the next step we insert $\hat{c}_{{\bf q}\uparrow}(s)$ into the expression for $\hat{d}_{\uparrow}(s)$ taken from Eq. 24a together with $\hat{c}_{{\bf k}\uparrow}(s)$ obtained from Eq. (24b). In result we have:

[TABLE]

We next repeat the procedure (24e-24), obtaining

[TABLE]

where

[TABLE]

From equations Eq. (25, 26) we obtain

[TABLE]

In order to find $\hat{d}_{\downarrow}(s)$ one should repeat the above described procedure for corresponding set of equations of motion for operators: $\hat{d}_{\downarrow}$, $\hat{d}_{\uparrow}^{\dagger}$, $\hat{c}_{{\bf k}\downarrow}$, $\hat{c}_{{\bf k}\uparrow}^{\dagger}$, $\hat{c}_{{\bf q}\downarrow}$, $\hat{c}_{-{\bf q}\uparrow}^{\dagger}$, $\hat{c}_{{\bf q}\uparrow}^{\dagger}$ and $\hat{c}_{-{\bf q}\downarrow}^{\dagger}(t)$, respectively. In result we have:

[TABLE]

Additionally, also the expression for $\hat{c}_{{\bf q}\sigma}$ required in calculations of the current flowing between the QD and the superconducting lead can be obtained from the set of Eqs (24a-24)

[TABLE]

where $\alpha=+(-)$ for $\sigma=\uparrow(\downarrow)$. The Laplace transforms of $\hat{d}^{\dagger}_{\sigma}(s)$ can be obtained taking the hermitian conjugation of the corresponding operator $\hat{d}_{\sigma}(s)$, Eqs. (31,32). Note, that in the wide-band-limit approximation the functions $M^{(\pm)}_{\sigma}(s)$ and $K(s)$ can be expressed in the superconducting atomic limit $\Delta=\infty$ in the forms $s\pm i\varepsilon_{\sigma}+\Gamma_{N}/2$ and $\Gamma_{S}/2$, respectively. Finally, as an example, we give here in the explicit form the Laplace transform of $\hat{d}_{\uparrow}(t)$:

[TABLE]

The expression for $\hat{d}_{\uparrow}^{\dagger}(s)$ can be obtained taking the hermitian conjugation of $\hat{d}_{\uparrow}(s)$ and making the replacement $(s_{3},s_{4})\leftrightarrow(s_{1},s_{2})$, where

[TABLE]

Note, that in the formula for $\hat{d}_{\sigma}(s)$ there is still present the superconduction energy gap parameter ($\Delta$) in all operator terms. The limit $\Delta\rightarrow\infty$ will be done later in calculations of the average values of the product of two corresponding operators, e.g. $\left<\hat{d}_{\sigma}^{\dagger}(t)\hat{d}_{\sigma}(t)\right>$ or $\left<\hat{d}_{\sigma}^{\dagger}(t)\hat{c}_{{\bf k}^{\prime}\sigma}(t)\right>$.

A.2 QD occupancy

We calculate the QD occupancy, $n_{\sigma}(t)$, according to the formula Eq. 6 taking the expectation value of the product of two corresponding inverse Laplace transforms. As at $t=0$ the QD is decoupled from the external reservoirs, the only nonvanishing expectation values would comprise the following averages $\left<\hat{d}_{\sigma}^{\dagger}(0)\hat{d}_{\sigma}(0)\right>=n_{\sigma}(0)$, $\left<\hat{d}_{\sigma}(0)\hat{d}_{\sigma}^{\dagger}(0)\right>=1-n_{\sigma}(0)$, $\left<\hat{c}_{{\bf k}\sigma}^{\dagger}(0)\hat{c}_{{\bf k}^{\prime}\sigma}(0)\right>=\delta_{{\bf k},{\bf k}^{\prime}}f_{N}(\varepsilon_{\bf k})$, and $\left<\hat{c}_{{\bf k}\sigma}(0)\hat{c}_{{\bf k}^{\prime}\sigma}^{\dagger}(0)\right>=\delta_{{\bf k},{\bf k}^{\prime}}\left[1-f_{N}(\varepsilon_{\bf k})\right]$, where $f_{N}(\varepsilon_{\bf k})$ is the Fermi-Dirac distribution of the normal lead. Other terms, corresponding to mobile electrons of the superconducting lead, can be neglected in the limit $\Delta=\infty$ (we return to this problem later), but they can be of course included when considering the finite energy gap. Finally, using Eq. (31), the QD occupancy can be expressed in the form:

[TABLE]

where for $\sigma=\downarrow$ the replacement $(s_{1},s_{2},s_{3},s_{4})\rightarrow(s_{3},s_{4},s_{1},s_{2})$ should be made. In the wide-band limit approximation we can recast the third and fourth terms to the form:

[TABLE]

The final explicit formula for $n_{\sigma}(t)$ is somewhat lengthy, therefore we present it here for the case $\varepsilon_{\sigma}=0$:

[TABLE]

Functions $F_{1}(\varepsilon,t)$ and $F_{2}(\varepsilon,t)$ have the following analytical forms:

[TABLE]

[TABLE]

and $A(\varepsilon)=\left(\frac{\Gamma_{N}^{2}}{4}+\varepsilon_{-}^{2}\right)\left(\frac{\Gamma_{N}^{2}}{4}+\varepsilon_{+}^{2}\right)$, $\varepsilon_{+/-}=\varepsilon\pm\frac{\Gamma_{S}}{2}$. It should be noted that for $\Gamma_{S}=0$ the formula (A.2) for $n_{\sigma}(t)$ reduces to the standard expression obtained by the non-equilibrium Green’s function technique, e.g. Jauho et al. (1994):

[TABLE]

Let us return to the discussion about the terms which appear in general formula, Eq. (6), but involve operators $\hat{c}_{{\bf q}\sigma}$. Let us analyze one of these terms e.g.

[TABLE]

which can be reduced to the form:

[TABLE]

where we have used the equality $\left<\hat{c}_{{\bf q}\uparrow}^{\dagger}(0)\hat{c}_{{\bf q}^{\prime}\uparrow}(0)\right>=\delta_{{\bf q}{\bf q}^{\prime}}f_{s}(\varepsilon_{{\bf q}})$. We have checked numerically (integrating the product of two corresponding inverse Laplace transforms) that this integral is smaller and smaller with increasing $\Delta$, so for $\Delta=\infty$ we put it equal to zero. Similarly, we have also checked all other terms involving operators $\hat{c}_{{\bf q}\sigma}(0)$ and found they can be neglected for $\Delta=\infty$.

A.3 QD pair correlation function

Using the formulas for $\hat{d}_{\sigma}(s)$, Eq. (31,32) and performing similar calculations as for the QD occupancy, the induced on-dot pairing, $\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$, can be written in the general form as:

[TABLE]

A.4 QD-superconducting lead current

Starting with the formula similar to the one given in Eq. (13) the QD-superconductor current can be obtained from

[TABLE]

where $\hat{d}_{\sigma}^{\dagger}(s)$ and $\hat{c}_{{\bf q}\sigma}(s)$ are given in Eqs. (31) and (33). Performing similar calculations as in previous subsections let us consider the non-vanishing term proportional to $n_{\uparrow}(0)$. It has the form

[TABLE]

The first part of this equation vanishes for $\Delta=\infty$ and the second part we calculate interchanging the summation over ${\bf q}$ with the Laplace transformation. In result we have:

[TABLE]

As $\sum_{\bf q}\frac{V_{\bf q}^{2}\Delta}{(s^{2}+\varepsilon_{\bf q}^{2}+\Delta^{2})}=\frac{\Gamma_{S}}{2}$ for $\Delta=\infty$, then finally the term proportional to $n_{\uparrow}(0)$ takes the form:

[TABLE]

In similar way we calculate the terms proportional to $\left<\hat{d}_{\uparrow}(0)\hat{d}_{\uparrow}^{\dagger}(0)\right>$, $\left<\hat{c}_{{\bf k}\uparrow}^{\dagger}(0)\hat{c}_{{\bf k}^{\prime}\uparrow}(0)\right>$ and $\left<\hat{c}_{{\bf k}\uparrow}(0)\hat{c}_{{\bf k}^{\prime}\uparrow}^{\dagger}(0)\right>$. All other terms containing the expectation values of two superconducting lead electron operators vanish (for $\Delta=\infty$), similarly as in the case for $n_{\uparrow}(t)$. Finally, we get the general equation for $j_{S\sigma}(t)$

[TABLE]

where

[TABLE]

For $\sigma=\downarrow$ the replacement $(s_{1},s_{2},s_{3},s_{4})\rightarrow(s_{3},s_{4},s_{1},s_{2})$ should be done. Note that for $\mu_{N}=0$ the formula for the superconducting current simplifies as $\mbox{\rm Re}\Phi_{\sigma}=0$. To show this property we assume the case $\varepsilon_{\sigma}=0$ and express $\Phi_{\sigma}$ in the following form:

[TABLE]

where

[TABLE]

[TABLE]

and $\varepsilon_{+/-}=\varepsilon\pm\frac{\Gamma_{S}}{2}$. For $\mu_{N}=0$ and zero temperature the function $\Phi_{\sigma}$ can be written in the form:

[TABLE]

and using the following properties of $A(\varepsilon)$ functions, $A_{1}(\varepsilon)=A_{1}(-\varepsilon)$, $A_{2}(\varepsilon)=-A_{2}(-\varepsilon)$, $\mbox{\rm Re}A_{3}(\varepsilon)=\mbox{\rm Re}A_{3}(-\varepsilon)$ and $A_{4}(-\varepsilon)=A_{4}^{*}(\varepsilon)$, we find that $\mbox{\rm Re}\Phi_{\sigma}=0$. This conclusion is also valid for $\varepsilon_{\sigma}\neq 0$.

Appendix B Mean field approximation

Let us consider the effective Hamiltonian of the proximitized QD coupled to the normal lead, treating correlations within the Hartree-Fock-Bogoliubov approximation

[TABLE]

In general, all parameters $\varepsilon_{\sigma}$, $\Gamma_{S}$, $V_{\bf k}$, $\varepsilon_{{\bf k}\sigma}$ can be time-dependent. We outline the algorithm for numerical computation of the QD charge $n_{\sigma}(t)$ and the induced pairing $\chi(t)=\left<\hat{d}_{\downarrow}(t)\hat{d}_{\uparrow}(t)\right>$. We have to solve numerically the following set of coupled equations of motion

[TABLE]

where $\bar{\varepsilon}_{\sigma}(t)=\varepsilon_{\sigma}+Un_{-\sigma}(t)$ and $\bar{\Delta}(t)=\frac{\Gamma_{S}}{2}+U\chi(t)$. Here we have used the wide-band-limit approximation and assumed $\varepsilon_{{\bf k}\sigma}$ to be time-independent. The new functions appearing in (55,56) can be determined solving corresponding equations of motion

[TABLE]

where $\alpha=+(-)$ for $\sigma=\uparrow(\downarrow)$ and $f_{N}(\varepsilon_{\bf k})$ is the Fermi-Dirac distribution of mobile electrons in the normal lead.

Appendix C Subgap Kondo effect

When the Coulomb potential $U$ is sufficiently large in comparison to $\Gamma_{S}$ the QD ground state evolves towards the spinful (doublet) configuration $\left|\sigma\right>$. Under such conditions the effective spin exchange between the correlated QD and mobile electrons of the metallic lead activate the subgap Kondo effect. It has been analyzed by many groups, using various techniques Martín-Rodero and Levy Yeyati (2011). In the present context we shall make use of basic facts, pointed out recently by R. Žitko et al Žitko et al. (2015) and independently by one of us Domański et al. (2016, 2017).

The exchange interaction $-\sum_{{\bf k},{\bf p}}J_{{\bf k,p}}\;\;\hat{{\bf S}}_{d}\cdot\hat{{\bf S}}_{{\bf k}{\bf p}}$ between the QD spin $\hat{{\bf S}}_{d}$ and spins $\hat{{\bf S}}_{{\bf k}{\bf p}}$ of the mobile electrons in normal lead can be determined by means of the generalizing canonical Schrieffer-Wolff transformation. Adopting it to the N-QD-S setup it has been found, that for the superconducting atomic limit the exchange coupling near the Fermi energy $J_{{\bf k}_{F},{\bf k}_{F}}$ is equal to Domański et al. (2016)

[TABLE]

For a spinful configuration the Kondo temperature can be estimated e.g. using the Bethe-Ansatz formula $T_{K}\propto\mbox{{\rm exp}}\left\{-1/\left[2\rho(\varepsilon_{F})J_{{\bf k}_{F}{\bf k}_{F}}\right]\right\}$, where $\rho(\varepsilon_{F}$) is the density of states of the normal lead at the Fermi level. We have compared such results with the unbiased NRG calculations and it has been found that the Kondo temperature is expressed by Domański et al. (2016)

[TABLE]

with $\eta\approx 0.6$. In particular, for the half-filled quantum dot ($\varepsilon_{\sigma}=-U/2$) the exchange couling (59) simplifies to

[TABLE]

where $J_{{\bf k}_{F},{\bf k}_{F}}^{(N)}$ stands for the normal case ($\Gamma_{S}\!=\!0$). Upon approaching a transition from the spinful doublet to the BCS-like (spinless) ground state the Kondo temperature is substantially enhanced Žitko et al. (2015); Domański et al. (2016)

[TABLE]

To get some insight into the transient phenomena related with the subgap Kondo regime we make use of the final conclusions inferred in Ref. Nordlander et al. (1999) from the time-dependent noncrossing approximation study. The characteristic time $\tau_{K}$ needed for the Abrikosov-Suhl peak to emerge at the Fermi level has been found to scale inversely with the Kondo temperatutre, i.e. $\tau_{K}\sim 1/T_{K}$. This information adopted to our N-QD-S setup implies the following relative ratio for the half-filled QD

[TABLE]

where $\tau_{K}^{(N)}$ stands for the normal state value ($\Gamma_{S}=0$). We plot this scaling in Fig. 12. Let us remark, that many-body screening (59) of the QD spin can be practically realized only in the doublet ground state (which for the half-filled QD occurs when $\Gamma_{S}<U$). By increasing the ratio $\Gamma_{S}/U$ the Andreev bound states tend to their crossing and simultaneously the Abrikosov-Suhl peak (62) quickly broadens Žitko et al. (2015); Domański et al. (2016). This explains why the characteristic time-scale $\tau_{K}$ strongly decreases with respect to $\Gamma_{S}/U$. More systematic analysis of this phenomenon is beyond a scope of the present paper.