Composite pulses for high fidelity population transfer in three-level systems

In this work, we propose a composite pulses scheme by modulating phases to achieve high fidelity population transfer in three-level systems. To circumvent the obstacle that not enough variables are exploited to eliminate the systematic errors in the transition probability, we put forward a cost function to find the optimal value. The cost function is independently constructed either in ensuring an accurate population of the target state, or in suppressing the population of the leakage state, or both of them. The results demonstrate that population transfer is implemented with high fidelity even when existing the deviations in the coupling coefficients. Furthermore, our composite pulses scheme can be extensible to arbitrarily long pulse sequences. As an example, we employ the composite pulses sequence for achieving the three-atom singlet state in an atom-cavity system with ultrahigh fidelity. The final singlet state shows robustness against deviations and is not seriously affected by waveform distortions. Also, the singlet state maintains a high fidelity under the decoherence environment.


Introduction
Implementation of high fidelity quantum coherent control is a top priority in quantum information processing (QIP) [1][2][3][4][5][6][7][8][9][10]. Many works [11][12][13], which design different kinds of pulse shapes, have been devoted to ensuring a remarkable quantum computing performance. Originally, the resonant pulse (RP) technique, where the frequency of radiation exactly matches the transition frequency, is regarded as a popular tool to accurately achieve coherent control, due to its fast operation and simple waveform [14]. For example, one can achieve complete population inversion through a π-pulse, or a maximum superposition state through a π/2-pulse [15]. However, the RP is extremely susceptible to external perturbations, such as fluctuations of control fields. To overcome this difficulty, the adiabatic passage (AP) technique has been developed [16,17]. The AP technique selects one eigenstate of the system Hamiltonian as an evolution path and propels the initial state to adiabatically evolve along this path. In this technique, the robustness against parameter perturbations is improved at the cost of the slow evolution rate and lower fidelity. To enjoy both the ultrahigh fidelity of RP and the robustness of AP, one can adopt the composite pulses (CPs) technique.
The CPs technique, comprised of a well-organized train of constant pulses with determined pulse areas and relative phases, was conceived in early nuclear magnetic resonance (NMR) [18][19][20][21]. To date, this technique has been used extensively in QIP [22][23][24][25][26][27][28][29][30][31]. One unique feature of CPs is that the pulse sequence is available to compensate for the systematic error in any physical parameters (e.g. pulse duration, pulse Figure 1. Infidelity F r (top panels) and population leakage P a e (bottom panels) vs the deviations 1 and 2 in the RPs, where the infidelity F r = |P a r − P r | quantifies the population drift in the desired population P r . Here, P r = 1/2 in the left column, while P r = 3/4 in the right column. Henceforth, the white-solid curves, the green-dashed curves, and the cyan-dotted curves correspond to F r (P a e ) = 0.001, 0.01, and 0.05, respectively. The rectangular region (labeled by the yellow-dot-dashed lines) are: | 1 | 0.2 and | 2 | 0.2. The results demonstrate that the population transfer suffers from the deviations in the RPs. Infidelity is particularly noticeable, since the region enclosed by white-solid curves is almost negligible.

Toy model and the general theory
Consider a Λ-type three-level quantum system interacting with two external fields, where the three-level system has two ground states |g and |r , and an excited state |e . The ground states |g and |r , acting as a qubit, cannot be immediately coupled to each other. Hence, we require an extra excited state |e to construct the indirect coupling between two ground states. More specifically, the transition |g ↔ |e (|r ↔ |e ) is driven by the external field with the coupling strength Ω 1 (Ω 2 ), the phase α (β), and the detuning Δ. In presence of deviations in the external fields, in the interaction picture, the system Hamiltonian is given by ( = 1) H = Δ| e e| + (1 + 1 )Ω 1 e iα |g e| + (1 + 2 )Ω 2 e iβ |r e| + h.c., ( 1 ) where the deviations 1 and 2 are random unknown constants, which would give rise to the systematic errors during quantum operations.
Here, we do not intend to specify a concrete physical system, because the three-level system studied here can be found in very different physics research fields ranging from high precision spectroscopy, to QIP, NMR, and metrology. For example, this physical model is quite familiar in the quantum system of a three-level atom interacting with two laser fields [14]. In the atomic system, the inhomogeneous distribution of the laser fields leads to different interaction coefficients between the atom and the laser fields. Also, the spatial position of the atom cannot be exactly determined due to its micro-vibration. As a result, the interaction coefficients in the system are susceptible to deviations. It is worth mentioning that similar energy-level structures can also be found in trapped ions [59], diamond nitrogen-vacancy centers [60], and superconducting circuits [61][62][63][64][65][66]. In addition, some complicated quantum systems, such as the electrons in semiconductors [67], the spin chain [68,69], neutral atoms interactions [70][71][72][73][74][75], and massive quantum particles in an optical Lieb lattice [76], can be reduced to the three-level physical model as well.
Assume that the system is initially in the ground state |g . If there are no deviations in the external fields (i.e. 1 = 0 and 2 = 0), one can easily achieve perfect population transfer for the qubit. To this end, we choose the evolution time T = 2π/ Δ 2 + 4(Ω 2 1 + Ω 2 2 ), then the population of the ground state |r becomes P r = 4Ω 2 1 Ω 2 2 cos 2 ΔT 4 (Ω 2 1 + Ω 2 2 ) 2 . (2) Thus, the value of P r can be modulated by altering either the detuning or the ratio of two coupling coefficients. However, in presence of deviations in the external fields, the actual population P a r would keep away from the desired value P r . This is verified in figures 1(a) and (c), which demonstrate that the infidelity F r = |P a r − P r | is high even for small deviations. Note that figures 1(a) and (c) also demonstrate that the value of infidelity becomes extremely low in some regions where the deviations are large, e.g. F r ≈ 8 × 10 −5 when 1 ≈ −0.46 and 2 ≈ 0.42 in figure 1(c). However, it is meaningless in practice, because the values of the deviations 1 and 2 are unknown in principle and they are usually independent of each other. On the other hand, population leakage from the ground states to the excited state invariably does happen in the presence of deviations, which can be observed in figures 1(b) and (d). Hence, the goal of this work is to precisely achieve the predefined population P r , and meanwhile, to dramatically suppress the population leakage by designing a composite N-pulse sequence.
First, we require to calculate the population of the ground state |r in absence of the deviations. Assume that the Hamiltonian of the nth pulse has the following form (n = 1, . . . , N) H n = Δ n |e e| + (1 + 1 )Ω 1n e iα n |g e| + (1 + 2 )Ω 2n e iβ n |r e| + h.c. (3) For the pulse duration T n = 2π/ Δ 2 n + 4(Ω 2 1n + Ω 2 2n ), the corresponding propagator, labelled as U(θ n , γ n ), can be written as (up to a global phase) Then, the total propagator of the composite N-pulse sequence can be written as The detailed derivation of the total propagator is given in appendix A. After some algebraic manipulations, the population of the ground state |r , labelled as P (0) r , becomes where the expression of ϑ N is also presented in appendix A. For simplicity, we set Δ n = 0 in the following, which means that the three-level system works in the resonant regime. Apparently, the following analysis can be easily generalized to the two-photon resonant regime. Furthermore, we choose θ n = π/4 (n = 1, . . . , N), which means that two coupling coefficients Ω 1 and Ω 2 remain unchanged in the N-pulse sequence. Note that it is also suitable for the other values of θ n . As a result, we only modulate the phases α n and β n , and the form of CPs now becomes In presence of deviations in the external fields, by the Taylor expansion, the actual population P a r of the ground state |r can be written as where C (l) N,k (k = 1, 2, . . .) are the coefficients of the lth-order term in the N-pulse sequence. To achieve the predefined population P r , one should design the phases α n and β n (n = 1, . . . , N) to guarantee the validity of the following equation That is, the zeroth-order term of the actual population P a r should be equal to the predefined population. This is the first condition that the phases α n and β n must satisfy. Furthermore, the global phase of the propagator is inessential for the system dynamics, and what really matters is the phase difference α mn = α m − α n and β mn = β m − β n (m, n = 1, . . . , N). Hence, one can set the phases α 1 and β 1 of the first pulse to be arbitrary. Nevertheless, they become extremely useful when the CPs are designed for single qubit gates, because the values of α 1 and β 1 can be used to modulate the relative phase between the ground states.
As a result, apart from equation (6), there still remain (2N − 3) free phases (variables) in the N-pulse sequence, which can be designed to compensate for the systematic errors caused by the deviations 1 and 2 .
Unlike the case of two-level systems [30][31][32][33][34][35][36][37][38], there exists a population leakage from the ground states to the excited state in the three-level system. Thus, we need to suppress population leakage (i.e. the population of the excited state |e ) as well. Similar to the treatment of P a r , the actual population P a e of the excited state |e can be expanded by the Taylor expansion, as follows where D (l) N,k (k = 1, 2, . . .) are the coefficients of the lth-order term in the N-pulse sequence. Contrary to equations (5) and (7) does not contain the first-order term.
Therefore, to obtain high fidelity population transfer in the qubit, the general method is to satisfy the following conditions: (i) eliminate the higher-order terms in P a r as much as possible. Namely, (ii) demand that the value of P a e be as small as possible, i.e. D (2) N,1 = D (2) N,2 = . . . = 0. Note that the design procedure of the phases is more complicated in three-level systems than that in two-level systems [30][31][32][33][34][35][36][37][38], since there are two coefficients [C (1) N, 1 and C (1) N,2 ] in the first-order term while there are six coefficients [C (2) N,k and D (2) N,k , k = 1, 2, 3] in the second-order term, and so forth. Particularly, when considering the small number of pulses, which is often the practical case, it does not have sufficient phases to eliminate the systematic errors up to the desired order.
To solve this issue, we can proceed as follows. First, it is worth noting that the influence of the higher-order term of systematic errors on dynamical evolution would gradually diminish. Therefore, the first-order terms have the most serious effect on dynamical evolution, and we should give priority to making these terms vanish by designing suitable phases, i.e. C (1) N,1 = C (1) N,2 = 0. Then, the remaining phases are designed to minimize the following cost function: where A l and B l are the weighting coefficients, satisfying 0 A l−1 A l and 0 B l−1 B l . Physically, the cost function F represents the trade-off between the accuracy of population transfer and the leakage to the excited state |e . The conditions 0 A l−1 A l and 0 B l−1 B l ensure that the coefficients of the low order terms have a high weight, and thus are preferentially eliminated. Remarkably, different weighting coefficients have different effects. Setting B l = 0, we have the cost function F = ∞ l=2 A l F (l) c aimed to keep the accuracy of the population transfer, while setting A l = 0 we have the cost function F = ∞ l=2 B l F (l) d aimed to reduce the leakage to the excited state |e . Note that the minimum value of the cost function is zero, corresponding to satisfy the equations:

Robust implementation of arbitrary population transfer by composite pulses
In this section, we present the design process of the phases α n and β n (n = 2, . . . , N) to achieve high fidelity population transfer by the N-pulse sequence. For the two-pulse sequence, the phases are derived analytically, while for the three-pulse sequence we combine both the analytical method and the numerical method to obtain the phases. For more than three pulses, we carry out the four(five)-pulse sequence by numerical calculations to demonstrate the feasibility of eliminating the higher-order terms of the systematic errors. It is worth mentioning that the cycle of the phases is 2π, thus we restrict the values of all phases in the interval [0, 2π) following from here.

Two pulses
The form of the two-pulse sequence is expressed by , γ 1 , and the zeroth-order coefficient in equation (5) becomes where α 12 = α 1 − α 2 and β 12 = β 1 − β 2 . As a result, the solution of equation (6) can be written as It is worth mentioning that the first-order coefficients C (1) 2,1 and C (1) 2,2 in equation (5) automatically vanish in the two-pulse sequence, i.e. C (1) 2,1 = C (1) 2,2 = 0. Therefore, we can further reduce the detrimental effect of the second-order terms, where the corresponding coefficients C (2) 2,k and D (2) 2,k (k = 1, 2, 3) are presented in appendix B. Note that there are only two controllable variables α 12 and β 12 , and β 12 is required to satisfy equation (12). Only single variable α 12 is left to eliminate the second-order terms of the systematic errors, and thus it is enough to consider the second-order coefficients in the cost function given by equation (8). In the following, we give out the optimal value of the variable α 12 for three types of cost functions (see appendix B for details): (a) The cost function is chosen as F = F (2) c . That is, we aim merely at eliminating the deviation in the actual population P a r . It is not hard to ascertain that F (2) c is minimal when where θ = arcsin 19π 2 / √ 2(16 + 19π 2 ) . (b) The cost function is chosen as F = F (2) d . That is, we aim merely at suppressing the population leakage P a e . Note that F (2) d is minimal when (c) The cost function is chosen as That is, we make an equal weight between the elimination of the deviation in the actual population P a r and the suppression of the population leakage P a e . Note that F = F (2) c + F (2) d is minimal when Here, Θ satisfies the following quartic equation: where the expressions for the coefficients a k (k = 1, . . . , 5) are given in appendix B. Figure 2 demonstrates the infidelity F r and the population leakage P a e as a function of the deviations 1 and 2 by different cost functions. Compared figures 2(a) and (b) with figure 1, both the robust behaviors and the population leakage are slightly improved when we choose the cost function F = F (2) c + F (2) d . An inspection of figures 2(c) and (d) demonstrates that population transfer would be more accurate when the phases are determined by equation (13), since the region enclosed by white-solid curves in figure 2(c) is much larger than that in figure 2(a). Nevertheless, there are plenty of population leakages in this situation. When adopting the phases determined by equation (14), the population leakage is suppressed at a very low level in a wide region, as shown in figure 2(f). These results show that the robust behaviors are quite different when we choose different cost functions.

Three pulses
The form of the three-pulse sequence is expressed by The results demonstrate that the robust behaviors are different by choosing different cost functions. and the zeroth-order coefficient in equation (5) becomes where α mn = α m − α n and β mn = β m − β n , m, n = 1, 2, 3. The first-order coefficients C (1) 3,1 and C (1) 3,2 in equation (5) are For the second-order coefficients, since the expressions are too complicated to present here, we give them in appendix C.
Three-pulse sequence Four-pulse sequence Five-pulse sequence where s = 0 when 0 < θ < π/4, while s = 1 when π/4 θ π/2. Note that there are four variables (α 12 , α 23 , β 12 , and β 23 ) in the three-pulse sequence, and the values of the variables β 12 and β 23 are given by equation (18). Therefore, only two variables α 12 and α 23 can be designed to eliminate the influence of the deviations. In the following, the cost function is chosen as d , that is, the impact of the accuracy of population transfer and the leakage to the excited state are equally weighted in the three-pulse sequence. Figure 3 shows the performance of the infidelity F r and the population leakage P a e in the three-pulse sequence. It can be seen that both the infidelity F r and the population leakage P a e keep a very small value in the presence of strong deviations. Compared with the case of the RPs (cf figure 1), the region F r (P a e ) 0.001 becomes much wider in the three-pulse sequence, and thus this sequence is more robust against the deviations 1 and 2 . Noting that it is hard to obtain the analytical expressions for the phases α n and β n in the three-pulse sequence, we present some numerical solutions for different populations (i.e. different θ) in table 1.

More than three pulses
The form of the four-pulse sequence reads , γ 1 , and the zeroth-order coefficient in equation (5) becomes where α mn = α m − α n and β mn = β m − β n , m, n = 1, 2, 3, 4. The first-order coefficients C (1) 4,1 and C (1) 4,2 in equation (5) are expressed by Thus, one solution of equations can be written as As a result, there are still four variables available to eliminate the influence of the deviations in the four-pulse sequence, and some numerical solutions for different θ are presented in table 1.
Similarly, the zeroth-order in the five-pulse sequence can be expressed by  Figure 4 shows the infidelity and the population leakage as a function of the deviations in the four(five)-pulse sequence. One can find from figure 4 that the infidelity and the population leakage are small in a wider region as the number of pulses increase. Apparently, the population transfer would be more accurate when more pulses are taken into account. Note that the cost function is chosen as figure 4, namely, the impact of the accuracy of the population transfer and the leakage to the excited state are equally weighted. We present in appendix E the detailed discussion on the performance of the accuracy and the leakage by choosing different forms of cost functions. Next, we demonstrate the evolutionary trajectory of the qubit (i.e. the ground states |g and |r ) on the Bloch sphere under the N-pulse sequence. Taking the four-pulse sequence as an example, suppose the initial state of this system is |g . Then, the final state of this system after the four-pulse sequence can be expressed in the basis {|g , |r , |e } by employing four angle variables: Since we focus on the evolutionary trajectory of the qubit, the phase ϕ is irrelevant and we can simply write the time evolution of qubit as: |Ψ = cos ϑ cos θ|g + sin θ exp(iφ )|r . When |cos ϑ | = 1, the evolutionary trajectory is on the Bloch sphere. While it is inside the Bloch sphere when |cos ϑ | = 1. Here, the state inside the Bloch sphere does not represent the mixed state, and means the population leakage from ground states to the excited state. Note that the smaller the value of |cos ϑ | is, the larger the population leakage will be. The center of the Bloch sphere, i.e. cos ϑ = 0, represents the excited state |e .
show the evolutionary trajectories of the qubit for different cost functions in the four-pulse sequence. By contrast, we plot in figure 5(a) the evolutionary trajectory for the RPs. We can observe from figure 5 that the population leakage always exists during the evolution process, since the trajectories are inside the Bloch sphere in some evolution stages. Nevertheless, the population leakage is strongly suppressed and the system state approaches the target state at the final time, as shown in figure 5(d). It is shown in figure 5(b) that the cost function F = F (2) c + F (3) c only guarantees the accuracy of P a r , while there exists leakage to the excited state. As a result, the final state deviates from the target state. A similar situation is also found in figure 5(c), because the cost function F = F (2) d + F (3) d + F (4) d only prevents the population leakage to the excited state. However, figure 5(d) demonstrates that the system is almost driven into the target state, because the cost function F = F (2) c + F (2) d is a tradeoff between the accuracy and the population leakage. Note that the corresponding robust region designed by F = F (2) c + F (2) d (e.g. the region enclosed by the white-curves in figure 4) is generally smaller than those designed by the other two cost functions. For details, one can see appendix E.
The initial state |g is labeled by , the target state 1/ √ 2[|g + exp(iφ )|r ] is labeled by , and the final state of the qubit after the four-pulse sequence (or the RPs) is labeled by . The trajectory inside the Bloch sphere means the population leakage from the qubit to the excited state.

Applications: robust preparation of three-atom singlet state by composite pulses
In above section, we have showed the design process of CPs in a general physical system, where the three-level structure can be found in atomic systems [14], trapped ions [59], diamond nitrogen-vacancy centers [60], and superconducting circuits [61][62][63][64][65][66], etc. In this section, we demonstrate that the CPs scheme can be also generalized to complicated systems for performing different quantum tasks, e.g. coherent conversion between two qubits [77] or preparation of entangled states [78], etc. The fundamental is to reduce complicated systems to the familiar three-level physical model. As an example, in presence of deviations, we next prepare the three-atom singlet state with ultrahigh fidelity in an atom-cavity system, and this approach is easily extended to other complicated systems.
As shown in figure 6(a), consider the atom-cavity system, where three identical four-level atoms are trapped in a bimodal cavity. The four-level atom has three ground states |1 , |2 , and |3 , and an excited state |e , as shown in figure 6(b). The transition |2 k ↔ |e k (|3 k ↔ |e k ) is resonantly coupled by the cavity mode a (b) with the coupling constant λ a k (λ b k ), where the subscript k represents the kth atom. The transition |1 k ↔ |e k is resonantly coupled by the laser field with the coupling strength Ω k and the phase α k . In the interaction picture, the Hamiltonian of the atom-cavity system reads ( = 1) whereâ andb are the annihilation operators of the cavity mode a and b, respectively. k and ζ a k (ζ b k ) are respectively the deviations of the laser fields and the cavity mode a (b) due to the inhomogeneities of spatial distribution.  1, 2, 3). (a) The three-pulse sequence scheme. (b) The RPs scheme shown in reference [78]. The fidelity is defined by F S = | Ψ S |Ψ | 2 , where |Ψ is the system state. The parameters of the three-pulse sequence scheme are λ k = 30Ω 1 , Apparently, the excited number is a conserved quantity in this system since [H,N e ] = 0, where the excited number operator is defined byN e = 3 k=1 (|e kk e| + |1 kk 1|) +â †â +b †b . When the condition Ω k λ a k (λ b k ) is satisfied, we can restrict the system dynamics into the single-excited subspace and the Hamiltonian is approximated by [78] where we set Ω 3 = Ω 2 , α 2 = α 3 = β 1 , 3 = 2 , and The state |lmn represents the first atom in |l , the second atom in |m , and the third atom in |n . Since the cavity mode a (b) is in vacuum state, we have ignored it in the Hamiltonian given by equation (24). It is easily found from equation (24) that the atom-cavity system can be reduced to the three-level physical model studied in section 2. Hence, if the initial state of this system is |Ψ 1 , by fixing the ratio of two coupling coefficients and modulating the phases according to the above CPs theory, one can obtain the following superposition state |Ψ = cos θ|Ψ 1 + e −iγ sin θ|Ψ 2 .
By setting γ = 0 and θ = arctan √ 2, the state |Ψ becomes which is actually the three-atom singlet state. We first explore how the deviations of two coupling constants λ a k and λ b k affect the fidelity F S of the singlet state. The numerical results simulated by the three-pulse sequence scheme are shown in figure 7(a), where we choose ζ a , ζ b ∈ [0, 1] to guarantee the validity of the effective Hamiltonian given by equation (24) (i.e. the condition Ω k λ a k , λ b k is well satisfied). As a comparison, we also plot in figure 7(b) the fidelity F S as a function of the deviations ζ a and ζ a by the RPs scheme shown in reference [78]. The results demonstrate that both schemes are robust against the deviations in two coupling constants, since the impact on the fidelity F S can be negligible no matter how large the deviations are, cf, F S 0.9978 in figure 7. The reason for this phenomenon is that both ζ a and ζ b are absent in the Hamiltonian given by equation (24). Besides, an interesting finding is that the fidelity F S tends to a higher value when the deviations ζ a and ζ b are much larger. A reasonable interpretation is given as follows. During the process of deriving the Hamiltonian given by equation (24), the strong coupling condition Ω k λ a k (λ b k ) is exploited. The coupling constant λ a k (λ b k ) increases with the raising of the deviation ζ a (ζ b ). As a result, the dynamics of this system is more and more restricted in the single-excited subspace so that the Hamiltonian given by equation (24) ensures a better validity for suppressing the leakage to other subspaces. Thus, a higher fidelity of the singlet state can be achieved when a larger deviation ζ a (ζ b ) appears. Next, we address the influence of the two deviations 1 and 2 on the fidelity F S . In figure 8(a), we display the fidelity F S versus the deviations 1 and 2 for the three-pulse sequence. For comparison, we also plot in figure 8(b) the performance of the fidelity for the RPs scheme shown in reference [78]. As shown by the area surrounded by the blue-solid (or pink-dotted) curve in figure 8(b), it is clear that the scheme shown in reference [78] is highly susceptible to the deviations 1 and 2 . However, the fidelity of the singlet state obtained by the three-pulse sequence remains in a high value (F S 0.99) even when a very large deviation occurs. This benefits from the fact that the second-order coefficients of systematic errors are restricted to an extremely low value by the specific set of phases. Furthermore, as the pulse number increases, the higher-order coefficients can be restricted to a low value that ensures a robust behavior against the deviations 1 and 2 . Thus, the CPs sequence can further improve an error-tolerant performance of the singlet state.
Another experimental issue, which inevitably impacts the performance of the final fidelity, is the waveform distortion. In the ideal condition, our input pulse shape is associated with a perfect square waveform. However, in practice, the input square wave always produces a tiny smooth rising and falling edge. Let us study this issue by taking the three-pulse sequence as an example. First of all, we design a group of functions in which the waveform distortion can be arbitrarily adjustable. The functions are given by , 0 t 3T/2, , 0 t 3T/2, where T = π/ Ω 2 1 + Ω 2 2 . Here, χ is a dimensionless parameter that defines the magnitude of the rising and falling edge of the square wave.
As χ increases, the functions (26) and (27) gradually become a square wave. Specifically from figures 9(a) and (b), we see that the functions are almost square waveforms for χ = 1000 (the blue-dotted curves). For χ → ∞, they describe perfect square waves. However, when χ is small, e.g. χ = 10, the duration of the rising and falling edge of the square wave becomes excessively long and the waveform has a serious shape change. Then, functions (26) and (27) turn into a smooth time-dependent pulse, as shown by the red-dashed curves in figures 9(a) and (b). Obviously, this kind of waveform does not satisfy the request of our CPs scheme. As illustrated in figure 9(c), even though the waveform suffers from serious distortion (the red-dashed curve), the fidelity F S can still maintain a high value (>0.99). Moreover, as revealed by the cyan-solid and pink-dashed-dotted curves in figure 9(c), a slight distortion in the waveform hardly makes an impact on the final fidelity. Hence, the CPs sequence has the ability to acquire a precise and robust singlet state even with a distorted waveform. Finally, we investigate the influence of the phase shift errors on the fidelity. In presence of the deviations in the phases, the actual phases can be rewritten as where α and β represent the phase shift errors. We plot in figure 9(d) the relation between the fidelity and the phase shift errors, which shows that the current CPs sequence is sensitive to phase shift errors [79]. The reason is that we do not consider the phase shift errors in the Hamiltonian given by equation (3). To tackle this problem, one requires to regard the phase shift errors as variables in equation (3), and then constructs a phase-error-corrected CPs sequence [80].
In a more realistic situation, one has to consider the influence of decoherence on the fidelity of the singlet state. In this system, the decoherence mainly comes from two aspects: (i) the atomic spontaneous emission from the excited state to the ground states. (ii) The decay of two cavity modes. Under the decoherence environment, the dynamics of this system are governed by the Lindblad master equation and the form can be written asρ whereσ le k = |l kk e| and the Lindblad super-operator is whereô is the standard Lindblad operator. γ le k and κ a (κ b ) are the dissipation rate from the excited state |e to the ground state |l and the decay rate of the cavity mode a (b), respectively. For simplicity, we set γ le k = γ/3 and κ a = κ b = κ. Figure 10 displays the fidelity F S as a function of the dissipation rate γ and the decay rate κ by both the three-pulse sequence scheme and the RPs scheme shown in reference [78]. We observe from figure 10(b) that the fidelity cannot maintain a relatively high value in the RPs scheme even if the system does not suffer from the decoherence, cf, F S < 0.956 in figure 10(b). However, in the three-pulse sequence scheme, the fidelity can reach a high value ( 0.99) when the dissipation rate is small, cf, the left region of the pink-dotted line in figure 10(a). Another intriguing aspect in figure 10 is that the fidelity does not sharply drop in the RPs scheme while the range of fidelity varies widely in the three-pulse sequence scheme. This is because the evolution time in the RPs scheme is much shorter than that in the three-pulse sequence scheme. As a result, the influence of the decoherence on the system is relatively small by the RPs scheme. From this point, we find that the longer pulse sequence is not the optimal choice in a severe decoherence environment, even though the longer pulse sequence can help to improve the robustness of the system. Furthermore, according to the pink-dotted line in figure 10(a), it is seen that the influence of the cavity decay is almost negligible, since F S 0.99 even when the decay rate κ increases to 0.1. This phenomenon can be explained by the fact that the strong coupling constant λ a k (λ b k ) ensures the decoupling of the single-photon state to minimize the loss of the cavity decay. From a physical point of view, the cavity is regarded as intermediary to realize indirect coupling between atoms, and it is almost in the vacuum state under the strong coupling condition.
On the other hand, for a given cavity decay rate, the fidelity F S decreases with the increase of the dissipation rate γ. Thus, the major influence on the preparation of single state is atomic spontaneous emission. The physical mechanism can be clarified as follows. Under the strong coupling condition Ω k λ a k (λ b k ), the dynamics of the system can be approximately described by the Hamiltonian H in equation (24). To obtain the singlet state, we need to utilize the intermediate state |Ψ 3 . Nevertheless, we can see that |Ψ 3 contains the excited states |e k (k = 1, 2, 3), which cannot be effectively eliminated during the evolution process. Hence, the dynamics of the system is sensitive to the spontaneous emission. An alternative approach to solve this problem is to adiabatically eliminate the excited state |e k [81,82], i.e. adopting the Raman transition process. To be specific, the one-photon detuning between the excited state and the ground states is large enough, while the system satisfies the two-photon resonance condition. Then, the population of the intermediate state is sharply suppressed, and the population can adiabatically transfer between ground states during the evolution. As a result, the influence of the spontaneous emission is effectively eliminated.

