Quench dynamics of Fano-like resonances in the presence of the on-dot superconducting pairing

Abstract

We explore the electron dynamics of a system composed of double quantum dot embedded between metallic and superconducting leads in a “T-shape” geometry. In nanoscopic systems, where electron transfer between electrodes can be realized via different paths, interference effects play an important role. For double quantum dot system in the chosen geometry, interference of electrons transferred between electrodes via the interfacial quantum dot and electrons scattered on the side dot gives rise to Fano-like interference. If such a system is additionally coupled to a superconducting electrode, together with the well-understood Fano resonance an additional resonance appears on the opposite side of the Fermi level. In the recent work (Barański et al. in Sci Rep 10:2881, 2020), we showed that this resonance occurs solely as a result of the local pairing of non-scattered electrons with scattered ones. In this work, considering the quench dynamics, we explore how much time is required for formation of each of these resonances. In particular, (i) we analyze the charge oscillations between subsystems; (ii) we estimate the time required for each resonance to achieve stable equilibrium upon an abrupt change of interdot connection; (iii) we discuss a typical energy and time scales for experiments on similar architectures.

Introduction

Understanding of electron dynamics in nanoscopic systems is important issue for development of modern information processing devices. To determine the upper limit of how fast such devices can perform logical operations, the system reaction time to a sudden change of parameters should be taken into account. Process in which one parameter of quantum system is changed much faster than time required for its thermalization is called a quantum quench¹.

Equilibrium conditions analysis of charge transport through mesoscopic heterostructures turns out to be invaluable in the context of understanding of quantum interference effects and their coexistence with other many-body effects. In particular, the Fano resonance originating from interference of electrons transmitted simultaneously through broad and narrow energy levels has been extensively explored both on theoretical as well as experimental ground. Interplay of the Fano-like interference with the Coulomb blockade^2,3, the Kondo effect^4,5,6, topological states^7,8,9,10, proximity-induced pairing^{11,12,13,14,15,16} or electron-boson interaction^17,18,19 are only a few examples. It was also established that sensitivity of Fano resonance for electron phase shift could be useful for probing of decoherences^20,21,22. Asymmetric Fano-line shapes have been also experimentally observed in various mesoscopic structures including double quantum dots⁶, Aharonov-Bohm rings^23,24 or quantum wires with a side-coupled quantum dot^25,26.

One particularly interesting experiment has been conducted by K. Kobayashi et al.²⁴. The authors realized controllable device in which the Fano effect can be electrostatically turned on and off. The architecture used by the authors is based on a ,,bridge” concept, where quantum dot is placed between the electrodes and an additional “bridge” electrode couples source and drain directly. Conceptually such realization is the same as double slit experiment for electrons with exception that one ,,slit” (namely QD) possess discrete energy levels, while bridge electrode has continuum of states. By electrostatic pinching on and off the bridge arm, the authors were able to tune the system from ballistic to interferometric regime. Although results presented in²⁴ are gathered in equilibrium conditions, the experiment seems to be a promising ansatz for analysis of time required for formation of the Fano-type interference. Implementation of similar experiment using modern techniques which provide an insight into electron time domain with picosecond resolution^27,28,29 would allow for an experimental inspection of the Fano resonance formation time. This perspective encourages to theoretical approaches concerning dynamics of interference effects in quantum dots.

Transient effects in quantum dot hybrids in various configurations have already been studied theoretically by a few authors^30,31,32,33. Dynamics of superconducting (SS) proximity effect has been inspected as well in Refs.^{34,35,36,37,38,39} including also topological phenomena⁴⁰. Recently, our group has also contributed to such studies⁴¹. Dynamics of Fano-type interference, however, has been presented only in Refs.^42,43. The aforementioned works focus on the case of two quantum dots coupled to non-superconducting/metallic electrodes in T-shape geometry. Steady state analysis of our group on Fano resonances in presence of superconductors shows that interplay of quantum interference with local pairing leads to unique effects such as anomalous Fano resonance¹⁴. In this work we will analyze such effects upon quench dynamics. To our knowledge time resolved analysis of interference effects in presence of local pairing has not been conducted yet.

One of convenient systems for theoretical analysis of such interferences in nanoscopic systems is composed of two quantum dots deposited in such a way that one quantum dot ($QD_1$) is strongly coupled to electrodes, while the other one ($QD_2$) is only side attached to the first quantum dot. This is so-called “T-shape” configuration. Realization of similar architecture have been done e.g. in²⁶. In this architecture lifetime of electron on $QD_1$ is relatively short and in accordance to Fermi golden rule its energy levels become broad. In contrary, due to lack of direct coupling of $QD_2$ to continuum of states, its energy levels remain quasidiscrete. Interference of electrons transmitted through broad level and those resonantly scattered on discrete level gives rise to well pronounced Fano-like resonances. Asymmetric Fano-like lineshapes can be observed both in density of states of interfacial dot ($QD_1$) near energy of discrete level, as well as in differential conductivity for gate voltages close to narrow level energy.

In the present work, we inspect the dynamics of the interference effects when such system is additionally connected to superconducting electrode. Thus, the full system considered here is composed of a double quantum dot deposited between two metallic (L and R) and one superconducting (S) electrodes as schematically depicted in Fig. 1. Our goal is to estimate how much time is required for formation of Fano-like resonances upon establishing abrupt connection of interdot coupling $t_{12}$.

In our recent work on a similar architecture system, we showed that, if the local pairing is present, along the well-understood Fano resonance another feature on the opposite side of the Fermi level is formed¹⁴. The shape of this feature and the nature of its origin reach beyond traditional descriptions of the Fano-like resonances. It was shown that this anomalous Fano (AF) feature arises as a result of pairing of non-scattered electrons with scattered ones. In particular, we considered a “toy model”, where we assumed spin-polarized interdot tunneling. Such an assumption ensures that only one spin (say $\uparrow$) electrons are directly scattered on the side level and eventual resonant features in the other spin component ($\downarrow$-spin) originate solely from pairing with scattered ones. For such a model, in the spectrum of directly scattered electrons ($\uparrow$-spin), we observe only the ordinary Fano resonance while for electrons that are not directly scattered ($\downarrow$-spin) only the anomalous resonance appears. This shows that information about scattering of a given electron on the side level is transferred to another electron via the local pairing and AF resonance arises as a response to scattering of opposite spin electron. In a case of non-polarized tunneling, both electrons experience simultaneously direct scattering (leading to the emergence of the ordinary Fano (OF) resonance near $\varepsilon _2$) and coupling to scattered (opposite spin) electrons leading to appearance of AF resonance near $-\varepsilon _2$. Consequently, in the non-polarized case, both resonances are present in density of states and differential conductivity.

