text:
# Weak Charge Quantization on Superconducting Islands

## I Introduction

Physics of interacting electronic systems in presence of disorder has been a subject of an intense study already for a few decades . Various theoretical approaches have been developed for the description of both metallic and insulator phases. The non–linear $`\sigma `$–model (NL$`\sigma `$M) in the replica  or dynamic (Keldysh)  formulation has proven to be the most powerful tool to deal with the weakly disordered (metallic) phase. It was shown that both versions of the $`\sigma `$–model may effectively treat the perturbation theory as well as the renormalization group (RG) formalism. Unlike the non–interacting models, where a whole specter of non–perturbative results is available , it was relatively little progress in the development of non–perturbative solutions of the interacting NL$`\sigma `$M. Our goal is to make a step in this direction, using the Coulomb blockade (CB) on a superconductive (SC) quantum dot as a prototypical example.

The Coulomb blockade on a quantum dot coupled to a certain number of leads proved to be an extremely rich and fascinating phenomena both theoretically and experimentally (see Refs.  for review). From the theoretical point of view it provides a model, where the Coulomb interactions, being spatially localized, may be treated in a non–perturbative way. The interactions strongly affect charge (and spin) fluctuations between the dot and the leads, that manifests in the peculiar transport and thermodynamic behavior of the coupled dot–lead system. In the case of practically isolated dot, such that the dimensionless conductance (measured in units of $`G_Qe^2/(2\pi \mathrm{})`$) of the dot–lead interface is small ($`G1`$), the fluctuations usually may be taken into account perturbatively. (The notable exceptions are provided by the vicinity of the charge degeneracy point and by the Kondo effect on the dot ). Here we concentrate on the opposite scenario of the dot strongly connected to the leads ($`G>1`$). In this case the Coulomb blockade is expected to be suppressed by the charge fluctuations and the overall effect of interactions on a single dot to be weak. It is, however, a challenging theoretical problem, lacking an obvious small parameter, to understand a remnants of the CB on a strongly connected dot. More importantly the weak CB may prove to be a strong phenomena in the granulated systems, with many (superconductive) dots incorporated in the conducting matrix.

For the normal dot the weak CB effect was relatively well understood from the various points of view. Matveev  gave a complete picture of the phenomena for the case of one or two conducting channels connecting the dot and the leads. For the multi–channel system the CB suppression was calculated with the exponential accuracy employing RG technique , the instanton calculus  and bosonization . The factor, suppressing the thermodynamic CB oscillations, was formulated by Nazarov  in terms of transmission coefficients, $`T_\alpha `$, of the dot–lead interface as $`_\alpha (1T_\alpha )^{1/2}`$, where $`\alpha =1\mathrm{}N`$ and $`N`$ is the number of spinless channels. The remarkable feature of this result is that the essentially many–body phenomena may be described via the single–particle (non–interacting) scattering matrix only. Yet the knowledge of the interface conductance $`G=_\alpha T_\alpha `$ alone is not sufficient to describe the CB oscillations. For example, in a dot coupled to a lead through the tunnelling barrier with $`T_\alpha 1`$, one obtains for the CB suppression factor $`\mathrm{exp}\{G_T/2\}`$. Another important case of a dot–lead interface is a coupling via a piece of diffusive metal. Employing Dorokhov statistical distribution  of the transmission coefficients $`P(T)=G_D/(2T\sqrt{1T})`$, where $`G_D`$ is the conductance of the diffusive area, one obtains for the typical CB suppression factor $`\mathrm{exp}\{\pi ^2G_D/8\}`$. Since the result is exponentially sensitive – the difference between e.g. $`G/2`$ and $`\pi ^2G/8`$ may be actually enormous.

The superconductive dot in contact with normal leads is even more challenging system. We consider the low temperatures $`T\mathrm{\Delta }`$, where $`\mathrm{\Delta }`$ is the SC gap, and therefore only very few (if at all) quasi–particle excitations are allowed to leave or enter the superconductor. The dominant mechanism of the charge transfer through the interface is thus the Andreev transport of the correlated pairs. The Andreev transmission of a given channel is known to be $`𝒯_\alpha =T_\alpha ^2/(2T_\alpha )^2`$ . For the tunnelling barrier setup ($`T_\alpha 1`$) this leads to an overall Andreev conductance that scales like $`G_AG_T^2/N`$ (pair tunnelling probability). In most cases this is a very small number. The presence of the diffusive normal metal adjacent to the tunnelling barrier increases the Andreev conductance to $`G_AG_T^2/G_D`$ (or, in the case $`G_T>G_D`$, to $`G_AG_D`$). The physical reason of this phenomena is multiple attempts of Andreev transmission due to the coherent back–scattering on the normal impurities. The natural question is whether it is $`G`$ or $`G_A`$ (or may be none of them) that determines the amplitude of the CB oscillations. In case $`G_A`$ is the relevant quantity, one may wonder whether the coherent back–scattering enhancement should be taken into account. The answers are not immediately obvious, since while the SC strongly prefers the pair transport – the Coulomb energy of the dot makes the entrance of two charges at once energetically costly. The interactions of the dot may also provide a dephasing mechanism which ruins coherent back–scattering. We show that parametrically it is indeed $`G_A`$ (including diffusive enhancement) which determines the CB amplitude. The coefficient in the exponent, however, is setup dependent, making the final answer sensitive to the ratio $`G_T/G_D`$ and not only to $`G_A`$. The effective charging energy turns out to be proportional to $`_\alpha (1𝒯_\alpha )^{1/2}`$, where $`𝒯_\alpha `$ are the Andreev transmissions defined above. This expression provides a remarkable analogy between the normal and Andreev results.

The other interesting phenomena brought by superconductivity is the parity effect . At very small temperature the pair tunnelling is the dominant mechanism and the period (in the gate voltage) of the CB oscillations is twice larger than in the normal dot. For a closed dot the parity phenomena is destroyed by the entropic effects at moderately small temperature $`T^{}=\mathrm{\Delta }/\mathrm{ln}(\mathrm{\Delta }/\delta )\mathrm{\Delta }`$, where $`\delta `$ is the mean single–particle level spacing on the dot . The physical reason for this temperature to be much less than $`\mathrm{\Delta }`$ is that it is enough to have a single excited quasiparticle to destroy the parity effect. At larger temperature the system exhibits the “normal” oscillation period. We show that for an open dot the transition temperature between the normal and doubled period is somewhat larger than for a closed dot and is given by $`T^{}=\mathrm{\Delta }/\mathrm{ln}(\mathrm{\Delta }/G_A\delta )`$, provided $`G_A\delta <\mathrm{\Delta }`$.

Technically we treat the problem from the two seemingly distinct perspectives. First we look for the spatially–dependent instanton solution of interacting NL$`\sigma `$M on SC dot in contact with the normal diffusive region. The finite action of such instanton configuration results in the exponential suppression of CB amplitude. Alternatively we can treat the problem employing the functional RG technique. The latter approach encodes the entire dot–lead interface in an (infinite) set of coefficients, whose values are subsequently renormalized by the quantum fluctuations of the phase on the dot. Renormalization is terminated at the cut–off energy scale where the conductance reaches unity. The renormalized (exponentially small) cut–off energy is the effective charging energy $`\stackrel{~}{E}_C`$, this dictates the amplitude of the CB oscillations. One of the messages of the present paper is that these two approaches lead to the identical results.

The paper is organized as follows: in section II we describe the setup and formulate the appropriate action in the interacting NL$`\sigma `$M language. Section III is devoted to the instanton treatment of the problem in the real space. In section IV we derive the proximity action functional and obtain CB amplitude for an arbitrary set of transmission eigenvalues (using instanton approach) at moderately low temperatures $`\stackrel{~}{E}_CTE_C`$. Then in section V the RG approach is employed to treat the CB in the same temperature range. It is demonstrated that the results coincide with those obtained by instanton techniques. Finally in section VI we discuss the physical results and their possible experimental signatures. Appendix contains derivation of the instanton action starting from the real-time Keldysh functional technique.

## II Problem Setup and Action

We consider a large diffusive (or chaotic) SC dot. The mean single–particle level spacing of the dot, $`\delta `$, is supposed to be the smallest energy scale in the problem. The SC gap, $`\mathrm{\Delta }`$, on the other hand, is the largest scale. The electrons on the dot interact via the capacitive interaction of the form

$$H_{int}=E_C(\widehat{N}q)^2,$$
(1)
where $`\widehat{N}`$ is the electron number operator and $`q`$ is the rescaled dimensionless gate voltage potential. The charging energy, $`E_C=e^2/(2C)`$, is assumed to satisfy the inequalities $`\delta <E_C<\mathrm{\Delta }`$.

The dot is separated from the normal diffusive metal by the tunnelling barrier with the conductance $`G_T`$, see Fig. 1. The piece of the quasi 1D or 2D metal having size $`L`$ and conductance $`G_D`$ is in turn connected to a clean bulk 3D lead. By the reasons which are explained later we shall assume that the Thouless energy of the diffusive region, $`E_{\mathrm{Th}}=D/L^2`$, where $`D`$ is the diffusion coefficient, is larger than the charging energy of the dot, $`E_{\mathrm{Th}}>E_C`$. The opposite limiting case requires a separate treatment and will be presented elsewhere.

We shall be interested in those thermodynamic characteristics of the dot, which depend in an oscillatory way on the gate voltage potential, $`q`$. Specifically, we look for the free energy $`F(q)=T\mathrm{ln}Z(q)`$, where $`Z(q)`$ is the partition function. Both of these quantities depend on a particular realization of disorder in the diffusive region and inside the dot. We shall therefore look for the disorder averaged quantity like $`F(q)`$. For an open dot $`F(q)`$ is an oscillatory function, with a phase sensitive to a disorder realization. It is therefore convenient sometimes to calculate also the correlation function $`F(q)F(q^{})`$ . The latter carries information about the typical $`F(q)`$, rather than the average one (that may spuriously vanish due to phase randomness).

The disorder averaging may be performed in two ways either by introducing $`p`$ replica of the system  and sending $`p0`$ in the end, or by dealing with the dynamical Keldysh formulation . We shall employ both of these approaches to demonstrates what they are consistent and interchangeable. In the replica formalism, the NL$`\sigma `$M is formulated in terms of the matrix field theory of the matrix field $`Q_{ij}^{ab}(𝐫,\tau ,\tau ^{})`$, where $`a,b=1,2\mathrm{}2p`$ are the replica indexes (one needs $`2p`$ replica to describe the correlation function $`Z^p(q)Z^p(q^{})`$) and $`i,j=1,2`$ are Gorkov–Nambu indexes. The correlation function of the free energies is given by $`F(q)F(q^{})=lim_{p0}p^2(Z^p(q)Z^p(q^{})1)`$. The $`Q`$–matrix obeys the constraint $`Q^2=1`$ and the fermionic anti–periodic boundary condition in both of its imaginary time arguments $`\tau ,\tau ^{}[0,\beta ]`$. The coordinate $`𝐫`$ runs over the volume of the SC dot and the normal diffusive region. The matrix $`Q`$–field describes dynamics of the electrons; it is coupled to a scalar bosonic vector field $`\mathrm{\Phi }_a(\tau )`$ which originates from the Coulomb interactions on the dot. Since we restrict ourselves by the simplest capacitive interaction Eq. (1), the field $`\mathrm{\Phi }_a(\tau )`$ is space independent throughout the dot and vanishes outside the dot. As a result, only spatially independent ($`𝐪=0`$) component of the superconductive $`Q`$–matrix of the dot appears to be coupled to the Coulomb field $`\mathrm{\Phi }`$. Therefore we assume that the $`Q`$–matrix is spatially constant inside the dot, while it may have a non–trivial spatial structure inside the normal diffusive region. The effective action for our geometry contains three terms

$$S=S_{dot}+S_T+S_D,$$
(2)
where $`S_{dot}`$ and $`S_D`$ are the effective actions of the isolated SC dot and normal diffusive region correspondingly, while $`S_T`$ is the tunnelling action which provides the coupling between the two. We shall examine separately these three contributions.

### A Action of the dot

According to our model – the dot is the only region, where electrons interact via the Coulomb interaction Eq. (1). As a result, the dot’s action contains the two coupled fields: $`\mathrm{\Phi }_a(\tau )`$ and $`Q_S^{ab}(\tau ,\tau ^{})`$. The first one originates from the Hubbard–Stratonovich decoupling of the Coulomb term , while the latter is the result of disorder averaging and integration of the fermionic degrees of freedom on the dot . Both of these fields are spatially independent, reflecting the fact that the Thouless energy of the dot is assumed to be large. The action takes the standard by now form:

$$S_{dot}[Q_S,\mathrm{\Phi }]=\underset{a=1}{\overset{2p}{}}\underset{0}{\overset{\beta }{}}𝑑\tau \left(\frac{\dot{\mathrm{\Phi }}_a^2}{4E_c}iq^a\dot{\mathrm{\Phi }}_a\right)\frac{\pi }{\delta }Tr\left\{\left(\sigma _3_\tau i\dot{\mathrm{\Phi }}+\widehat{\mathrm{\Delta }}\right)Q_S\right\}.$$
(3)
Here the replica vector of gate voltages, $`q^a`$, is defined in a way that $`q^a=q`$ for $`a[1,p]`$ while $`q^a=q^{}`$ for $`a[p+1,2p]`$. The SC order parameter $`\widehat{\mathrm{\Delta }}=(\mathrm{\Delta }_a\sigma _+\mathrm{\Delta }_a^{}\sigma _{})\delta ^{ab}`$, is the same in all replica and is written as a matrix in the Gorkov–Nambu space, where $`\sigma _\pm =(\sigma _1\pm i\sigma _2)/2`$. Notice, that the gate voltage dependence of the partition function originates entirely from $`\mathrm{exp}\{i\pi _aq^aW_a\}`$, where $`\pi W_a𝑑\tau \dot{\mathrm{\Phi }}_a(\tau )=\mathrm{\Phi }_a(\beta )\mathrm{\Phi }_a(0)`$ is the zeroth Matsubara component of the $`\dot{\mathrm{\Phi }}_a(\tau )`$ field.

Since all the energy scales we consider are larger than $`\delta `$, one may evaluate the second term in the dot’s action (3) in the saddle point approximation over the field $`Q_S`$. In so doing, one disregards the mesoscopic conductance fluctuations of the dot–lead interface . The saddle point value of the $`Q_S`$ field is given by the Gorkov Green function gauged by the phase $`\mathrm{\Phi }(\tau )`$:

$$Q_S^{ab}(\tau ,\tau ^{})=e^{i\sigma _3\mathrm{\Phi }_a(\tau )}\mathrm{\Lambda }_S^{ab}(\tau ,\tau ^{})e^{i\sigma _3\mathrm{\Phi }_b(\tau ^{})}.$$
(4)
The Gorkov Green function, $`\mathrm{\Lambda }_S`$, has the standard form, which is more familiar in the Matsubara basis (we assume the phase of the SC dot without the Coulomb interactions to be zero)

$$\mathrm{\Lambda }_S^{ab}(n,m)=\delta ^{ab}\delta _{nm}\left(\begin{array}{cc}\mathrm{cos}\theta _n& \mathrm{sin}\theta _n\\ \mathrm{sin}\theta _n& \mathrm{cos}\theta _n\end{array}\right),$$
(5)
where

$$\mathrm{cos}\theta _n\frac{ϵ_n}{\sqrt{ϵ_n^2+|\mathrm{\Delta }|^2}};\mathrm{sin}\theta _n\frac{|\mathrm{\Delta }|}{\sqrt{ϵ_n^2+|\mathrm{\Delta }|^2}}$$
(6)
and $`ϵ_n=\pi T(2n+1)`$ is the fermionic Matsubara frequency. The phase rotation (4) preserves the fermionic anti–periodic boundary conditions if all $`W_a`$ are even integers. At small temperatures, $`T\mathrm{\Delta }`$, however, one has $`\mathrm{cos}\theta _n0`$ and $`\mathrm{sin}\theta _n1`$ and therefore the Gorkov matrix, $`\mathrm{\Lambda }_S`$ is (almost) off–diagonal in the Nambu space and local in time. As a result, only $`2\mathrm{\Phi }(\tau )`$ participate in the phase rotation, Eq. (4). Therefore the odd–integer $`W_a`$ preserve the fermionic boundary conditions, as well. The odd–integers $`W_a`$ result in the doubling of the period of the $`F(q)`$ function with respect to the normal case. This reflects the fact that the pair transfer is the dominant mechanism of the charge exchange.

In fact, one has to be more careful and recall that, according to Eq. (3), one has to perform integrations over all Matsubara components of $`\dot{\mathrm{\Phi }}_a(\tau )`$ fields . All non–zero Matsubara components may be eliminated by the gauge transformation, Eq. (4). As a result, they have no effect on the thermodynamic of an isolated dot at all (they are of major importance, of course, once the dot is coupled to the leads). The remaining (usual) integral over the zeroth Matsubara component, $`W_a`$, must be performed explicitly. To this end one notices that $`\pi TW_a`$ enters the action as an imaginary chemical potential in the replica $`a`$. The (replicated) free energy of an isolated dot is thus a periodic function of each of $`W_a`$ with the period $`2`$ (indeed the chemical potential always enters as $`\mathrm{exp}\{\mu \widehat{N}/T\}`$, since the number operator, $`\widehat{N}`$, has only integer eigenvalues – the periodicity is apparent). The free energy possesses deep minima at even integer values of $`W_a`$ with the quadratic behavior in their vicinity $`F_{dot}(W)\pi ^2T^2W^2/(2\delta )`$, where $`\delta =(^2F_{dot}/\mu ^2)^1`$ is the mean level spacing of the dot. For a sufficiently large dot, where $`\delta <T`$, the integrals over $`W_a`$ may be performed in the saddle point approximation, which results in the even–integer quantization of $`W_a`$. In a SC dot there are additional minima at odd–integer values of $`W_a`$. Consequently the integration over $`W_a`$ is substituted by the summation over all integers. At $`T>0`$ the additional minima of the free energy at odd integer $`W_a`$ are not as deep as at the even integers, reflecting the fact that the addition of an odd number of electrons is possible by creating a quasi–particle. We shall evaluate now the action cost of the odd minima with respect to the even ones.

For the even values, say $`W_a=0`$, one substitute the saddle point solution, Eq. (4), into the action, Eq. (3), and obtains (we disregard for a moment the first term in Eq. (3))

$$S_{dot}(W_a=0)=\frac{\pi }{\delta }Tr\{(\sigma _3_\tau +\widehat{\mathrm{\Delta }})\mathrm{\Lambda }_S^{aa}\}=\frac{2\pi }{\delta }\underset{n}{}\sqrt{ϵ_n^2+|\mathrm{\Delta }|^2}.$$
(7)
This sum is divergent. However, it is only the difference of the action between different replica that has a physical significance. The later quantity is convergent as we shall see momentarily. For the odd integers, say $`W_b=1`$, there is an imaginary component $`i\pi T`$ of the chemical potential in Eq. (3). It may be eliminated by the gauge transformation which converts the antiperiodic boundary conditions for the fermions into the periodic one. As a result, one arrives to the same expression as Eq. (7) with the fermionic Matsubara frequency $`ϵ_n=\pi T(2n+1)`$ substituted by the bosonic one $`\omega _n=2\pi Tn`$:

$$S_{dot}\left(W_b=1\right)=\frac{2\pi }{\delta }\underset{n}{}\sqrt{\omega _n^2+|\mathrm{\Delta }|^2}.$$
(8)
Employing the Poisson summation formula, one obtains for the difference

$$𝒜(T)S_{dot}(1)S_{dot}(0)=\frac{8|\mathrm{\Delta }|}{\delta }\underset{l=0}{\overset{\mathrm{}}{}}\frac{1}{2l+1}K_1\left((2l+1)\frac{|\mathrm{\Delta }|}{T}\right)\frac{4\sqrt{2\pi |\mathrm{\Delta }|T}}{\delta }e^{|\mathrm{\Delta }|/T},$$
(9)
where the last equality assumes $`T\mathrm{\Delta }`$. In the opposite limit of the normal dot, $`T>\mathrm{\Delta }`$, one finds $`𝒜(T)\pi ^2T/\delta `$, reflecting the fact that in a normal dot the odd–integer minima are absent. Although these minima persist up to $`T\mathrm{\Delta }`$, their contribution is exponentially suppressed at $`T>T^{}`$, where $`T^{}`$ is determined from the condition $`S_{dot}(1)S_{dot}(0)1`$

$$T^{}\frac{|\mathrm{\Delta }|}{\mathrm{ln}|\mathrm{\Delta }|/\delta }|\mathrm{\Delta }|.$$
(10)
As a result, there is an important temperature dependence associated with the odd–integer winding numbers at the scale $`TT^{}|\mathrm{\Delta }|`$. We shall see below that for the dot strongly coupled to the leads, the corresponding temperature scale is slightly different.

We summarize now our discussion of the action of a large SC dot with the Coulomb interaction at $`T\mathrm{\Delta }`$. The scalar potential in each replica obeys the boundary condition

$$\mathrm{\Phi }_a(\beta )\mathrm{\Phi }_a(0)=\pi W_a,$$
(11)
where $`W_a`$ is an integer winding number. The corresponding saddle point value of the $`Q_S`$–matrix field is given by Eq. (4). The action of the dot takes the form:

$$S_{dot}[\mathrm{\Phi }]=\underset{a=1}{\overset{2p}{}}\left[i\pi q^aW_a+𝒜(T)\delta _{W_amod2,1}+\underset{0}{\overset{\beta }{}}𝑑\tau \frac{\dot{\mathrm{\Phi }}_a^2}{4E_C}\right],$$
(12)
where $`𝒜(T)`$ is given by Eq. (9).

### B Tunnelling barrier and diffusive region action

The tunnelling action couples the $`Q_S`$–field on the dot with the $`Q(𝐫=0)`$–field at the point $`𝐫=0`$ adjacent to the tunnelling barrier from the normal metal side. It has the standard form

$$S_T=\frac{G_T}{8}Tr\{Q_SQ(0)\},$$
(13)
where $`G_T`$ is the tunnelling conductance.

The action of the normal diffusive region also has the standard form

$$S_D=\frac{\pi \nu }{4}\underset{0}{\overset{L}{}}𝑑𝐫DTr\{(Q(𝐫))^2\},$$
(14)
where $`D`$ is the diffusion constant and $`\nu `$ the density of states (per spin) of the normal region. The total conductance of the normal region is given by $`G_D=4\pi \nu D/L`$ for the quasi 1D geometry of the normal region and $`G_D=8\pi ^2\nu D/\mathrm{ln}(L/d)`$ for the 2D geometry. Here $`L`$ is the length (radius) of the 1D (2D) region and $`d`$ is the radius of the SC dot. We have omitted the frequency term $`\nu Tr\{ϵQ\}`$ on the r.h.s. of Eq. (14), because of the assumption that $`E_{\mathrm{Th}}>E_C>T`$. At the point where the diffusive region is attached to the normal bulk lead one has to impose the boundary condition

$$Q(𝐫=L)=\mathrm{\Lambda }_N,$$
(15)
where $`\mathrm{\Lambda }_N^{ab}(n,m)=\sigma _3\delta ^{ab}\delta _{nm}\text{sign}(ϵ_n)`$ is the appropriate $`Q`$–matrix of the normal bulk lead.

Alternatively, one may imagine integrating out the $`Q(𝐫)`$–field of the normal region, subject to the boundary condition Eq. (15) and weighted by the action $`S_T+S_D`$. This procedure (we shall describe it in details in section IV) leads to the effective action $`S_{TD}`$ of the interface plus diffusive region written in terms of $`Q_S`$ and $`\mathrm{\Lambda }_N`$ only. If all the relevant energy scales are less than the Thouless energy of the diffusive metal – the general form such action may take is

$$S_{TD}=\frac{G_D}{8}\underset{l=1}{\overset{\mathrm{}}{}}\gamma _lTr\{(Q_S\mathrm{\Lambda }_N)^l\},$$
(16)
where $`\gamma _l`$ are coefficients which depend on the details of the interface (in our case the ratio $`G_T/G_D`$). The largeness of the Thouless energy is necessary to disregard the retardation effects and thus possible time non–local coupling between $`Q_S`$ and $`\mathrm{\Lambda }_N`$. Under such condition, the proximity action, Eq. (16), is completely equivalent to those given by Eqs. (13) and (14) upon the proper choice of the set, $`\gamma _l`$ .

## III The Instanton Approach

