{'fraction_non_alphanumeric': 0.0433650048040995, 'fraction_numerical': 0.01146578413579588, 'mean_word_length': 4.549289099526066, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Dynamics and phonon-induced decoherence of Andreev level qubit

We present detailed theory for Andreev level qubit, the system consisting of a highly transmissive quantum point contact embedded in a superconducting loop. The two-level Hamiltonian for Andreev levels interacting with quantum phase fluctuations is derived by using the path integral method. We also derive kinetic equation describing qubit decoherence due to interaction of the Andreev levels with acoustic phonons. The collision terms are non-linear due to fermionic nature of the Andreev states, leading to slow non-exponential relaxation and dephasing of the qubit at temperature smaller than the qubit level spacing.


I. INTRODUCTION
The possibility to employ Andreev bound levels in superconducting contacts for quantum computation has been suggested in Refs. 1, 2. The proposed Andreev level qubit (ALQ) consists of a highly transmissive, with reflectivity R ≪ 1, quantum point contact (QPC) embedded in a low inductance superconducting loop. In the ALQ, quantum information is stored in the microscopic two-level system of Andreev bound states in the contact. Hybridization of the clockwise and counterclockwise persistent current states in the ALQ loop is produced by the microscopic processes of electronic back scattering in the QPC. This is different from the macroscopic superconducting flux qubits [3][4][5] and charge-phase qubit 6 , where the hybridization is provided by charge fluctuations on the tunnel junction capacitors. Thus requirement of the large charging energy or large loop inductance is not critical for the ALQ. A single ALQ consists of a pair of the Andreev bound levels belonging to the same normal conducting mode in the QPC; a multimode QPC will form a qubit cluster. The way of the ALQ operation is similar to the one of the experimentally tested flux qubits 3-5 -the Andreev levels can be excited by driving biasing magnetic flux through the qubit loop. 7,8,1 The read out method is also similar to the flux and charge-phase qubits: the quantum state of the Andreev levels determines the magnitude and direction of the persistent current circulating in the loop, and also the magnitude of the induced flux. Since the quantum information is stored in the microscopic system, Andreev levels, while the access for manipulation and readout is provided by macroscopic persistent currents, ALQ occupies an intermediate place between the microscopic solid state qubits (like localized spins on impurities 9 or quantum dots 10 ) and macroscopic superconducting qubits.
During 90-s, the Josephson transport in superconducting QPCs has been intensively investigated, and a number of remarkable experiments has been performed on atomic size metallic QPCs using controllable break junction technique 11,12 , as well as on gated quantum constrictions in 2D electron gas confined between superconductors. 13 The critical Josephson current and current-voltage characteristics have been thoroughly examined in these experiments by applying current or voltage bias. [13][14][15][16] There was, however, one experiment, particularly important in the qubit context, where the flux bias has been implemented: Koops et al. 17 have inserted metallic QPC in a SQUID, and evaluated the Josephson current-phase dependence by measuring the induced flux with an inductively coupled SQUID magnetometer. The measurements have been only reported for the equilibrium state. Unfortunately, no experimental attempts to drive the QPC out of equilibrium to some coherent or incoherent excited state have been performed so far.
The purpose of this paper is to present a theory for the Andreev level qubit in more detail than that was possible in short communication. 2 We will also consider the electron-phonon interaction as a potential "intrinsic" source of the qubit decoherence, and derive a kinetic equation for the qubit density matrix. The ALQ Hamiltonian and the kinetic equation are derived by using a path integral method. [18][19][20] The central technical difficulty here is to extend the method to contacts with large transparency. This difficulty is overcome by incorporating exact boundary condition in the QPC action. Another important point discussed in the paper is the role of the charge electro-neutrality in the junction electrodes, which affects the qubit Hamiltonian.
The fermionic nature of the Andreev levels does not affect the qubit operation and qubit-qubit coupling but it plays important role for the qubit decoherence. We find that the electron-phonon collision terms in the kinetic equation for the ALQ differ qualitatively from the Bloch-Redfield master equation 21 commonly applied to study decoherence of the macroscopic superconducting qubits. 22 This results in long phonon-induced decoherence time for the ALQ, namely, at the temperature smaller than the qubit level spacing, both the relaxation and dephasing processes are governed by a power rather than exponential law.
The structure of the paper is the following. We discuss the model description of the QPC in Section II, and then explain in detail, in Section III, the path integral approach for transmissive QPC: we consider the action for the contact and derive the effective Hamiltonian for Andreev levels interacting with the quantum phase fluctuations; the single-particle density matrix, and effective current operator are also discussed in this section. In Section IV we discuss averaging over fast phase fluctuations and derive an effective Hamiltonian for the qubit, this procedure is extended in Section V to the two inductively coupled qubits to derive direct qubit-qubit interaction. Section VI is devoted to the electron-phonon interaction: we derive effective action for the Andreev level-phonon interaction and derive corresponding collision terms in the equation for the qubit density matrix; then we present solution of the kinetic equation, and evaluate the decoherence rate.

II. CONTACT HAMILTONIAN
Let us consider a superconducting quantum point contact with bulk 3D electrodes. We model the contact with a smooth on the atomic scale (adiabatic) constriction and assume a local scatterer situated in the neck of constriction (see Fig. 1) producing weak electronic back scattering with the reflectivity R ≪ 1. We further assume that the constriction supports a single conducting mode. We adopt the mean field approximation for the electrons in the contact, which are described with the Hamiltonian, where the first term is the BCS Hamiltonian for the bulk superconducting electrons, Ψ(r, t) being the two-component Nambu field operator, and the second term describes the charging energy of the contact capacitor C 18-20 . The single-particle Hamiltonian h in Eq.
(1) has the form, where ∆(r, t) and χ(r, t) are, respectively, the modulus and phase of the superconducting order parameter, the potential U (r) accounts for the confinement of electrons within the contact as well as the electron scattering, while ϕ(r, t) and A(r, t) are electromagnetic potentials. The voltage drop at the contact, V (t), in Eq. (1) is related to the discontinuity of the electric potential at the contact, V (t) = ϕ(−0, t)−ϕ(+0, t). To investigate the decoherence effects, we allow the electrons in the electrodes to interact with acoustic phonons, and include corresponding electron-phonon interaction and phonon terms in the total Hamiltonian of the contact, It is convenient to introduce quasiclassical (Andreev) approximation for the superconducting electrons. Following standard procedure, we eliminate rapidly varying in space potential, U (r), by introducing quasiclassical wave functions of the single conducting mode in the left and right electrodes, and couple these wave functions by means of a normal electron scattering matrix. In Eq. (4), ψ (σ) (x) are slowly varying 1D envelopes for the longitudinal electron motion (σ = ± indicates the direction of the motion), ψ ⊥ (r ⊥ , x) is a normalized wave function of the transverse motion with energy is the quasiclassical longitudinal electronic momentum. The coupling of the quasiclassical envelopes in the left (L) and right (R) electrodes, ψ L,R , is conveniently described by the transfer matrix, 7,23 Here d and r are the energy independent transmission and reflection amplitudes, respectively. Since any observable quantity is expressed through a bilinear combination of the envelopes with the same σ, the energy independent scattering phases can be eliminated from the boundary condition (5), hence, without loss of generality, the scattering amplitudes will be further assumed to be real, d = √ D, r = √ R, where D and R are the transmission and reflection coefficients of the contact, respectively.
The electromagnetic potentials, ϕ(r, t), A(r, t), and the complex order parameter in Eq.
(2) are to be found from the Maxwell equations and the self-consistency equation. It is convenient to present the Hamiltonian in a gauge invariant form by extracting the phase of the order parameter using a local gauge transformation, Ψ(r, t) → e iσz χ(r,t)/2 Ψ(r, t) .
Then the superfluid momentum p s =h∇χ/2 − (e/c)A, and the gauge-invariant electric potentialφ =hχ/2 + eϕ, appear in the quasiclassical Hamiltonian of the electrode, while the phase difference φ(t) = χ R (0, t) − χ L (0, t), appears in the boundary condition, In the bulk metallic electrodes with good screening, and at small frequencies relevant for the problem, the gauge-invariant fields,φ(r, t), p s (r, t), are to be found from the electro-neutrality condition, and the current conservation, [24][25][26] δn(r, t) = 0, ∇j(r, t) = 0, where n(r, t) is the electronic density. In the electrodes, the charge imbalance relaxation yields the equilibrium relation, ϕ = (∂n/∂µ) −1 δn, on the distance exceeding the electric field penetration length. Furthermore, in the absence of normal dissipative current, p s is proportional to the total current density, j = (en/m)p s , which is negligibly small far from the contact due to rapid spreading out of the current in the point contact geometry. Thus the conditions (9) yield complete cancellation of the electromagnetic potentials in the electrodes, Taking into account that the modulus of the order parameter far from the contact is equal to the equilibrium value, 27 ∆ = const, we conclude that the Hamiltonian (7) approaches the equilibrium form. Potentialφ can be expressed through the electric field E,φ =hχ/2 + e dx E, by introducing the gauge invariant phaseχ = χ − (2e/ch) dx A = (2/h) dx p s . The spatial distribution of these quantities is illustrated on Fig.  2. Sinceφ vanishes far from the contact, andχ is constant according to Eq. (10), then the gauge-invariant phase difference across the contactφ =χ (Josephson relation). In the SQUID, this relation is equivalent to the phase versus flux relation,φ = 2eΦ/hc, since the voltage drop across the contact is generated by time variation of the magnetic flux Φ threading the SQUID. Finally, we notice that the gauge-invariant phase difference,φ, rather than φ enters the boundary condition (8), which can be explicitly seen by extracting the Aharonov-Bohm phase from the transfer matrix. Thus the gauge-invariant phase difference remains the only free collective variable whose dynamics is determined by the electrodynamic environment of the contact. Below we will not distinguish betweenφ and φ, because the difference is negligibly small in the QPC. Proceeding with discussion of the interaction of electrons with phonons, we consider only longitudinal acoustic phonons and describe the interaction within the deformation potential approximation, where γ is the deformation potential constant, u(r, t) is the phonon field operator, s is the sound velocity, and ρ is the crystal mass density. The Hamiltonian of free phonons has a standard form, Our strategy now will be to derive effective Hamiltonian for the Andreev levels including interaction with phonons. If the phase difference would be a classical variable, this derivation can be in principle done by direct truncation of the Hamiltonian (3). However, in the presence of quantum phase fluctuations it is convenient to apply the path integral technique.

III. CONTACT EFFECTIVE ACTION
Let us consider the whole system, the QPC and superconducting loop (see Fig. 3), and introduce the path integral representation for the propagator, The Lagrangian of the system L tot consists of the contact part L c , and the part describing the circulating current in the loop. The latter is conveniently combined with the charge term from the electronic Hamiltonian (1) giving the Lagrangian of the loop oscillator L osc , where φ e corresponds to the bias magnetic flux, and L is the loop inductance. The remaining part of the contact Lagrangian consists, similar to Eq. (3), of the electronic part, phonon part, and electron-phonon interaction, In the quasiclassical approximation, the electronic Lagrangian splits into two parts, L α , α = L, R, corresponding to the left and right electrodes, and the third part L BC discussed in detailed in the next section, which accounts for the boundary condition. Noting that relaxation processes are caused by phonons with small wave vectors compared to the Fermi wave vector, q ≪ k, transitions between the states ψ (+) and ψ (−) are forbidden, and the electron-phonon Lagrangian can be written on the form, The phonon Lagrangian is given by equation,

A. Boundary condition
The boundary condition (5) is valid for any contact transparency. To include this boundary condition in the path integral formulation, we introduce an additional term in the Lagrangian 2 , where η is an auxiliary fermionic Nambu field playing the role of a Lagrange multiplier. Correspondingly, the integration over η is to be included in the propagator in Eq. (15), giving the following form for the electronic part of the propagator, Let us prove that such a Lagrangian indeed generates the boundary condition (5), (8). To this end, it is convenient to make a step back and restore a non-quasiclassical form for the fermionic field, in the bulk part of the Lagrangian, The dynamic equations and the boundary condition result from the zero variation of the action S e = dt (L L + L R + L BC ) with respect toψ L,R andη, or in the explicit form, and Integrating Eqs. (26) over x in a small vicinity of x = 0 yields the relations, Then introducing again the quasiclassical envelopes and combining Eqs. (27), (28) with the quasiclassical relation, we get the boundary condition equivalent to Eqs. (5), (8).

B. Effective action for Andreev levels
Now we are prepared to derive an effective action for the Andreev levels. Following the procedure of Ref. 18, we integrate out fast electronic fields ψ α in Eq. (22), Then we are left with the effective action S 0 η , which contains only variables, η(t) and φ(t), where g(t) is given by the Fourier component, Connection between the effective action (31) and the Andreev levels can be established by considering the case of time independent phase, φ = const. Indeed, by writing the effective action in the Fourier representation, we identify the spectrum of the system,hω = ±E a (φ), by calculating eigenvalues of the matrix in the brackets, this equation coincides with well known Andreev level spectrum 28,29 (see Fig. 4). Thus we conclude that the fermionic field η represents the Andreev levels.
Proceeding to a time-dependent φ(t), we restrict its time variation speed to small values,h∂ t φ/4 ≪ ∆. Furthermore, the dynamics of the Andreev levels is also to be slow, E(φ) ≪ ∆, which implies that the contact reflectivity must be small, R ≪ 1, ALQ must be biased at φ e ≈ π, and the amplitude of the quantum phase fluctuationsφ(t) = φ(t) − φ e , must be sufficiently small,φ ≪ φ e . We emphasize that the constraint on the amplitude of the phase fluctuations is actually provided in our case due to the loop geometry of the electrodes having sufficiently small inductance; the constraint is important to suppress the Landau-Zener transitions to the continuum states. Under the imposed conditions, the non-local in time kernel g(ω), Eq. (32), can be replaced by a constant value (adiabatic approximation), , leading to the equation, We eliminate the term with the phase time derivative in Eq. (35) by transforming η, Then we finally arrive at the effective action, where presents the effective single-particle Hamiltonian for the two-level Andreev system. 2 The Hamiltonian in Eq. (38) differs by the exponential pre-factor from the two-level Hamiltonian derived in Refs. 30, 31 and further employed in Refs. 32,33. This factor appears in the present derivation after the electric potential has been included into consideration in Eq.
(1) to provide the electro-neutrality in the electrodes [see text after Eq. (8)]. Both Hamiltonians are equivalent under stationary conditions, ∂ t φ = 0, and the difference is not important for the adiabatic dynamics. However, in general, the pre-factor is important, e.g. for derivation of correct equation for the current operator in Eq. (50). It is instructive to compare the case of the transparent contact considered here with the case of a tunnel contact extensively studied in the MQC theory. [18][19][20] Physical difference between the two cases is that the Andreev level system in transparent contacts is a slow one, while in tunnel contacts it is a fast one because the Andreev level energy in tunnel contacts is close to ∆. Within the present formalism, the integration over ψ fields is similar in both cases. However, the next step, the adiabatic approximation in Eqs. (31), (32) is not allowed in the tunnel limit; instead one should perform also the integration over η in Eq. (22), and make an expansion over small D. The result of this calculation, presented in Appendix A, coincides with the results of Refs. 18, 19, 34.

C. Andreev levels density matrix
The Hamiltonian (38) governs the free evolution of the Andreev variable η, The same Hamiltonian governs the free evolution of the single-particle density matrix of the Andreev levels ρ a (t) defined via statistical average, whereη denotes the fermionic operator corresponding to the Grassmann field η, and the average is taken over all electronic states. The density matrix ρ a is a 2 × 2 matrix in the Nambu space, which satisfies the normalization condition Tr ρ a = 1. The statistical average in Eq. (40) is represented, after the averaging over ψ L,R , by a path integral, Then calculating the time derivative, using Eq. (39) and the conjugated equation, we get, Thus the Andreev level dynamics is described by a Liouville equation similar to ordinary quantum mechanical two-level systems.

D. Andreev level current
We conclude this section with derivation of an effective current operator for the Andreev levels. A common quantum mechanical expression for the current applied to the fields ψ L,R at x = ∓0 gives, or using Eq. (28), [the current is continuous, I L (t) = I R (t), by virtue of Eq. (27)]. The same equation can be obtained by varying the electronic part of the action, S e , with respect to the phase difference, Using this relation the statistically averaged current, can be written as or, after tracing out the fields, ψ L,R , This leads to the following equation for the current in terms of the Andreev variable, in which the kernel I corresponds to the effective single particle current operator of Andreev levels, Apparently the current operator does not commute with the Andreev level Hamiltonian, [h a , I] = 0, which is the consequence of a normal electron reflection at the QPC. Hence the Andreev level eigen states consist of superpositions of the current eigen states, unless R = 0 (see Fig. 4). Correspondingly, the Andreev level current, defined as an expectation value of the current operator over the Andreev state,

IV. AVERAGING OVER PHASE FLUCTUATIONS
Equation (42) describes dynamics of the Andreev levels for a given realization of the time dependent phase across the QPC. However, the phase dynamics is strongly coupled to the Andreev levels. Intrinsic dynamics of the phase is governed by the quantum Hamiltonian of the loop oscillator [cf. Eq. (16)], where ω p = 8E L E C /h 2 is a plasma frequency of the oscillator determined by the contact charging energy, E C = e 2 /2C, and the loop inductive energy, E L = (hc/2e) 2 (1/L). Thus the whole system is generally a multilevel one. On the other hand, the phase dynamics and the coupling between the loop oscillator and the Andreev levels is a vital part of the ALQ operation: The Andreev levels cannot be manipulated if the phase dynamics was frozen because any manipulation requires variation of the current, and second, the read out of the Andreev levels can be only performed via measuring the quantum state of the loop oscillator. An obvious way to solve this problem and preserve the qubit property of the coupled system is to "enslave" the loop oscillator by choosing the oscillator level spacinghω p much larger than the Andreev level spacing,hω p ≫ 2E a . Then the Andreev state evolution will not excite the oscillator, which will remain in the ground state and adiabatically follow the evolution of the Andreev levels. This implies that the phase is a fast variable which should be averaged out, leading to an effective qubit Hamiltonian. To facilitate the averaging procedure we take advantage of the small amplitude of the phase fluctuations,φ = φ − φ e ≪ φ e , which was already assumed when proceeding from Eq. (33) to Eq. (35). This assumption is justified by the large inductive energy, E L ≫ E J ≈ ∆, and it allows us to expand the Andreev level Hamiltonian (38) over small φ; then proceeding to the current eigen basis, η → e −iσx √ Rφe/4 η, we get, Averaging over fast phase fluctuations can be done directly in Eq. (15) by performing explicit integration over phase. However, there is a simpler way to get the same result. Equation (42) holds for a fluctuating phase provided the oscillator is not excited. It can be viewed as the diagonal with respect to φ part of a more general equation for a full density matrix, ρ(σφ, σ ′ φ ′ ), for the Andreev levels plus oscillator system, The interaction between the Andreev levels and the oscillator given by the second term in Eq. (54) displaces the oscillator steady state from the origin by ±I/2M ω 2 p depending on the direction of the current in the junction or, equivalently, the state of the Andreev levels (see Fig. 5). We eliminate this term in (54) by applying the transformation, and then average the resulting Hamiltonian over the oscillator ground state, taking into account the relation, e iA p 0 = exp −A 2 p 2 0 /2 [ ... 0 indicates averaging over the oscillator ground state]. As the result we get an effective Hamiltonian, where the bare contact reflectivity, R, is renormalized by the phase fluctuations, This renormalization effect can be understood as the effect of inertia of the loop oscillator, which hinders the current variations, i.e. it works against the effect of the electronic back scattering at the contact responsible for the hybridization of the current states [see discussion in the end of the Section IIID]. The renormalization effect becomes increasingly strong in the limit of a classical oscillator with large "mass". Because of renormalization of the contact reflectivity, the Andreev level spectrum is modified, and the qubit frequency reduces. This might be important for practical applications, because it would allow one to reduce the qubit frequency by choosing the circuit parameters rather than by tuning the contact reflectivity. Eq. (57) gives effective Hamiltonian h q for the ALQ in the absence of interaction with phonons. Similarly, averaged over phase fluctuation density matrix ρ 0 gives density matrix for the ALQ. Keeping the same notation, ρ, for the qubit density matrix, we finally arrive at the equation for the free evolution of the ALQ, This equation is sufficient to describe the ALQ manipulation (by driving the biasing flux) and read out (by measuring the induced current) 1,2 , and also the qubit-qubit interaction, which is discussed in the next section.

V. INDUCTIVE QUBIT-QUBIT COUPLING
Our treatment of the interaction of the Andreev levels with phase fluctuations can be easily extended on the case of several inductively coupled SQUIDs to describe direct qubit-qubit coupling. Let us consider as an example two SQUIDs with different QPC reflectivities, R 1 = R 2 , and identical circuit parameters, C and L, mutual inductance is M. Then the inductance terms in the Lagrangian (16) written for the two qubits will take the form, whereφ 1,2 are fluctuating phases in the first and second qubits. By introducing the normal modes for the LCoscillators, φ → e −iτyπ/4 φ, we present the two-qubit Hamiltonian on the form similar to Eqs. (54), (56), where H osc,i describes the normal oscillator with frequency ω pi = c/ C (L ± M), and the index i = 1, 2 refers to the first and second qubit. Then we apply the similar transformations as in the previous section to Eq. (62), namely, we eliminate the linear inφ i terms by means of a canonical transformation, , and then average the transformed Hamiltonian over the ground state of the normal oscillators. As the result, we obtain an effective two-qubit Hamiltonian including direct qubit-qubit coupling, The renormalized contact reflectivities are now given by equation, ]. It is worth mentioning that the two-qubit configuration may be also realized with a single SQUID containing QPC with two conducting modes. In this case, we have the two Andreev level Hamiltonians coupled to the same loop oscillator, Averaging over the phase fluctuations, we arrive at the same interaction Hamiltonian as in Eq. (63) but with the loop inductance L substituting for −M. Thus we come to an interesting conclusion that the effect of the phase fluctuations not only reduces the bare contact reflectivity but also introduces effective coupling of the Andreev levels of different conducting modes in multimode QPC.

