Nearly-optimal quantum control: analytical approach

TL;DR

This paper introduces an analytical approach to nearly-optimal quantum control, enabling efficient state switching in quantum systems without complex computations.

Contribution

It presents a new analytical framework for quantum control optimization within simple protocol families, reducing complexity and enhancing physical transparency.

Findings

Designed control protocols for switching between orthogonal states

Achieved near-optimal control with reduced computational complexity

Demonstrated effectiveness on a biased two-level quantum system

Abstract

We propose nearly-optimal control strategies for changing states of a quantum system. We argue that quantum control optimization can be studied analytically within some protocol families that depend on a small set of parameters for optimization. This optimization strategy can be preferred in practice because it is physically transparent and does not lead to combinatorial complexity in multistate problems. For demonstration, we design optimized control protocols that achieve switching between orthogonal states of a naturally biased quantum two-level system.

Figures

Equations

$$i\frac{d\Psi }{dt}=\hat{H}(t)\Psi ,$$

$${\hat{H}}^{\mathrm{free}}=-\epsilon {\hat{\sigma }}_z,$$

$$\hat{H}=-\epsilon {\hat{\sigma }}_z+B_x^{\mathrm{ex}}(t){\hat{\sigma }}_x+B_y^{\mathrm{ex}}(t){\hat{\sigma }}_y+B_z^{\mathrm{ex}}(t){\hat{\sigma }}_z,$$

$$C={\int }_{-\infty }^{0}dt\left\{A⟨\Psi (t)|{\hat{\sigma }}_x|\Psi (t){⟩}^2+(B_x^{\mathrm{ex}}(t){)}^2+(B_y^{\mathrm{ex}}(t){)}^2+(B_z^{\mathrm{ex}}(t){)}^2\right\}.$$

$$C={\int }_{-\infty }^0\mathcal{L}dt,$$

$$⟨\Psi (t)|{\hat{\sigma }}_x|\Psi (t){⟩}^2={\sin }^2(2B_y^{\mathrm{ex}}t).$$

$$T_u=\pi /(2B_y^{\mathrm{ext}}),$$

$$C_{\mathrm{u}}={\int }_{-T_u}^0\mathcal{L}dt=\frac{{\pi }^2}{4T_u}+\left({\epsilon }^2+\frac{A}{2}\right)T_u.$$

$$C_u^{\min }=\pi \sqrt{{\epsilon }^2+\frac{A}{2}}.$$

$$2{\int }_{-T_0}^0{B}_y^{\mathrm{ex}}(t)dt=\pi .$$

$$⟨\Psi (t)|{\hat{\sigma }}_x|\Psi (t){⟩}^2={\sin }^2\left({\int }_{-T_0}^{t}2B_y^{\mathrm{ex}}(t^{\prime })d{t}^{\prime }\right).$$

$$\varphi (t)\equiv -\pi /2+2{\int }_{-T_0}^t{B}_y^{\mathrm{ex}}(t^{\prime })d{t}^{\prime },$$

$$C_{sG}[\varphi (t)]={\int }_{-T_0}^0\left\{\frac{{\dot{\varphi }}^2}{4}+{\epsilon }^2+A{\cos }^2(\varphi )\right\}dt,$$

$$\varphi (-T_0)=-\pi /2,\varphi (0)=\pi /2,$$

$$\frac{1}{4}{\dot{\varphi }}^2-{\epsilon }^2-A{\cos }^2\varphi =e=\mathrm{const}.$$

$$\dot{\varphi }=2\sqrt{{\epsilon }^2+A{\cos }^2\varphi +e}.$$

$$\varphi =\arcsin \left[\mathrm{sn}\left(2\sqrt{{\epsilon }^2+A+e}(t+T_0)|\frac{A}{{\epsilon }^2+A+e}\right)\right],$$

$$T_0(e)=\frac{1}{\sqrt{{\epsilon }^2+A+e}}K\left(\frac{A}{{\epsilon }^2+A+e}\right),$$