We are interested in the limit of strong coupling between the dot and the leads, meaning $`G_D,G_T>1`$ (the weak coupling limit may be treated in the spirit of Refs. ). For $`G_D1`$ the fluctuations of the $`Q`$–field around its optimal value are suppressed. One may employ therefore the stationary phase treatment of the NL$`\sigma `$M for the dot–lead interface . Taking the variation of the action Eqs. (13), (14) under the condition $`Q^2=1`$ one obtains the Usadel equation

$$2\pi \nu (DQQ)\delta (𝐫)G_T[Q,Q_S]=0.$$
(17)
This equation is to be solved for a fixed $`Q_S=Q_S[\mathrm{\Phi }]`$ given by Eq. (4) and with the boundary conditions Eq. (15). The solution $`Q=Q[\mathrm{\Phi }]`$ after being substituted back into the action Eqs. (13), (14) – results in the semiclassical phase action $`S[\mathrm{\Phi }]`$. The later may then be investigated using the instanton approach applied to the $`\mathrm{\Phi }(\tau )`$–field.

### A Zero winding number

As a warm–up exercise, consider the zero winding number sector of the theory, $`W_a=0`$ for $`a=1,\mathrm{}2p`$. Obviously it does not produce an oscillatory dependence of $`F(q)`$, cf. Eq. (12), and therefore serves only an axillary purpose. The lowest energy configuration in zero winding number sector is $`\mathrm{\Phi }_a(\tau )=0`$ and therefore $`Q_S`$ on the dot is simply given by the BCS $`\mathrm{\Lambda }_S`$, Eq. (5). Solution of the Usadel equation (17) may be written as

$$Q(𝐫)=\mathrm{\Lambda }_N\mathrm{exp}\left\{iu(𝐫)\theta \sigma _2\right\},$$
(18)
where $`u(𝐫)`$ is the normalized “voltage drop” inside the normal region: $`u(r)=(Lr)/L`$ in 1D and $`u(r)=\mathrm{ln}(L/r)/\mathrm{ln}(L/d)`$ in 2D. The absolute value of the voltage drop, $`\theta `$, is coordinate independent diagonal (in replica and Matsubara space) matrix. It has a physical meaning of the SC rotation angle of the normal metal in the direct proximity to the dot. One may substitute the solution back into the action to obtain (for a single replica and Matsubara component)

$$S_0=\frac{1}{8}\left[G_D\theta ^22G_T\mathrm{cos}(\theta _n\theta )\right],$$
(19)
where the subscript “$`0`$” stresses that we work with zero winding numbers and $`\theta _n`$ is defined by Eq. (6). For small energy, $`ϵ_n\mathrm{\Delta }`$, one has $`\theta _n\pi /2`$ and therefore the corresponding action takes the form

$$S_0=\frac{G_D}{8}\left[\theta ^22t\mathrm{sin}\theta \right],$$
(20)
where $`tG_T/G_D`$. This action is minimized when $`\theta =\theta (t)`$ satisfies the equation

$$\theta =t\mathrm{cos}\theta .$$
(21)
The lowest energy solution of this equation smoothly interpolates between $`\theta =0`$ for $`t1`$ and $`\theta =\pi /2`$ for $`t1`$. Finally, this solution has to be substituted into Eq. (20) to find the action cost, $`S_0=S_0(t)`$, for the zero winding number configuration.

### B Non–zero winding numbers

To calculate the oscillatory component of the free energy, $`F(q)`$, one has to consider the $`\mathrm{\Phi }`$–field configurations with non–zero winding numbers, cf. Eq. (12). Consider, thus, the simplest even configuration of winding numbers with $`W_1=2`$ and all others $`W_a=0`$ in the remaining $`2p1`$ replica. With the exponential accuracy it is sufficient to consider the “straight” windings: $`\mathrm{\Phi }_1(\tau )=2\pi T\tau `$. The saddle point of the SC $`Q`$–field on the dot is given by Eq. (4) and takes the form

$$Q_S^{11}(n,m)=\left(\begin{array}{cc}\delta _{n,m}\mathrm{cos}\theta _{m1}\hfill & \delta _{n,m+2}\mathrm{sin}\theta _{m+1}\hfill \\ \delta _{n,m2}\mathrm{sin}\theta _{m1}\hfill & \delta _{n,m}\mathrm{cos}\theta _{m+1}\hfill \end{array}\right),$$
(22)
where the $`2\times 2`$ structure is the Nambu space. In all other replica except of $`a=1`$ the $`Q_S`$–matrix has the form Eq. (5). One may check that the Usadel equation is solved by exactly the same $`O(2)`$ rotation as in $`W=0`$ case, Eq. (18), performed in each of the following $`2\times 2`$ Nambu blocks

$$\left(\begin{array}{cc}Q_{11}(n+1,n+1)\hfill & Q_{12}(n+1,n1)\hfill \\ Q_{21}(n1,n+1)\hfill & Q_{22}(n1,n1)\hfill \end{array}\right)=\left(\begin{array}{cc}\mathrm{cos}\theta _n& \mathrm{sin}\theta _n\\ \mathrm{sin}\theta _n& \mathrm{cos}\theta _n\end{array}\right);$$
(23)
the Nambu indexes are explicitly stated on the l.h.s. In this $`2\times 2`$ block $`\mathrm{\Lambda }_N=\sigma _3`$ and therefore the matrix from the r.h.s of Eq. (23) may indeed be rotated into $`\mathrm{\Lambda }_N`$ as in Eq. (18). The only exception are the two lowest Matsubara components $`n=0`$ and $`n=1`$, which form $`4\times 4`$ block of the form

$$\left(\begin{array}{cccc}Q_{11}(1,1)\hfill & & Q_{12}(1,1)\hfill & \\ & Q_{11}(2,2)\hfill & & Q_{12}(2,0)\hfill \\ Q_{21}(1,1)\hfill & & Q_{22}(1,1)\hfill & \\ & Q_{21}(0,2)\hfill & & Q_{22}(0,0)\hfill \end{array}\right)=\left(\begin{array}{cccc}\mathrm{cos}\theta _0& & \mathrm{sin}\theta _0& \\ & \mathrm{cos}\theta _1& & \mathrm{sin}\theta _1\\ \mathrm{sin}\theta _0& & \mathrm{cos}\theta _0& \\ & \mathrm{sin}\theta _1& & \mathrm{cos}\theta _1\end{array}\right).$$
(24)
The corresponding $`4\times 4`$ block of the $`\mathrm{\Lambda }_N`$–matrix on the normal lead is the unit matrix. Obviously the unit matrix cannot be unitary rotated into the block Eq. (24) and therefore solution of the Usadel equation in this block is impossible. The difficulty originates from the fact that, due to the random phase of the CB oscillations, the average free energy, $`F(q)`$, is not an oscillatory function. We need therefore to consider winding number configuration of the form $`W_1=2`$ and $`W_{p+1}=2`$, while all other $`W_a=0`$. The contribution to the correlation function from such a configuration is proportional to $`\mathrm{exp}\{2\pi i(qq^{})\}`$, that is the lowest normal harmonic of the correlation function $`F(q)F(q^{})`$. The “dangerous” $`4\times 4`$ block in the $`p+1`$–st replica is given by

$$\left(\begin{array}{cccc}Q_{11}(1,1)\hfill & & Q_{12}(1,1)\hfill & \\ & Q_{11}(0,0)\hfill & & Q_{12}(0,2)\hfill \\ Q_{21}(1,1)\hfill & & Q_{22}(1,1)\hfill & \\ & Q_{21}(2,0)\hfill & & Q_{22}(2,2)\hfill \end{array}\right)=\left(\begin{array}{cccc}\mathrm{cos}\theta _0& & \mathrm{sin}\theta _0& \\ & \mathrm{cos}\theta _1& & \mathrm{sin}\theta _1\\ \mathrm{sin}\theta _0& & \mathrm{cos}\theta _0& \\ & \mathrm{sin}\theta _1& & \mathrm{cos}\theta _1\end{array}\right),$$
(25)
while the corresponding $`4\times 4`$ block of the $`\mathrm{\Lambda }_N`$–matrix is minus one times the unit matrix. Here as well the unitary rotation between the points $`𝐫=L`$ and $`𝐫=0`$ is impossible. However, if one combines $`a=1`$ and $`a=p+1`$ replica and allows rotation between them – then the unitary rotation may be found . Indeed, combining both “dangerous” blocks into the single $`8\times 8`$ block (i.e. combine Eq. (24) with Eq. (25) on the dot and unit matrix with the minus unit matrix on the normal lead), one readily see that they may be unitary connected (since they possess the same set of eigenvalues).

The calculations are simplified in the low–temperature case, $`T\mathrm{\Delta }`$, where $`\theta _n\pi /2`$. In this case the $`8\times 8`$ block takes the form

$$Q_S=\left(\begin{array}{cc}1\hfill & 0\hfill \\ 0\hfill & 1\hfill \end{array}\right)\left(\begin{array}{cc}\sigma _1\hfill & 0\hfill \\ 0\hfill & \sigma _1\hfill \end{array}\right)$$
(26)
on the dot and

$$\mathrm{\Lambda }_N=\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)\left(\begin{array}{cc}\sigma _0& 0\\ 0& \sigma _0\end{array}\right)$$
(27)
on the normal lead ($`\sigma _0`$ is the unit matrix in the Nambu space). Here the outer $`2\times 2`$ unit matrix represents the space of $`n=0`$ and $`n=1`$ Matsubara components, while the inner one – the replica space of $`a=1`$ and $`a=p+1`$ (Pauli matrices act in the Nambu space). We seek thus for the solution of the Usadel equation in the $`8\times 8`$ sub–space, having in mind that the solution for all other Matsubara components and replica is exactly the same as in $`W=0`$ sector, Eqs. (18)–(21). The general solution for $`𝐫0`$ may be written as

$$Q(𝐫)=\mathrm{\Lambda }_N\mathrm{exp}\left\{iu(𝐫)\left(\begin{array}{cc}0\hfill & B\hfill \\ B^{}\hfill & 0\hfill \end{array}\right)\right\},$$
(28)
where $`u(𝐫)`$ function was defined after Eq. (18) and $`B`$ is a coordinate–independent Nambu matrix. Employing the singular value decomposition, it may be written as

$$B=U^1\left(\begin{array}{cc}\mathrm{\Theta }_1\hfill & 0\hfill \\ 0\hfill & \mathrm{\Theta }_2\hfill \end{array}\right)V,$$
(29)
where $`U=e^{iu_3\sigma _3}e^{iu_2\sigma _2}e^{iu_1\sigma _1}`$ and $`V=e^{iv_3\sigma _3}e^{iv_2\sigma _2}e^{iv_1\sigma _1}`$ are $`SU(2)`$ matrices and $`\mathrm{\Theta }_1\mathrm{\Theta }_2`$ are real singular values. Substituting the solution back into the action, one obtains (for each of the two involved Matsubara frequencies)

$$2S_{\pm 2}=\frac{1}{4}\left[G_D(\mathrm{\Theta }_1^2+\mathrm{\Theta }_2^2)G_T(\mathrm{sin}2u_2\mathrm{sin}2v_2)(\mathrm{cos}\mathrm{\Theta }_1\mathrm{cos}\mathrm{\Theta }_2)\right],$$
(30)
where the subscript $`\pm 2`$ stays for the corresponding winding numbers and the coefficient two on the l.h.s. reminds that two replica were involved. The next step is to minimize the action over $`\mathrm{\Theta }_{1,2}`$, $`u_2`$ and $`v_2`$ angles. Three of the four equations for the minima have only trivial (parameter independent) solutions: $`\mathrm{\Theta }_1=0`$; $`u_2=\pi /4`$ and $`v_2=\pi /4`$. The action in terms of the single non–trivial angle $`\mathrm{\Theta }_2`$ finally takes the form

$$2S_{\pm 2}=\frac{G_D}{4}\left[\mathrm{\Theta }_2^2+2t(\mathrm{cos}\mathrm{\Theta }_21)\right],$$
(31)
where, as above, $`tG_T/G_D`$. The corresponding saddle point equation for $`\mathrm{\Theta }_2=\mathrm{\Theta }_2(t)`$ is

$$\mathrm{\Theta }_2=t\mathrm{sin}\mathrm{\Theta }_2.$$
(32)
For $`t1`$ the only solution of this equation is $`\mathrm{\Theta }_2=0`$, while for $`t1`$ the angle $`\mathrm{\Theta }_2(t)`$ interpolates between zero (for $`t=1`$) and $`\pi `$ (for $`t1`$).

### C CB suppression