In the context of quench dynamics, an important question arises: Is there a difference in the timing of the formation of these two resonances? On one hand, if one foresees the appearance of the AF features in the spectrum of a given spin electrons, the electron of opposite spin should be scattered on the side level firstly, and only then information about this event should be passed to the first electron giving rise to the AF resonance. From this perspective, one could expect that the feature of direct scattering will be formed faster than its superconducting response. On the other hand, due to the local pairing, two electrons of opposite spin form one quantum object. If the superconducting-like correlations on the QD’s are instantaneous, than any process that involve one electron should have its instantaneous effect on the other. Therefore, if information transfer between paired electrons is infinitely fast, there should not be any delay between formation of the OF and the AF resonances. In this regard, our analysis can be understood as answering the question whether the information transfer through local pairing have instantaneous character?

In the present work, we explore the electron dynamics of a system composed of double quantum dot embedded between metallic and superconducting leads in a “T-shape” geometry. Considering the quench dynamics, we explore how much time is required for formation of both the ordinary Fano resonance (associated with direct scattering of electrons) and the anomalous Fano resonant features (associated with indirect scattering of electrons due to the superconducting pairing). In particular, we analyze the charge oscillations between subsystems in the absence of normal electrodes. The main part of the paper is devoted to estimation of the time required for both resonances to achieve stable equilibrium upon an abrupt change of interdot connection. Moreover, we discuss typical energy and time scales for experiments on similar architectures.

Figure 1

The schematic illustration of the analyzed system. It consists of the double quantum dot (i.e., two quantum dots: $QD_1$ and $QD_2$) embedded between two metallic (L, R) electrodes and one superconducting electrode (S).

Model

A setup investigated in this work and schematically depicted in Fig. 1 can be described by the following Hamiltonian:

$$\begin{aligned} \hat{H}= \sum _{\beta =L,R,S} (\hat{H_{\beta }} + \hat{H}_{T \beta }) + \sum _{i=1,2} \hat{H}_{QDi} + \hat{H}_{t}, \end{aligned}$$

(1)

where $\hat{H}_{\beta }$ correspond to three external electrodes (two metallic $\beta =L,R$ and one superconducting $\beta =S$). Two single-level quantum dots and the tunnel coupling between them are represented by $\hat{H}_{QDi}= \sum _{ \sigma } \varepsilon _{i} \hat{d}^{\dagger }_{i \sigma } \hat{d}_{i \sigma } + \sum _{ \sigma } U_{i} \hat{n}_{i \sigma } \hat{n}_{i \bar{\sigma }}$ and $\hat{H}_{t}=\sum _{\sigma } t_{12 \sigma } \left( \hat{d}^{\dagger }_{1 \sigma } \hat{d}_{2 \sigma } + \hat{d}^{\dagger }_{2 \sigma } \hat{d}_{1 \sigma }\right)$, respectively ($\sigma = \uparrow , \downarrow$). $\hat{H}_{T \beta }$ stands for hybridization of each electrode to interfacial quantum dot ($QD_1$). We treat metallic reservoirs as free electron gas, thus for $\beta =L,R$ we have $\hat{H}_{\beta } = \sum _{k,{\beta }, \sigma } \xi _{k {\beta } } \hat{c}^{\dagger }_{k {\beta } \sigma } \hat{c}_{k {\beta } \sigma }$, where electron energies $\xi _{k {\beta } }$ are measured with respect to chemical potentials $\mu _\beta$. Assuming isotropic energy gap ($\Delta _S$) of superconducting electrode, $\hat{H}_{\beta }$ for $\beta =S$ reads as $\hat{H}_{S}=\sum _{k,S, \sigma } \xi _{k S} \hat{c}^{\dagger }_{k S \sigma } \hat{c}_{k S \sigma } - \sum _k ( \Delta _S \hat{c}^{\dagger } _{k S \uparrow } \hat{c}^{\dagger }_{-k S \downarrow }+h.c.)$. The hybridization part is described by $\hat{H}_{T \beta }= \sum _{k, \sigma } V_{k \beta } \left( \hat{c}_{k \beta \sigma } \hat{d}^{\dagger }_{1 \sigma } + h.c.\right)$ and represents the coupling of $QD_1$ to $\beta$ electrode.

Restricting our analysis to the energies deep inside the superconductor energy gap (i.e., $\Delta _S \rightarrow \infty$), we can neglect the existence of single-electron energy levels in the SC electrode. This, so-called superconducting atomic limit approximation, allows to reduce the influence of SC electrode on $QD_1$ by a static terms representing induced local pairing ($\hat{d}_{1 \uparrow } \hat{d}_{1 \downarrow } + h.c.$)^44,45,46. Studies on similar architecture in static conditions^14,16,20 show that correlations have a marignal effect on appearance of interference features apart from two specific cases. First one, when energy level of a side dot is close to 0 and interference patterns coincide with Kondo resonance leading to the emergence of so-called Fano–Kondo resonance⁵. The second case is when Coulomb satellite state of $QD_2$ (with energy $\varepsilon _2+U_{2}$) leads to appearance of additional interference structures²⁰. Transient effects in the Kondo regime were also studied in single-quantum dot systems, e.g., Refs.^47,48. We note that in experimental realizations correlations are usually much stronger than effective pairing and thus satellite state lies far beyond the considered energy scale. Let us underline that transient effects of laterally coupled quantum dots between normal and superconducting electrodes in both uncorrelated and correlated regimes was studied in Ref.³⁶. As it was pointed there, for the hybridization of quantum dots to the SC electrode is stronger than correlations on the dots, the system is dominated by the Andreev scattering. In the mentioned work, even for correlations strength of $U = 0.5 \Gamma _{S}$ ($U_1=U_2=U$, the Hartree-Fock-Bogoliubov decoupling scheme used), its influence on conductance was rather marginal and it was only limited to the influence on the final shape of transmittance and not the time of its evolution. In the case of strongly correlated system (i.e., $U=1.5 \Gamma _{S}$), two Andreev states form a single Lorentzian, however, the evolution in the time domain of all features remains the same. Thus, one expects that the Coulomb repulsion, although it rearranges the static background transmittance, has rather a negligible effect on evolution of the resonant features in time. Detailed analysis of this issue is out of the scope of the present work.