$$C_{sG}=-e{T}_0(e)+2\sqrt{{\epsilon }^2+A+e}E\left(\frac{A}{{\epsilon }^2+A+e}\right),$$

$$\varphi (t)=\arcsin \left[\mathrm{sn}\left(2\sqrt{{\epsilon }^2+A}(t+T_0)|\frac{A}{{\epsilon }^2+A}\right)\right],$$

$$\pi /2\ge \varphi \ge -\pi /2,$$

$$T_0=\frac{1}{\sqrt{{\epsilon }^2+A}}K\left(\frac{A}{{\epsilon }^2+A}\right),$$

$$C_{sG}^{\min }=2\sqrt{{\epsilon }^2+A}E\left(\frac{A}{{\epsilon }^2+A}\right).$$

$${\hat{H}}_{ct}=\frac{i}{2}\left(|{\partial }_{t}u(t)⟩⟨u(t)|-|u(t)⟩⟨{\partial }_{t}u(t)|\right).$$

$${\hat{H}}_0(t)=B_z(t){\hat{\sigma }}_z+B_x(t){\hat{\sigma }}_x,$$

$${\hat{H}}_{ct}(t)=\frac{1}{2}\frac{d}{dt}\left(\arctan \frac{B_x}{B_z}\right){\hat{\sigma }}_y=-\frac{{\dot{B}}_z{B}_x-B_z{\dot{B}}_x}{2\left[(B_z{)}^2+(B_x{)}^2\right]}{\hat{\sigma }}_y,$$

$${\hat{H}}_{ct}=-\frac{1}{2}\dot{\varphi }{\hat{\sigma }}_y,$$

$$\varphi (-\infty )=-\frac{\pi }{2},\varphi (0)=\frac{\pi }{2}.$$

$$\mathcal{L}=\frac{1}{4}{\dot{\varphi }}^2+B^2+2B\epsilon \sin \varphi +{\epsilon }^2+A{\cos }^2\varphi .$$

$$\frac{\partial \mathcal{L}}{\partial B}=2B+2\epsilon \sin \varphi =0,$$

$$B=-\epsilon \sin \varphi .$$

Full text

I Introduction

The future miniaturization technologies will face with the challenge to design elementary information-processing units that employ quantum effects and operate with simultaneously high speed and low dissipation. Such requirements are usually contradictory to each other. Therefore, they raise the need for quantum optimal control, i.e., the design of protocols for application of time-dependent fields that lead to the desired behavior of a quantum system while minimizing some cost function grigorenko-book .

There are two major complications with control optimization at the quantum level. First, such a control has to be nonlinear. This follows already from the fact that physical characteristics of a quantum system are nonlinear functions of dynamic variables, which are the state vector components. Moreover, control fields multiply the state vector in the Schrödinger equation, so this equation is nonlinear if we treat both the state vector and control fields as dynamic variables. Even for a spin-1/2 system, conventional optimal control already becomes complex, requiring involved mathematical treatment and sometimes numerical calculations to find the optimal protocols Boscain1 ; Dalessandro ; Boscain2 ; Garon . The second complication is the exponentially fast growth of the number of variables that parametrize the state vector. For example, $N$ qubits are described quantum mechanically by a state vector with $2^{N}$ components. So, beyond the most elementary cases, even numerical solution of a quantum optimization problem by conventional methods of the control theory is hard to achieve.

In this article, we argue that complexity of the optimal quantum control problem can be considerably reduced in many applications. We propose to search for the desired optimization only within some families of control protocols that allow analytical solution of the nonstationary Schödinger equation. If such a family is sufficiently broad, one can expect that optimization within this family will lead to a reasonable cost while the problem will be tractable.

Finding large classes of exactly solvable models with desired properties may look impossible at first view. Indeed, there is even no known general analytical solution for the dynamics of a spin-1/2 system in an arbitrary time-dependent magnetic field. Higher dimensional cases are even more complex because they require solutions of higher than 2nd order differential equations with time-dependent parameters, which remain poorly understood.