We are now on the position to discuss the suppression of the CB. Consider first the component of the correlation function $`Z^p(q)Z^p(q^{})`$ which is proportional to $`\mathrm{cos}2\pi (qq^{})`$. As was explained above the relevant field configurations are those having a single replica with $`W_a=\pm 2`$, where $`a[1,p]`$ and a single replica with $`W_a=2`$, where $`a[p+1,2p]`$. The corresponding action is given by $`4𝒮2S_{\pm 2}+(2p2)2S_0`$; (here $`S_0`$ is multiplied by the number of replica with zero winding number – $`2p2`$, and by factor of two, because two Matsubara components are different between $`W=\pm 2`$ and $`W=0`$). Taking the replica limit $`p0`$ and employing Eqs. (20), (31), one finds

$$2𝒮(t)=S_{\pm 2}2S_0=\frac{G_D}{8}\left[\mathrm{\Theta }_2^2(t)2\theta ^2(t)+2t(\mathrm{cos}\mathrm{\Theta }_2(t)+2\mathrm{sin}\theta (t)1)\right],$$
(33)
where $`\theta (t)`$ and $`\mathrm{\Theta }_2(t)`$ are the solutions of Eqs. (21) and (32) correspondingly. As a result, the contribution to the correlation function with the unit period has the form: $`\mathrm{exp}\{4𝒮\}\mathrm{cos}2\pi (qq^{})`$, where the factor of $`4`$ stays for the fact that two replica are involved in the correlation function and the winding number is two. We shall discuss this result in more details in section VI, after deriving it using other methods. For the later reference we need to calculate $`t^2_t(2𝒮/t)`$; employing the saddle point equations (21) and (32), we obtain

$$t^2\frac{}{t}\left(\frac{2𝒮}{t}\right)=\frac{G_D}{8}\left[\mathrm{\Theta }_2^22\theta ^2\right].$$
(34)
As was explained in section II, the parity effects sets in at small temperature $`T<T^{}\mathrm{\Delta }`$. Technically it manifests itself in the appearance of the odd winding numbers. The calculations for the field configuration with $`W_1=1`$, and $`W_{p+1}=1`$ are exactly parallel to the one for the even winding numbers. The only difference is that there is the single “dangerous” Matsubara frequency, $`n=0`$, consequently there is no outer $`2\times 2`$ Matsubara structure as in Eqs. (26) and (27) (the replica and Nambu structures are exactly the same). As a result, the action is exactly twice smaller than for the oscillations with the unit period. We thus obtain for the component of the correlation function with doubled gate voltage period: $`\mathrm{exp}\{2(𝒮+𝒜(T))\}\mathrm{cos}\pi (qq^{})`$, where the factor of two in the exponent originates from the two replica involved.

Finally, let us discuss the technicalities of the replica limit, $`p0`$ . There are altogether $`p^2`$ distinct configurations of the winding numbers which contribute to the correlation function $`Z^p(q)Z^p(q^{})`$ with the same oscillatory factor, say $`\mathrm{exp}\{2\pi i(qq^{})\}`$. They correspond to the $`p`$ possible choices for $`W_a=2`$ and $`p`$ independent possible choices of $`W_a=2`$. There is thus $`p^2`$ combinatorial factor in the $`Z^p(q)Z^p(q^{})`$ that is cancelled against the same factor in the denominator for $`F(q)F(q^{})=lim_{p0}p^2(Z^p(q)Z^p(q^{})1)`$.

## IV Proximity action

In this Section we discuss the derivation and simple applications of the proximity action approach to treat the interface between the SC dot and normal lead. This method was recently suggested  for the analysis of superconductive-normal metal transition in a 2D proximity-coupled array . Notice that our definition of the running parameter $`\zeta `$ and charges $`\{\gamma _l\}`$ differ from those used in Refs. , see Appendix A for details.

### A Derivation of the proximity action

Our immediate goal is to derive the proximity action, Eq. (16), starting from Eqs. (13) and (14). To this end, let us imagine adding the diffusive metal step by step with infinitesimal shells having conductance $`\delta GG_D>1`$. At some intermediate step of this procedure one has the tunneling barrier with the part of the normal metal incorporated into the set of coefficients $`\gamma _l(\zeta )`$. Here $`\zeta G_D/G_D(\zeta )[0,1]`$ is the running parameter such that $`\gamma _l(0)=\delta _{l,1}G_T/G_D`$, and $`S(0)=S_T`$ is the bare tunneling barrier action; while $`\gamma _l(1)=\gamma _l`$ and $`S(1)=S_{TD}`$ is the full proximity action of the barrier “dressed” with the diffusive metal. At each step one adds another thin shell of the diffusive metal. The action of the entire system takes the form

$$S(\zeta +\delta \zeta )=\frac{G_D}{8}\underset{l=1}{\overset{\mathrm{}}{}}\gamma _l(\zeta )Tr\{(Q_SQ)^l\}+\frac{\delta G}{16}Tr\{(Q\mathrm{\Lambda }_N)^2\},$$
(35)
where $`Q`$ is the field right on the boundary between the new shell and the already integrated region. The action of the newly added shell has the simple form, $`(\delta G/16)Tr(Q\mathrm{\Lambda }_N)^2`$, since $`\delta GG_D>1`$ and therefore the matrix $`Q`$ has to be rather close to $`\mathrm{\Lambda }_N`$ (in the leading order in $`(\delta G)^1`$ one may disregard all higher order terms $`(Q\mathrm{\Lambda }_N)^l`$ with $`l>2`$). The next step is to integrate out the $`Q`$ field and obtain the new set $`\gamma _l(\zeta +\delta \zeta )`$. We perform this procedure in the saddle point approximation. Taking variation of the action, Eq. (35), under the condition $`Q^2=1`$ one finds the saddle point (Usadel) equation

$$G_D\underset{l=1}{\overset{\mathrm{}}{}}l\gamma _l(\zeta )\left[(Q_SQ)^l(QQ_S)^l\right]\delta G\left[Q\mathrm{\Lambda }_N\mathrm{\Lambda }_NQ\right]=0.$$
(36)
In the leading order in $`G_D/\delta G1`$ one obtains for the saddle point solution

$$Q=\mathrm{\Lambda }_N+\frac{G_D}{2\delta G}\underset{l=1}{\overset{\mathrm{}}{}}l\gamma _l\left[(Q_S\mathrm{\Lambda }_N)^l\mathrm{\Lambda }_N\mathrm{\Lambda }_N(Q_S\mathrm{\Lambda }_N)^l\right].$$
(37)
Putting it back into the action, one finally finds

$$S(\zeta +\delta \zeta )=\frac{G_D}{8}\left[\underset{l=1}{\overset{\mathrm{}}{}}\gamma _lTr\left\{(Q_S\mathrm{\Lambda }_N)^l\right\}+\frac{\delta \zeta }{2}\underset{l,k=1}{\overset{\mathrm{}}{}}lk\gamma _l\gamma _kTr\left\{(Q_S\mathrm{\Lambda }_N)^{l+k}(Q_S\mathrm{\Lambda }_N)^{|lk|}\right\}\right],$$
(38)
where we employed the fact that $`\delta G=G_D/\delta \zeta `$. One ends up with the following evolution equations for $`\gamma _l(\zeta )`$:

$$\frac{d\gamma _l}{d\zeta }=\frac{1}{2}\underset{k=1}{\overset{\mathrm{}}{}}k(l+k)\gamma _k\gamma _{l+k}\frac{1}{4}\underset{k=1}{\overset{l1}{}}k(lk)\gamma _k\gamma _{lk}.$$
(39)
To solve this set of equations it is convenient to define the function $`u(\zeta ,\mathrm{\Theta })\underset{l=1}{\overset{\mathrm{}}{}}l\gamma _l(\zeta )\mathrm{sin}(l\mathrm{\Theta })`$. With it’s help Eq. (39) takes the form

$$u_\zeta =uu_\mathrm{\Theta },$$
(40)
supplemented with the initial condition $`u(0,\mathrm{\Theta })=t\mathrm{sin}\mathrm{\Theta }`$, where, as before, $`tG_T/G_D`$. Employing the method of characteristics, one may write the solution of the latter problem in the implicit form:

$$u(\zeta ,\mathrm{\Theta })=t\mathrm{sin}(\mathrm{\Theta }\zeta u(\zeta ,\mathrm{\Theta })).$$
(41)
As a result, the function $`u_0(\mathrm{\Theta })u(1,\mathrm{\Theta })`$, that has $`l\gamma _l`$ as it’s Fourier coefficients, may be found as the solution of the following differential or algebraic equations

$$t^2\frac{}{t}\left(\frac{u_0}{t}\right)+u_0_\mathrm{\Theta }u_0=0;u_0(\mathrm{\Theta })=t\mathrm{sin}(\mathrm{\Theta }u_0(\mathrm{\Theta })).$$
(42)
The last expression solves the problem of writing the action of the barrier and the diffusive region in the form of Eq. (16). Note, that to derive this result we have not used a specific form of $`Q_S`$ and $`\mathrm{\Lambda }_N`$ (other than the fact that $`Q_S^2=\mathrm{\Lambda }_N^2=1`$). Therefore the action and the coefficients $`\gamma _l`$ are equally applicable to any other setup. The only two conditions that were employed are $`G_D1`$ to disregard localization effects and perform the saddle point calculations and $`E_{\mathrm{Th}}E_C`$ to disregard the frequency term in the action and justify the structure of Eq. (16). The condition $`E_{\mathrm{Th}}E_C`$ also allows to neglect Cooper-channel interaction in the normal conductor, while deriving proximity action (cf.  for the discussion of general situation).

For the later use we define $`\theta =\theta (t)u_0(\pi /2)`$ and $`\mathrm{\Theta }_2=\mathrm{\Theta }_2(t)u_0(\pi )`$. According to Eq. (42), these two angles satisfy $`\theta =t\mathrm{cos}\theta `$ and $`\mathrm{\Theta }_2=t\mathrm{sin}\mathrm{\Theta }_2`$ correspondingly. Comparing these relations with Eqs. (21) and (32), we conclude that the angles $`\theta `$ and $`\mathrm{\Theta }_2`$ introduced here coincide with the instanton angles introduced in section III.

### B Physical observables

In this subsection we repeat briefly derivation presented originally in .

Normal state conductance. Consider for a moment the dot in the normal state with a small applied bias $`V(\tau )=V\mathrm{cos}\omega _m\tau `$. The dot’s Green function takes the form $`Q_S=\mathrm{exp}\{i\mathrm{\Phi }(\tau )\}\mathrm{\Lambda }_N(\tau \tau ^{})\mathrm{exp}\{i\mathrm{\Phi }(\tau ^{})\}`$, where $`\mathrm{\Phi }=𝑑\tau V(\tau )`$. Substituting such $`Q_S`$ in the action and taking the second variation with respect to $`V(\tau )`$ at $`V=0`$, one finds for the linear normal state conductance

$$G_N=G_D\underset{l=1}{\overset{\mathrm{}}{}}l^2\gamma _l=G_D_\mathrm{\Theta }u_0(0).$$
(43)
For small $`\mathrm{\Theta }`$ one has $`u_0(\mathrm{\Theta })\mathrm{\Theta }_\mathrm{\Theta }u_0(0)`$; employing Eq. (42), one finds $`_\mathrm{\Theta }u_0(0)=t/(1+t)`$. As expected, one obtains the resistance addition rule

$$G_N^1=G_D^1+G_T^1.$$
(44)
Andreev conductance. We shall concentrate now on the temperature region $`\mathrm{\Delta }T\delta `$, where the SC $`Q`$–matrix may be written in the form Eq. (4) with $`\theta _n=\pi /2`$, that corresponds to infinite $`\mathrm{\Delta }`$ (the leading effect of the finite $`\mathrm{\Delta }`$ is to renormalize $`E_C^1E_C^1+G/2\pi \mathrm{\Delta }`$ ). In this case the $`Q_S`$ matrix is off–diagonal in the Nambu space and time–local. Since the $`\mathrm{\Lambda }_N`$–matrix is diagonal in the Nambu space, only even powers of $`(Q_S\mathrm{\Lambda }_N)`$ may contribute to the action. As a result, one finds for the proximity action