In this work, we focus solely on building up of interference effects in case when ordinary and anomalous Fano resonances are well separated (this is $\varepsilon _2 \ne 0$). Therefore, in our calculations, we will omit onsite Coulomb interactions. Consequently, constituents of the Hamiltonian representing the interfacial quantum dot, the SC electrode, and their mutual interactions can be rewritten in the following form:

$$\begin{aligned} \hat{H}_{QD1}\!+\!\hat{H}_{S}\!+\!\hat{H}_{T S}\!=\!\sum _{\sigma }\! \varepsilon _{1} \hat{d}_{1 \sigma }^{\dagger }\hat{d}_{1 \sigma }\!-\!\Delta _d \left( \hat{d}_{1 \downarrow } \hat{d}_{1 \uparrow } \!+\!\text{ h.c. } \right) \!, \ \end{aligned}$$

(2)

where $\Delta _d$ is dependent solely on the $QD_1-SC$ electrode coupling constant: $\Delta _d=\Gamma _S/2=\pi \sum _{k} |V_{k S}|^2 \delta (\omega -\xi _{k S})$. In our work, we use the total coupling to conduction electrodes ($\Gamma _N$) as energy unit: $\Gamma _N=\Gamma _L + \Gamma _R= 2 \pi \sum _{k\beta } |V_{k,\beta }|^2 \delta (\omega -\xi _{k \beta })$. Our time unit (t.u.) is given by the inverse of this energy scale $1\ \text {t.u.} = \hbar /\Gamma _{N}$. More detailed discussion on the time and energy scales is given in Section Quantified values of time.

Methods

Dynamics of coupled complex systems with discrete energy levels can be described in both time and energy scales. One of the scales is related to charge oscillations between subsystems. The period of these oscillations is usually associated with the energy of coupling between subsystems. In our case, apart from interdot coupling $t_{12\sigma }$, oscillations originating from the coupling of the system to superconducting electrode play an important role. All this oscillations do not thermalize unless the system is hybridized to continuum of states. Such continuum acts as damping force for classical oscillations^41,49. Therefore, if one is looking for a time of achieving steady-state solutions, the time scale is predominantly defined by energy of the coupling to metallic lead. In our approach we first calculate charge oscillations on both quantum dots.

Using the Laplace transformation, we conduct analytic calculations of the exactly solvable model with the time-independent parameters, where coupling to continuum is assumed to be negligible, i.e., $\Gamma _N \rightarrow 0^{+}$. The main goal of the work, i.e., analysis of the above described double quantum dot system for arbitrary $\Gamma _N \ne 0$ and with the spin-independent interdot coupling ($t_{12\uparrow }(t) = t_{12\downarrow }(t) \equiv t_{12}(t)$), will be achieved by numerical calculations using the 4th order Runge-Kutta (4RK) method applied to equations of motion obtained in the Heisenberg representation.

Time dependent occupancies and charge oscillations

In order to derive the time dependent occupation number $n_{i \sigma } (t) = \langle {\hat{d}}^{\dagger }_{i \sigma }(t){\hat{d}}_{i \sigma }(t)\rangle$, we start with equation of motion in Heisenberg representation for a given operator $\hat{O}$:

$$\begin{aligned} \textbf{i} \hbar \frac{d}{d t} \hat{O} (t) = [\hat{O}(t), \hat{H}]. \end{aligned}$$

(3)

We calculate the set of such differential equations for creation ${\hat{d}}^{\dagger }_{1 \sigma }, {\hat{d}}^{\dagger }_{2 \sigma }, {\hat{c}}^{\dagger }_{k \sigma }$ and annihilation ${\hat{d}}_{1 \sigma }, {\hat{d}}_{2 \sigma }, {\hat{c}}_{k \sigma }$ operators. Next, for each fermion operator we conduct a Laplace transform ${{\mathscr {L}}}[\hat{O}(t)](s)$ defined as

$$\begin{aligned} \hat{O}(s) = \int _0^\infty dt \exp {\left( -st \right) } \hat{O}(t) \equiv {{\mathscr {L}}}[\hat{O}(t)](s), \end{aligned}$$

(4)

where s is a complex variable. Such a procedure transforms a set of differential equations into a set of linear equations dependent on s which are solvable analytically. We get analytical expression for each operator $\hat{O}(s)$. For example, in the case of time-independent $t_{12\sigma }$’s, for the annihilation operators ${\hat{d}}_{1 \uparrow }$ and ${\hat{d}}_{1 \downarrow }$, the s-dependent form looks as follows

$$\begin{aligned} {\hat{d}}_{1 \uparrow }(s) = - \frac{ \textbf{i} \left[ {\hat{A}}^{\dagger }_{\downarrow }(0)\Gamma _{S} + \textbf{i} {\hat{A}}_{\uparrow }(0) M_{\downarrow }^{*}(s) \right] }{\Gamma _{S}^2 + M_{\downarrow }^{*}(s) M_{\uparrow }(s)}, \qquad {\hat{d}}_{1 \downarrow }(s) = - \frac{ \textbf{i} \left[ {\hat{A}}^{\dagger }_{\uparrow }(0)\Gamma _{S} + \textbf{i} {\hat{A}}_{\downarrow }(0) M_{\uparrow }^{*}(s) \right] }{\Gamma _{S}^2 + M_{\uparrow }^{*}(s) M_{\downarrow }(s)}, \end{aligned}$$

(5)

where

$$\begin{aligned} {\hat{A}}_{\sigma }(0)= & {} -\textbf{i} \sum _{k,\beta } \frac{V_{k \beta }}{s + \textbf{i} \xi _{k \beta }} {\hat{c}}_{k \beta \sigma }(0) - \frac{ \textbf{i} t_{12\sigma }}{s + \textbf{i} \varepsilon _2} {\hat{d}}_{2 \sigma }(0) + {\hat{d}}_{1 \sigma }(0),\nonumber \\ M_{\sigma }(s)= & {} s + \textbf{i} \varepsilon _1 + \sum _{k,\beta } \frac{|V_{k \beta }|^2}{s + \textbf{i}\xi _{k\beta }} + \frac{t_{12\sigma }^2}{s + \textbf{i} \varepsilon _2}. \end{aligned}$$