Contrary to such expectations, we would like to point out that it is actually not hard to generate families of time-dependent control protocols whose effects on quantum systems can be understood analytically. One possibility to do this is to apply control fields that simply compensate for terms in the Hamiltonian that are responsible for complex behavior. For example, if there are controlled couplings between some qubits, we can simply set such couplings to zero in order to make individual qubit dynamics easy to control by simple local fields for a while. It is then possible to determine and optimize costs within the family of such protocols analytically.

Another possibility to design a family of solvable Shrödinger equations is to solve an inverse problem. Namely, in equation

[TABLE]

we can prescribe the functional form for some desired types of dynamics of $\Psi(t)$ and then treat (1) as the linear equation for unknown elements of the matrix Hamiltonian $\hat{H}(t)$. Since the number of components of the state vector scales linearly with the size of the phase space $N$ and the number of components of the matrix $\hat{H}$ scales as $N^{2}$, this equation can generally be solved with large redundancy. Many families of such obtained exactly solvable models are known garanin ; sarma ; berry . Models generated this way describe somewhat unnatural dynamics but this is not a drawback for quantum control purposes. After parametrizing the class of solutions of such an inverse problem, we can explore the cost function within the resulting parameter space.

In the rest of the article, we demonstrate these two strategies using the model of a spin-1/2 system control by a time-dependent magnetic field. Instead of straightforward global optimization, we are going to use the fact that there are classes of protocols for which analytical solutions of the Schrodinger equation are known. Within such classes the optimization problem strongly simplifies, but generally the result will not coincide with a fully optimized solution. In this sense, our result is “nearly-optimal”.

II The model

Consider a two level system shown in Fig. 1(a). We assume that the system is initially in the ground state that is separated from the excited state by a natural energy splitting bias $2\varepsilon$. So, the free Hamiltonian of the system is

[TABLE]

where $\hat{\sigma}_{z}$ is the Pauli operator. Assume that the goal is to design pulses of some control field that will lead to the definite transition from the ground state to the excited state. It is convenient to think about a two-level system in terms of a spin-1/2 in a time-dependent magnetic field. So, the general Hamiltonan is

[TABLE]

where $B_{x}^{\rm ex}$, $B_{y}^{\rm ex}$, and $B_{z}^{\rm ex}$ are control parameters, which we will call components of the external magnetic field. There are numerous known ways to design pulses of a magnetic field in order to flip the spin. However, practical situation may impose conditions that favor some protocols over the others. Specifically, we will assume that following four conditions must be met:

(i) During the spin-flip, nonzero average spin polarization along $y$-axis should be avoided by all means possible. Thus we require that $\langle\Psi(t)|\hat{\sigma}_{y}|\Psi(t)\rangle=0$ at any time.

(ii) We want to minimize time that the spin spends not being strictly polarized along $z$-axis. In other words, we want to flip the spin as fast as possible.

(iii) The size of the external control field that we apply should be minimized.

(iv) We are looking for protocols whose properties can be written in terms of some known special functions.

Here we would like to stress that our choice of conditions (i)-(iv) is not motivated by some immediate contemporary experimental research. We chose them for purely illustrative reasons because they provide an example that is sufficiently nontrivial mathematically to illustrate advantages and possible weaknesses of our approach.

For definiteness, we will assume that the control protocol starts at some negative time moment so that spin appears in the state “up” by time $t=0$, after which all external fields are switched off. Condition (i) may be justified in real situation by requirement to avoid unwanted coupling with ambient qubits that are placed along the $y$-axis. We assume that this condition is not negotiable so it merely restricts the allowed phase space for dynamics. In particular, this condition forbids application of popular Rabi-like pulses that induce circulation of spin polarization in the $xy$-plane.