$`S_{TD}[\mathrm{\Phi }]`$ $`=`$ $`{\displaystyle \frac{G_D}{8}}{\displaystyle \underset{l=1}{\overset{\mathrm{}}{}}}\gamma _{2l}Tr\left\{\left(Q_S\mathrm{\Lambda }_N\right)^{2l}\right\}`$ (45)
$`=`$ $`{\displaystyle \frac{G_D}{4}}{\displaystyle \underset{l=1}{\overset{\mathrm{}}{}}}(1)^l\gamma _{2l}{\displaystyle \underset{a=1}{\overset{2p}{}}}\underset{2l}{\underset{}{{\displaystyle \underset{0}{\overset{\beta }{}}}𝑑\tau _1\mathrm{}{\displaystyle \underset{0}{\overset{\beta }{}}}𝑑\tau _{2l}}}e^{2i\mathrm{\Phi }_a(\tau _1)}G(\tau _1\tau _2)e^{2i\mathrm{\Phi }_a(\tau _2)}G(\tau _2\tau _3)\mathrm{}e^{2i\mathrm{\Phi }_a(\tau _{2l})}G(\tau _{2l}\tau _1),`$ (46)
where $`G(\tau \tau ^{})=iT/\mathrm{sin}\pi T(\tau \tau ^{})`$ is the Green function of the normal lead (the Fourier transform of $`\mathrm{\Lambda }_N(ϵ_n)`$ to the time domain is $`\mathrm{\Lambda }_N(\tau ,\tau ^{})=G(\tau \tau ^{})\sigma _3`$). If the external voltage $`V(\tau )`$ is applied to the SC dot, the phase is to be understood as $`\mathrm{\Phi }=𝑑\tau V(\tau )`$. Calculating the second variation of the above action with respect to $`V(\tau )`$, one finds for the linear (Andreev) conductance

$$G_A=G_D\underset{l=1}{\overset{\mathrm{}}{}}(1)^l(2l)^2\gamma _{2l}=G_D_\mathrm{\Theta }u_0(\pi /2).$$
(47)
Employing Eq. (42) for $`\mathrm{\Theta }\pi /2`$ one finds

$$G_A=\frac{G_T\mathrm{sin}\theta }{1+t\mathrm{sin}\theta },$$
(48)
where $`\theta =t\mathrm{cos}\theta =u_0(\pi /2)`$. One therefore obtains $`G_A=G_T^2/G_D+O(t)`$ for $`t1`$ and $`G_A=G_D+O(1/t)`$ for $`t1`$.

### C Instanton treatment of the proximity action

To calculate the proximity action on the simplest instanton trajectory $`\mathrm{\Phi }_1(\tau )=\pi T\tau W`$ and $`\mathrm{\Phi }_a(\tau )=0`$ for $`a[2,2p]`$, one may employ e.g. the Matsubara basis (see Appendix A for calculations in the time domain). Employing Eq. (22) (modified in the obvious way for arbitrary $`W`$), one finds

$$(Q_S\mathrm{\Lambda }_N)^{2l}(n,m)=(1)^l\delta _{n,m}\left(\begin{array}{cc}\text{sign}(ϵ_{nW})\text{sign}(ϵ_n)& 0\\ 0& \text{sign}(ϵ_{n+W})\text{sign}(ϵ_n)\end{array}\right)^l;$$
(49)
in the replica $`a=1`$, while $`(Q_S\mathrm{\Lambda }_N)^{2l}=(1)^l\sigma _0`$ in all other replica. As a result, for all even $`l=2k`$ the action in $`a=1`$ replica is not different from that in the other $`2p1`$ replica and therefore does not contribute to the total action in the replica limit $`p0`$. On the other hand, for odd $`l=2k1`$, there are $`|W|`$ Matsubara components, where $`a=1`$ and $`a1`$ replica come with the opposite signs. Employing Eq. (45), one finds for action on the instanton trajectory (in the replica limit)

$$𝒮_W(t)=|W|\frac{G_D}{2}\underset{k=1}{\overset{\mathrm{}}{}}\gamma _{4k2}=|W|\frac{G_D}{8}\left[\underset{0}{\overset{\pi /2}{}}𝑑\mathrm{\Theta }u_0(\mathrm{\Theta })\underset{\pi /2}{\overset{\pi }{}}𝑑\mathrm{\Theta }u_0(\mathrm{\Theta })\right],$$
(50)
where the last equality is a direct consequence of the definition $`u_0(\mathrm{\Theta })=_{l=1}^{\mathrm{}}l\gamma _l\mathrm{sin}(l\mathrm{\Theta })`$. Employing the differential equation (42), one finds

$$t^2\frac{}{t}\left(\frac{𝒮_2}{t}\right)=\frac{G_D}{8}\left[\underset{0}{\overset{\pi /2}{}}𝑑\mathrm{\Theta }_\mathrm{\Theta }u_0^2\underset{\pi /2}{\overset{\pi }{}}𝑑\mathrm{\Theta }_\mathrm{\Theta }u_0^2\right]=\frac{G_D}{8}\left[2\theta ^2\mathrm{\Theta }_2^2\right].$$
(51)
Comparing with Eq. (34), we find an exact coincidence of the real–space result, Eq. (33), and the proximity action, Eq. (50), $`𝒮(t)=𝒮_1(t)`$. Notice that, since the (random) diffusive region was integrated out upon derivation of the proximity action, the phase of the CB oscillations is not random. In a sense the proximity action represents a typical diffusive region (mesoscopic fluctuations are omitted in the derivation given above). One may therefore calculate directly typical $`Z(q)`$, without resorting to the correlation function. From the two smallest winding numbers $`W=1`$ and $`W=2`$ one finds $`Z(q)\mathrm{exp}\{(𝒮_1+𝒜(T))\}\mathrm{cos}\pi q+\mathrm{exp}\{𝒮_2\}\mathrm{cos}2\pi q+\mathrm{}`$, in agreement with the result of the previous section.

It is worth mentioning here that the general expression for the instanton action in the case of the normal island, derived in Ref. , may be presented in the form similar to Eqs. (33) and (50):

$$𝒮_W^N=|W|\frac{G_D}{4}\underset{k=1}{\overset{\mathrm{}}{}}\gamma _{2k1}=|W|\frac{G_D}{8}\underset{0}{\overset{\pi }{}}𝑑\mathrm{\Theta }u_0(\mathrm{\Theta })=|W|\frac{G_D}{16}\left[\mathrm{\Theta }_2^2+2t(\mathrm{cos}\mathrm{\Theta }_2+1)\right],$$
(52)
and only even $`W`$ are allowed for the normal case. It is easy to see, using Eq. (32), that in the entire interval $`G_D>G_T`$ the value of the action $`𝒮_2^N`$ is the same as in the tunnelling limit: $`𝒮_2^N=G_T/2`$. As $`t=G_T/G_D`$ approaches unity, a kind of phase transition occurs , which reveals itself as a non-analytic behavior of $`𝒮_2^N(t)`$ at $`t=1`$.

### D Landauer approach

Employing Nazarov’s techniques , one may show that the function $`u_0(\mathrm{\Theta })`$ is directly related to the generating function of the transmission coefficients

$$\frac{1}{2}G_Du_0(\mathrm{\Theta })=F(\mathrm{\Theta })\mathrm{sin}\mathrm{\Theta }\underset{\alpha =1}{\overset{N}{}}\frac{T_\alpha \mathrm{sin}\mathrm{\Theta }}{1T_\alpha \mathrm{sin}^2(\mathrm{\Theta }/2)},$$
(53)
where $`T_\alpha `$ are eigenvalues of the $`N\times N`$ matrix $`t^{}t`$ describing transmission between the dot and the lead. Factor $`1/2`$ on the l.h.s. takes into account that $`G_D`$ is defined for the two spin components. Performing the Fourier transform one finds

$$\frac{1}{2}G_D\gamma _l=\frac{4(1)^{l+1}}{l}\underset{\alpha =1}{\overset{N}{}}e^{2\mu _\alpha l},$$
(54)
where $`T_\alpha =\mathrm{cosh}^2\mu _\alpha `$. As a result, the proximity action takes the compact form:

$$S_{TD}=\underset{\alpha =1}{\overset{N}{}}Tr\mathrm{ln}\left(1+e^{2\mu _\alpha }Q_S\mathrm{\Lambda }_N\right)=\frac{1}{2}\underset{\alpha =1}{\overset{N}{}}Tr\mathrm{ln}\left(1\frac{T_\alpha }{4}(Q_S\mathrm{\Lambda }_N)^2\right).$$
(55)
In the last expression we have omitted a constant which goes to zero in the replica limit. Exactly the same interface action is known in the supersymmetric formalism  (notice that here we deal with spin $`1/2`$ electrons). One may calculate now the action on the instanton trajectory to find the renormalized charging energy. Employing Eq. (54) and recalling that $`𝒮_1=(G_D/2)_k\gamma _{4k2}`$, one finds

$$𝒮_1=\frac{1}{2}\underset{\alpha =1}{\overset{N}{}}\mathrm{ln}\left(1\frac{T_\alpha ^2}{(2T_\alpha )^2}\right).$$
(56)
This expression is to be compared with the one for the Andreev conductance obtained from Eqs. (47) and (54) :

$$G_A=4\underset{\alpha =1}{\overset{N}{}}\frac{T_\alpha ^2}{(2T_\alpha )^2};$$
(57)
(c.f. with the normal conductance $`G_N=G_D_ll^2\gamma _l=2_\alpha T_\alpha `$, where $`2`$ stays for spin). Factor $`4`$ on the r.h.s. of Eq. (57) stems from the fact that Andreev particles carry charge $`2e`$ and no spin. One may thus identify the combination $`𝒯_\alpha =T_\alpha ^2/(2T_\alpha )^2`$ as the Andreev conductance in channel $`\alpha `$. The renormalisation of the charging energy is therefore given by

$$e^{𝒮_1}=\underset{\alpha =1}{\overset{N}{}}(1𝒯_\alpha )^{1/2}.$$
(58)
This may be compared with the corresponding factor for the normal spinless particles : $`_\alpha (1T_\alpha )^{1/2}`$. Indeed, in our notations the normal instanton action is $`𝒮_2^N=(G_D/2)_l\gamma _{2l1}=_\alpha \mathrm{ln}(1T_\alpha )`$; for spinless particles one has $`\mathrm{exp}\{𝒮_2^N/2\}`$ for the CB suppression. There is thus a perfect analogy between Andreev and normal spinless charge fluctuations mechanisms, provided that the normal transmission coefficients, $`T_\alpha `$, are substituted by the Andreev ones, $`𝒯_\alpha `$.

From Eq. (58) one may extract some general results. In the tunnelling limit: $`𝒯_\alpha 1`$ in all channels, one finds: $`𝒮_1=G_A/8`$. In the diffusive metal interface the transmission eigenvalues are distributed according to the Dorokhov distribution : $`P(T)=G_D/(4T\sqrt{1T})`$. The corresponding distribution for the Andreev transmissions is given by $`P(𝒯)=G_D/(8𝒯\sqrt{1𝒯})`$, as a result $`G_A=G_D`$. For the typical action of the diffusive interface one finds therefore $`𝒮_1=\pi ^2G_A/32`$.

### E Fluctuation determinant and summation over many-instanton configurations

We now expand the action, Eq. (45), to the second order in deviations from the instanton trajectory $`\mathrm{\Phi }_a(\tau )=\pi T\tau W_a+\delta \mathrm{\Phi }^a(\tau )`$, where $`\delta \mathrm{\Phi }^a(0)=\delta \mathrm{\Phi }^a(\beta )`$. After diagonalization of the resulting quadratic form one finds the following spectrum of the small fluctuations :

$$\lambda _n^{(W_a)}=\frac{\pi ^2T}{E_C}n^2+\frac{G_A}{4}\left(|nW_a|+|n+W_a|2|W_a|\right),$$
(59)
where $`n[\mathrm{},\mathrm{}]`$. For $`TE_C`$ there are $`2|W_a|+1`$ (almost) zero modes labelled by $`n|W_a|`$. One of them $`n=0`$ is the trivial shift of $`\mathrm{\Phi }_a(\tau )`$ by a constant, whereas the remaining $`2|W_a|`$ zero modes are associated with the deformation of the instantons. The general solution of the saddle point equations (for $`TE_C`$) known as the Korshunov instanton  may be written as