(6)

Note that equations above are obtained for a general case of spin-dependent $t_{12\sigma }$ ($t_{12\uparrow } \ne t_{12\downarrow }$) and arbitrary $\Gamma _{N}$, $\Gamma _{S}$, $\varepsilon _1$, and $\varepsilon _2$.

Now, we can get expressions for the occupation numbers in s-dependent form (i.e., expressions for $n_{i \sigma } (s) = \langle {\hat{d}}^{\dagger }_{i \sigma }(s){\hat{d}}_{i \sigma }(s)\rangle$) in the limit of $\Gamma _N \rightarrow 0$. The general equations are rather complex, therefore, we present here the results obtained in the conditions: $\varepsilon _1=\varepsilon _2=0$ and $n_{2 \sigma }(0)=0$. After removing pairs of operators that give zero and applying the inverse Laplace transform for a given operator ${{\mathscr {L}}}^{-1}[\hat{O}(s)](t)$, for the spin-independent inter dot coupling, i.e., $t_{12\uparrow } = t_{12\downarrow }=t_{12}$, we formulate the time dependence of expectation values $n_{1 \sigma }(t)$ and $n_{2 \sigma }(t)$:

$$\begin{aligned} n_{1 \sigma } (t)= & {} - [1 - {n}_{1 \bar{\sigma }}(0)] L_1^2 + {n}_{1 \sigma }(0) (L_2 + L_3)^2, \end{aligned}$$

(7)

$$\begin{aligned} n_{2 \sigma } (t)= & {} [1 - {n}_{1 \bar{\sigma }}(0)] L_4^2 - {n}_{1 \sigma }(0) (L_5 + L_6)^2, \end{aligned}$$

(8)

where $n_{i \sigma } (0) = \langle \hat{d}^{\dagger }_{i \sigma }(0)\hat{d}_{i \sigma }(0)\rangle$ and

$$\begin{aligned} L_1= & {} {{\mathscr {L}}}^{-1} \left[ \frac{\textbf{i} \Gamma _s}{\Gamma _s^2 + \left( s + \frac{t_{12}^2}{s} \right) ^2} \right] , \qquad \ L_2 = {{\mathscr {L}}}^{-1} \left[ \frac{s}{\Gamma _s^2 + \left( s + \frac{t_{12}^2}{s}\right) ^2} \right] , \quad L_3 = {{\mathscr {L}}}^{-1} \left[ \frac{t_{12}^2}{s \left[ \Gamma _s^2 + \left( s + \frac{t_{12}^2}{s} \right) ^2\right] } \right] , \qquad \\ L_4= & {} {{\mathscr {L}}}^{-1} \left[ \frac{- t_{12} \Gamma _s}{s \left[ \Gamma _s^2 + \left( s + \frac{t_{12}^2}{s} \right) ^2 \right] } \right] , \quad L_5 = {{\mathscr {L}}}^{-1} \left[ \frac{\textbf{i} t_{12}}{\Gamma _s^2 + \left( s + \frac{t_{12}^2}{s} \right) ^2} \right] , \quad L_6 = {{\mathscr {L}}}^{-1} \left[ \frac{ \textbf{i} t_{12}^3}{s^2 \left[ \Gamma _s^2 + \left( s + \frac{t_{12}^2}{s} \right) ^2 \right] } \right] . \nonumber \end{aligned}$$

(9)

Figure 2

Charge oscillations in the limit of $\Gamma _N \rightarrow 0$. Occupation numbers $n_{1 \sigma }(t)$ (blue lines) and $n_{2 \sigma }(t)$ (orange lines) obtained for $t_{12}=0.3\Delta _d$ (left) and $t_{12}=10\Delta _d$ (right) with the initial conditions: $n_{1 \sigma }(0)=1$ and $n_{2 \sigma }(0)=0$ [cf. Eqs. (10) and (11)]. Green dashed lines indicate charge transferred (per spin) into the region of the SC electrode.

Here, we will find the explicit equations for the time-dependence of occupation number $n_{i \sigma }(t) = \langle \hat{d}^\dagger _{i \sigma }(t) \hat{d}_{i \sigma } (t) \rangle$ in the limit of $\Gamma _N \rightarrow 0$ and $\varepsilon _1=\varepsilon _2=0$, time-independent $t_{12\sigma } \equiv t_{12}$ for $\sigma =\uparrow ,\downarrow$. We also assume that at the initial moment ($t=0$), the side attached quantum dot is empty, i.e., $n_{2 \sigma }(0)=0$ (for both $\sigma$’s). Note that, in such a case, time dependence of $n_{i\sigma }$ is spin-independent ($n_{i\uparrow }(t)=n_{i\downarrow }(t)$). Taking into account this simplifications, we obtain from (7) and (9) (and $\Delta _d=\Gamma _S/2$)

$$\begin{aligned} n_{1 \sigma }(t)= & {} \frac{n_{1 \sigma }(0)\left[ \omega _{+}\cos (\omega _{+} t) + \omega _{-}\cos (\omega _{-} t)\right] ^2}{\Delta _d^2+4t_{12}^2} + \frac{h_{1 \bar{\sigma }}(0) \left[ \omega _{+}\sin (\omega _{+} t)-\omega _{-}\sin (\omega _{-} t)\right] ^2}{\Delta _d^2+4t_{12}^2}, \end{aligned}$$

(10)

where $h_{1 \sigma }(0)=1-n_{1 \sigma }(0)$ and $\bar{\sigma }$ denotes spin direction opposite to $\sigma$. Charge oscillations are defined by two frequencies

$$\begin{aligned} \omega _{\pm } = \sqrt{\tfrac{1}{2}\left( \Delta _d^2+2t_{12}^2 \pm \Delta _d\sqrt{\Delta _d^2+4t_{12}^2} \right) }. \end{aligned}$$