Condition (iv) expresses the desire to work with pulse shapes with easily characterizable properties. This means that we are not looking for a complete solution of the optimization problem in terms of some nonlinear differential equation that can be solved only numerically. Instead, we restrict our search to a specific class of protocols with known analytical solutions. This class is quite large so we hope that the deviation of the cost of our protocol is not essentially different from the cost of a numerically exact solution of the optimization problem.

Conditions (ii) and (iii) are conflicting with each other. In order to minimize time of spin flip, we should apply strong control fields but such fields induce unwanted dissipation that we also want to minimize. In order to resolve this conflict, we should quantify conditions (ii) and (iii) by introducing some cost functional $C$. We will assume the simplest form of $C$ that is consistent with symmetry of the problem and condition that this functional should involve not higher than 2nd order powers of expectations of Pauli matrices and control field components:

[TABLE]

Here the first term inside the integral penalizes all values of the state vector $|\Psi(t)\rangle$ that have nonzero expectation value of $\hat{\sigma}_{x}$ operator. Since we are interested in protocols that asymptotically have $\langle\Psi(t)|\hat{\sigma}_{x}|\Psi(t)\rangle^{2}=0$ at $t\rightarrow-\infty,0$, the time integral of this expression has the meaning of the effective time of spin flip. If the external field were the true magnetic field, the physical meaning of the rest of the cost would be merely the energy of the pulse. So, naturally, this cost favors application of smaller magnetic field values. Parameter $A$ describes the relative importance of the two cost contributions. Bigger values of $A$ favor faster protocols that require larger external field amplitudes. We will not consider other possible restrictions that have been encountered in the literature other .

We can now quantify the problem. Our goal is to find the time-dependent external field that induces dynamics of the state vector that minimizes the functional of the form

[TABLE]

under condition that during the evolution we have $\langle\Psi(t)|\hat{\sigma}_{y}|\Psi(t)\rangle=0$ and the spin is fully polarized in opposite directions along $z$-axis at the beginning and the end of the protocol.

III Optimal square pulse protocol

Without the natural bias (i.e., at $\varepsilon=0$), we would be able to rotate the spin by merely applying the external field along the $y$-axis. This problem is exactly solvable and its optimization was solved in Ref. grigorenko . At $\varepsilon\neq 0$, the formal approach to optimization, however, would involve introduction of Lagrange multipliers that account for evolution of the state vector and constraints. In contrast, as our first strategy, we consider a straightforward way to rotate the spin while keeping its polarization within $xz$-plane: we just apply an external field in the $z$-direction, i.e., momentarily set $B_{z}^{\rm ex}=\varepsilon$, in order to compensate for the intrinsic field bias. Simultaneously, we apply a constant field along the $y$-direction while keeping $B_{x}^{\rm ex}=0$. As a result, the spin will precess around $y$-axis. Let $\phi$ be the angle that our spin makes with the $x$-axis, as shown in Fig. 1(b). For precession in this constant magnetic field we have

[TABLE]

We should switch the external field off when the spin completes rotation by an angle $\pi$. Obviously, for uniform rotation, our protocol should have a time duration

[TABLE]

where the subscript “u” stands for “uniform”, and the cost is

[TABLE]

Using (7), this cost is minimized at $B_{y}^{\rm ex}=\sqrt{\varepsilon^{2}+\frac{A}{2}}$:

[TABLE]

The square-pulse protocol is simple but it is expected to be strongly suboptimal because minimization is performed only within a family of field pulses parametrized by a single parameter $T_{u}$. In what follows, we are going to explore more complex strategies. First, we can improve the above strategy by allowing a more complex time-dependence of $B_{y}^{\rm ex}(t)$. Second, we can consider a different family of protocols that have nonzero values of $B_{x}^{\rm ex}(t)$.

IV Protocols with time-dependent $B_{y}^{\rm ex}$

Let us now assume that we again apply the field with $B^{\rm ex}_{z}=\varepsilon$ and $B^{\rm ex}_{x}=0$ during some time $T_{0}$. Now, we treat $T_{0}$ as a free parameter and we also allow arbitrary time-dependence of $B_{y}^{\rm ex}(t)$. The only constraint on the transverse field is that

