Optimization of STIRAP-based state transfer under dissipation

Using a perturbative treatment, we quantify the influence of non-adiabatic leakage and system dissipation on the transfer fidelity of a stimulated Raman adiabatic passage (STIRAP) process. We find that, optimizing transfer time rather than coupling profiles, leads to a significant improvement of the transfer fidelity. The upper bound of the fidelity has been found as a simple analytical function of system cooperativities. We also provide a systematic approach to reach this upper bound efficiently.

Introduction.-State transfer, where an arbitrary quantum state is transmitted from a source to a target system, is a fundamental task in quantum state engineering. While sometimes it is possible to directly couple quantum systems (even of different nature), in many circumstances the transfer must be mediated by a third party (a quantum bus), whose coherence properties play a crucial role for the transfer process. To mitigate the effect of coupler decoherence, a particularly interesting class of indirect transfer protocols is based on the idea of stimulated Raman adiabatic passage (STIRAP), as it allows to perform the state transfer with negligible excitation on the quantum bus.
While the optimization of STIRAP for atomic population transfer was performed in a decoherence-free subsystem [31], with the objective of reducing the non-adiabatic leakage and of minimizing the pulse area or duration, the interplay between decoherence and non-adiabatic transitions is a crucial aspect for many current applications. Different from atomic systems, most solid-state quantum devices suffer significant dissipation. Hence, the prolonged operation time of STIRAP (required by adiabaticity) becomes a severe drawback: even if the source and target are relatively weakly coupled to the environment, the accumulation of errors over a long transfer time could still significantly degrade the transfer fidelity. Thus, in the presence of dissipation, optimization of the coupling profiles is essential to reach a trade-off between the non-adiabatic leakage and system dissipation, and to guarantee a high-fidelity state transfer.
To the best of our knowledge, a general understanding of this trade-off, as well as the resulting fundamental limitations, is missing in the literature. Some previous works considering the effect of dissipation on STIRAP have only included the decoherence of the intermediate level [32][33][34], often with a phenomenological approach [32,35]. In other studies, the dissipation is numerically simulated (see, e.g., [28,[36][37][38]). In order to provide physical insight and determine the power of STIRAP in the most realistic scenario -and especially to clarify the crucial interplay by the two types of dissipation mechanisms -we pursue here an analytical treatment, with full consideration of the system dissipation. We solve the master equation with a perturbative approach which is accurate in the desired high-fidelity regime. Although it is usually believed that the coupling profiles play an important role in optimizing the protocol, we quantify their influence and find it is relatively minor. Instead, the upper bound of the STIRAP fidelity critically depends on the transfer time and it is decided by a simple function of the cooperativities. We also provide a systematic optimization procedure to reach this upper bound efficiently. Our results can be applied to a variety of physical implementations such as optomechanics, circuit QED, and hybrid systems [39,40].
System and dynamics.-We consider the setup schematically illustrated in Fig. 1(a), where two qubits interact resonantly with a common bosonic mode. Using the rotating wave approximation, the Hamiltonian in the interaction picture is with G i (t) the tunable coupling strengths, a the bosonic annihilation operator of the bus, and σ ± = σ x ± iσ y . Such a Jaynes-Cummings Hamiltonian has been realized in cavity QED and various circuit QED architectures. Considering dissipation, the bus and qubits decay to their ground states with rates γ and κ i , respectively. The system dynamics is described by the master equation Our dis-cussion can also be easily extended by including pure dephasing terms 1 2 [41]. This description of the dynamics is commonly adopted in cavity QED, superconducting qubits, optomechanical systems, and might capture the main features of other types of qubits as well, although their decoherence dynamics in some cases can be more complex (see, e.g., [42]).
The general goal is to transfer an arbitrary state |ψ = c g |g + c e |e from qubit 1 to qubit 2. Considering a sufficiently low environmental temperature (to allow at most one excitation), the state transfer protocol is confined in a 4-level subspace formed by |1 = e (1) , 0, g (2) , |2 = g (1) , 1, g (2) , |3 = g (1) , 0, e (2) , |4 = g (1) , 0, g (2) . It is clear that in this low-excitation limit, each party can be either a qubit or a bosonic mode, i.e., the discussion below is also applicable to the case of two qubits coupled through a third qubit [13,43] or a solid-state qubit coupled to an optical cavity through a mechanical oscillator [44,45] (as illustrated in Fig. 1(b)).
An ideal transfer corresponds to c g |4 + c e |1 → c g |4 + c e |3 . At zero temperature, state |4 is stable and all the loss of fidelity is due to the transfer of the excited state. Hence in the following we only discuss c g = 0, which is the most demanding case and it is equivalent to a pure STIRAP process (marked by the blue shadow in Fig. 1(c)). The idea of STIRAP is to adiabatically tune G 1 (t) from zero to a finite value, while G 2 (t) is tuned from a finite value to zero such that the system evolves from |1 to |3 . In the whole process, the system adiabatically follows the instantaneous eigenstate (dark state): with tan θ(t) = G 1 (t)/G 2 (t) and where t f is the final operation time and it determines the overall speed of the transfer. The density matrixρ = U † (t) ρ I U (t) (where U (t) = k |k (t) k| and {|k(t) } are the instantaneous eigenstates of H I [41]) satisfies: The transformed dissipator is defined byLρ = U † L(UρU † )U and its matrix form can be found in a straightforward way [41]. The last two terms in Eq. (3) corrupt the desired transfer process: the first one represents the non-adiabatic leakage out of 2 (t) , which dominates for a fast-changing pulse; while the second term, i.e., the environment dissipation, dominates for a slow-changing pulse. How to reach an optimal trade-off between the two effects will be the central issue in the following sections. Before moving to that discussion, it is worth mentioning that the unitary part of Eq. (3) describes a fast spin precession around M (t) = G(t)ê z +θ(t)ê y , using a spin-1 language [41]. A perfect STIRAP state transfer can be realized when G, θ are constant and the transfer time is a multiple of the precession period (see Fig. 1(d)). In the absence of dissipation, there are also various other strategies to realize a perfect state transfer [41].
Perturbative treatment.-A successful transfer requires that the adiabatic dynamics plays a dominant role. Hence the effect of the non-adiabatic leakage and dissipation in Eq. (3) can be treated perturbatively. To do this, the density matrix is expanded asρ(t) =ρ (0) (t)+ρ (1) (t)+ ρ (2) (t) + · · · , which yields a corresponding expansion for 22 (t f ) (using the boundary condition Eq. (2)). The lowest-order result is F (0) = 1 and the higher-order contributions are obtained by iterative solution of the equation forρ (k) (t). Calculation details and the expressions including dephasing can be found in Ref. [41] and we report here only the simplified final results without dephasing. The first-order correction reads: which describes the loss of fidelity via qubit decay. The integrands has a rather transparent physical meaning: it simply reflects the decay of the time-dependent dark state Eq. (1). The next order contribution is: which is due to the non-adiabatic correction (i.e., the leakage out of the instantaneous eigenstate |2 (t) ). It vanishes when the system approaches the deep adiabatic limitθ(t) G(t). The last term has an oscillating dependence with respect to t f due to the spin-1 precession mentioned earlier. Finally, the dissipation of the quantum bus only enters the 3rd order contribution as the quantum bus can only be populated through second-order non-adiabatic leakage [see Fig. 1(c)]. This dissipation effect can be suppressed by long operation time and, not surprisingly, shows the same type of oscillating behavior of the non-adiabatic correction.
Optimizing the STIRAP operation.-Based on the perturbative treatment, we investigate the upper bound of the transfer fidelity and the best optimization strategy. We start with a special choice of the coupling profiles known as parallel adiabatic passage (PAP) [46], where G 1 (t) = G 0 sin θ(t) and G 2 (t) = G 0 cos θ(t). PAP is characterized by a constant energy splitting G(t) = G 0 and is commonly adopted by STIRAP protocols, as it allows one to suppress leakage errors by avoiding anticrossing points (see, for example, Refs. [31,47]). Due to the equal maximum couplings, PAP is a natural choice in the case of identical or similar qubits, which motivates us to take κ 1 = κ 2 = κ. The extension to the asymmetric case will be discussed later.
For PAP, F (1) = −κt f and the dependence on θ(t) vanishes.θ(t) is the only function to optimize and it can be written as a Fourier series: where α 0 = π/(2t f ) is fixed by the boundary condition Eq. (2) and α n (n > 0) represent a set of optimization parameters. At a certain transfer time t f , optimizingθ(t) requires ∂F/∂α n = 0, which can be easily solved since the fidelity is a quadratic form of the α n . The solution is: i.e., to obtain the optimal fidelity, the coefficients of the even terms are all equal; while the odd-n Fourier components vanish, which is a consequence of the symmetric setup (κ 1 = κ 2 ). In Eq. (8), N even is the total number  Fig. 2(a). Including more optimization parameters yields a higher fidelity. In the limit N even → ∞, we get α opt n∈even −π/ (2t f N even ) and the optimized fidelity over coupling profiles is F opt (t f ) = 1 − κt f − π 2 γ/ 4G 2 0 t f , which shows the competition between qubit decay and non-adiabatic transition with respect to the transfer time. If we further optimize over t f , the maximum transfer fidelity is: with the corresponding optimal transfer time t opt f = π/(κ √ C), and C = 4G 2 0 /γκ is the system cooperativity. This expression provides the largest attainable fidelity when an arbitrary state is transferred with a STIRAPlike pulse in a symmetric configuration and for a given set of parameters G 0 , κ and γ. Defining the loss of fidelity ∆F = 1 − F , we show in Fig. 2(b) that Eq. (9) is in excellent agreement with the numerically optimized fidelity, even with γ ∼ G 0 . Thus, Eq. (9) shows that the ultimate power of STIRAP is limited by the cooperativity whereas the transfer by sequential swapping requires the strong coupling condition (the fidelity to swap |1 → |2 followed by |2 → |3 is F ∼ 1 − (κ + γ)/G 0 , thus requires γ, κ G 0 ).
The maximum fidelity for a finite number of optimization parameters can also be found as F max (N even ) = 1 − 2π ξ (N even ) /C with ξ(N even ) ≡ (N even +3/2)/(N even + 1) [41]. Thus, for a given setup, the goal of the fidelity (i.e., the acceptable deviation form the upper bound) determines N even , which subsequently determines the coupling profiles and operation time. Systematic improvements of the pulse shape can be gained by progressively increasing the number N even of optimization parameters.
It is interesting to notice that, even taking N even = 0, i.e., adopting the simple dependence θ(t) = πt/(2t f ), the transfer fidelity is F max (N , which is just sightly smaller than F max . This observation suggests that, for a STIRAP-based state transfer under dissipation, optimizing the operation time (i.e., the operation speed) is far more efficient than introducing a complicated pulse. This is demonstrated in Fig. 2(c), where the maximum fidelity grows quickly with N even and moderate values are sufficient to achieve a small deviation from the upper bound.
We now go beyond PAP by considering generic coupling profiles. Notice that (see Eq. (5-6)) the key variable in F (2) and F (3) isθ(t)/G(t), which can be similarly decomposed in a Fourier series: The main difference is that β 0 = t −1 f t f 0 dτθ(τ )/G(τ ) cannot be fixed by the boundary condition, but it is constrained by the maximum coupling strengths G 1,max and G 2,max [41], so we can still perform the minimization with respect to β n . The procedure closely patterns the PAP case, for example, β opt n∈even is given by an expression similar to Eq. (8) except that the pre-factor π/(2t f ) is replaced by β 0 , and cos G 0 t f by cos( t f 0 G(τ )dτ ). The optimization over the pulse shape yields (N even → ∞) [41]: Further optimizing over t f leads to: where we indicate with C i = 4G 2 i,max /γκ i the cooperativity of side i = 1, 2. In the limit C 1 = C 2 = C, this result allows for a slightly larger fidelity than Eq. (9) (4 √ 2/(2π) 0.9). This can be attributed to the fact that here G(t) is not kept constant. Equation (12) also shows that in the limit when C 1 and C 2 are very different, the less coherent system (with lower C i ) dominates the fidelity loss. In Fig. 3(a) we compare our analytic expression to the numerical optimization, and show that Eq. (12) is indeed an accurate characterization of the maximum fidelity. Finally, we consider a large difference in the qubit coherence (e.g., κ 1 κ 2 ). Since F (1) also depends on θ(t) (see Eq. (4)), the previous analytic approach is not easily applicable. Still, one can understand the general parametric dependence by replacing the second term in Eq. (11) with F (1) ∼ −(κ 1 + κ 2 )t f /2. Note that large deviations from this estimate occur when θ(t) spends a significant fraction of the transfer time close to the initial or final angles θ = 0, π/2. However, this situation is essentially equivalent to a shorter t f . The argument indicates that the relevant figure of merit is the smallest cooperativity C min = √ γκ max /G 2 min , where κ max = max[κ 1 , κ 2 ] and G min = min[G 1,max , G 2,max ]. We find that this conclusion is generally in agreement with the numerics. For example, Fig. 3(a) shows that when G 1,max = G 2,max , the maximum fidelity is well approximated by: which is in agreement with Eq. (12) when C 1 = C 2 . When C 1 is very different from C 2 , Eqs. (12) and (13) differ in the numerical prefactors but confirm the general argument that ∆F ∼ 1/ √ C min , for an optimized transfer time t opt f ∼ 1/(κ max √ C min ). Conclusion and Remarks.-We have analyzed the STIRAP-based state transfer between two qubits interacting with a common harmonic mode. The trade-off between the non-adiabatic transitions and dissipation needs careful optimization of the time-dependent couplings. Instead of a numerical approach, we have developed an analytical treatment on the general parameter dependence, which has allowed us to reach a physical understanding of the optimal transfer time and the upper bound of the fidelity. We also showed how to reach this bound efficiently. Our results provide a useful guideline for implementations of such protocol in a variety of physical systems, as well as generalizations to alternative setups and more complicated adiabatic transfer schemes.
YDW and SC acknowledge the support from Chinese Youth 1000 Talents  Supplemental material for "Optimization of STIRAP-based state transfer under dissipation"

UNITARY DYNAMICS AND TRANSFER FIDELITY
We discuss briefly here a simple optimization of the fidelity based on an exact solution of the unitary dynamics, and compare the STIRAP scheme with other types of state-transfer protocols.
First we consider Eq. (3) of the main text and note that the first two terms on the right side yield an effective HamiltonianH where we map |1 , |2 , |3 to the spin-1 angular momentum states (| + 1 , |0 , | − 1 , respectively). As seen, the time evoultion is equivalent to the rotation of a spin in a time-dependent magnetic field (see Fig. 1(d) in the main text). For a coupling profile of the simple sin-cos form: the magnetic field is constant and tilted from the z direction by an angle φ = arctan[π/(2G 0 t f )]. Since the initial state is an eigenstate of J z (i.e., is not aligned to M ), an oscillatory dynamics follows, with angular frequency given by | M | = G 2 0 + (π/2t f ) 2 . The exact solution for the fidelity is: Besides recovering F → 1 in the adiabatic limit t f 1/G 0 , we also notice that, even for a much shorter operation time t f ∼ 1/G 0 , a perfect transfer is possible at the discrete values G 0 t f = 2π n 2 − (1/4) 2 (n = 1, 2, . . .). Such a condition can be matched by tuning the transfer time t f or the coupling strength G 0 . Physically, these times correspond to the periodic return of the state to the z direction.
It is certainly not surprising that an ideal state transfer can be realized in the absence of dissipation. In fact, more efficient and natural ways exist. Among these, we would like to mention the well-known sequential swapping, where a constant value of G 1 (with G 2 = 0) transfers |1 → |2 and then a second swap leads to |2 → |3 . The total transfer time is π/G 0 (if both coupling strengths are equal: G 1 = G 2 = G 0 ), which is slightly shorter than the minimum transfer time 2π 1 − (π/4) 2 /G 0 1.2π/G 0 implied by Eq. (S3). Another way which is operationally simplest is to set both couplings to a constant G 0 and let the system evolve for a time π/ √ 2G 0 , which in the spin-1 language corresponds to a π-rotation around the x axis.
Nevertheless, the choice Eq. (S2) is interesting in the context of the present discussion. On one hand, it shows that a STIRAP-like time-dependence of the couplings still allows for a perfect state transfer deep in the non-adiabatic regime (i.e., with t f ∼ 1/G 0 ). On the other hand, this exact solution can serve as a useful reference for our discussion of the general case, including dissipation and more complicated choice for θ(t). In fact, we can use the solution Eq. (S3) to estimate the fidelity loss due to the decay of the intermediate system: Without dissipation, the population of the excited quantum bus due to the unitary time evolution is: (S4) The population loss through the bus is approximately: As discussed in the main text, this result is in agreement with the perturbative calculation. Combining this result with the cavity damping κt f , one can optimize the transfer time and find the fidelity F max (N even = 0) = 1 − 1.22 2π √ C . Thus, this simple example allows one to understand the relevant physics and typical time scales characterizing the state transfer process, while further optimization of θ(t) leads to a modest improvement of the fidelity.
In the reference frame defined by |k(t) , we write the density matrix asρ(t) =ρ (0) (t) +ρ (1) (t) +ρ (2) (t) + · · · , where the lowest-order satisfies: The solution is simplyρ where the non-adiabatic term in matrix form is (for simplicity, we omit the explicit time dependence, i.e.ρ =ρ(t)) The full expression of the dissipator is too cumbersome to present here. We only show the contribution from the quantum bus:L which is independent onρ 22 . Since at zero orderρ Eq. (S10) shows thatρ (1) (t) has no contribution proportional to γ. To evaluate the first-order correction to the fidelity, the full expression of Lρ 22 is necessary: Equation (S8) can be solved iteratively, to yield the perturbative expansion of the transfer fidelity F =ρ 22 (t f ). In particular, integrating Eq. (S8) gives: where ξ (k−1) (τ ) ≡ [µ,ρ (k−1) (τ )], L (k−1) ab (τ ) ≡ (Lρ (k−1) (τ )) ab , and ∆E ab (t) = E a (t) − E b (t), with E a (t) = G(t)(δ a,1 − δ a,3 ). Equation (S12), together with Eq. (S11), yields the 1st order correction: This represents a generalization of Eq. (4) of the main text (where we assumedγ ϕ = 0). At 2nd-order the disspation of the bus does not contribute to the fidelity either, because [L mρ ] 22 is identically zero [see Eq. (S10)]. This observation, together with the independence ofL mρ onρ 22 , reflects the fact that |2 does not involve excitations of the bus (it is a dark state). Thus, the populationρ 22 is not directly affected by a finite γ and ρ 22 has no influence onL mρ . The second order correction can be written as: which takes into account the corruption of fidelity due to purely non-adiabatic leakage. While a full calculation of the 2nd-order result should take into account the qubits dissipation, these corrections are neglected here. This is due to the fact that the qubit dissipation already enters the 1st-order result and higher order terms involving κ i , γ (i) ϕ should be much smaller.
We then consider the 3rd order correction, which finally yields a contribution proportional to γ. As illustrated schematically by Fig. 1(c) of the main text, the loss of fidelity via the mechanical damping is through the non-adiabatic leakage. This is also clear from Eqs. (S9) and (S10). Consider, for example, that at 1st order the non-adiabatic leakage leads to a finite valueρ (1) 12 ∝θ. At 2nd order, this correction gives a contribution toρ (2) 12 ∝θγ (cf. Eq. (S10)). Finally, the non-adiabatic terms leads to a correction toρ 22 fromρ  (S15) The above formulas for F (2) and F (3) can be simplified in the relevant case of a sufficiently large G(τ ). In fact, they are given by integrals of the form t 0 dt f (t ) exp ±i t 0 dτ G(τ ) , where f (t ) is a relatively smooth function while exp ±i t 0 dτ G(τ ) is a fast oscillating factor. By performing multiple integrations by parts, a systematic expansion of such integrals in powers of G(t) −1 can be derived. A straightforward but tedious calculation yields the leading-order results, cited in Eqs. (5) and (6) of the main text.

MAXIMUM FIDELITY WITH UNEQUAL COUPLINGS
In the PAP case, the discussion is based on equal maximum couplings and equal damping rates for both qubits.
Here we relax the constrain of equal maximum couplings and investigate again the upper bound of the fidelity.
To do this, we writeθ(t)/G(t) into a Fourier expansions with β n the new set of optimization parameters (see Eq. (10) in the main text). This makes the optimization very similar to the PAP case, except two difficulties. The first one is that not onlyθ(t)/G(t), but also cos t f 0 G(τ )dτ is affected by the coupling profiles (while in the PAP case this gives a simple constant, cos G 0 t f ). The second complication is thatθ(t)/G(t) must satisfy more involved constrains thanθ(t). In fact, at each value of θ(t) the coupling G has a maximum value G max (θ) [see Eq. (S23) below]. It is not immediately clear how the condition G(t) ≤ G max (θ(t)) (together with the old one, t f 0θ (τ )dτ = π/2) can be simply written in terms of the β n . To avoid these difficulties, we first perform the maximization by considering the β n>0 and cos t f 0 G(τ )dτ as arbitrary parameters. This approach is useful because removing these constrains gives an upper bound to the fidelity, at fixed β 0 . Furthermore, in performing the calculation, we will also see how the above two points can be resolved.