A. Effective action
In this Section we include electron-phonon interaction in the consideration. We start with the derivation of an effective action for the Andreev level-phonon interaction. To this end we repeat the calculation of the previous section adding the Lagrangian L e−ph , Eq. (19), to Eq. (30). By keeping in the electrode Green functions only the first-order correction in small interaction constant γ, we arrive at the following action, where the Green functions G L,R read, In this equation, the quantities g (σ) L,R refer to the different parts, g (σ) (x < 0, x ′ < 0; t) and g (σ) (x > 0, x ′ > 0; t), respectively, of the translation-invariant free electron Green function, g (σ) (x − x ′ , t), which obeys the equation The Green functions g L,R (x, 0; t) = −Θ(∓x) where ζ t is given by the Fourier component, ζ ω = ∆ 2 − (hω) 2 . Proceeding to the adiabatic approximation discussed in the previous section (ζ ω → ζ e ), and performing the transformation (36), we arrive at the following effective action for the Andreev level-phonon interaction, where and Taking into account the zero-order dynamic equation with respect to γ, Eq. (39), and putting φ = φ e , we obtain for the quantity n(x, t), Eq. (70), the following expression, where κ = ζ e /hv. Finally, integrating over r in the action (69), we present the effective interaction on the form, where γ q is the effective constant of the Andreev level-phonon coupling, We notice that the Andreev level-phonon interaction in Eq. (74) has purely transverse origin, i.e. while inducing interlevel transitions and hence the relaxation it does not produce any additional dephasing to the one associated with the relaxation. It is important to mention that the effective coupling constant is proportional to √ R, and it turns to zero in the case of perfectly transparent constriction. This results from the already mentioned fact that the relevant phonons have small wave vectors and are not able to provide large momentum transfer (∼ 2hk) during scattering with electrons.
The effective action in Eqs. (74) and (75) was derived for a given realization of the time dependent phase. To take into account the effect of phase fluctuation, we have to apply the transformation, Eq. (56), to the action (74); the Lagrange form of the transformation reads, Then the integration over the phase adds the factor e −λ to the action, cf. Eq. (58), which simply implies a renormalization of the coupling constant, This result can be expected: since the Andreev level-phonon coupling is transverse in the current basis and depends on the contact reflectivity, the renormalized reflectivity R * rather than the bare one has to enter the coupling constant, γ * q = γ q (R * ); this is consistent with Eqs. (77), (75).