Conclusion
In this work, we achieved a high fidelity population transfer for the three-level system by the CPs sequence. We derived the transition probability from the total propagator and used the Taylor expansion to disassemble it into a series of derivatives. For a small number of pulses, since not enough variables (phases) can be employed to eliminate anticipated deviations, we constructed a cost function to control the fidelity of the result. Then, we showed that high fidelity of the target state, or suppressing the population of the leakage state, or both two objectives are robustly achievable by choosing different forms of cost functions. Particularly, the method of the short pulse sequence can be further extended to a longer one if necessary.
As an example, we applied the CPs sequence to implement the three-atom singlet state with ultrahigh fidelity in the atom-cavity system, where we only needed to individually control the phase of the inputting pulses. The numerical results indicate that the final singlet state is robust against the deviations in two coupling coefficients. Moreover, we found that significant distortion in the pulse sequence has a very small impact on the fidelity of the singlet state. As the simulated results further reveal, the fidelity of the final singlet state keeps a relatively high value in the decoherence environment. Our CPs scheme provides a selectable way for robust controlling in the error-prone environment and shows an enormous potential application in various physical models, e.g. implementation of coherently convertible dual-type qubits with high fidelity [77].

Data availability statement
All data that support the findings of this study are included within the article (and any supplementary files).