$$e^{2i\mathrm{\Phi }_a(\tau )}=\underset{k=1}{\overset{|W_a|}{}}\left[\frac{e^{2\pi iT\tau }z_k}{1e^{2\pi iT\tau }z_k^{}}\right]^{sign(W_a)},$$
(60)
where $`z_k`$ is a set of $`|W_a|`$ complex numbers (instanton coordinates) which parametrize the $`2|W_a|`$ dimensional zero–mode manifold. The instanton coordinates, $`z_k`$, are complex numbers from inside of the unit circle. One may thus characterize any deviation $`\delta \mathrm{\Phi }_a(\tau )`$ by $`2|W_a|`$ zero–mode coordinates, $`z_k`$, and remaining transversal modes having non–zero masses, Eq. (59). The Jacobian of such transformation may be found in the standard way  and is given by

$$J=\frac{1}{|W_a|!}\text{det}\frac{1}{1z_lz_k^{}},$$
(61)
where $`l,k[1,|W_a|]`$. The fact that the Jacobian vanishes if two of the coordinates coincide signals the repulsive interaction between the instantons.

The Gaussian integration over the massive degrees of freedom normalized by those in the zero winding number sector results in

$$\frac{\underset{n=1}{\overset{\mathrm{}}{}}(\lambda _n^{(0)}/\pi )}{\underset{n=|W_a|+1}{\overset{\mathrm{}}{}}(\lambda _n^{(W_a)}/\pi )}=\left(\frac{G_A}{2\pi }\right)^{|W_a|}\mathrm{exp}\left\{|W_a|\mathrm{ln}\left(\frac{G_AE_C}{2\pi ^2T}\right)\right\},$$
(62)
up to the terms of the order $`T/(G_AE_C)1`$. Combining together all the factors, one finds for the one instanton, $`W_a=1`$, trajectory

$$\frac{Z_1}{Z_0}=e^{i\pi q}\underset{|z|1}{}\frac{d^2z}{1|z|^2}\left(\frac{G_A}{2\pi }\right)e^{\{𝒮_1\mathrm{ln}(G_AE_C/2\pi ^2T)\}}=e^{i\pi q}\frac{\stackrel{~}{E}_C}{2T}\mathrm{ln}\frac{E_C}{T},$$
(63)
The cutoff of the logarithmically divergent $`z`$–integral originates from the following consideration. In presence of the charging energy term the Korshunov instantons, Eq. (60), are not true zero–modes of the action. The $`z`$–dependence of the action brings the factor $`\mathrm{exp}\{(T/2E_C)|z|^2/(1|z|^2)\}`$ to the integral. As a result the effective integration range shrinks to $`|z|^2<1T/E_C<1`$. The renormalized charging energy $`\stackrel{~}{E}_C`$ is defined as

$$\stackrel{~}{E}_C=\frac{E_CG_A^2}{2\pi ^2}e^{𝒮_1}.$$
(64)
It is also instructive to look at the $`W_a=2`$ contribution to the partition function

$`{\displaystyle \frac{Z_2}{Z_0}}`$ $`=`$ $`e^{2i\pi q}{\displaystyle \frac{1}{2!}}{\displaystyle d^2z_1d^2z_2\left[\frac{1}{(1|z|_1^2)(1|z|_2^2)}\frac{1}{(1z_1z_2^{})(1z_1^{}z_2)}\right]\left(\frac{G_A}{2\pi }\right)^2e^{\{𝒮_22\mathrm{ln}(G_AE_C/2\pi ^2T)\}}}`$ (65)
$`=`$ $`{\displaystyle \frac{1}{2!}}e^{2i\pi q}\left({\displaystyle \frac{\stackrel{~}{E}_C}{2T}}\right)^2\left[\mathrm{ln}^2{\displaystyle \frac{E_C}{T}}{\displaystyle \frac{\pi ^2}{6}}\right].`$ (66)
The leading power of the large logarithm comes from the diagonal term of the Jacobian, Eq. (61). This term corresponds to the non–interacting instantons approximation, which is justified therefore by the large parameter $`\mathrm{ln}(E_C/T)1`$. Notice, that the (repulsive) interaction correction (the off–diagonal part of the Jacobian, Eq. (61)) comes with $`\mathrm{ln}^0`$ and not $`\mathrm{ln}^1`$, as may be naively expected.

Due to the large logarithmic factor, originating from the integration over sizes of each individual instanton (i.e. neglecting their interactions), one may treat positive and negative instantons (anti–instantons) on the equal footing. Strictly speaking, a combination of instantons and anti–instantons is not a saddle point solution. However, relative weakness of instanton interactions (no log-factor in the second term in Eq. (66)) makes it possible to consider an ideal gas of instantons and anti–instantons. The contribution to the partition function with a given winding number $`W`$ is given by all configurations having $`m+|W|`$ instantons and $`m`$ anti–instantons in an arbitrary order. Taking only the terms with the leading power of the logarithms (diagonal term of the Jacobian) one finds

$$Z(q)=Z_0\underset{W=\mathrm{}}{\overset{\mathrm{}}{}}e^{\pi iWq}\underset{m=0}{\overset{\mathrm{}}{}}\frac{1}{m!(m+|W|)!}\left(\frac{\stackrel{~}{E}_C}{2T}\mathrm{ln}\frac{E_C}{T}\right)^{2m+|W|},$$
(67)
This ugly looking series may be summed up into an unexpectedly simple expression:

$$Z(q)=Z_0\mathrm{exp}\left\{\frac{\stackrel{~}{E}_C}{T}\mathrm{ln}\left(\frac{E_C}{T}\right)\mathrm{cos}(\pi q)\right\}.$$
(68)
The simplest way to check it is to expand back Eq. (68) in a double series in powers of $`\stackrel{~}{E}_C`$ and $`e^{i\pi q}`$. The role of the anti–instantons, therefore, is to convert $`\mathrm{cos}(W\pi q)`$, which may be expected from the instantons only, into $`\mathrm{cos}^W(\pi q)`$.

Under the condition $`\stackrel{~}{E}_CTE_C`$ the interaction between instantons can be treated perturbatively; for the lowest-order correction we use the last term in Eq. (66) to find for the gate voltage - dependent free energy

$$F(q)=\stackrel{~}{E}_C\left[\mathrm{ln}\left(\frac{E_C}{T}\right)\mathrm{cos}(\pi q)\frac{\pi ^2}{24}\frac{\stackrel{~}{E}_C}{T}\mathrm{cos}(2\pi q)+\mathrm{}\right].$$
(69)
We retain the second term in (69) since it will be important at temperatures comparable with the parity-effect temperature $`T^{}`$, see Sec. VI. It also demonstrates that the instanton–instanton interaction corrections are small only if $`\stackrel{~}{E}_C<T`$.

## V RG treatment of the quantum fluctuations

We next approach the problem of the CB from the different perspective. Instead of calculating the action on the instanton field configurations, we shall perform the RG analysis of the coefficients $`\gamma _l`$ upon integrating out fast fluctuations of $`\mathrm{\Phi }(\tau )`$. To this end we write $`\mathrm{\Phi }(\tau )=\mathrm{\Phi }^s(\tau )+\mathrm{\Phi }^f(\tau )`$, where superscripts $`s`$ and $`f`$ stay for slow ($`\omega _m<\mathrm{\Omega }`$) and fast ($`\omega _m>\mathrm{\Omega }`$) component of the field correspondingly. The running cutoff $`\mathrm{\Omega }`$ runs from $`\mathrm{\Omega }G_AE_C`$ down to the lowest Matsubara frequency $`\mathrm{\Omega }T`$. We next substitute the $`\mathrm{\Phi }`$–field into the SC proximity action, Eq. (45), expand to the second order in $`\mathrm{\Phi }^f`$ and integrate out the fast field fluctuations. The bare propagator of the $`\mathrm{\Phi }^f`$–fields may be read out from the action Eq. (45):

$$\mathrm{\Phi }_a^f(\omega _m)\mathrm{\Phi }_b^f(\omega _m)=\frac{\delta _{ab}}{|\omega _m|G_A(\zeta _\mathrm{\Omega })},$$
(70)
where $`G_A(\zeta _\mathrm{\Omega })=G_D_{l=1}^{\mathrm{}}(1)^l(2l)^2\gamma _{2l}(\zeta _\mathrm{\Omega })`$ is the running value of the Andreev conductance, and $`\zeta _\mathrm{\Omega }\mathrm{ln}(G_AE_C)/\mathrm{\Omega }`$. (Note that this definition of $`\zeta `$ differs from that used to derive the proximity action in Section IV.) The further calculation is essentially similar to that in the Keldysh formalism . One obtains the following set of the evolution equations for $`\gamma _l(\zeta _\mathrm{\Omega })`$:

$$\frac{d\gamma _{2l}}{d\zeta }=\frac{8}{G_A(\zeta )}\left(l\gamma _{2l}+2\underset{k=1}{\overset{\mathrm{}}{}}(1)^k(l+k)\gamma _{2l+2k}\right),$$
(71)
where we have omitted the subscript $`\mathrm{\Omega }`$. Since only the even coefficients are involved it is convenient to define the function $`\stackrel{~}{u}(\zeta ,\mathrm{\Theta })=[u(\zeta ,\mathrm{\Theta })u(\zeta ,\pi \mathrm{\Theta })]/2`$ which contains only even harmonics. In terms of this function the RG equations (71) take the form

$$\stackrel{~}{u}_\zeta =\frac{4}{G_A(\zeta )}\left[\stackrel{~}{u}(\mathrm{\Theta })\mathrm{tan}\mathrm{\Theta }\right]_\mathrm{\Theta },$$
(72)
where the Andreev conductance is $`G_A(\zeta )=G_D\stackrel{~}{u}_\mathrm{\Theta }(\zeta ,\pi /2)`$. The initial condition for Eq. (72) is $`\stackrel{~}{u}(0,\mathrm{\Theta })=\stackrel{~}{u}_0(\mathrm{\Theta })`$, with the $`u_0(\mathrm{\Theta })`$ function evaluated above. Let us mention for completeness, that in the normal case the corresponding RG equation takes the form $`u_\zeta =(2/G_N(\zeta ))[u(\mathrm{\Theta })\mathrm{cot}(\mathrm{\Theta }/2)]_\mathrm{\Theta }`$, where the running normal conductance according to Eq. (43) is given by $`G_N(\zeta )=G_D_{l=1}^{\mathrm{}}l^2\gamma _l(\zeta )=G_Du_\mathrm{\Theta }(\zeta ,0)`$.

The key parameter that determines the strength of phase fluctuations is the Andreev conductance $`G_A(\zeta )`$. We consider the large bare value, $`G_A(0)1`$. Upon integration over fluctuations of $`\mathrm{\Phi }(\tau )`$, the effective value of $`G_A(\zeta )`$ decreases (the physical reason of this phenomena is a partial loss of the phase coherence between multiple Andreev reflections), and becomes comparable to unity at some time scale $`\mathrm{\Omega }_c^1=(G_AE_C)^1e^{\zeta _c}`$. At longer times coupling between the dot and the lead is weak, thus the phase fluctuations grow and autocorrelation function $`C(\tau )=e^{i(\mathrm{\Phi }(\tau )\mathrm{\Phi }(0))}`$ decays fast. It is natural to associate $`\mathrm{\Omega }_c`$ with an effective Coulomb energy of the island connected to the wire, $`\stackrel{~}{E_C}`$. We will show now that the estimate for $`\mathrm{\Omega }_c`$ that follows from the RG equation coincides (within exponential accuracy) with the results of the instanton analysis.

Equation (72) may be solved with the method of characteristics. We define the function $`v(\mathrm{\Theta })=\stackrel{~}{u}(\mathrm{\Theta })\mathrm{tan}\mathrm{\Theta }`$ and then introduce new “space” and “time” variables: $`\xi =\mathrm{ln}\mathrm{sin}\mathrm{\Theta }`$ and $`\tau `$ which is defined via the relation

$$\frac{4}{G_A(\zeta )}d\zeta =d\tau .$$
(73)
In terms of these new variables the RG equation (72) takes the simple form: $`v_\tau =v_\xi `$, yielding $`v(\xi ,\tau )=v_0(\xi +\tau )`$. As a result,

$$\stackrel{~}{u}(\zeta (\tau ),\mathrm{\Theta })=\frac{e^\tau \mathrm{cos}\mathrm{\Theta }}{\sqrt{1e^{2\tau }\mathrm{sin}^2\mathrm{\Theta }}}\stackrel{~}{u}_0(\mathrm{arcsin}(e^\tau \mathrm{sin}\mathrm{\Theta })),$$
(74)
The dependence $`\zeta (\tau )`$ can be now found from Eq. (73) with the Andreev conductance given by $`G_A(\zeta )=G_D\stackrel{~}{u}_\mathrm{\Theta }(\zeta ,\pi /2)`$.