B. Kinetic equation
Qubit decoherence is usually described with collision terms in the Liouville equation for the qubit density matrix, which take into account interaction with an environment. Our goal now will be to derive the phonon-induced collision terms in the equation (60) for ALQ density matrix, and to evaluate the decoherence of the ALQ.
While the description of free qubit evolution was possible in the terms of the single-particle density matrix, evaluation of the collision terms goes beyond the single-particle approximation and requires the knowledge of electronic twoparticle correlation functions. This is because the Andreev levels are not a rigorously isolated system but rather belong to a large fermionic system of the superconducting electrons in the contact electrodes. Thus to derive the collision terms, we apply the many-body Keldysh-Green function technique 36 combined with the path integral approach. The method described below automatically takes into account the many-body effects, namely the Pauli exclusion principle, leading to a non-linear form of the collision terms and eventually to suppression of the decoherence.
The starting point of the derivation is Eqs. (15), (22) for the propagator, in which the integration over the fast fermionic fields, ψ L,R , and phase, φ, has been performed while integration over the phonons and Andreev states remains, The time evolution in the action follows now along the Keldysh contour 37 C K , S = CK dt L, which goes from −∞ to +∞, and then backwards. 20,38,39 The interaction is supposed to be switched on and off adiabatically at the remote past, t = −∞, and the phonon bath is supposed to be in the thermal equilibrium. The Lagrangian L has the form, where L ph is given by Eq. (20). To reduce the time integration along the Keldysh contour to an ordinary time integral, we distinguish forward and backward branches of the contour by labeling them with index s = 1, 2, and introduce the two-component fields η s and X s q . Since the action is local in time, it can be rewritten as S = . The first step of the derivation is to integrate out the phonon fields, which will give rise to an effective self-interaction for the field η, with the kernelĎ(t − t ′ ) given by whereĎ q (t − t ′ ) is the equilibrium Keldysh Green function of the phonons represented with the 2 × 2 matrix in the Keldysh space. In these equations, T C is the time-ordering operator on the Keldysh contour, and τ z is the Pauli matrix operating in the Keldysh space; summation over repeated indices is implied. The next step is to take advantage of the weak electron-phonon interaction, and to decouple the four-fermionic interaction term in Eq. (80) by introducing the Hubbard-Stratonovich field G ss ′ αβ (t, t ′ ), which is a matrix in the Keldysh-Nambu-time space. Before doing this, it is convenient to explicitly extract the small parameter λ ph , which determines the electron-phonon coupling strength, from the kernelĎ in Eq. (80) by redefining the kernel,Ď → λ phĎ (λ ph ∼ (v/s)(ω/ω D ) 2 ≪ 1, ω D is the Debye frequency). As the result, we get, where L(t, t ′ ) = (ih∂ t − E * a σ z ) δ(t − t ′ ) is diagonal in the Keldysh space, and Equation (82) describes dynamics of the field η interacting with the effective fieldǦ. In terms of the Keldysh-Green function for the field η, this evolution is described by the Dyson equation, where the self-energyΣ depends on the effective fieldǦ. A closed equation forǦ can be derived by integrating out the field η in Eq. (82). This procedure leads to the equation, where Tr denotes both the matrix trace in the Keldysh-Nambu space and the integration over the time variables.
Noticing that the action is large, S[Ǧ ] ∼ λ −1 ph , we evaluate the integral in Eq. (87) within the saddle-point approximation (cf. Ref. 38). The corresponding saddle-point equation is derived by varying the action with respect toǦ, which yields, Comparing Eqs. (88) and (86), we obtain the relation,Ǧ = λ phǦ . Written in the terms ofǦ, where G R(A) = G 11 − G 12 (21) is the retarded (advanced) Green function, and G K = G 11 + G 22 = G 12 + G 21 is the Keldysh component. Similar relations also hold for the self-energy. Then Eq. (89) takes the form, Kinetic equation is obtained by considering the difference between Eq. (91) and its Hermitian conjugate for the Keldysh component, 36 The right hand side of Eq. (92) describes the qubit decoherence as well as dynamic corrections due to the phonons. 40 In the absence of the coupling to phonons, the solution of Eq. (92) has the form, where F is a time-independent matrix determined by the initial state of the qubit. When a weak interaction with the phonons is switched on, an asymptotical solution to the equation (92) can be written on the form, where the matrix F is a slowly evolving on the time scale,h/E * a , function of the global time (t + t ′ )/2, andF is a rapidly oscillating small-amplitude correction (see Appendix B). Equations for the matrix elements of F are derived in Appendix B, Eqs. (B14), (B15), and for the diagonal matrix elements they read, while for the off-diagonal matrix element, F 12 = F * 21 , the equation has the form, where ν is the phonon-induced transition rate between the qubit levels, N = e 2βE * a − 1 −1 is the phonon distribution function at the frequency equal to the qubit level spacing, β = 1/kT , and δ and δ 0 are small dynamic corrections defined in Eq. (B16). For equal times, t = t ′ , G K (t, t) is related to the qubit density matrix, Eq. (40), as follows: , and therefore Eqs. (95), (96) in fact give kinetic equation for the qubit density matrix in the interaction picture,ρ = e iσz E * a t/h ρ e −iσzE * a t/h , F σ = 2ρ σσ − 1, F 12 = 2ρ 12 . It is instructive to write Eq. (95) in terms of the qubit occupation numbers n σ = 1 − ρ σσ , The right hand side of this equation has the standard form of the electron-phonon collision term, yielding the Fermi distribution for the equilibrium occupation numbers, This conclusion is consistent with the well known fact that the Fermi distribution of the Andreev levels gives correct magnitude for the equilibrium Josephson current 28,29 . Furthermore, it follows from Eqs. (98) and (99) that Tr ρ(t) = 1. Then equations for the two independent components of the density matrix, ρ z = (ρ 11 − ρ 22 )/2, andρ 12 , omitting the dynamic corrections, read, These nonlinear equations are drastically different from the linear Bloch-Redfield equation describing decoherence of the macroscopic superconducting qubits, 22 and they have qualitatively different solutions, as illustrated in Fig. 6. Exact solutions for Eqs. (100), (101) read, where δρ z (t) = ρ 0 z − ρ z (t) is the deviation from the equilibrium, ρ 0 z = (1/2) tanh(βE * a /2). The evolutions of the diagonal (relaxation) and off-diagonal (dephasing) parts of the density matrix are qualitatively similar. One may distinguish the linear regime, δρ z (0) sinh(βE * a ) ≪ 1, when the decoherence is determined by the exponential law, δρ z (t) = δρ z (0) e −Γt ,ρ 12 (t) = ρ 12 (0) e −Γt/2 .
However, the decoherence rate, Γ, becomes exponentially small at temperature smaller than the qubit level spacing, βE * a ≫ 1. At this temperature, the most interesting is the opposite, nonlinear regime, δρ z (0) sinh(βE * a ) ≫ 1. In this case, there is a wide time interval, t ≪ sinh(βE * a )/ν, where both the relaxation and dephasing follow the power law (see Fig. 6), and only at very large times, t ≫ sinh(βE * a )/ν, the decoherence undergoes a crossover to an exponential regime similar to Eq. (104). We note that the exponentially small relaxation rate in the linear regime is well known for the quasiparticle recombination in bulk superconductors.