[TABLE]

We have then

[TABLE]

To simplify this expression, we will introduce the phase variable

[TABLE]

in terms of which the cost functional now reads

[TABLE]

where the dot denotes differentiation with respect to $t$ and the subscript “$sG$” stands for “sine-Gordon” because Eq. (IV) coincides with a Langrangian of a classical mechanical motion of a particle along coordinate $\phi$ in the sine-Gordon potential $-A\cos^{2}(\phi)$. The term with $\varepsilon^{2}$ in (IV) should be kept because it penalizes protocols with larger duration $T_{0}$.

Varying the cost functional over $\phi(t)$ one can find the equation of motion. Instead of solving this equation, we invoke the classical mechanical analogy and say that the energy $e$ of the particle motion is conserved:

[TABLE]

Hence,

[TABLE]

Formally, duration of the protocol $T_{0}$ should be also optimized. However, varying cost over this parameter we find the value of energy that cannot be used to satisfy the boundary conditions. Instead, we will treat $e$ as an additional parameter for optimization. Integrating over time we find (definitions of special functions are provided in appendix A):

[TABLE]

and using the boundary conditions, $\phi(t=-T_{0})=-\pi/2$ and $\phi(t=0)=\pi/2$, we obtain

[TABLE]

where $K(m)$ is the complete elliptic integral of the first kind (appendix A).

Substituting (13)-(15) into (IV) we find

[TABLE]

where $E(m)$ is the complete elliptic integral of the second kind. Finally, we should minimize the cost over the remaining free parameter $e$. Differentiating (16) with respect to $e$, we find that the minimum is achieved at $e=0$ so,

[TABLE]

where $T_{0}$ is given by

[TABLE]

and the minimal cost within this family of protocols is given by

[TABLE]

V Optimal shortcut to adiabaticity

A bigger class of analytically tractable protocols can be obtained by solving the inverse problem to the nonstationary Schrödinger equation. Many techniques to do this have been developed for spin-1/2 systems garanin ; sarma ; berry . Here we will focus on the class of such solutions called shortcuts to adiabaticity berry ; Demirplak_2008 ; Deffner_2014 ; Bason_2011 ; Malossi_2013 ; Zheng_2016 . This class covers not only spin-1/2 but also many practically interesting multistate situations delCampo .

Let $\hat{H}_{0}(t)$ be some time-dependent Hamiltonian and $|u(t)\rangle$ be one of the instantaneous eigenstates of this Hamiltonian, which continuously depends on $t$. Note that $|u(t)\rangle$ is generally not a solution of the nonstationary Schrödinger equation with $\hat{H}_{0}(t)$. The Hamiltonian of a shortcut to adiabaticity is obtained by adding a counterterm $\hat{H}_{ct}(t)$ such that the evolution of the state vector with the Hamiltonian $\hat{H}(t)=\hat{H}_{0}(t)+\hat{H}_{ct}(t)$ does follow the path of eigenstates $|u(t)\rangle$. It turns out that there is a formal expression for such a counterterm berry ; Bason_2011 ; Zheng_2016 :

[TABLE]

In other words, the instantaneous eigenvector $|u(t)\rangle$ of $\hat{H}_{0}(t)$ becomes the exact solution of the nonstationary Schrödinger equation (1) with the Hamiltonian $\hat{H}(t)=\hat{H}_{0}(t)+\hat{H}_{ct}(t)$, where the counterterm is given by (20).

For the following discussion, it will be convenient to introduce a new vector ${\mathbf{B}}$ with absolute value $B$ and components $B_{x}\equiv B\cos(\phi)=B_{x}^{\rm ext};\,\,\,B_{y}=0;\,\,\,B_{z}\equiv B\sin(\phi)=B_{z}^{\rm ext}-\varepsilon.$ If now we choose

[TABLE]