We note that, for the vanishing interdot coupling ($t_{12}\rightarrow 0$), frequency $\omega _{-}$ is equal to 0, while $\omega _{+} \rightarrow \Delta _d$. This brings the Eq. (10) to the form of $n_{1 \sigma }(t)=n_{1 \sigma }(0)\cos ^2(\Delta _d t)+h_{1 \bar{\sigma }}(0)\sin ^2(\Delta _d t)$, which reproduces the result obtained for a single QD placed between metallic and superconducting leads for vanishing $\Gamma _N$ (see Ref.⁴⁹). In the case of the most interest to us, when the interdot coupling is significantly smaller then the effective pairing (but not negligible), the resultant oscillations have low- and high-frequency modes. High-frequency mode is predominantly governed by the coupling of the dot to the superconducting reservoir $\omega _{+} \sim \Delta _d$. Low frequency ($\omega _{-}$), in turn, is related to the interdot coupling $t_{12}$. It is also interesting to investigate the charge fluctuations that occur on $QD_2$. Oscillations of $n_{1 \sigma }(t)$ are associated to a charge flow from $QD_1$ to $QD_2$ (and vice versa) and a charge flow between $QD_1$ and the SC electrode. In our model, $QD_2$ is tunnel coupled only to $QD_1$, therefore changes of $n_{2 \sigma }(t)$ originates solely from the interdot charge flow. Considering the same conditions as previously (i.e., $\varepsilon _{1}=\varepsilon _{2}=0$, and $n_{2 \sigma }(0)=0$), occupation $n_{2 \sigma }(t)$ can be calculated from (8) by performing the inverse Laplace transformation of the Eq. (9) given previously. One obtains

$$\begin{aligned} n_{2 \sigma }(t)\simeq & {} h_{1 \bar{\sigma }}(0) t_{12}^2\frac{\left[ \cos \left( \omega _{+} t\right) - \cos \left( \omega _{-} t \right) \right] ^2}{\Delta _d^2+4t_{12}^2} + n_{1 \sigma }(0)t_{12}^2\frac{\left[ \sin (\omega _{+} t)-\sin (\omega _{-} t)\right] ^2}{\Delta _d^2+4t_{12}^2}. \ \end{aligned}$$

(11)

The occupancies of both dots and the dot-lead charge flow (per spin) calculated from Eqs. (10) and (11) are presented in Fig. 2. Assuming that at the initial conditions, $QD_1$ was full and $QD_2$ was empty (i.e., $n_{1 \sigma }(0)=n_{1 \bar{\sigma }}(0)=1$, $n_{2 \sigma }(0)=n_{2 \bar{\sigma }}(0)=0$), the charge transferred to the SC electrode can be expressed by deficiency of charge in both QD’s, i.e., $1-n_{1 \sigma }(t)-n_{2 \sigma }(t)$.

For weak interdot coupling (i.e., $t_{12} < \Delta _d$), the amplitude of charge oscillations on $QD_{2}$ is strongly reduced (compared to the amplitude of $n_1(t)$), however, frequencies of oscillations for both quantum dots are very close (cf. the left panel of Fig. 2). As the occupation number of $QD_1$ oscillates in a full range (from 0 to 1) while the amplitude of oscillations in $QD_2$ is significantly reduced, we conclude that only part of the charge flowing from SC reservoir to $QD_1$ is transmitted to $QD_2$.

Although in our work we mainly focus on the limit of the weak (but not negligible) interdot coupling, it is very interesting to investigate charge oscillations for $t_{12} \gg \Delta _d$ (cf. the right panel of Fig. 2). In this so-called molecular regime, two quantum dots form a single molecule placed on the top of the superconductor. Consequently, charge oscillates are very fast between dots (cf. oscillations of blue and orange line in the figure), while the weak coupling of the molecule to the SC reservoir causes slow charge flow from the dots to SC electrode and back to the double-quantum-dot (DQD) molecule (cf. green line on the figure). Charge flow from the DQD to the SC is possible only if a charge on $QD_1$ dot is nonzero. Conversely, charge flow from the SC to the DQD is possible only if the $QD_1$ is at least partially empty. Therefore, in the case of $t_{12}\gg \Delta _d$ (where charge between dots is swapped rapidly), the charge transfer to and from the SC electrode appears sequentially after each short-term cycle between the dots. Such behavior can be inspected in Fig. 2 (right), where charge transfer between the DQD molecule and the SC (marked by the green line) occurs sequentially after each interdot charge transfer. In discussed here limit (i.e., $t_{12} \gg \Delta _d$), frequencies $\omega _{+}$ and $\omega _{-}$ become very close to each other ($\omega _{+} \approx \omega _{-}$). Introducing $\frac{1}{2}(\omega _{+}+\omega _{-}) \equiv \omega _m$ and $\frac{1}{2}(\omega _{+}-\omega _{-}) \equiv {\delta _\omega }$ (with ${\delta _\omega } \ll \omega _m$), we obtain, in this limit, the equations for charge oscillations on both dots in the form of a square of classic beats $n_{1 \sigma }(t)\simeq n_{1 \sigma }(0)\left[ \cos (\omega _{m} t)\cos (\delta _{\omega } t)\right] ^2+ h_{1 \bar{\sigma }}(0) \left[ {\cos } (\omega _{m} t)\sin (\delta _{\omega } t)\right] ^2$, $n_{2 \sigma }(t)\simeq n_{1 \sigma }(0)\left[ \cos (\omega _{m} t)\sin (\delta _{\omega } t)\right] ^2+ h_{1 \bar{\sigma }}(0) \left[ \sin (\omega _{m} t)\sin (\delta _{\omega } t)\right] ^2$. High-frequency beats (given by $\omega _m$) refer to interdot oscillations, while the low-frequency mode (given by $\delta _\omega$) is related to a charge transfer between the molecule and the SC electrode. We note that, at the limit of $t_{12} \rightarrow \infty$, low-frequency mode $\delta _{\omega }$ is equal to $\Delta _d/2$, which reproduces the case of a single molecule placed on top of superconductor, but with the twice weaker coupling⁴⁹.

Numerical approach

Asymmetric line-shapes being a symptom of electron scattering on a side level can be practically investigated by inspection of time-dependent charge current induced by a bias voltage. Our aim is to calculate the charge current $I_{L}(t,eV)$ flowing through one of the metallic electrodes upon bias voltage V applied to another one (e denotes the elementary charge), i.e., the source-drain voltage between the normal electrodes. We assume that chemical potential of the SC electrode and the metallic R electrode are equal and we measure the energy with respect to these potentials ($\mu _{SC}=\mu _{R}=0$). The pronounced resonant features can be inspected in differential conductivity $G(t,eV)=\frac{d }{d V}I(t,eV)$. The charge flowing through L electrode is given by average change of the electron number in the L lead

$$\begin{aligned} I_{L \sigma } (t, V)= -e \left\langle \frac{d N_{L \sigma }}{dt} \right\rangle , \end{aligned}$$