Appendix A. The detailed derivation of the total propagator without deviations in the composite N-pulse sequence
In this appendix, we first calculate the general expression of the propagator for each pulse, and then give out the total propagator of the composite N-pulse sequence.
(A2) Therefore, the total propagator of the composite N-pulse sequence becomes where ϕ N , φ N , and ϑ N satisfy the recursive relations given by equation (A2). The population of the ground state |r reads
As a result, it is not hard to calculate that the cost function F = F (2) c is minimal when where θ = arcsin 19π 2 / √ 2(16 + 19π 2 ) . Second, the expression of the cost function F = F (2) d is Obviously, F (2) d is minimal when Finally, the cost function F = F (2) c + F (2) d is minimal when Here, Θ is the solution of the following quartic equation: where the coefficients a k (k = 1, . . . , 5) are:

Appendix D. The first-order coefficients in the five-pulse sequence
In this appendix, we provide the expressions of the first-order coefficients in the five-pulse sequence, which are . Infidelity F r (top panels) and population leakage P a e (bottom panels) vs the deviations 1 and 2 in the four-pulse sequence. The parameters are designed according to the cost function F = F (2) c + F (2) d , and the values are presented in table 2. The results demonstrate that both the infidelity and the population leakage maintain a small value ( 0.001, cf the white-curves) in a very wide region, but the corresponding regions are not larger than those in figure 11 or figure 12. Table 2. The parameters to achieve the predefined population P r = sin 2 θ in the four-pulse sequence, where α 1 = β 1 = 0.   where α mn = α m − α n and β mn = β m − β n (m, n = 1, 2, 3, 4, 5).

Appendix E. Detailed discussion on the performance of the accuracy and leakage in the four-pulse sequence
In this appendix, taking the four-pulse sequence as an example, we discuss in detail the influence of different cost functions on the performance of the accuracy and the leakage. As stated in section 2, the cost functions are mainly divided into three types: (i) F = ∞ l=2 A l F (l) c is meant to eliminate the deviations in P a r ; (ii) F = ∞ l=2 B l F (l) d is aimed to eliminate the population leakage P a e ; (iii) F = ∞ l=2 A l F (l) c + B l F (l) d is a combination of both, intended to eliminate both the deviations in P a r and the population leakage P a e . Figures 11-13 are the plotted infidelity and population leakage as a function of deviations through the cost functions F = F (2) c + F (3) c , F = F (2) d + F (3) d + F (4) d , and F = F (2) c + F (2) d , respectively. The results demonstrate that P a r in figure 11 is much more accurate than those in figures 12 and 13, while the population leakage P a e in figure 12 is much smaller than those in figures 11 and 13. In figure 13, there is a tradeoff between the accuracy and the leakage. Therefore, different cost functions have different effects, and one can choose a specific cost function according to the particular task in QIP.