C. Evaluation of transition rate
We conclude our study with the evaluation of the phonon-induced transition rate, ν in Eq. (97). To evaluate the transition rate, one needs to specify geometry of the junction in more detail. Let us suppose that our adiabatic constriction, Eq. (4), is formed by a hard-wall potential and has an axial symmetry. Under these assumptions, the Fourier component of transverse wave function in Eq. (75) has the form, where r ⊥ (x) is the radius of the constriction cross section. The magnitude of the relaxation rate essentially depends on the parameter r ⊥ (0)Q, where r ⊥ (0) is the radius in the neck of the constriction, and Q = 2E * a /hs is the wave vector of phonons responsible for the interlevel transitions; for atomic-size constrictions, this parameter is small, r ⊥ (0)Q ≪ 1. Let us assume that the qubit level spacing is not too small, E * a /∆ ≫ s/v, then the phonon wave vector Q is large compared to the inverse penetration length of the Andreev level wave function, Q ≫ κ.
Let us for a moment assume that the Andreev level wave function does not spread out in the electrodes, but remains confined in the transverse direction, r ⊥ (x) = const = r ⊥ (0) (see Fig. 7); then the Fourier component in Eq. (106) is close to unity, and the interaction region in Eq. (75) is limited by the penetration length of the Andreev state, x ∼ 1/κ, restricting relevant phonon longitudinal wave vectors to small values, q x ∼ κ ≪ Q. The transition rate in this case reads, where τ −1 ph (E * a ) ∼ E * 3 a /hΘ 2 D is a bulk electron-phonon relaxation rate at the Andreev level energy (Θ D is the Debye temperature). This result has been derived in Ref. 32 (although neglecting the renormalization effect), and it can be qualitatively applied to long constrictions, whose length exceeds the coherence length. For short constrictions considered here, L ≪ 1/κ, the effect of spreading out of the Andreev level wave function is essential, and the approximation, r ⊥ = const, is not appropriate. Let us adopt the following model for the constriction shape, r ⊥ (x) = r ⊥ (0) (1 − |x|/L) −a , a > 1. It is easy to see then that the integral in Eq. (75) will have a cut off at |x| = L ≪ 1/κ. Function J 1 (z a )/z a in Eq. (106) can with good accuracy be approximated with the function, (1/2) θ(z 0 − z), where z 0 ≈ 1, and the integral in Eq. (75) is easily evaluated, giving, i.e., the transition rate in short constrictions significantly reduces. This is the effect due to small spatial region available for the Andreev level-phonon interaction in short QPC.