(12)

where $N_{L \sigma } = \Sigma _k \hat{c}^{\dagger }_{k L \sigma } \hat{c}_{k L \sigma }$.

In order to find time-dependent statistical averages we derive the closed set of ordinary differential equations of motion in terms of the Heisenberg notation $\frac{d}{dt}\hat{O}(t)=\frac{\textbf{i}}{\hbar }[\hat{H},\hat{O}]$ and apply the 4th order Runge–Kutta (RK4) method for numerical calculations of time evolution for each average. Details of this procedure were introduced in the previous work on dynamics of Majorana-QD hybrid⁴¹. The equations of motion in the Heisenberg representation have been derived using sneg library created by R. Žitko⁵⁰.

Quench protocol

We assume that initially ($t \le 0$) all parts of the system (both QD’s and three electrodes) are separated. Such initial condition can be met by setting the averages comprising operators referring to different parts of the system like $\langle \hat{d}^{(\dagger )}_{i \sigma }(0) \hat{c}^{(\dagger )}_{k \beta \sigma '}(0) \rangle$, $\langle \hat{d}^{(\dagger )}_{i \sigma }(0) \hat{d}^{(\dagger )}_{j \sigma '}(0) \rangle$ to be equal 0. The assumption that the quantum dots are separated from the SC electrode for $t<0$ also imposes that there was no pairing potential in the region of the QDs. This requirement is met by $\langle \hat{d}^{\dagger }_{i \sigma }(0) \hat{d}^{\dagger }_{i \bar{\sigma }}(0) \rangle = \langle \hat{d}_{i \bar{\sigma }}(0) \hat{d}_{i \bar{\sigma }}(0) \rangle =0$. The average numbers of electrons in L and R leads are given by the Fermi distribution function $\langle \hat{c}^{\dagger }_{k \beta \sigma }(0) \hat{c}_{k \beta \sigma }(0) \rangle = \lbrace 1 + \exp \left[ (\xi _{k \beta }-\mu _{\beta }) / (k_B T) \right] \rbrace ^{-1}$. In order to do not confuse the effect of formation of quasiparticles with interference effects we developed two step procedure. First, we assume that at time $t=0$ $QD_{1}$ is connected only to the external electrodes (i.e., L, R and S) keeping $t_{12}=0$. Static results for single QD in such heterostructure show that evolution should lead to the formation of two quasiparticle Andreev states located near $\pm \sqrt{\varepsilon _{1}^{2} +\Delta _{d}^{2} }$ (cf., e.g., Refs.^44,51). When these states achieve its static values and time fluctuations will vanish (let say at time $t=t_0$), we abruptly connect the second quantum dot [with constant interdot coupling term, $t_{12\sigma }(t) \equiv t_{12}(t)=t^{0}_{12} \theta ( t-t_0 )$ for $\sigma = \uparrow ,\downarrow$]. From this moment, scattering of electrons on the side level becomes possible and we can observe the evolution of interference patterns emerging for voltages (eV) close to $\pm \varepsilon _2$. In all our calculations we assume the quench time to be $t_0=20 \ \hbar /\Gamma _N$, which, in the considered energy scale, is much larger than the relaxation time after connecting the dot to the external electrodes.

Figure 3

Differential conductivity G $[2e^2/h]$ as a function of bias voltage eV $[\Gamma _{N}]$ and time t $[\hbar /\Gamma _{N}]$. The results are obtained for $\Delta _{d}=2\Gamma _N$, $\varepsilon _{2}=\Gamma _{N}$ and $t^{0}_{12}=0.3{\Gamma _{N}}$. The red dashed line indicates the moment ($t_{0}=20$ $\hbar /\Gamma _{N}$) at which tunneling between dots is turned on. The blue dashed line shows the temporary shape of the evolving resonant states at $t=60$ $\hbar /\Gamma _{N}$.

Static calculations

In order to verify obtained results we can calculate also the static (time-independent) conductivity. Conductivity calculated using the RK4 method for ($t \rightarrow \infty$) should reproduce these results. In three terminal heterostructure comprising two metallic and one superconducting electrode, charge transport is provided by three types of processes: (i) ballistic single electron transfer (ET) from L to R electrode, (ii) direct Andreev reflection (DAR), where single electron from L electrode is converted to a Cooper pair propagating in the SC electrode with simultaneous reflection of hole back to L electrode, and (iii) cross Andreev reflection (CAR), in which the hole is reflected to second metallic electrode (R). For the static case these three processes can be evaluated using the following Landauer-like formulas^52,53

$$\begin{aligned} J^{ET}_{L}= & {} \frac{2e}{h} \int d \omega \Gamma _{L} \Gamma _{R} |G^{r}_{11}|^2 (f_{L} - f_{R}), \end{aligned}$$

(13)

$$\begin{aligned} J^{DAR}_{L}= & {} \frac{2e}{h} \int d \omega \Gamma _L^{2} |G^{r}_{12}|^{2} (f_{L} - \tilde{f}_{L}), \end{aligned}$$

(14)

$$\begin{aligned} J^{CAR}_{L}= & {} \frac{2e}{h} \int d \omega \Gamma _{L} \Gamma _{R} |G^{r}_{12}|^{2} (f_{L} - \tilde{f}_{R}), \end{aligned}$$

(15)

where $G^r_{ij}$ are matrix elements of the retarded Green functions of $QD_1$ in the Nambu representation (given, e.g., in Ref.¹⁴), $f_{\beta } =\left\{ 1+\exp {\left[ (\omega -\mu _{\beta })/(k_B T)\right] }\right\} ^{-1}$ and $\tilde{f}_{\beta }=\left\{ 1+\exp { \left[ (\omega +\mu _{\beta })/(k_B T) \right] } \right\} ^{-1}$ are the Fermi distributions of electrons and holes, respectively. In our calculations we assume that voltage (V) is applied to L electrode while the chemical potential of S and R electrodes are equal and energies are measured with respect to them ($\mu _{R}=\mu _{S}=0$). Time-dependent current calculated using the RK4 method accounts for all these three processes together.

Results and discussion

Figure 4

The comparison of the resonant features near $\varepsilon _{2}$ (red lines) and $-\varepsilon _{2}$ (orange lines) obtained for several different times ($\tilde{t} =t-t_{0}$) with the static ($t \rightarrow \infty$) conductance (the blue lines). The results are presented for $\tilde{t}= 10, \ 40,\ 70,\ 130$ $[\hbar /\Gamma _N]$, respectively (as labelled). V denotes the bias voltage. Other parameters used here are the same as in Fig. 3.