Our goal is to find a scale $`\zeta _c`$, where the Andreev conductance $`G_A`$ becomes small. To find $`\zeta _c`$ we notice that in terms of the RG time $`\tau `$ the smallness of $`G_A`$ means the limit $`\tau \mathrm{}`$. Then, using definition of $`\tau `$ given by Eq. (73), we find

$$\zeta _c=\frac{G_D}{4}\underset{0}{\overset{\mathrm{}}{}}𝑑\tau \frac{\stackrel{~}{u}_0(\mathrm{arcsin}e^\tau )}{\sqrt{e^{2\tau }1}}=\frac{G_D}{4}\underset{0}{\overset{\pi /2}{}}𝑑\mathrm{\Theta }\stackrel{~}{u}_0(\mathrm{\Theta }).$$
(75)
Recalling the definition of $`\stackrel{~}{u}_0(\mathrm{\Theta })`$ and comparing this result with Eq. (50) we conclude that $`\zeta _c=𝒮_1`$. As a result, the amplitude of the lowest oscillatory component ($`\mathrm{cos}\pi q`$) of the free energy, from the RG point of view, is given by the precisely the same exponential factor as in the framework of the instanton calculation.

Yet another instructive way to approach the problem is to follow the renormalization of the instanton action. According to Eq. (50), $`𝒮_W(\zeta )=|W|(G_D/2)_k\gamma _{4k2}(\zeta )=|W|(G_D/4)_0^{\pi /2}𝑑\mathrm{\Theta }\stackrel{~}{u}_0(\mathrm{\Theta })`$. Employing Eq. (72), one finds

$$\frac{d𝒮_W}{d\zeta }=\frac{|W|G_D}{G_A(\zeta )}\left(\stackrel{~}{u}(\mathrm{\Theta })\mathrm{tan}\mathrm{\Theta }\right)|_{\mathrm{\Theta }\frac{\pi }{2}}=|W|,$$
(76)
where in the last equality we have used that $`G_A(\zeta )=G_D\stackrel{~}{u}_\mathrm{\Theta }(\zeta ,\pi /2)`$. (In the normal case one finds $`dS_W^N/d\zeta =|W|/2`$; only even $`W`$ are allowed.) From the dimensional analysis one expects $`F(q)\mathrm{\Omega }\mathrm{exp}\{𝒮_1\}\mathrm{cos}\pi q+\mathrm{}`$, where $`\mathrm{\Omega }E_CG_A`$ is the cutoff energy. We may now renormalize down the cutoff energy $`\mathrm{\Omega }\mathrm{\Omega }(\zeta )\mathrm{\Omega }e^\zeta `$, simultaneously following renormalization of the action, $`𝒮_1𝒮_1(\zeta )`$. According to Eq. (76), $`𝒮_1(\zeta )=𝒮_1\zeta `$. As a result,

$$\mathrm{\Omega }(\zeta )e^{𝒮_1(\zeta )}=\text{const}.$$
(77)
This means that the lowest oscillatory component of the free energy is expected to be temperature independent, in agreement with Eq. (69) (up to log–factor, originating from the zero modes). The instanton calculation is done at $`\zeta =0`$, alternatively the RG calculation presented above looks for $`\zeta _c`$, such that $`𝒮_W(\zeta _c)0`$ (the scale where the instantons do not cost anything). In view of Eq. (77) it is not surprising that they yield the same result. Moreover, one can renormalize down to any intermediate scale $`\zeta `$ (such that $`G_A(\zeta )>1`$) and then calculate the instanton action with the current parameters $`\{\gamma _l(\zeta )\}`$ – the result is still guaranteed to be the same. Notice also that the renormalized action $`S_W(\zeta )=S_W|W|\mathrm{ln}(G_AE_C/T)`$ coincides with the result that comes from the instanton action together with the massive Gaussian fluctuations, c.f. Eqs. (62)–(66). The effect of the zero modes ($`\mathrm{ln}E_C/T`$ in the prefactor) is missed in the one–loop RG calculation.

Finally we comment upon the comparison of the approach used in this Section with Ref.  where the same kind of RG approach was employed in the case of 2-dimensional normal conductor with very low Thouless energy, $`E_{\mathrm{Th}}E_C`$. In the last case the RG procedures involves simultaneous integration over phase fluctuations $`\mathrm{\Phi }(\tau )`$ and diffuson/Cooperon modes in the normal conductor. As a result, the full functional RG equation contains three terms in the r.h.s.: one coming from Eq. (40), the second is similar to Eq. (72), and the third contribution accounts for the Cooper-channel repulsion in the normal conductor (cf. ). The resultant RG equation cannot be reduced to the set of even harmonics $`\stackrel{~}{u}`$ only; we are not aware of any method to solve it apart from “brut-force” numerical integration. Considerable simplification of the problem treated here stems from the fact that, due to the condition $`E_{\mathrm{Th}}E_C`$ all diffuson/Cooperon modes may be integrated out before (means at larger energies than) the phase fluctuations $`\mathrm{\Phi }(\tau )`$ become relevant. This procedure leads to the action functional $`S_{TD}`$ that depends upon $`\mathrm{\Phi }(\tau )`$ trajectory only, as defined in (45). Thus the results presented in this Section can be considered as “minimal generalization” of those obtained in Ref. .

## VI Discussion of the results

The main result of our study, given by Eq. (69), applies in the temperature range

$$\stackrel{~}{E}_CT\mathrm{min}\{E_C,T^{}\},$$
(78)
where $`T^{}`$ is the parity-effect temperature given by Eq. (10). The condition $`TT^{}`$ ensures that both even and odd winding numbers contribute to the partition function, since $`𝒜(T)1`$, c.f. Eq. (9). The two lowest oscillatory components of the free energy are given then by Eq. (69). At higher temperature $`TT^{}`$ an addition of odd number of electrons to the dot becomes possible, changing the relative amplitude of the harmonics. On the level of the partition function the components with odd winding numbers acquire the temperature dependent factor $`\mathrm{exp}\{𝒜(T)\}`$. As a result, one finds for the gate voltage dependent free energy:

$$F(q)=\stackrel{~}{E}_C\left[e^{𝒜(T)}\mathrm{ln}\left(\frac{E_C}{T}\right)\mathrm{cos}(\pi q)+\left((1e^{2𝒜(T)})\mathrm{ln}^2\left(\frac{E_C}{T}\right)\frac{\pi ^2}{6}\right)\frac{\stackrel{~}{E}_C}{4T}\mathrm{cos}(2\pi q)+\mathrm{}\right].$$
(79)
Finally at $`T>T^{}`$ the parity effect disappears and only the normal (unit) oscillation period remains in the free energy

$$F(q)=\stackrel{~}{E}_C\frac{\stackrel{~}{E}_C\left(\mathrm{ln}^2\left(\frac{E_C}{T}\right)\frac{\pi ^2}{6}\right)}{4T}\mathrm{cos}(2\pi q)+\mathrm{}$$
(80)
Notice that if $`\mathrm{ln}E_C/T^{}>\pi /\sqrt{6}`$ there is the sign change of the $`\mathrm{cos}(2\pi q)`$ component at $`TT^{}`$. Moreover, unlike the low temperature case, where the amplitude of the dominant harmonics was weakly (logarithmically) temperature dependent, Eq. (69), at larger temperature there is the stronger dependence ($`T^1\mathrm{ln}^2(E_C/T)`$) of the amplitude. The characteristic crossover temperature is determined from the relation $`𝒜(T^{})\mathrm{ln}(T^{}/\stackrel{~}{E}_C)`$ and is given by

$$T^{}=\frac{\mathrm{\Delta }}{\mathrm{ln}\frac{\mathrm{\Delta }}{G_A\delta }},$$
(81)
where we have used the fact that $`𝒮_1G_A`$.

We emphasize an interesting feature of the weak Coulomb blockade compared to the usual one: although the role of effective Coulomb energy is taken by the effective energy $`\stackrel{~}{E}_CE_C`$, the oscillation amplitude depends relatively weakly upon temperature within the range given by Eq. (78). This is to be contrasted with the usual case of an island connected by highly resistive contacts: at temperatures $`T>E_C`$ oscillation amplitude vanishes exponentially fast with $`T/4E_C`$.

Next we shall discuss the quantitative value of the renormalized charging energy $`\stackrel{~}{E}_C`$, which determines both the magnitude of the CB and the region of applicability of our results. To this end we need to evaluate $`𝒮_1(t)`$ (where $`t=G_T/G_D`$), which is given by one of the two equivalent expressions, Eqs. (33) or (50). In the two limiting cases $`t1`$ (the tunneling barrier limit) and $`t1`$ (the diffusive metal limit) the answer may be found analytically. One obtains

