Completely Positive Map for Noisy Driven Quantum Systems Derived by Keldysh Expansion

Accurate modeling of decoherence errors in quantum processors is crucial for analyzing and improving gate fidelities. To increase the accuracy beyond that of the Lindblad dynamical map, several generalizations have been proposed, and the exploration of simpler and more systematic frameworks is still ongoing. In this paper, we introduce a decoherence model based on the Keldysh formalism. This formalism allows us to include non-periodic drives and correlated quantum noise in our model. In addition to its wide range of applications, our method is also numerically simple, and yields a CPTP map. These features allow us to integrate the Keldysh map with quantum-optimal-control techniques. We demonstrate that this strategy generates pulses that mitigate correlated quantum noise in qubit state-transfer and gate operations.


Introduction
The ubiquity of decoherence errors in current quantum computing platforms poses a bottleneck for performing error-correctable quantum computation [1].Further reducing these errors relies on the accurate modeling of them, which is challenging due to the presence of complicated drive and noise background.For example, if the quantum system is strongly driven, or the noise (quantum and classical) is correlated, the widely-used Lindblad master equation is generally not applicable [2][3][4][5][6][7][8].
To obtain more accurate predictions, recent research is exploring generalizations of the Lindblad master equation, where the constant damping operators and rates in the original form are replaced by time-dependent ones [2][3][4][9][10][11][12][13][14][15].Impressively, after a more careful treatment of the bath degrees of freedom than that in the Lindblad formalism, the master equation not only becomes compatible with drives and correlated noise, but also maintain the property of generating completely positive and trace-preserving (CPTP) maps [2,3,[12][13][14][15].In addition to this route, a formalism based on filter functions has been developed to model errors caused by correlated noise.This formalism can predict the sensitivity of the driven system to noise at different frequencies [6,7,[16][17][18][19][20][21].As a comparison, the filter-function method usually requires fewer integrals, and Ziwen Huang: zhuang@fnal.govhas a clearer physical picture of how noise at difference frequencies contributes differently.However, this method mostly focuses on classical (or dephasing) noise, and does not always guarantee the CPTP character of the map.
In this paper, we present a decoherence model which combines the advantages of the two routes mentioned above, and is tailored for optimizing gate operations.Our method belongs to the filter-function category, while the Keldysh technique [9,22] used here extends the scope of the formalism in Ref. [6,17] to quantum noise.Furthermore, the map derived by this method is guaranteed to be CPTP, after a special secular approximation to the filter functions.Such approximation also significantly simplifies the calculation, which further allows us to explore error-mitigation strategies by integrating our method with the technique of quantum optimal control [7,23].Using a few examples, we show that such a combination can generate pulses that suppress decoherence errors induced by correlated quantum noise.
The paper is structured as follows.In Sec. 2, we outline the derivation of the Keldysh map.The main results are summarized in Eqs.(19) and (20).In Sec. 3, we apply our method to a variety of quantum systems, which not only reproduces some familiar results, but also extends the prediction of decoherence errors to several less familiar situations.In Sec. 4, we integrate the Keldysh method with the quantum-optimal-control technique, and demonstrate improvement of gate and state-transfer fidelities via the optimization of drive pulses.

Formal Keldysh expansion
We start by deriving the formal map for the qubit density matrix using the Keldysh expansion.The Hamiltonian of the full system is Ĥ () = Ĥ () + Ĥ +  Ĥ , where Ĥ (), Ĥ , Ĥ denote the Hamiltonians for a driven quantum system, bath, and interaction.The system Hamiltonian Ĥ () = Ĥ0 + Ĥ () comprises the static Hamiltonian Ĥ0 and the drive operator Ĥ ().For simplicity, we specify the system-bath interaction by Ĥ = x η, where x and η are the system and bath operators, respectively1.The small dimensionless parameter  is used to keep track of the order.In Eq. ( 1), we assume that the interaction is weak and can be treated perturbatively.To conveniently perform the perturbative calculation, we first move to the interaction picture with the unperturbed propagator Û0 () = Û () ⊗ Û (), where the partial propagator for the sys- Here, the interaction-picture propagator is given by Ũ () = T exp[− ∫  0  H ()] and H () = Û † 0 () Ĥ Û0 () is the system-bath coupling term in the interaction picture.We use ρ (0) and ρ (0) to denote the initial partial density matrices for the system and bath.Note that in this work, we assume that there is no entanglement between the system and bath initially, and the bath is prepared in its thermal equilibrium ρ,eq .
To evaluate this formal expression, we expand Ũ () =  Ũ ()  () as a Dyson series, where the th term Ũ ()  () is given by Inserting this into Eq.( 2 To simplify this expression, we define the th-order map and the sum map as which casts Eq. ( 4) into ρ () = () ρ (0).
1In the main text, we focus on time-independent x and η for a concise description, but the derivation in Sec. 2 can be straightforwardly generalized to cases where the system operator is time-dependent.Such a situation becomes important in several experiments, e.g., Refs.[24,25] Although we can in principle use Eq. ( 5) to calculate the map () to arbitrary order, it is usually not the most convenient quantity to extract physical measurables from, according to the discussion in Refs.[9,26,27].Instead, we follow the Keldysh theory and define the self-energy () ≡ ln[()], (7) where redundant terms in higher-order expansions can be conveniently identified, and the derivation of quantities such as relaxation rates is easier [9,26,27].Below, we will focus on ().Similar to (), the self-energy () can be expanded in powers of  by () =   () (), where  () () can be derived from  () () by a Taylor expansion.For example, the lowest two orders are related by  (1) () =  (1) (),  (2) () =  (2) () − 1 2  (1) ()
These relations can be further simplified, if we assume that the noise has a zero mean, Tr  { η() ρ (0)} = 0.In that case, the first-order map  (1) () vanishes, resulting in the following simplified relations

Second-order truncation
For most experiments involving gate operations and state transfer, it is usually sufficient to estimate the decoherence error up to leading order.In the following, we focus on the leading-order self-energy  (2) ().Explicitly, the secondorder self-energy takes the form With the knowledge of the noise spectrum, we can further expand the right-hand side of Eq. (9).For example, the first term can be expressed as and the other two can be derived similarly.These integrals can be conveniently summarized by four Keldysh diagrams, which we show and explain in Appendix A. In the equation above, is the noise spectrum, and the interaction-picture operators are derived as x() = Û †  () x Û () and η() = Û †  () η Û ().
Given the information of the noise spectrum   () and system propagator Û (), we can use Eqs.( 9) and (10) to calculate the approximated dynamical map [2,[13][14][15], i.e.2, This calculation is reminiscent of the filter-function method shown in Refs.[6,7,[16][17][18].For example, the double integral is closely related to the filter functions studied there.(See Appendix B for a detailed discussion of the connection to these theories.)For comparison, the time ordering of the two coupling operators η( 1 ) and η( 2 ) in the Keldysh expansion resolves the asymmetric noise spectrum for a non-classical noise source.However, evaluating the triple integral in Eq. ( 10) accurately is usually not numerically efficient.In the next section, we show that appropriate approximations can simplify the calculation.

Fourier expansion and secular approximation
Our strategy for solving the aforementioned problem is based on the Fourier expansion of where we define the fundamental frequency   = 2/.Inserting it into Eqs.( 10) and ( 9), we find Here,  , ′ () is the filter function defined by 2For a noise bath initiated as  ,eq =    | ⟩ ⟨  | (  > 0), the formal map (11) without any approximation or truncation is already formally CPTP.(The proof of this is analogous to that in [2] for a static system.)Note that sufficiently accurate calculation of the integrals (13) may be required toward a numerical CPTP map.(The analytical evaluation of this integral is discussed in Appendix A.) For  =  ′ and  ≠  ′ , the filter functions  , ′ () behave differently, which is worth careful inspection.
Diagonal filter functions.-For =  ′ , the filter functions  , ′ () can be cast into the following form where the real and imaginary parts are given by: We illustrate the behavior of   () (solid line) and   () (dashed line) in Fig. 1.Panel (a) shows the real and imaginary parts of  0,0 (), respectively.Visibly, the function   () has a predominant peak located at  = 0.According to Eq. ( 16), the width of this peak is ∼ 2/.Different from the real part, Im{ 0,0 ()} =   () flips its sign at  = 0, showing both a peak and a valley.Compared to   () ≲ 2|| All these properties can be observed in Fig. 1  Based on these three features, we arrive at the following conclusion: if variations of   () are insignificant over the frequency scale of a few   , the off-diagonal elements of with  −  ′ ≠ 0 have negligible amplitude.We justify this claim in three steps.First, for large | −  ′ |, the amplitudes of  , ′ () are small, rendering a negligible  , ′ .Second, for terms with small but nonzero | −  ′ |, the slow-varying   () allows us to treat it as quasi-constant.Third, using the property of net-zero area in Eq. ( 17), the integral  , ′ vanishes for small but non-zero | −  ′ |.Since  , ′ are the coefficients of terms shown in Eq. ( 13), we conclude that all off-diagonal terms are less important than the diagonal ones in that expansion, if the spectrum is sufficiently smooth at the resolution determined by   .
After neglecting the terms with off-diagonal filter functions, we simplify Eq. ( 13) to  (2) () ρ (0) ≈  (2)  CP () where  (2)  CP () is the simplified second-order self-energy.This step resembles the secular approximation performed in the derivation of the Lindblad master equation -in both cases, small off-diagonal terms are neglected.Since the coefficient ∫ ∞ −∞ (/2)  ()2  ( −   ) is strictly positive, the self-energy  (2)  CP () has the form of a Lindbladian (up to an extra time dimension) [28].Then, according to Ref. [28], the exponential of  (2)  CP () yields a CPTP map In the following, we refer to x = [ ∫  0  x() exp(  )]/ and   =   as the filter operator and its corresponding filter frequency, respectively.
According to the justification of the secular approximation, the map (20) tends to be more accurate for smoother noise spectrum spectra.However, we observe that even for an   () that exhibits strong peaks, the magnitude of the terms in Eq. ( 13) with diagonal filter functions can still dominate those with the off-diagonal ones.As a result, the secular map (20) is found to qualitatively agree with the full map for many common noise spectra, including those showing strong peaks.This is illustrated in Sec.3.3 for the example of 1/  noise.Therefore, although we find it challenging to quantify the magnitude of the approximation error for arbitrary noise spectra, we still adopt Eq. ( 20) for an estimation if an extreme spectrum is considered; then, the agreement between the secular and full maps can be checked for validation.
We append two remarks to compare our method with several existing ones.First, although our method is not a differential equation, the derived map is reminiscent of the dynamical map generated by the Lindblad equation [28,29].Specifically, the first and second lines in Eq. ( 19) resemble the damping terms and Lamb shift in the master equation, respectively.A difference is that, our map considers noise contributions from the frequency set {  | ∈ Z}, while the Lindblad master equation only includes noise at system transition frequencies.For a more intuitive comparison, we illustrate the decoherence channels for an undriven and driven qubit in Fig. 2 (a) and (b), respectively (see a detailed description of the qubit in the caption).While three damping operators fully describe the decoherence processes in an undriven qubit according to the Lindblad master equation, more damping terms are relevant for a driven one according to Eq. (19).In fact, we show in Sec.3.1 that the Lindblad map is a special case of the dynamical maps derived by our method.
Second, the form of Eq. ( 19) is also reminiscent of the coarse-grained master equation [2,3].We understand this similarity as follows: the second-order secular Keldysh expansion used here is comparable to the coarse-graining step detailed in Refs.[2,3].Differently, our framework focuses on a single map of the reduced density matrix from  = 0 to , rather than its full evolution during  ∈ [0, ].The calculation of the Fourier series of x() can also be easily carried out numerically, e.g., using the method of fast Fourier transformation.Therefore, our model tends to require less computational resource due to fewer compu-  , which describe the excitation, decay, and pure-dephasing channels of the qubit, respectively [29].In (b), because an arbitrary drive is applied on the qubit, the previous simple decomposition of x() no longer holds; instead, one can use Eq. ( 12) to perform the frequency decomposition over the time period  ∈ [0, ].As a result, the decoherence channels for a driven qubit are given by the first line of the self-energy (19).tational steps, if the map is only needed for one final time .

Total decoherence error
Using the map (20), we can conveniently derive the process infidelity for a noisy quantum processor.Following Ref. [30], this error is expressed as Above,   is the dimension of the system Hilbert space, V tg ≡ Ûtg ⊗ Û † tg denotes the target superoperator, where Ûtg is the target unitary, and V  () is the closed-system map defined by V  () ≡ Û () ⊗ Û †  ().If we only focus on the decoherence contribution to  gate , we can neglect the possible coherent errors by setting  tg =   (), and reduce Eq. ( 21) to Above, we have used the leading-order approximation () ≈ I  +  (2)  CP ().The trace Tr{•} is for the superoperator for the system density matrix, while Tr  {•} denotes the usual trace for regular system operators3.This approximation evaluates  dh to the order  2 .Up to this order, only the real part of  , contributes.
We interpret the sum in the second line of Eq. ( 22) as follows.The total decoherence error  dh is a sum of contributions by noise from frequency bands indexed by .The th band has the approximate bandwidth ∼   and is centered at   [see the filter function in Fig. 1 (a)].The total noise amplitude over this bandwidth is given by the integral 2Re{ , }.The driven qubit, however, is not equally sensitive to noise from all these frequency bands -according to Eq. ( 22), we can quantify the sensitivity by the filter strength which satisfies the conservation rule The conservation rule implies that, if only white noise is present, i.e., the noise spectrum   () =  is a constant over frequency, the decoherence error is which increases with  but is independent of the shape of drive applied during  ∈ [0, ] up to order  2 .For non-Markovian noise [  () is not a constant], however, Eq. ( 24) does not hold.In this case, different drives generally result in different magnitudes of  dh .To reduce decoherence errors, one should design pulses such that   is suppressed where the integrated noise amplitude Re{2 , } is large.In the following sections, most strategies discussed for reducing decoherence are centered around this strategy.

Applications
In this section, we demonstrate the power of our framework by a few examples from a wide range of applications.Our framework not only reproduces some of the established conclusions, but also extends the prediction to situations that have not been carefully studied by previous theories.

Static quantum systems
We first apply our method to an undriven system, with the main purpose of reproducing the dynamical map derived by the Lindblad master equation.In this case, we set the drive Ĥ () = 0 in Eq. ( 1), which makes the system Hamiltonian Ĥ () = Ĥ0 time-independent.Following the procedure described in Sec. 2, we first derive the interaction-picture coupling operator x(), and use it to find the filter operators x needed in the derivation of the self-energy (19).For the undriven system, the propagator is given by which yields the coupling operator in the interaction picture Above, the frequencies   associated with different terms are contained in the frequency set For this expression,   is the dimension of the system Hilbert space, and   is the eigenenergy of the th eigenstate for Ĥ .We refer to the elements in F as the transition frequencies, and the corresponding x(  ) as the damping operator.(Note that   = 0 is also included in this transition-frequency set F.) Using Eq. ( 12), we obtain the th filter operator where we define (  , ) Inserting the expansion (27) into Eq.( 19), we obtain the self-energy for the undriven system as where we define the damping operator Inserting it into Eq.( 20), we obtain the dynamical map for the undriven system.In deriving this map, we only specify a time-independent Hamiltonian, but do not make further assumptions such as those usually required by the Lindblad master equation.
If we do enforce these assumptions, then our framework reproduces the Lindblad dynamical map.In detail, these conditions are: ① the difference between the transition frequencies is much larger than   , i.e., |  −  ′  | ≫   for   ,  ′  ∈ F,   ≠  ′  .② the system evolution time  is sufficiently long such that the spectral variation in   () is negligible over the small frequency scale   = 2/; These two conditions can be translated to the following two more familiar statements: ① the system's characteristic time is much shorter than the system evolution time  of interest, which usually has a similar timescale as the system relaxation time   ; ② the bath correlation time   is also much shorter than the evolution time  ∼   (see Appendix C for a more detailed explanation).Under ① and ②, the fundamental frequency   is by far the smallest frequency scale.This allows us to consider the limit   → 0, and perform the following three approximations.
First, we approximate the function   () ≈ () and   () ≈ − P (1/), where () is the Dirac delta function and P denotes the Cauchy principal value.Using the approximated  / (), we simplify the integral  , and obtain where we define Second, the infinitesimal   also justifies the replacement of the summation over  by an integral over .This step transforms the self-energy in Eq. ( 19) to an integral  (2)  CP ρ where we define x ≡ x⌊ /  ⌋ / √   .Third, using the expansion (26) for the undriven system and the definition of x , we find the approximation A similar delta-function approximation holds for x †  ⊗ x in D [ x𝜔 ].Inserting these approximations into the integral Eq. ( 31) and carrying out the integral over the delta functions, we finally arrive at the self-energy The secular CPTP map (20) generated by the self-energy above is identical to that predicted by the Lindblad master equation.
As a minimal example, we apply the map above to a static qubit that is described by the Hamiltonian where   is the qubit frequency.The transition-frequency set for this qubit is F = {0, ±  }.If it is transversely coupled to a bath through the operator x = σ , the expansion ( 26) is given by Then, following the steps described above, we find the self-energy under conditions ① and ② approximated as If the coupling operator is x = σ , we instead have x() = σ , which yields the approximated self-energy

Weakly driven systems
Although the derivation of the Lindblad dynamical map [29] assumes no drive on the system, in the literature, weak drives are sometimes naively added to the master equation with the damping rates and operators unaffected.Different from the Lindblad method, our framework rigorously includes the drive in the derivation.Below, we use the Keldysh framework to investigate the change of the decoherence map (20) if a weak drive is added.Because the specific noise spectrum may vary in different experiments, here we choose to focus on the filter operators x , which determine the map (20) up to the specific noise spectrum.We start by considering a general system, which is described by the Hamiltonian Ĥ () = Ĥ0 +  Ĥ ().This Hamiltonian consists of the static part Ĥ0 and a sufficiently weak driving term  Ĥ () ( is a small dimensionless parameter).Due to the small amplitude of the latter term, we can perturbatively calculate Û () using the Magnus expansion [31]: (38) where Û0 () is the propagator for the undriven system given in Eq. ( 25), and Ω() is the Magnus exponent.To leading order of , this exponent is approximately By inspecting Eq. ( 38) and the definition of x(), we note that: if Ω() is small, x() can be approximated as (40) which is only slightly modified from Eq. (26).In the limit | Ω()| → 0, the filter operator x can still be approximated by the undriven expansion (27).Therefore, the noise channels and the resulting dynamical map (20) should also approach those obtained in the undriven case.By contrast, if the condition of negligible Ω() is not satisfied, such approximation may be invalid.
In the following, we will use a concrete example to concretely demonstrate both scenarios.We consider a qubit driven by a sinusoidal tone.The Hamiltonian of this driven system is given by The coupling operator for this qubit is taken to be x = σ , which corresponds to the transverse coupling between the qubit and the noise bath.The drive strength  is assumed to be weak, i.e.,  ≪   .In that case, the Magnus expansion of the qubit propagator is applicable, with the exponent given by where the detuning is defined by   ≡   −   .From Eq. ( 42), we find that the Magnus exponent has a negligible magnitude, i.e., | Ω()| ≲ 4|/  |, if the drive is off-resonant ( ≪ |  |).For such off-resonant drive, Eq. ( 40) predicts that the expansion of x can be approximated by Eq. (35).To verify this, we numerically [32] calculate x and the resulting filter strength   for both undriven (red coloring) and off-resonantly driven (blue coloring) qubits.The resulting filter strengths   versus filter frequencies  =   are shown in Fig. 3 (a).[We only focus on the frequency range  ≈   as an example, where noise induces energy decay in the qubit.For qubit excitation, the discussion is analogous.]As shown in the plot, the filter strengths for the driven qubit only insignificantly differ from those for the undriven qubit.The filter operator associated with the most prominent filter strength is approximately σ− , which is the decay operator for the undriven qubit (red peak).
For the resonantly-driven qubit, however, the exponent Ω() can grow significantly, even in the limit  ≪   .Choosing   = 0 as an example, we find that the exponent grows linearly with time.This exponent leads to the expression for the coupling operator Certainly, Eq. ( 44) cannot be approximated by Eq. ( 26), rendering the previous approximation of x by Eq. ( 27) invalid.The difference in x between the driven and undriven cases causes distinctive behaviors of   's, as shown in Fig. 3 (b).Compared to the plot of   for the undriven qubit (red coloring), the plot for the driven qubit (blue coloring) exhibits two additional peaks located at frequencies  =   ± .These extra peaks imply additional decoherence channels, rendering the qubit sensitive to noise at frequencies   ±  in addition to its transition frequency   .These additional damping channels are missed by the standard Lindblad method [Eq.(36)], but can be crosschecked by a Golden-rule type of calculation in the rotating frame [33], or the Floquet theory [5].(In Appendix D, we explain the appearance of the side peaks in the framework of the Floquet master equation.)In summary, the filter operators x for the driven system in general differ from those in the undriven case.This further results in different dynamical maps (20).The approximation of one set of x by the other is possible if the drive is off-resonant such that | Ω()| is sufficiently small4.
Finally, we use a third example to demonstrate the predictive power of our method in more complicated situations involving non-periodic drives.For example, we calculate x and   for a qubit driven by a pulse with a hyperbolic envelope [see inset of Fig. 3 (c)].Compared to the plot of filter strengths in (b), the ramping up and down of the drive result in two wider side peaks that are not centered at the maximal driving strength, as shown in (c).Such distinctive feature implies the difference in the decoherence processes between systems with periodic and non-periodic drives.The latter case is thought to go beyond the description by the rotating-frame analysis or the Floquet theory, but is conveniently captured by the Keldysh method.

Ramsey, echo and 1/ 𝑓 noise
Along with the introduction of Eq. ( 19) in Sec.2.3, we claim that the secular approximation is applicable even for spectra that show strong variation within the frequency scale characterized by   .Here, as a supporting example, we study the state evolution of a qubit which is coupled to a 1/  noise source, whose spectrum is strongly peaked at  ≈ 0. Particularly, we compare the prediction by the secular CPTP map (20) and the full-wave version (13) for this example.
We consider a qubit longitudinally coupled to a noise bath and subject to a transverse drive.The Hamiltonian for this qubit is Ĥ () =   σ /2 +  () σ , and the coupling operator is x = σ .The 1/  noise spectrum is given by   () = 2A 2  /||, where we set an infrared cutoff frequency  ir to regularize the singularity at  = 0.
We first calculate the map (20) for the simple case of an undriven qubit [ () = 0], which is relevant for a Ramsey experiment [34].Different from the discussion in Sec.3.1, the presence of the strong peak in the noise spectrum violates condition ②, rendering the Lindblad prediction (37) invalid.This difficulty, instead, can be overcome by our Keldysh framework, which takes advantage of the filter functions [6].To perform the Keldysh calculation, we first derive the filter operator for the undriven qubit x =  ,0 σ , meaning that x=0 is the only non-vanishing filter operator.Such decomposition enables the analytical evaluation of both Eqs. ( 19) and (13), which predicts the 4If only white noise is present, however, the total decoherence errors in the driven and undriven cases equal.See discussion around Eq. ( 24).
The log-quadratic scaling of the self-energy implies a wellknown sub-Gaussian dephasing profile [34,35], which differs from the exponential one predicted by the Lindblad map (37).
As the secular and full-wave maps agree well for the undriven qubit, we next check whether that agreement extends to the driven case.We focus on a qubit undergoing a spin echo.This protocol uses a -pulse to help refocus the phase of the qubit and, as a result, mitigate the qubit dephasing due to 1/  noise.
For this simulation, we choose the drive such that a finitewidth  pulse is applied at the middle of the whole echo  20) and Eq. ( 11), respectively.The sketch of the echo pulse used for this simulation is plotted by the gray dashed curve.The two insets show the filter strengths   as functions of the filter frequencies   , at the middle and end of the echo duration.Along with the filter strengths, the 1/  noise spectrum used for simulation is also sketched.
duration  (gray curve in Fig. 4).Due to the application of the pulse, the filter strengths for different  ∈ [0, ] differ characteristically.For example, we calculate x for  = /2 and  = , and show the corresponding filter strength   in the insets of Fig. 4. For  = , the filter strength  =0 vanishes, and the most prominent peaks of   are located at  = ±1.By contrast, for  = /2, the only prominent filter strength is  =0 , indicating strong sensitivity to noise from  ≈ 0. Using the filter operators and the 1/  noise spectrum, we further calculate the decoherence maps predicted by both Eq. ( 11) and Eq. ( 20), with  varied over  ∈ [0, ].Note that since this framework is not based on a differential equation, each ρ () with  ∈ [0, ] is calculated separately rather than recursively.Then, we evaluate the off-diagonal matrix element ρ (), which is plotted in Fig. 4. The magnitude of this matrix element indicates the phase coherence of the qubit [36], if the qubit is initially prepared in an equal superposition state.
Visibly, the two calculations show qualitative agreement, while a small deviation exists as the consequence of neglecting the off-diagonal filter functions in Eq. (19).The comparison suggests that, even for a non-trivial decomposition x and a highly structured noise spectrum, the secular approximation can still be applicable5.
Besides the agreement between the two sets of simulations, we also observe the interesting rebound of the matrix element | ρ ()|.Specifically, | ρ ()| decreases during the first half of the echo period and then increases after the pulse is applied.[For the first half period, the evolution of 5For numerical simplicity, we can always choose to truncate the expansion in Eq. ( 13) by only keeping terms with |  −  ′ | smaller than a certain integer, since those with a larger |  −  ′ | are less important due to their diminishing amplitude according to Eq. (55).But one should be alerted that such calculation has the risk of generating a non-CPTP map.ρ () indeed shows the sub-Gaussian dephasing behavior predicted by Eq. ( 45).]Such "inverse dephasing" of the qubit implies a negative decoherence rate, which is described in more detail by a time-local master equation introduced in Ref. [4].This behavior also sends another useful message: although the map () ≈ exp[ (2)  CP ()] is guaranteed to be CPTP, the intermediate map () [( ′ )] −1 (0 <  ′ < ) is not necessarily so.Noticing this, one may naturally ask whether it is possible to follow the procedure in Sec. 2 to derive a CPTP map from  =  ′ to  (0 <  ′ < ), if the system and bath are initialized at  = 0; however, we point out that the basis for such derivation may not hold -for  =  ′ > 0, the two subsystems may already be entangled, while the derivation starting from Eq. ( 2) requires isolation between them.(See discussion of the relation between initial entanglement and the CPTP character of the map in Ref. [37].)

Floquet qubits
We designate this final subsection to test the Keldysh method in studying the Floquet qubit [5,24,[38][39][40].This type of qubit uses the Floquet states of a periodically driven system to store and manipulate quantum information, which can offer advantages such as increased coherence times and more convenient gate operations than the static qubits.Although the open-system Floquet theory is developed to calculate the decoherence rates in such systems, its applicability is limited to the idle Floquet qubit.In the following, we show that the Keldysh method not only reproduces some results by such theory, but also explores situations that are beyond its application.
We start by studying an idle Floquet qubit, where the drive Ĥ () in Eq. ( 1) is periodic, i.e., Ĥ ( +   ) = Ĥ () (  = 2/  is the drive period).In this case, the closedsystem propagator can be expressed as where |  ()⟩ and   are the th independent Floquet state and its corresponding quasi-energy.In the interaction picture, the coupling operator is transformed as where the set F = {  −   ′ +   | 0 < ,  ≤   ,  ∈ Z} contains all possible quasi-energy differences; the operators x(  ) are damping operators in the basis of Floquet states (time-independent in the interaction picture), rather than the eigenstates of the undriven qubit discussed in Sec.

3.1.
With these preparations, we next show how our Keldysh method reproduces the prediction of the decoherence process via the Floquet theory in Ref. [5].In that instance, a Floquet qubit is coupled to both the 1/  flux-noise bath and dielectric noise bath.The spectrum of the former has a strong peak at  ≈ 0, and that of the latter has a smoother spectrum.
Similar to the steps in that reference, we first disregard the peak at  = 0 in the spectrum, and then correct the resulting dynamical map with a more careful consideration of the peak.For the first step, both conditions ① and ② are satisfied if we assume a sufficiently large evolution time .If so, the map takes on the same form as Eq. ( 33), while the operators x(  ) are updated according to Eq. (47).Such map reproduces that generated by the Markovian Floquet master equation [29].Then, to address the strong peak at  = 0 due to the 1/  spectrum, we carefully evaluate the coefficient Re{2 0,0 } = ∫ ∞ −∞ (/2)2  ()  (), which is found to be approximately 2A 2    2 ln |2 ir |.This coefficient should replace   (  = 0) multiplying the damping term D[ x(  = 0)] in Eq. (33).After these two steps, the decoherence map derived in Ref. [5] is reproduced exactly.
The discussion above focuses on an idle Floquet qubit.For gate operations and readout on a Floquet qubit, nonperiodic control is required, and the Floquet theory is no longer applicable.Remarkably, the Keldysh method is still useful in predicting the decoherence map, since the knowledge of the Floquet states and their quasi-energies are not prerequisites for our numerical calculation of the Fourier expansion of x().Below, we show one such example, where a Floquet qubit undergoes an adiabatic evolution from a dynamical sweet spot to an unprotected static working point.(This previously has been used for readout of Floquet qubits [24,39,41].) For a concrete simulation, we reuse the qubit model and parameters in Ref. [5] [Fig.3 (b) and (c) of that reference].Specifically, we consider a fluxonium qubit with its external flux   biased slightly away from the half-flux-quantum sweet spot.Under a periodic drive that is carefully tuned, the qubit can be operated at a so-called dynamical sweet spot, where the derivative  01 /  vanishes.This leads to the first-order insensitivity of the qubit to the 1/  flux noise.(See more details in the caption of Fig. 5 in the current paper.) As references, we first calculate the filter operators x for the dynamical sweet spot (initial) and the static point (final), and plot the resulting   in (a) and (b), respectively.At the dynamical sweet spot, the qubit is to the first order insensitive to 1/  noise, as shown by the vanishing filter strength  =0 with the corresponding filter frequency   = 0.In turn, it is sensitive to noise at other Floquet transition frequencies contained in F [the locations of these frequencies are pointed to by pink arrows in Fig. 5  (a)].By contrast, the qubit at the static working point is strongly sensitive to 1/  noise, as indicated by the large filter strength  =0 .In addition, the qubit is also sensitive to noise at the qubit frequencies ±  .These frequencies are marked by the blue arrows in Fig. 5 (b).We note that the plots for   in both (a) and (b) are reminiscent of the plot for the filter weights in Ref. [5] [Fig.3 (b) and (c) of that reference].In fact, the resulting dynamical maps re- where the parameters are chosen as   /Δ = 1.17 and /Δ = 1.37.The static qubit frequency is   = √ Δ 2 +  2 .The evolution time is set as  = 20 • 2/  for all three plots.For (a), we set /Δ = 2.27.The filter strength  =0 vanishes, which corresponds to a dynamical sweet spot.For (b), the qubit is undriven ( = 0).In this case, the filter strength  =0 is predominant, implying strong sensitivity to 1/  noise.For (c), the drive used in (a) is continuously switched off according to a hyperbolic envelope F () = [1 + tanh[( mid − )/ ramp ]]/2, as shown in the inset ( mid and  ramp are used to tune the starting time and duration of the switch-off, respectively).The resulting plot of   differs from those in both (a) and (b).In (a)-(c), the noise spectra assumed in Ref. [5] is also plotted in gray.
produce those obtained in Ref. [5], as long as the evolution time  is taken to be sufficiently large.
Finally, we calculate   for the adiabatic process connecting the two working points, and show the results in (c).For this case, the plot of   differs from those from both (a) and (b).Interestingly, the locations of the peaks in (c) overlap with those from both (a) and (b).[We mark the peak locations with arrows with different coloring to indicate their apparent origin.]This feature suggests that, the qubit undergoing the adiabatic evolution is subjected to a combination of decoherence channels from both the Floquet and static regimes.Such a feature cannot be predicted by the open-system Floquet theory.
The examples used in this and the previous subsections are limited to qubits with only two energy levels, while the theory developed in Sec. 2 applies to general quantum systems.In Appendix E, we apply this framework to derive the dynamical map for an arbitrarily driven harmonic oscillator as an example.

Quantum optimal control
Above, we have introduced the secular CPTP map and applied it in studying a variety of driven systems.Besides predicting the decoherence maps, the framework can also be used to design drive pulses that mitigate decoherence errors, once it is combined with the technique of quantum optimal control [7,17,18,[42][43][44][45][46][47][48][49].
Optimal-control techniques utilize computer-aided optimization of pulses to minimize state-transfer or gate infidelities.For open-system optimization, the widely used decoherence model is the Lindblad master equation [43,44].Besides this model, Refs.[7,17,18,[46][47][48][49] have used the filter-function approach to reduce the error in quantum operations induced by correlated classical noise.
In the following, we implement our Keldysh decoherence model and extend the filter-function-aided optimization to the mitigation of quantum noise.The optimization with our method is also relatively simpler.In fact, once the total evolution time  is set, the integrals for evaluating the coefficients  , can be precalculated, which are independent of the form of the drives.The CPTP maps (20) also avoid the risk of optimization over unphysical maps.In the following, we use two examples to showcase the capability of the Keldysh-assisted quantum optimal control.

State transfer in the presence of Ohmic noise
We first study the state transfer in a qubit coupled to a typical quantum-noise bath, an Ohmic noise.The spectrum of such noise is   () = A   Θ() [assuming zero temperature; see inset of Fig. 6 (a) for the spectrum], where A  denotes the noise strength and Θ() is the Heaviside function.The qubit coupled to this bath is described by the Hamiltonian Ĥ () =   σ /2 +  () σ , where  () denotes the drive field.The coupling operator is set as x = σ .
In this setup, the spectrum exhibits clear asymmetry between the amplitudes of noise at positive and negative frequencies, which implies distinctive excitation and decay rates in the idle qubit [see Eq. (36)].Because of this, if the decay rate overwhelms the excitation rate, one can leverage the natural system-bath interaction to realize the |⟩ → |⟩ transfer; for the reverse transfer, however, such decay should instead be carefully mitigated.We note that the usual filter-function formalism is not designed to resolve such asymmetry in the noise spectrum, which is a difficulty we can overcome using our method.
For our example, we consider the more difficult |⟩ → |⟩ transfer.To mitigate the error caused by the energy decay, we program the optimizer to minimize the infidelity Above, Û () is the closed-system propagator, () is the CPTP map (20), and the two density matrices are ρ() ≡ |()⟩⟨()|.The implemented optimization algorithm is Gradient Ascent Pulse Engineering, which is commonly used in quantum optimal control [23,50].
For comparison, we first optimize the pulse assuming that the noise is absent.In this case, the optimizer chooses a pulse reminiscent of a typical resonant Rabi drive [red curve in Fig. 6 (a)], whose amplitude is almost constant.(The step-like structure of the pulse is a result of our requirement of the piecewise-constant drive.)For a closedsystem simulation, the pulse induces a smooth increase of the excited-state population, resulting in a negligible statetransfer error (< 10 −6 ) at the end of the pulse.However, for the open-system simulation including the noise bath, the calculation by Eq. ( 20) predicts a much higher error ( st = 8.8 × 10 −2 ), which is caused by the interaction between the qubit and the Ohmic noise bath.Especially, the smooth increase in population renders the qubit prone to the energy loss for the whole state-transfer duration [red dashed curve in (b)].
We next optimize the pulse with the noise included in Eq. (48).The optimized pulse is shown by the blue curve in Fig. 6 (a).Different from the closed-system version, the amplitude of the open-system optimized drive is held close to zero until the latter half of the duration, where the amplitude is ramped up rapidly.In this way, the qubit stays in the excited state for a shorter time than in the previous version.Such behavior of the excited-state population reduces the decoherence error [see the comparison between blue and red curves in Fig. 6 (b)], yielding a 4.4× reduction in the state-transfer infidelity ( st = 2.0 × 10 −2 ).
We note that a similar result is obtained by Ref. [44], where the optimization is based on quantum trajectories.For comparison, that optimization presumed the knowledge of the damping rates and operators, which is derived for an idle qubit using the Lindblad equation rather than a driven system [see discussion in Sec.3.2].From this aspect, our optimization method tends to be more accurate, since it avoids the potential inaccuracy in the damping rates and operators.

Avoiding two-level-system losses in gate operations
Besides state transfer, the Keldysh-assisted optimizer can also help improve gate fidelities.Especially, if the fidelities of certain intuitive gates are limited by one or several resonance peaks in the noise spectrum, our optimizer can offer solutions that reduce the system sensitivity to noise associated with those peaks.
To demonstrate this, we consider a noise bath that consists of one Ohmic bath and a few two-level systems (TLSs) [51,52] [see gray curve for the spectrum of the bath in  The spectrum of the Ohmic noise is given in the inset.In (b), we simulate the evolution of the excited-state populations during the |⟩ → |⟩ transfer using the two pulses, respectively.The solid curves show the open-system evolution of the population, and the dashed one shows the closed-system evolution.For this simulation, we choose the amplitude of the Ohmic noise as A  = 0.001.Fig. 7 (b)].These discrete-level defects have been widely believed to limit the coherence times of many solid-state qubits [52].Therefore, their mitigation is currently an indispensable task.
For the concrete simulation, we choose the same Hamiltonian and coupling operator as in the previous case.We consider a situation where the qubit has a fixed frequency but is accidentally in close resonance with the TLSs.For this setup, the gate operations enabled by idling or weakly driving the qubit should suffer significantly from the TLS loss, because the locations of the peaks of   overlap with those of the resonance peaks (see Fig. 3).This situation motivates us to explore pulses that can mitigate the TLS loss.
Toward this goal, we use the Keldysh-assisted optimalcontrol technique to optimize the pulse  (), with the cost function set as the gate infidelity given in Eq. ( 21).In the following, we focus on optimizing the identity gate as an example, which is the key for the quantum-information storage and multi-qubit operations [53].For the simple identity gates enabled by idling the qubit, the interaction between the qubit and the TLSs limits the fidelity of such operation to  gate = 7.0 × 10 −2 over the duration  = 40 • 2/  [see the overlap of the filter-strength peak of the idle qubit (red) and noise resonance peaks in Fig. 7 (b)].The Keldysh-assisted optimizer, by contrast, proposes to Figure 7: Optimization of the identity gates for a qubit coupled to a few near-resonant TLSs.(a) shows the pulse obtained by the Keldysh-assisted optimal control (blue) and the resulting evolution of the excited-state population (purple) with the qubit initialized in the excited state.(b) plots the filter strengths for the free-induced identity gate (red) and Keldysh-optimized (blue) version.The noise spectrum used for the simulation is also plotted in gray, which is the sum of an Ohmic noise spectrum (A  = 3 × 10 −4 ) and those of three transversely coupled TLSs.(See the form of the spectral function of TLSs in Ref. [51].)The frequencies of the TLSs are chosen to be close to that of the qubit (0.980, 0.995, and 1.020   ), and their relaxation times are chosen as 10 • 2/  .(c) compares the fidelities of the two identity operations after they are each repeated multiple times.
drive the qubit strongly [ ()/  ∼ 0.25] according to the solid curve shown in Fig. 7 (a).As a result of the application of this drive, the populations in the two qubit states oscillate over the whole gate period, and return to the original values at the end of the pulse [see dashed curve in (a) for the excited-state population].These oscillations lead to the appearance of multiple peaks in   located away from qubit frequency   [blue plot in (b)], while the values of   at the TLS resonance frequencies are suppressed.As a result, the error in the identity operation is reduced by 3.2× to  gt = 2.2 × 10 −2 .
For a clear visual comparison between the two schemes, in (c) we show the identity fidelities after multiple such operations are applied.One can observe that the optimized gate yields much higher fidelity for such repeated applica-tions, which corresponds to a longer effective coherence time in the qubit.We also perform optimization for  and phase gates, and find improvement of a similar magnitude.

Conclusion and outlook
In conclusion, we introduce a decoherence model for evaluating errors in a noisy driven system subjected to correlated quantum noise.The second-order Keldysh expansion and the secular approximation lead to a simple CPTP map (20) for the system density matrix.Using this map, we study decoherence errors in a variety of quantum systems with both periodic and non-periodic drives.The clear physical picture of the noise sensitivity described by   provides useful information for developing noise-mitigation strategies, especially if the noise spectrum is only qualitatively understood but cannot be accurately measured.The simplicity of the map after the secular approximation makes this decoherence model suitable to be integrated with the quantum-optimal-control technique.Using the examples of both state-transfer and single-qubit gate, we show that the combination can help mitigate non-classical and correlated noise in state transfers and gate operations.
In the future, one may consider using the technique developed in Ref. [54] for calculating Û () for an even more numerically efficient optimization, since that technique is also based on a Dyson series (also the basis for our Keldysh calculation).For capturing higher-order decoherence effects, it is also useful to explore a higher-order CPTP map [27].Finally, our analytically simple map (20) can provide hints for studying decoherence processes for more complicated systems, e.g., nonlinear oscillators [55] and composite system in the ultrastrong-coupling regime [56].
B Connection to Ref. [17] We use this appendix to connect our theory to previous filter-function research.Particularly, we show that if only classical noise is present, the map Eq. ( 11) can reproduce some of the formulas used in Ref. [17].

C Bath correlation time and spectral variation
In this appendix, we investigate the relation between the variation of   () and the bath correlation time.Such a relation is useful for interpreting condition ② as a comparison between the bath correlation time and the system evolution time.The variation of   () can be roughly quantified by the second-order derivative of the spectrum.Specifically, we define the spectral roughness () by To relate this quantity to the correlation time, we insert the inverse Fourier transform of   () into the expression () and express it as The right-hand side of Eq. ( 59) appears to be related to the bath correlation time.To understand this expression more clearly, we consider a two-level-system defect as an example.The time-domain correlation function of such a two-level system (at zero temperature) has the form [51]  (, 0) =  (0, 0) −    −/  , where   and   are its resonance frequency and coherence time (usually considered as its correlation time), respectively.For this function, the right-hand side of Eq. ( 59 which indicates that more significant variations in   () correspond to longer bath correlation times.

D Floquet master equation for a weakly driven qubit
In this appendix, we use the Markovian Floquet master equation to explain the appearance of the side peaks in Fig.

(b).
The two Floquet states for the weakly driven qubit are and their quasi-energies are  ± = ±/2 −   /2.Inserting them into Eq.( 47), we can again find the expression of the rotated coupling operator Eq. (44).Following the derivation of the Floquet master equation, we extract the transition frequencies   and their corresponding damping operators x(  ) from Eq. (44) The damping terms present in this map indeed capture the noise channels predicted in Fig. 3 (b).

E Driven harmonic oscillator
In the main text, the examples we present are limited to qubits with only two levels.Here, we demonstrate the applicability of our framework in a quantum harmonic oscillator, which has infinite levels.
The system Hamiltonian of this oscillator is specified by and the coupling operator for the oscillator is x = â + â † .For this linear system, the closed-system propagator can be analytically derived as [57] Û () =  ( ) â † −  * ( ) â  −   â † â  −Φ( ) , where the displacement  is calculated by  () =   () +  (), and the additional phase acquired is given by In the interaction picture, the coupling operator is transformed as x() = [ â −    + ()] + H.c., which leads to only two frequency components, x(  ) = â +  and x(−  ) =  † + ᾱ * .We note that the c-numbers (),  * () in x(±  ) only contribute to the Lamb-shift terms to the map (20), which we choose to omit due to the weak noise strength.
Then, if we again assume that the two conditions ① and ② in Sec.3.1 hold, the self-energy is given by  (2)  CP () =    (  )D identical to the prediction by the Lindblad master equation.The message sent by such analysis is that, the Lindblad map is a good approximation for the harmonic oscillator under arbitrary linear drives, as long as   and   ∼ 2/  are much smaller compared to  ∼   .Such a conclusion is in clear contrast to those in Sec.3.2 and 3.4 for the driven qubits.Note that this conclusion may be invalid if nonlinearity in the cavity is induced by its coupling to qubits [55].This conclusion may also be invalid if the drive also affects the noise bath [58], which results in a varying   () during the drive time.

Figure 2 :
Figure 2: A cartoon comparing the decoherence channels of an undriven and driven qubit.The qubit Hamiltonian is given by Ĥ () =   σ /2 + Ĥ (), where Ĥ () describes an arbitrary drive on the qubit.This qubit is coupled to the bath via the operator x =   σ +   σ , and the bath spectrum is denoted by   ().Panels (a) and (b) illustrate the qubit decoherence channels if the qubit is undriven or driven, respectively.In (a), the three terms in the decomposition x() = ±   σ± exp(±  ) +   σ result in the damping operators D[ σ± ] and D[ σ ], which describe the excitation, decay, and pure-dephasing channels of the qubit, respectively[29].In (b), because an arbitrary drive is applied on the qubit, the previous simple decomposition of x() no longer holds; instead, one can use Eq.(12) to perform the frequency decomposition over the time period  ∈ [0, ].As a result, the decoherence channels for a driven qubit are given by the first line of the self-energy(19).

Figure 3 :
Figure 3: Filter strength   for a driven qubit described by Hamiltonian Eq. (41).In (a)-(c), the horizontal axis shows the filter frequency   =   , and the vertical shows the filter strength   [Eq.(23)].The widths of the columns are given by the fundamental frequency   = 2/.The parameters are chosen as follows.The drive amplitudes for all three simulations are chosen as /  = 0.02, and the frequency   used for each plot is given in each figure.The duration is set as  = 200•2/  for all three simulations.For (a) and (b), we choose sinusoidal drives with a constant amplitude; the results of   for the driven qubit are shown in blue, while those for an undriven qubit are shown in red for reference.For (c), the sinusoidal drive used for (b) is multiplied by a hyperbolic envelope F () = [1 + tanh[( −  mid1 )/ ramp ]] [1 + tanh[( mid2 − )/ ramp ]]/4, as shown in the inset.( mid1,2 and  ramp are used to tune the location and duration of the ramps.)

Figure 4 :
Figure 4: Evolution of ρ () for a qubit during a spin echo experiment.The blue solid and red dashed curves correspond to the predictions by Eq. (20) and Eq.(11), respectively.The sketch of the echo pulse used for this simulation is plotted by the gray dashed curve.The two insets show the filter strengths   as functions of the filter frequencies   , at the middle and end of the echo duration.Along with the filter strengths, the 1/  noise spectrum used for simulation is also sketched.

Figure 6 :
Figure 6: Optimization of the state-transfer fidelity for a qubit coupled to an Ohmic noise bath.(a) shows pulses from both the closed-system (red) and open-system (blue) optimizations.The spectrum of the Ohmic noise is given in the inset.In (b), we simulate the evolution of the excited-state populations during the |⟩ → |⟩ transfer using the two pulses, respectively.The solid curves show the open-system evolution of the population, and the dashed one shows the closed-system evolution.For this simulation, we choose the amplitude of the Ohmic noise as A  = 0.001.