In Fig. 3 we present a time evolution of the total conductance G (in $2e^2/h$ units) versus applied bias (source-drain) voltage eV (in units of $\Gamma _{N}$). In accordance to the assumed quench protocol (cf. Section Quench protocol), at the beginning $QD_{1}$ is connected only to the external electrodes and two Andreev states build up over time. We noticed that, after approximately 10–15 $\hbar /\Gamma _{N}$, these states saturate and all fluctuations are suppressed. To be sure that the process of building up of the Andreev states does not affect formation of the interference features, we set the quench time safely later, i.e., at $t_{0}=20$ $\hbar /\Gamma _{N}$. At this moment, the interdot connection is established with coupling strength $t^{0}_{12}=0.3\Gamma _{N}$ and the resonant characteristics start to evolve. The moment of abrupt establishing of the interdot coupling is underlined by the red dashed line in Fig. 3. The blue dashed line highlights the shape of both resonances obtained at $t=60$ $\hbar /\Gamma _{N}$.

We note that a well-pronounced AF feature (i.e., this near $-\varepsilon _{2}$) starts to develop almost instantaneously after the abrupt connection is established, while only small fluctuations appear near $eV \simeq \varepsilon _{2}$. It is surprising that even though the direct scattering feature is not developed yet, the feature announcing its superconducting response is building up so vigorously. The non-equilibrium response for scattering indicates that whenever electron being a component of a local pair is involved in a given process even by a small fluctuations, the second spin component reacts instantaneously and robustly. Considering time required for achieving equilibrium conditions, however, one need to have in mind that, in the static solutions (i.e., for $t \rightarrow \infty$), the maximum of the feature announcing the SC response for scattering, i.e., the AF feature, is magnitude of order higher than the maximum of the ordinary Fano resonance feature. Measures of development for each resonance should be thus considered relative to its final amplitude rather than one to another.

In this regard, in Fig. 4, we overlay the resonant features obtained for several time parameters on plots of static ($t \rightarrow \infty$) conductance obtained using equations (13)–(15). We note that, after $\tilde{t}=10\ \hbar /\Gamma _{N}$ from the quenching of QD$_2$, the feature announcing direct resonance (red line) is given only by tiny fluctuations, while for the indirect resonance, a well-pronounced peak is already visible. However, after $\tilde{t}=40\ \hbar /\Gamma _{N}$, development of both resonances (relative to their static amplitude) appears to be comparable. At $\tilde{t}=70\ \hbar /\Gamma _{N}$, the Fano-like resonance is close to the stable solution and only small adjustment of its shape is contributed between $\tilde{t}=70\ \hbar /\Gamma _{N}$ and $\tilde{t}=130\ \hbar /\Gamma _{N}$. In contrary, the amplitude of the AF resonance is still building up even after $\tilde{t}=130\ \hbar /\Gamma _{N}$. From the picture above, it seems that although short after quench the AF resonance builds up very fast (comparing to the ordinary Fano feature), the stable solution for the ordinary Fano resonance is achieved slightly faster than for the AF resonance. To deliver quantified data for development of both these resonances, we will inspect closely the time required for development of a stable peak for each resonance.

Figure 5

Amplitudes of the local maxima near the ordinary Fano resonance (blue solid line) and the anomalous Fano feature (black solid line) as a function of time $\tilde{t}=t-t_0$ [$\hbar /\Gamma _{N}$]. The corresponding fitted exponential curves are shown by dashed red and orange lines. Results obtained for $t^{0}_{12}=0.4\Gamma _{N}$. Other parameters are the same as in Fig. 3.

Quantified analysis of time evolution of the resonant features

Table 1 Relaxation factors $\tau _F$ and $\tau _{AF}$ obtained for several couplings $t^{0}_{12}$ and corresponding equilibrium times $t^{F}_{eq}$, $t^{AF}_{eq}$.

We noticed that evolution of the OF and the AF peaks in time resembles the exponential growth (cf. Fig. 5 for $t^{0}_{12}=0.4\Gamma _{N}$). In order to quantify the time required for formation of the stable resonant features, we fit the exponential function describing decay of a difference between the actual and the final value of maxima. This function is defined as

$$\begin{aligned} \bar{G}(t) = \bar{G}(\infty ) -[\bar{G}(\infty )- \bar{G}(t_0)] \exp { \left[ (t_0-t )/ \tau \right] }, \end{aligned}$$

(16)

where $\bar{G}(\infty )$ is the peak amplitude in the static conditions (i.e., at $t\rightarrow + \infty$), $\bar{G}(t_0)$ is a value of the local maxima at the initial moment of the quench (at $t_0$), and the fitted parameter $\tau$ represents the interval, in which a mismatch between the initial conductivity and the equilibrium conductance diminishes $e \approx 2.71$ times. The calculations have been performed for a few amplitudes of the interdot coupling strengths ranging from $t^{0}_{12}=0.3\Gamma _{N}$ to $t^{0}_{12}=0.6 \Gamma _{N}$. We find that the characteristic time parameter $\tau _{AF}$ obtained for the anomalous Fano is larger than its counterpart $\tau _{F}$ obtained for the ordinary Fano resonance. The difference $\tau _{AF}-\tau _{F}$ is larger for weak interdot coupling $t^0_{12}$. Data of obtained $\tau$ for both resonances are collected in the second and third column of Table 1. The dependencies of both $\tau$’s as a function of the interdot coupling strength are presented in Fig. 6.

Figure 6

Dependencies of the equilibrium times on the interdot coupling. Times $\tau _F$ and $\tau _{AF}$ (in the units of $\hbar /\Gamma _{N}$) as a function of the interdot coupling $t^{0}_{12}$ obtained for the ordinary Fano resonance (blue line) and the anomalous Fano feature (red line), respectively.

Time of achieving the equilibrium conditions $t_{eq}$ could be considered as time after which the difference between the initial and the final amplitudes of a given resonant feature decrease by $95\%$. In order to deliver information on a typical timescale of the resonances in tangible units (i.e., nanoseconds), in the table, we present also $t_{eq}$ assuming that the dot-lead coupling $\Gamma _{N}$ is equal to 50 $\mu$eV. The energy scale of this magnitude is typical for experiments considering nanoscale objects coupled with superconductors^54,55.