VII. CONCLUSION
Let us summarize the outlined theory for the Andreev level qubit (ALQ). ALQ belongs to the family of superconducting flux qubits, but it differs from the macroscopic flux qubits 3-6 in several important respects. First, the quantum hybridization of the flux (and persistent current) states in the ALQ loop is produced by electronic back scattering in the quantum point contact (QPC) rather than by charge fluctuations on the junction capacitors in the case of the macroscopic qubits. Thus, in principle, neither small junction capacitance nor large loop inductance is critical for the ALQ operation. Second, the ALQ is based on a QPC with large, almost full, transparency, in contrast to classical tunnel junctions employed in the macroscopic superconducting qubits. Large contact transparency is required for placing Andreev levels deep within the superconducting energy gap to achieve good decoupling from the continuum electronic states. To guarantee good separation of the qubit levels from the continuum, the amplitude of the phase fluctuations around the biasing point, φ ≈ π, must be restricted to the small values, which implies small inductance of the qubit loop.
In the tunnel junctions of the macroscopic qubits, Andreev levels are fast variables whose effect, after averaging, reduces to the Josephson potential energy added to the Hamiltonian of the loop oscillator. In the transparent junction of the ALQ, Andreev levels are slow variables which cannot be averaged out, and the full description includes the Andreev two-level Hamiltonian strongly coupled to the quantum loop oscillator. Derivation of the effective two-level Hamiltonian goes beyond the tunnel model approximation, which is done by incorporating exact boundary condition into the action of the contact.
The qubit read out is achieved via measuring fluctuating persistent current or induced flux in the qubit loop. To simplify the read out, the loop plasma frequency is chosen to be large compared to the qubit frequency so that the loop oscillator is "enslaved" by the Andreev levels. Other regimes, e.g. resonance between the Andreev levels and loop oscillator, could be also considered, however they remained outside the scope of the paper. Typical circuit parameters for the ALQ could be chosen as follows: E J ≈ ∆ for the single open mode, R ≤ 0.01 with I c ∼ 400nA for Nb; L ∼0.1nH, and C ∼ 0.1pF, giving ω p ∼ 10 11 sec −1 , and E L /h ∼ 10 13 sec −1 ∼ 10∆ N b /h, while E * a /h ∼ 10 10 sec −1 . Full description of the ALQ dynamics and qubit-qubit coupling, including qubit manipulation and read out is given by the single-particle density matrix.
Even if a good decoupling of the qubit states from the continuum electronic states is provided, there are still soft microscopic modes in the junction, which could couple to the Andreev levels. These modes present potential source of "intrinsic" decoherence, in addition to the commonly considered external decoherence, e.g., due to fluctuating biasing and read out circuits. We considered such an intrinsic decoherence of the ALQ related to acoustic phonons. It turns out that the collision terms in the kinetic equation are nonlinear, in contrast to the linear master equation for the macroscopic superconducting qubits. 22 This reflects the fermionic nature of the Andreev states and leads to considerable enhancement of the decoherence time at small temperature. One can understand this effect in the following way. Andreev levels belong to a many-body system of superconducting electrons. Although the macroscopic behaviour of the ALQ can be expressed in terms of the single-particle density matrix, the microscopic interaction with phonons involves two-particle correlation functions, which are sensitive to the fermionic nature of the Andreev states and obey the Pauli exclusion principle. This leads to the reduced probability of the phonon-induced interlevel transitions and hence to a slower decoherence.
Furthermore, the rate of the phonon-induced transitions between the Andreev levels is significantly reduced compared to the bulk transition rate. The reason is that both the Andreev levels belong to the same normal electronic mode; this together with a rapid spreading out of the Andreev level wave function in the contact electrode strongly reduces the relevant phonon phase space.
In the paper, only the case of a single-mode QPC was considered for clarity. However, the approach might be also relevant for macroscopic Josephson qubits with tunnel junctions. In junctions with disordered tunnel barriers, the open conducting modes with large transmissivity are present. 42,43 This introduces low-energy Andreev levels, which implies that quantum phase fluctuations become coupled to these Andreev levels, and the system must be described with the effective ALQ-type Hamiltonian. Another kind of qubits where the effective ALQ Hamiltonian might be appropriate are d-wave qubits: for the most interesting junction geometries, low-energy Andreev levels, midgap states, build up in the junction, which also makes the effective ALQ Hamiltonian to be relevant.