$$𝒮_1=\mathrm{ln}\frac{2\pi ^2\stackrel{~}{E}_C}{G_A^2E_C}=\frac{G_A}{8}\{\begin{array}{cc}1;\hfill & \text{in the tunneling limit,}\hfill \\ \pi ^2/4;\hfill & \text{in the diffusive limit,}\hfill \end{array}$$
(82)
where the Andreev conductance is given by Eq. (48). Recall that for a normal dot the corresponding exponent is given by $`G_N/2`$ in the tunnelling limit and $`G_N\pi ^2/8`$ in the diffusive limit. We observe that upon the same dissipative conductance the SC dot exhibits factor of four (in the two limiting cases) smaller action. Therefore one may observe a sizeable CB in the superconductive state, while in the analogous normal dot the CB is practically suppressed. The crossover behavior for the SC action $`𝒮_1(t)`$, the normal action $`𝒮_2^N(t)`$ and the Andreev conductance $`G_A(t)`$ as function of $`t=G_T/G_D`$ is depicted in Fig. 2. Dependence of the same quantities on the resistance $`G_D^1`$ of the normal region at fixed $`G_T`$ is demonstrated in Fig. 3.

In the case $`\stackrel{~}{E}_CT^{}`$ the crossover between $`T<T^{}`$ and $`T>T^{}`$ regimes leads to the sharp drop of the oscillation amplitude by the factor $`\frac{\stackrel{~}{E}_C}{4T^{}}\mathrm{ln}\frac{E_C}{T^{}}1`$ The same drop in the residual Coulomb blockade may occur as function of magnetic field, due to the suppression of the $`\mathrm{\Delta }`$ and thus of $`T^{}`$ by the magnetic field. This may result in a strong negative magnetoresistance at low temperatures of a granular media made out of small superconductive grains (cf. e.g.  for the discussion of relevant experiments). This effect can occur if the low-temperature oscillation amplitude $`\stackrel{~}{E}_C`$ is sufficiently large to destroy Josephson coupling between grains, so transport of Cooper pairs between grains is blocked by $`\stackrel{~}{E}_C`$. However, quantitative theory of such an effect is still to be developed.

Another open question concerns the behavior of $`F(q)`$ in the zero temperature limit, $`T<\stackrel{~}{E}_C`$. As temperature decreases below $`\stackrel{~}{E}_C`$, all approximate methods, we used, run out of their applicability range: the renormalized conductance $`G_A`$ drops below unity, fluctuations become strong and there is no apriori reason to treat them within Gaussian approximation. Whereas the overall scale of the oscillation amplitude is probably given by its first harmonic, Eq. (69) at $`T\stackrel{~}{E}_C`$, and is of the order $`\stackrel{~}{E}_C\mathrm{ln}\frac{E_C}{\stackrel{~}{E}_C}`$, the precise shape of the oscillations is still to be determined. It is possible that in the $`T0`$ limit the function $`F(q)`$ becomes nonanalytic at the degeneracy point $`q=1/2`$. An extreme case of such a nonanalytic behavior – finite steps in $`dF/dq`$ at half-integer $`q`$ – was found in Ref.  where the Andreev conductance was completely neglected.

###### Acknowledgements.

We have greatly benefitted from the numerous discussions with A. V. Andreev, I. S. Beloborodov, L. I. Glazman, and K. A. Matveev. MVF and MAS were supported by the SCOPES program of Switzerland, Dutch Organization for Fundamental Research (NWO), Russian Foundation for Basic Research under grant 01-02-17759, the program “Quantum Macrophysics” of the Russian Academy of Sciences, the Russian Ministry of Science, and the Swiss National Foundation. AK was partially supported by the BSF grant N 9800338. AIL was partially supported by the NSF grant N 0120702. MAS was partially supported by the Russian science support foundation.

## A Instantons in the Imaginary-Time Proximity Action

### 1 Transition to imaginary time

The method of the multicharge proximity action was initially developed  in the Keldysh real-time representation. Here we show how the Matsubara proximity action, Eq. (16), may be obtained from its Keldysh analog by analytic continuation to imaginary time. Since this Appendix serves illustrative purposes we will consider the simplest case $`T=0`$. As discussed in the bulk of the paper, the form of $`F(q)`$ is unknown at $`T=0`$ since at $`T<\stackrel{~}{E}_C`$ the dilute instanton gas approximation fails and one has to consider an interacting instanton liquid. Nevertheless, the overall scale of oscillations in $`F(q)`$ may be inferred from the single instanton action which is temperature independent and may be calculated at $`T=0`$.

To proceed we need to establish a correspondence between our notations and those adopted in Refs. . The latter notations will be designated by a prime. Firstly, our conductance quantum $`G_Q=e^2/2\pi \mathrm{}`$ is different from $`G_Q^{}=e^2/\mathrm{}`$ used in those papers. Secondly, the running parameter $`\zeta `$ and charges $`\{\gamma _n\}`$ are related by

$$\frac{\zeta }{\zeta ^{}}=\frac{\gamma _n^{}}{\gamma _n}=\frac{G_D}{8\pi ^2g^{}}.$$
(A1)
Thus, the Keldysh multicharge proximity action may be written in our notations as

$$S_{\mathrm{prox}}[\mathrm{\Phi }(t)]=i\frac{G_D}{8}\underset{n=1}{\overset{\mathrm{}}{}}\gamma _n(\zeta )Tr(\stackrel{ˇ}{Q}_S\stackrel{ˇ}{\mathrm{\Lambda }}_N)^n,$$
(A2)
where the trace is taken over time indexes, Nambu and Keldysh spaces. Here $`\stackrel{ˇ}{Q}_S`$ is a matrix in the superconducting dot acting in the Nambu space as

$$\stackrel{ˇ}{Q}_S(t)=\left(\begin{array}{cc}0& e^{2i\stackrel{}{\mathrm{\Phi }}(t)}\\ e^{2i\stackrel{}{\mathrm{\Phi }}(t)}& 0\end{array}\right),$$
(A3)
where $`\stackrel{}{\mathrm{\Phi }}=diag(\mathrm{\Phi }_>,\mathrm{\Phi }_<)`$ is a matrix in the Keldysh space, with $`\mathrm{\Phi }_>(t)`$ and $`\mathrm{\Phi }_<(t)`$ being the fields residing on the forward and backward branches of the Keldysh contour. In the normal lead, $`\stackrel{ˇ}{\mathrm{\Lambda }}_N(E)=\stackrel{ˇ}{\mathrm{\Lambda }}_0(E)\sigma _3`$, where $`\sigma _3`$ is the Pauli matrix in the Nambu space, and the matrix $`\stackrel{ˇ}{\mathrm{\Lambda }}_N(E)`$ acts in the Keldysh space:

$$\stackrel{ˇ}{\mathrm{\Lambda }}_0(E)=\left(\begin{array}{cc}12f& 2f\\ 2(f1)& 2f1\end{array}\right),$$
(A4)
$`f(E)`$ being the distribution function.

Tracing over the Nambu space reduces the action (A2) to the form:

$$S_{\mathrm{prox}}[\mathrm{\Phi }(t)]=i\frac{G_D}{4}\underset{n=1}{\overset{\mathrm{}}{}}(1)^n\gamma _{2n}\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑t_1\mathrm{}𝑑t_{2n}\mathrm{tr}_Ke^{2i\stackrel{}{\mathrm{\Phi }}(t_1)}\stackrel{ˇ}{\mathrm{\Lambda }}_0(t_1t_2)e^{2i\stackrel{}{\mathrm{\Phi }}(t_2)}\stackrel{ˇ}{\mathrm{\Lambda }}_0(t_2t_3)\mathrm{}e^{2i\stackrel{}{\mathrm{\Phi }}(t_{2n})}\stackrel{ˇ}{\mathrm{\Lambda }}_0(t_{2n}t_1),$$
(A5)
where the trace is taken only over the Keldysh space.

Transition to the imaginary time is achieved by the deformation of the Keldysh contour $`𝒞`$ which initially run over the time axis from $`\mathrm{}`$ to $`\mathrm{}`$ and then in the backward direction. The desired deformation introduces a vertical segment of length $`\beta `$ at some time $`t=t_0`$ so that the contour originates at $`\mathrm{}+i\beta /2`$, runs through the points $`t_0+i\beta /2`$, $`t_0`$ to $`\mathrm{}`$ and then back through $`t_0`$, $`t_0i\beta /2`$ to $`\mathrm{}i\beta /2`$, see Fig. 4. For the purpose of the evaluation of the instanton action the choice of the point $`t_0`$ is somewhat arbitrary. Note, however, that were we to consider, e.g., the correlation function of the form $`\stackrel{~}{N}(q)\stackrel{~}{N}(q^{})`$, it would be necessary to allow for adiabatic evolution of the gate voltage $`q=q(t)`$, $`q^{}=q(t^{})`$ and to introduce two vertical segments at times $`t`$ and $`t^{}`$.

The analysis is simplified in the zero-temperature case. Choosing the position of the vertical segment at $`t_0=0`$, we deform the forward segment $`(\mathrm{},0)`$ of the Keldysh contour to the upper imaginary half-axis and the backward segment $`(0,\mathrm{})`$ to the lower imaginary half-axis. The remaining horizontal appendix of the Keldysh contour running from 0 to $`\mathrm{}`$ and then back to $`t=0`$ can be neglected as it does not contribute to thermodynamic quantities.

In order to continue the function $`\stackrel{ˇ}{\mathrm{\Lambda }}_0(t)`$ to the complex plane in time one should employ its analytic properties . In deforming the contour to the imaginary axis the element $`\stackrel{ˇ}{\mathrm{\Lambda }}_0^{12}(\tau )`$ enters only with $`\tau >0`$, and the element $`\stackrel{ˇ}{\mathrm{\Lambda }}_0^{21}(\tau )`$ enters only with $`\tau <0`$. Under this condition the function $`\stackrel{ˇ}{\mathrm{\Lambda }}_0(\tau )`$ can be substituted by

$$\stackrel{ˇ}{\mathrm{\Lambda }}_0(\tau )F(\tau )\left(\begin{array}{cc}1& 1\\ 1& 1\end{array}\right),F(\tau )=\frac{1}{\pi \tau }.$$
(A6)
Eq. (A6) solves the problem of analytic continuation of the Green function $`\stackrel{ˇ}{\mathrm{\Lambda }}_0(t)`$.

The trivial matrix structure of Eq. (A6) ensures that each term in Eq. (A5) can be written as a multiple integral of the single function

$$\mathrm{\Phi }(\tau )=\{\begin{array}{cc}\mathrm{\Phi }_>(\tau ),\hfill & \text{for }\tau >0\text{,}\hfill \\ \mathrm{\Phi }_<(\tau ),\hfill & \text{for }\tau <0\text{,}\hfill \end{array}$$
(A7)
defined on the whole imaginary axis. The infinitesimal element $`dt`$ is transformed to

$$dt\{\begin{array}{cc}id\tau ,\hfill & \text{for }\tau >0\text{,}\hfill \\ id\tau ,\hfill & \text{for }\tau <0\text{.}\hfill \end{array}$$
(A8)
Finally, we obtain

$$S_{\mathrm{prox}}[\mathrm{\Phi }(\tau )]=\frac{G_D}{4}\underset{n=1}{\overset{\mathrm{}}{}}\gamma _{2n}_{\mathrm{}}^{\mathrm{}}𝑑\tau _1\mathrm{}𝑑\tau _{2n}e^{2i\mathrm{\Phi }(\tau _1)2i\mathrm{\Phi }(\tau _2)+\mathrm{}2i\mathrm{\Phi }(\tau _{2n})}F(\tau _1\tau _2)F(\tau _2\tau _3)\mathrm{}F(\tau _{2n}\tau _1),$$
(A9)
that coincides with Eq. (45) in the zero-temperature limit.

### 2 Instantons in time domain

Here we calculate the action on the trajectories for which $`e^{2i\mathrm{\Phi }(\tau )}`$ is an analytic function in the upper half-plane. For such a solution, the instanton’s winding number $`\pi W=\mathrm{\Delta }\mathrm{\Phi }=𝑑\tau (\mathrm{\Phi }/\tau )`$ is positive. Integration over $`t_1`$, $`t_3`$, …, $`t_{2n1}`$ is performed with the help of

$$\frac{d\tau _1}{\pi }e^{2i\mathrm{\Phi }(\tau _1)}\frac{1}{\tau _0\tau _1}\frac{1}{\tau _1\tau _2}=i\frac{e^{2i\mathrm{\Phi }(\tau _0)}e^{2i\mathrm{\Phi }(\tau _2)}}{\tau _0\tau _2}\pi e^{2i\mathrm{\Phi }(\tau _0)}\delta (\tau _0\tau _2),$$
(A10)
that is obtained by regularizing $`1/t=[1/(ti0)+1/(t+i0)]/2`$ and making use of analyticity of $`e^{2i\mathrm{\Phi }(\tau )}`$ in the upper half-plane. As a result, the integral in Eq. (A9) is transformed to

$$_{\mathrm{}}^{\mathrm{}}\underset{k=1}{\overset{n}{}}\left(\frac{d\tau _{2k}}{\pi }e^{2i\mathrm{\Phi }(\tau _{2k})}\right)\underset{k=1}{\overset{n}{}}\left[i\frac{e^{2i\mathrm{\Phi }(\tau _{2(k1)})}e^{2i\mathrm{\Phi }(\tau _{2(k+1)})}}{\tau _{2(k1)}\tau _{2(k+1)}}\pi e^{2i\mathrm{\Phi }(\tau _{2k})}\delta (\tau _{2(k1)}\tau _{2(k+1)})\right],$$
(A11)
where $`\tau _0\tau _{2n}`$ and $`\tau _{2(n+1)}\tau _2`$. Expanding the last product and combining similar terms we rewrite it as

$$(1)^n\delta (0)+\underset{m=1}{\overset{n}{}}C_n^m(1)^{nm}K_m,$$
(A12)
where

$$K_m=(i)^m_{\mathrm{}}^{\mathrm{}}\underset{k=1}{\overset{m}{}}\left(\frac{d\tau _{2k}}{\pi }e^{2i\mathrm{\Phi }(\tau _{2k})}\right)\underset{k=1}{\overset{m}{}}\frac{e^{2i\mathrm{\Phi }(\tau _{2(k1)})}e^{2i\mathrm{\Phi }(\tau _{2(k+1)})}}{\tau _{2(k1)}\tau _{2(k+1)}}.$$
(A13)
The integrals for $`K_m`$ are calculated recursively with the help of the relation

$$\frac{d\tau _2}{\pi }e^{2i\mathrm{\Phi }(\tau _1)}\frac{e^{2i\mathrm{\Phi }(\tau _1)}e^{2i\mathrm{\Phi }(\tau _2)}}{\tau _1\tau _2}\frac{e^{2i\mathrm{\Phi }(\tau _2)}e^{2i\mathrm{\Phi }(\tau _3)}}{\tau _2\tau _3}=2i\frac{e^{2i\mathrm{\Phi }(\tau _1)}e^{2i\mathrm{\Phi }(\tau _3)}}{\tau _1\tau _3}.$$
(A14)
Thereby we get $`K_m=2^mW`$. For anti–instantons with $`W<0`$ the same analysis yields $`K_m=2^m|W|`$.

Now summation in Eq. (A12) becomes trivial and we obtain for the instanton action:

$$S[\mathrm{\Phi }(\tau )]=\frac{G_D}{4}\underset{n=1}{\overset{\mathrm{}}{}}\gamma _{2n}\left\{(1)^n\delta (0)+\right|W|[1(1)^n]\}.$$
(A15)
Eq. (A15) formally contains the divergent $`\delta `$-function of zero argument. This part of the answer gives the action of the non-instanton configuration $`\mathrm{\Phi }=\mathrm{\Phi }_0=\mathrm{const}`$ and thus drops from the difference $`𝒮_W=S[\mathrm{\Phi }(\tau )]S[\mathrm{\Phi }_0]`$ which coincides with Eq. (50) obtained in the frequency domain.