then, as desired, we will have evolution of the eigenstate $|u(t)\rangle$ such that spin polarization will be always directed along the field ${\mathbf{B}}$, which is restricted to $xz$-plane. For such a spin-1/2 Hamiltonian, the counterterm is Bason_2011 ; Malossi_2013 :

[TABLE]

which can also be expressed as

[TABLE]

i.e., it is induced by the magnetic field $B_{y,ct}=-\dot{\phi}/2$ directed along the $y$-axis. The boundary conditions in terms of the angle $\phi$ read:

[TABLE]

We will also assume that the spin rotates counterclockwise. The cost function (5) now has the Lagrangian

[TABLE]

Since the Lagrangian (25) does not depend explicitly on $\dot{B}$, variation over $B$ produces an algebraic constraint:

[TABLE]

from which we find

[TABLE]

At this point, we encounter a complication. By definition, $B$ is the absolute value of the field vector, so it is always non-negative. Hence, constraint (27) can be satisfied only for $\sin\phi\leq 0$. Since negative values of $B$ are impossible and larger positive values of $B$ are unfavored at $\sin\phi>0$, we have to construct our protocol so that at $\phi\leq 0$ condition (27) is satisfied and at $\phi>0$ we just set $B=0$. Let us now consider these two stages of the protocol separately.

Case $\bm{-\pi/2\leq\phi<0}$: In this region, the field is given by Eq. (27). Substituting (27) into (25), we obtain

[TABLE]

This is a Lagrangian of the classical motion of a particle with mass $m=2$ in a potential $V(\phi)=-\varepsilon^{2}\cos^{2}\phi-A\cos^{2}\phi$. Corresponding energy conservation equation reads:

[TABLE]

Energy $e$ is found here from the observation that at $t\rightarrow-\infty$ we must have $\dot{\phi}\rightarrow 0$, and consequently $e=0$. So,

[TABLE]

From which we find

[TABLE]

At this point, we treat $T$ as a free parameter, whose value we will obtain self-consistently later.

The field components $B_{x}$, $B_{z}$, and $B_{y,ct}$ are found as

[TABLE]

We plot them in Fig. 2(b). We can also calculate the cost of rotating the spin from $\phi=-\pi/2$ to $\phi=0$ along this protocol:

[TABLE]

Case $\bm{\pi/2\geq\phi>0}$: For such values of the rotation angle we have $B=0$, and the Lagrangian (25) reads

[TABLE]

This Lagrangian is the same as the one that we studied in Sec. IV. The only difference is that starting boundary conditions are now $\phi(t=-T)=0$. Using results of previous section and the symmetry of evolution with Lagrangian (36), we find that

[TABLE]

[TABLE]

and the cost of the 2nd part of the protocol:

[TABLE]

Finally, we combine two costs to find

[TABLE]

Figure 2 summarizes dependence of the angle $\phi$ on time as given by Eqs. (V) and (38), as well as the components of the field that induce this spin rotation.

Let us now compare costs of optimal protocols in Eq. (9) for uniform rotation, Eq. (19) for fully compensated bias $B_{z}^{\rm ex}=\varepsilon$ but flexible time-dependence of $B_{y}(t)$, and the cost (40) of the protocol that used shortcuts to adiabaticity. In terms of the dimensionless ratio $\gamma=\varepsilon^{2}/A$, these costs can be written as

[TABLE]

We find that ratios of these costs depend only on $\gamma$. Figure 3 presents the ratios $C^{\rm min}/C_{u}^{\rm min}$ and $C^{\rm min}/C_{sG}^{\rm min}$ vs. $\gamma$. It shows that these two ratios are smaller than 1 for any $\gamma$. Thus, for any choices of parameters $A$ and $\varepsilon$, the protocol using shortcuts to adiabaticity has a lower cost than the protocols based on instantaneous bias compensation. For small and large $\gamma$, the limits of the two ratios are

[TABLE]