APPENDIX A: TUNNEL LIMIT FOR QPC
In this Appendix, we consider a low transparency QPC, D ≪ 1, and connect our method to the results of Refs. 18, 19 for tunnel Josephson junctions. In this limit, the Andreev levels lie very close to the edge of the superconducting gap, E a ≈ ∆, and therefore the field η is a fast variable, which has to be integrated out along with all the other electronic degrees of freedom. This will result in the effective action for the phase difference alone. After the integration over electronic fields, the propagator in Eq. (22) will take the form, Here g is the Green function defined in Eq. (32), and the matrix product in Eq. (A2) includes also the time convolutions. Taking advantage of the small D, and expanding the action (A2), the lowest order term reads: where the trace refers to the Nambu space. After taking the trace, the action can be written on the following form, with the kernels α and β given by where H (1) 0,1 are the Hankel functions of the first kind. Analytical continuation to the imaginary time in Eqs. (A6) and (A7) using relation, K ν (t) = (π/2)i 1+ν H (1) ν (it), leads to the same expressions for the kernels, α and β, as derived in Refs. 18,19,34, for tunnel Josephson junctions [with normal resistance of the tunnel junction being replaced with the normal resistance of the single-channel QPC: R N = (De 2 /πh) −1 ].
For a slow varying phase, φ, on the time scale of the kernel variations,h/∆, the both cosine terms in Eq. (A5) can be expanded over the relative time coordinate, τ = t 1 − t 2 , up to the second order, and the effective action takes a simpler form, where I c = 2e +∞ −∞ dτ β(τ ) = eD∆/2h is the critical Josephson current of the single-channel tunnel point contact, and is the correction to the contact capacitance due to the quasiparticle tunneling. 34 whereǦ ω (t) = dτ e iωτǦ (t + τ /2, t − τ /2) , and the products include the convolutions defined by the equation, Neglecting effect of the phonons on the retarded and advanced Green functions of the qubit, G R,A , we replace them in the right hand side of Eq. (B1) by the free Green functions, G R,A ω = (hω − E * a σ z ± i0) −1 , which are time-independent in the mixed representation. Thus, the time-dependence in the self-energy comes only from the Keldysh component, G K . The self-energy components take the form, where N Ω = e βhΩ − 1 −1 is the equilibrium phonon distribution function, and D Ω is the spectral weight function of the phonon bath, D Ω = 2π sgn(Ω) q |γ * q | 2 δ(|Ω| − Ω q ) .