The Fano-type interference peaks are characterized by asymmetric line-shapes originating from the close coexistence of resonant enhancement and resonant suppression of transmission. The function resembling the Fano-like shapes reads as

$$\begin{aligned} F(\omega )=\alpha (\omega ) \frac{(q \Gamma _{K}/2+\omega -\omega _{res})^2}{(\Gamma _{K}/2)^2+(\omega -\omega _{res})^2} \end{aligned}$$

(17)

where $\alpha (\omega )$ is the background transmission (which, in our case, is the undisturbed Andreev feature and can be obtained from Eqs. (13–15) at $t=t_0$, i.e., before the quench), $\omega _{res}$ is the energy of the resonant (quasi-discrete) level, $\Gamma _K$ is the broadening of the resonant level and q is the so-called asymmetry parameter (cf. also Ref.¹⁴ as well as original works of the resonance on a flat background^56,57). One of the crucial parameters describing the Fano shape is an asymmetry parameter q. If $q=0$, the Fano function resembles symmetric deep, whereas, for $q \rightarrow \infty$, the Fano function develops into the Breit-Wigner (Lorentz) distribution. It is interesting to investigate the time evolution of the asymmetry parameter q. In order to study time dependence of q parameter, the Fano function (17) is fitted to the features evolving around $eV=\varepsilon _{2}$ after the quench. From such a fit, the values of q are extracted. For this purpose, we use the procedure developed and described in earlier work¹⁴. Such fitting makes sense if the original shape of the transmission has well-pronounced both local minima and maxima. At the moment of the quench, the shape of conductivity near the resonant level (i.e., at $eV\approx \Gamma _{N}$) is given by a smooth monotonous function, which does not resemble the Fano shape. Thus, short after the quench only small fluctuations appear near the resonant level and it is not possible to properly fit the curve from formula (17), cf. Figure 4 for $\tilde{t} = 10$ $\hbar /\Gamma _{N}$. Therefore, we inspect the evolution of the asymmetry parameter from the moment where local maxima and minima are well distinguished from small disturbances. For the set of parameters used in Fig. 7 ($t^{0}_{12} = 0.4 \Gamma _{N}$), this time was estimated for about $7\ \hbar / \Gamma _{N}$ after the quench.

Figure 7

The evolution of the ordinary Fano resonant feature in time. The four panels on the left show the comparison of the resonant feature (solid red lines) with the best fit of the Fano shape curve (blue dashed lines) obtained for several different times ($\tilde{t} =t-t_{0}$, as labeled). The right panel presents the time evolution of the parameter q of the fitted Fano curve to the data (solid green line). The orange dashed line shows the corresponding fitted exponential curve with $\tau _q$ parameter. Results obtained for $t^{0}_{12}=0.4\Gamma _N$. Other parameters are the same as in Fig. 3.

On four left panels of Fig. 7, we present the actual shape of conductivity near $eV \simeq \varepsilon _2$ (solid red lines) obtained for the same set of parameters as in Fig. 5 for times $\tilde{t} = 30,\ 50,\ 80,\ 130$ $\hbar /\Gamma _N$ after quenching (as labelled) and the fitted Fano functions according to Eq. (17) (dashed blue lines). The right panel of Fig. 7 shows the evolution of asymmetry parameter q. The relation between q and the value of $F_{max} \equiv F(\omega _0)$ at the maximum of $F(\omega )$ at $\omega =\omega _0$ is given by

$$\begin{aligned} q =\pm \frac{\sqrt{F_{max}-\alpha _{max}}}{\sqrt{\alpha _{max}}}, \end{aligned}$$

(18)

where we assumed that $\alpha (\omega )$ is almost not dependent on $\omega$ in the neighbourhood of $\omega _0$ and $\alpha (\omega ) \approx \alpha _{max} \equiv \alpha (\omega _0)$ near $\omega _0$ (cf. also equation (10) from Ref.¹⁴). As we showed in Fig. 5 with fitting of Eq. (16), the local maxima of the ordinary Fano resonant features evolve exponentially. Consequently, the asymmetry parameter grows exponentially as well, but with different (reduced) characteristic time $\tau _{q}$. For the same set of parameters as used for the analyses presented in Fig. 5, the corresponding time for the asymmetry parameter was estimated to $\tau _{q}=13.76\ \hbar /\Gamma _N$ (with the fitting function analogous to that given in (16)). It differs from the characteristic time for the Fano feature maxima found as $\tau _{F}=17.06\ \hbar /\Gamma _N$, but there is no simple relation between $\tau _{q}$ and $\tau _{F}$ because of the formula (18). A small decrease of q, for $\tilde{t}> 90\ \hbar /\Gamma _{N}$, is caused by minor oscillations of the local maxima.

Although switching on and off the interdot coupling in considered time scales would be difficult to realize experimentally, effectively the quench protocol presented in this work could be realized electrostatically by gate voltage applied to the side dot (cf., e.g., Refs.^35,41). In such a case, the quench protocol might assume to firstly set the gate voltage, such that energy level of $QD_2$ lies beyond considered energy scale. At a given moment the gate potential should be changed to a desired energy (within considered energy scale) allowing electron scattering. Our calculations on this type of quench protocol indicate that the time scale for formation of both resonances remains the same as establishing an abrupt connection between dots.

Conclusions

In this work, we estimated the time required for formation of ordinary Fano resonance and its superconducting response on quantum dot ($QD_{1}$) region upon abrupt connection of additional quantum dot ($QD_2$). We found that upon abrupt interdot connection nonequilibrium SC response for scattering has its instantaneous effect pronounced by a high-magnitude asymmetric peak observables in differential conductivity. Careful inspection of their amplitudes relatively to static results reveals that stable solution for the anomalous feature develops longer than the feature representing direct scattering. The difference in time of this “saturation” decreases with increasing of the coupling between dots. A comparison of a time scale to typical values of the dot-lead coupling in experiments on quantum dot – superconductor hybrids shows that time of reaching equilibrium ranges from few to few hundreds of nanoseconds depending on strength of coupling to the Fermi sea and interdot connection. Comparing our results to those obtained when $QD_2$ was substituted by one end of topological chain hosting Majorana particles⁴¹ we found that establishing of both ordinary and anomalous resonances takes considerably more time than it is required for Majorana mode to leak into the region of QD.