Between the two limits, $C^{\rm min}/C^{\rm min}_{u}$ has a maximum at $\gamma=0.0504$, with $C^{\rm min}/C^{\rm min}_{u}=0.905$. Thus, for any $A$ and $\varepsilon$ the protocol that uses shortcuts to adiabaticity leads up to $18.2$ percent reduction of the cost compared to the protocols based on instantaneous compensation of the natural bias term in the Hamiltonian.

VI Conclusions

We demonstrated the possibility of a mostly analytical approach to quantum control problems. Instead of pursuing for numerically exact but complex and nontransparent solutions we explored the possibility to restrict the optimization problem to some sufficiently broad classes of analytically tractable protocols. By construction, such protocols are generally suboptimal. The precise cost function is, however, unknown in most practical situations. There is no reason then to know the numerically exact solution of an optimization problem. It is much more desirable for practically useful protocols to be physically meaningful and analytically tractable. Such a control can be relatively easily planned and adjusted upon receiving additional information about performance. This simplicity is what our approach can provide.

By restricting the search for optimal protocols to a class of solutions of some exactly solvable model, we automatically resolve the problem of combinatorial complexity of quantum systems because state vector components no longer have to be treated as independent parameters. As an additional bonus, such a restriction leads to control field pulses whose shapes can be written in terms of known special functions, and consequently can be easily characterized. Numerical solution may be needed only for a small number of saddle point equations for control parameters,

Using the spin-1/2 example, we demonstrated two basic possibilities for choosing the family of solutions for optimization. One is based on simple compensation of terms in the Hamiltonian that prevented an analytical solution to exist in the first place. Another approach is based on using families of solutions obtained by solving an inverse problem, i.e., assuming some functional form of the state vector evolution (a shortcut to adiabaticity in our case) and then finding the family of Hamiltonians that produces this evolution. By restricting the search for optimal control among such classes of analytically solvable Schödinger equations, we were able to express the cost functional only in terms of control parameters.

The worked out example revealed not only the simplicity of our approach but also potential problems that may become important in more complicated situations. For example, protocols that are based on straightforward removal of unwanted couplings are simple and generally need only control of the local degrees of freedom. We showed, however, that more complex protocols can be substantially less costly. We also found that, interestingly, it was impossible to build the optimal protocol fully within the family of shortcuts to adiabaticity because equations of motion drove one of the parameters beyond its range of definition. We resolved this problem by combining different types of protocols in order to describe different time intervals. Resulting protocol, however, had the best performance among all strategies that we explored. It is expected that similar properties of optimized control protocols will be generally found in multistate problems too. So, it seems that more complex models will still require such manual resolution of complications.

Finally, we note that we left many characteristics of our approach unstudied. For example, in addition to extending it to more complex quantum problems, it should be useful to explore the robustness of our protocols with respect to uncertainty of parameters, e.g., using the integral nullification method introduced in Daems . It should be also interesting to compare our protocols with ones based on more complex, than shortcuts to adiabaticity, classes of solvable models that can be derived by inverting the Schrödinger equation garanin .

Acknowledgements

The work was carried out under the auspices of the National Nuclear Security Administration of the U.S. Department of Energy at Los Alamos National Laboratory under Contract No. DE-AC52-06NA25396. Authors also thank the support from the LDRD program at LANL.

Appendix: Elliptic integrals and Jacobi elliptic functions

Here we list special functions that we have used in this text, namely, the elliptic integrals and Jacobi elliptic functions.

The incomplete elliptic integral of the first kind, $F(\phi|m)$, with modulus $m$ is defined as

[TABLE]

The complete elliptic integral of the first kind, $K(m)$, is defined as

[TABLE]

The incomplete elliptic integral of the second kind, $E(\phi|m)$, is defined as

[TABLE]

and the complete elliptic integral of the second kind, $E(m)$, is defined as

[TABLE]

The Jacobi elliptic functions are the inverses of the incomplete elliptic integrals of the first kind wiki-Jef . If $u=F(\phi|m)$, then the Jacobi elliptic function $\operatorname{sn}$(x—m) is defined as wiki-Jef :

[TABLE]