# Optimal control and ultimate bounds of 1:2 nonlinear quantum systems

Jing-jun Zhu  Laboratoire Interdisciplinaire Carnot de Bourgogne, CNRS UMR
6303, Université de Bourgogne, BP 47870, 21078 Dijon, France  International
Center of Quantum Artificial Intelligence for Science and Technology
(QuArtist) and Department of Physics, Shanghai University, Shanghai 200444,
China  Kaipeng Liu  International Center of Quantum Artificial Intelligence
for Science and Technology (QuArtist) and Department of Physics, Shanghai
University, Shanghai 200444, China  Xi Chen  Department of Physical Chemistry,
University of the Basque Country UPV/EHU, Apartado 644, 48080 Bilbao, Spain
EHU Quantum Center, University of the Basque Country UPV/EHU, Barrio Sarriena,
s/n, 48940 Leioa, Spain  Stéphane Guérin<EMAIL_ADDRESS>Laboratoire
Interdisciplinaire Carnot de Bourgogne, CNRS UMR 6303, Université de
Bourgogne, BP 47870, 21078 Dijon, France

###### Abstract

Using optimal control, we establish and link the ultimate bounds in time
(referred to as quantum speed limit) and energy of two- and three-level
quantum nonlinear systems which feature 1:2 resonance. Despite the unreachable
complete inversion, by using the Pontryagin maximum principle, we determine
the optimal time, pulse area, or energy, for a given arbitrary accuracy. We
show that the third-order Kerr terms can be absorbed in the detuning in order
to lock the dynamics to the resonance. In the two-level problem, we determine
the non-linear counterpart of the optimal $\pi$-pulse inversion for a given
accuracy. In the three-level problem, we obtain an intuitive pulse sequence
similar to the linear counterpart but with different shapes. We prove the
(slow) logarithmic increasing of the optimal time as a function of the
accuracy.

##  I Introduction

The accurate control of quantum dynamics is at the core of the quantum world.
Quantum control protocols have been developed in order to design specific
shaped pulses including composite CPlevitt ; Wimperis ; CPprl106 ; Jones ;
CPprl129 ; CPnjp24 , adiabatic STIRAP ; STIRAP2 ; UltrafastSTIRAP , shortcut
to adiabatic STA ; robustNJP ; STA2 ; pra103  and single-shot shaped pulse
SSprl111 ; Hybrid ; Laforgue  techniques. However, these protocols, even when
accelerated compared to standard adiabatic passage, do not specifically
control the time of operation, which can lead to severe obstructions to
experimental implementation. In this context, optimal control theory (OCT) OCT
has emerged as a powerful tool to mitigate intensities of pulses allowing one
to attain the ultimate time bound in the system, which is also interpreted as
quantum speed limit (QSL) prl103_Tomaso ; prl111 ; Frey ; pra98 . Besides
numerical implementation of OCT, such as monotonically convergent iteration
algorithm OPalgorRuths134 ; OPalgorzhu108 ; OPalgorzhu110 ; nonOppra78 ,
global Krotov method Krotov , gradient ascent algorithms (GRAPE) GRAPE , one
can highlight Pontryagin maximum principle (PMP) PMP ; Boscain ; Extended ;
PMPprx , which, transforming the initial infinite-dimension control problem
into a finite dimension problem, allows analytic derivation of the optimal
controls (typically with respect to time or energy). One can also mention
recent geometric approaches Barnes ; Dridi ; Dridi2  treating simultaneously
robust and optimal control.

The extension of quantum control techniques to non-linear quantum systems
relevant to describe BEC, e.g. when one considers the conversion from atomic
to molecular BEC leading to a so-called 1:2 Fermi resonance pra65_Drummond ,
is a non-trivial issue. The system has to be reinterpreted and analyzed with
tools from classical mechanics, where the concept of integrability, without
counterpart in the standard linear quantum physics, plays an important role
classical_H ; Henrard . When the system is integrable, adiabatic passage
techniques can be formulated with trajectories formed by the instantaneous
(stable) elliptic fixed points defined at each value of the adiabatic
parameters and continuously connected to the initial condition. Obstructions
to classical adiabatic passage are given by the crossing of a separatrix
prl99_Itin ; pra_SG ; epl_stephane . In addition, for a two-level problem with
a 1:2 resonance, the north pole of the generalized Bloch sphere (associated to
the upper state and thus corresponding to a complete population transfer from
the ground state) is unstable since it is associated to an hyperbolic fixed
point in the classical phase space representation pra102_Jingjun . This
prevents adiabatic passage to be robust when it approaches the north pole. The
system is not controllable at this point in the sense that the nonlinearity
prevents reaching the upper state exactly pra_SG . However, one can approach
it as closely as required, and inverse-engineering techniques prl119_Stephane
; STAOPnon  have been developed for that purpose.

Ultimate bounds, e.g. quantum speed limit prl111 , can be defined via the
minimization of a given cost (such as time, pulse area or energy) determined
from optimal control; their extension to nonlinear systems is an open
question. The purpose of this work is to establish and link these ultimate
bounds in terms of time and energy using optimal control. We present a
complete study of optimal control via the PMP for the two- and three-level
systems featuring a 1:2 resonance, considering the cost as time or energy.
Since the complete inversion from the ground state is unreachable, we define
the target with a given (arbitrary) accuracy. This work completes the analysis
made in pra94xc  for the nonlinear two-level system. Among the results, we
extend the notion of (optimal) $\pi$-pulse Boscain  to nonlinear systems (for
a given accuracy) and show the asymptotic logarithmic increasing of the
optimal time as a function of the accuracy instead of Rabi oscillations. We
also show the similarity with the linear case: the nonlinear dynamics is shown
to be identical for time or energy optimum with a constant pulse (or constant
generalized pulse in the case of the three-level system).

Section II and III are devoted to two- and three-states problems,
respectively. We conclude in Section IV.

##  II 1:2 nonlinear two-level model

###  II.1 The model

The two-level model including second-order (with a 1:2 resonance) and third-
order Kerr nonlinearities is characterized by following equations of motion
pra_SG

$\displaystyle
i\dot{\psi}_{1}=\left[-\frac{\Delta}{3}+\Lambda_{11}|\psi_{1}|^{2}+\Lambda_{12}|\psi_{2}|^{2}\right]\psi_{1}+\frac{\Omega}{\sqrt{2}}{\psi}^{*}_{1}\psi_{2},$
(1a) $\displaystyle
i\dot{\psi}_{2}=\left[\frac{\Delta}{3}+\Lambda_{21}|\psi_{1}|^{2}+\Lambda_{22}|\psi_{2}|^{2}\right]\psi_{2}+\frac{\Omega}{2\sqrt{2}}{\psi}^{2}_{1},$
(1b)

where $\psi_{1}$ and $\psi_{2}$ are the state probability amplitude,
satisfying $|\psi_{1}|^{2}+2|\psi_{2}|^{2}=1$ that can vary in the respective
ranges $|\psi_{1}|^{2}\in[0,1]$, $|\psi_{2}|^{2}\in[0,1/2]$. The controls are
time-dependent: $\Delta\equiv\Delta(t)$ and $\Omega\equiv\Omega(t)$,
representing the detuning and Rabi-frequency, respectively. Here
$\Lambda_{ij}(i,j=1,2)$ denote the third-order nonlinearities (in units of
angular frequency) and $\Lambda_{21}=\Lambda_{12}$. Second-order
nonlinearities appear in the form of the coupling. We will use units such that
$\hbar=1$.

We can describe the dynamics on a generalized nonlinear Bloch sphere (see,
e.g., pra102_Jingjun ; Efstathiou ) by introducing the nonlinear coherences
and the population inversion, respectively:

$\displaystyle\eta_{1}$
$\displaystyle=\sqrt{2}\,\text{Re}\bigl{(}\psi^{2}_{1}\bar{\psi}_{2}\bigr{)},\quad\eta_{2}=\sqrt{2}\,\text{Im}\bigl{(}\psi^{2}_{1}\bar{\psi_{2}}\bigr{)},$
(2a) $\displaystyle\eta_{3}$
$\displaystyle=|\psi_{2}|^{2}-\frac{1}{2}|\psi_{1}|^{2},\quad\eta_{3}\in\Bigl{[}-\frac{1}{2},\frac{1}{2}\Bigr{]},$
(2b)

leading to

$\displaystyle|\psi_{1}|^{2}=\frac{1}{2}(1-2\eta_{3}),\quad|\psi_{2}|^{2}=\frac{1}{4}(1+2\eta_{3}).$
(3)

The generalized 1:2 nonlinear Bloch sphere is characterized by the following
surface equation

$\displaystyle\eta^{2}_{1}+\eta^{2}_{2}-\Bigl{(}\frac{1}{2}-\eta_{3}\Bigr{)}^{2}\Bigl{(}\frac{1}{2}+\eta_{3}\Bigr{)}=0.$
(4)

The south and north poles correspond to $|\psi_{1}|^{2}=1$,
$|\psi_{2}|^{2}=0$, i.e. $\eta_{3}=-1/2$ and $|\psi_{1}|^{2}=0$,
$|\psi_{2}|^{2}=1/2$, i.e. $\eta_{3}=1/2$, respectively. It has been proved
that $\eta_{3}=1/2$ is an unreachable (unstable) target with or without Kerr
terms in pra_SG ; epl_stephane . Using Eq. (1), we have

$\displaystyle\dot{\eta}_{1}=\bigl{(}-\Delta+\Lambda_{a}-2\Lambda_{s}|\psi_{2}|^{2}\bigr{)}\eta_{2},$
(5a)
$\displaystyle\dot{\eta}_{2}=\frac{\Omega}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}+\bigl{(}\Delta-\Lambda_{a}+2\Lambda_{s}|\psi_{2}|^{2}\bigr{)}\eta_{1},$
(5b) $\displaystyle\dot{\eta}_{3}=\Omega\eta_{2}$ (5c)

with the effective third-order nonlinearities

$\Lambda_{s}=2\Lambda_{11}+\Lambda_{22}/2-2\Lambda_{12},\quad\Lambda_{a}=2\Lambda_{11}-\Lambda_{21}.$
(6)

It can be seen that $\Lambda_{a}$ can be trivially compensated by a static
shift of the detuning, while the term proportional to $\Lambda_{s}$ depends on
the dynamical variable $|\psi_{2}|^{2}$. However it has been shown in
prl119_Stephane  that one can lock the resonance using the freedom in the
choice of the time-dependance of $\Delta$ by incorporating the term
$2\Lambda_{s}|\psi_{2}|^{2}$. Hence, we define the effective detuning (which
includes a change of sign of $\Delta$ for convenience):

$\tilde{\Delta}=-\Delta+\Lambda_{a}-2\Lambda_{s}|\psi_{2}|^{2}=-\Delta+\Lambda_{a}-\Lambda_{s}\Bigl{(}\frac{1}{2}+\eta_{3}\Bigr{)},$
(7)

such that the set (5) of differential equations only features the second-oder
nonlinearity:

$\displaystyle\dot{\eta}_{1}=\tilde{\Delta}\eta_{2},$ (8a)
$\displaystyle\dot{\eta}_{2}=\frac{\Omega}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}-\tilde{\Delta}\eta_{1},$
(8b) $\displaystyle\dot{\eta}_{3}=\Omega\eta_{2}.$ (8c)

This generalizes the standard Bloch equations to the 1:2 nonlinear system with
the detuning $\tilde{\Delta}$ (7) featuring an explicit (linear) dependance on
$\eta_{3}$ via the effective third-order nonlinear term $\Lambda_{s}$. To
determine the expression of optimal trajectories in this system, we apply the
PMP taking for the cost the time (time optimal control) or the energy (energy
optimal control). The 1:2 nonlinear two-state time optimal control has been
solved in pra94xc , and we closely follow its derivation in order to apply it
for the energy optimal control and in next section for the three-level
problem.

To achieve a given transfer from a population
$\eta_{\text{3i}}\equiv\eta_{3}(t_{i})$ at the initial time $t_{i}$, we define
the targeted final population
$\eta_{\text{3f}}\equiv\eta_{3}(t_{f})=\frac{1}{2}-\epsilon$ at the final time
$t_{f}$, i.e. $|\psi_{2}(t_{f})|^{2}=\frac{1}{2}(1-\epsilon)$. When one
targets the upper state, the deviation $\epsilon$ will be taken small but
different from zero since the second-order nonlinearity prevents to reach it
exactly.

###  II.2 Time and area optimal control

Various situations have been considered in pra94xc . Here, we particularly
focus on the (almost) complete transfer from the ground state as a function of
$\epsilon$, and will show that the minimum time (or pulse area) increases in a
(slow) logarithmic way with respect to small deviations $\epsilon$, or
reciprocally that the deviation $\epsilon$ decreases exponentially with
respect to the minimum time (or pulse area).

We consider the time minimizing functional

$\displaystyle J=\int_{t_{i}}^{t_{f}}dt.$ (9)

The corresponding control (or pseudo) Hamiltonian from the set (8) of
differential equations is (where we have added a constant $p_{0}$):

$\begin{split}&h_{c}=\tilde{\Delta}(\lambda_{1}\eta_{2}-\lambda_{2}\eta_{1})+\Omega\Bigl{[}\frac{\lambda_{2}}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}+\lambda_{3}\eta_{2}\Bigr{]}\end{split}$
(10)

with Hamiltonian’s equation for the (dimensionless) costate
$\Lambda=[\lambda_{1},\lambda_{2},\lambda_{3}]^{\mathsf{T}}$ gathering the
conjugate momenta of $\eta_{1}$, $\eta_{2}$, and $\eta_{3}$, respectively:

$\displaystyle\dot{\lambda}_{1}=-\frac{\partial
h_{c}}{\partial\eta_{1}}=\lambda_{2}\tilde{\Delta},$ (11a)
$\displaystyle\dot{\lambda}_{2}=-\frac{\partial
h_{c}}{\partial\eta_{2}}=-\lambda_{1}\tilde{\Delta}-\lambda_{3}\Omega,$ (11b)
$\displaystyle\dot{\lambda}_{3}=-\frac{\partial
h_{c}}{\partial\eta_{3}}=-\frac{\lambda_{2}}{2}\Omega(6\eta_{3}-1)-\frac{\partial\tilde{\Delta}}{\partial\eta_{3}}(\lambda_{1}\eta_{2}-\lambda_{2}\eta_{1}).$
(11c)

In order to prevent arbitrary large field amplitude detrimental for
experimental implementation, we impose a boundary on the field
$\Omega\leq\Omega_{0}$ as a constraint. The maximization of $h_{c}$ according
to the PMP corresponds thus to the necessary condition

$\displaystyle\frac{\partial h_{c}}{\partial\tilde{\Delta}}$
$\displaystyle=0,$ (12)

i. e.

$\lambda_{1}\eta_{2}-\lambda_{2}\eta_{1}=0.$ (13)

This leads to the shape of the external field $\Omega$:

$\Omega=\frac{2}{\lambda_{2}\left(3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\right)+2\lambda_{3}\eta_{2}}$
(14)

from (10,13) and the fact that the system (10) is autonomous, i.e.
$h_{c}=const.$ In Eq. (14), we have renormalized
$\lambda_{2}/h_{c}\to\lambda_{2}$ and $\lambda_{3}/h_{c}\to\lambda_{3}$
without loss of generality. Differentiating Eqs. (13) and (14) gives

$\displaystyle\frac{\lambda_{1}}{2}\left(3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\right)+\lambda_{3}\eta_{1}=0,$
(15) $\displaystyle\dot{\Omega}=0,$ (16)

from which we conclude that $\Omega$ is constant, taken at its maximum
$\Omega=\Omega_{0}$. Multiplying Eq. (14) by $\eta_{1}$ and Eq. (15) by
$\eta_{2}$ also using (13), we obtain a linear system of equations for the
variables $\lambda_{2}$ and $\lambda_{3}$:

$\displaystyle\lambda_{2}\frac{\eta_{1}}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}+\lambda_{3}\eta_{1}\eta_{2}=\frac{\eta_{1}}{\Omega},$
(17a)
$\displaystyle\lambda_{2}\frac{\eta_{1}}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}+\lambda_{3}\eta_{1}\eta_{2}=0$
(17b)

of determinant zero, which can give a solution when the inhomogeneous terms
are zero, i.e. $\eta_{1}=0$. This implies $\lambda_{1}=0$ from Eq. (15) (for a
non-constant $\eta_{3}$) and $\dot{\eta}_{1}=0$ in Eq. (8a) gives
$\tilde{\Delta}=0$, i.e. for the original detuning $\Delta$ from (7):

$\Delta=\Lambda_{a}-\Lambda_{s}\Bigl{(}\frac{1}{2}+\eta_{3}\Bigr{)}.$ (18)

This leads to an optimal trajectory along the meridian connecting the south to
the target near the north pole (of distance $\epsilon$ from it) in the
$(\eta_{2},\eta_{3})$ plane. The dynamics can be solved exactly from (8). For
instance, when we consider a population transfer from the ground state (south
pole), i.e. $\eta_{\text{3i}}=-1/2$, we obtain (taking $t_{i}=0$)

$\eta_{3}(t)=\tanh^{2}\Bigl{(}\frac{1}{2}\Omega_{0}t\Bigr{)}-\frac{1}{2}.$
(19)

Figure 1: Populations (lower frame), coherence and population inversion (upper
frame) history, given by (19) and integration of (8b), governed by the time
optimal (constant) pulse $\Omega_{0}$ for $\epsilon=0.002$, giving
$\mathcal{A}\approx 7.60$, i.e. $T_{\text{min}}\approx 7.60/\Omega_{0}$.

From Eq. (14), taken at initial time, we obtain
$\Omega_{0}=2/\lambda_{2\text{i}}$ (independently of the initial value of
$\lambda_{3,\text{i}}$). To obtain the explicit expression of the minimum time
$T_{\min}=\min(t_{f}-t_{i})$ for a given $\Omega_{0}$, we calculate the
corresponding minimum pulse area by integrating (8c) using (4), which is fully
determined by the (given) initial and final boundaries of $\eta_{3}$,

$\mathcal{A}_{\min}\equiv\Omega_{0}T_{\min}=\pm\int_{\eta_{\text{3i}}}^{\eta_{\text{3f}}}\frac{d\eta_{3}}{\sqrt{\left(\frac{1}{2}-\eta_{3}\right)^{2}\left(\frac{1}{2}+\eta_{3}\right)}}.$
(20)

The $\pm$ sign ensures a non-negative pulse area, i.e. the sign $+$ ($-$)
corresponds to $\eta_{\text{3i}}<\eta_{\text{3f}}$
($\eta_{\text{3i}}>\eta_{\text{3f}}$). Hereby, we consider
$\eta_{3\text{i}}<\eta_{3\text{f}}=\frac{1}{2}-\epsilon$ and finally get the
minimum time for given $\epsilon$ and $\Omega_{0}$,

$\displaystyle
T_{\text{min}}=\frac{2}{\Omega_{0}}\left|\text{atanh}\sqrt{\frac{1}{2}+\eta_{\text{3f}}}-\text{atanh}\sqrt{\frac{1}{2}+\eta_{\text{3i}}}\right|.$
(21)

This gives for the non-linear final transfer probability from the ground state
in optimal time $T_{\min}$:

$p=2|\psi_{2}(T_{\min})|^{2}=\tanh^{2}\Bigl{(}\frac{1}{2}\Omega_{0}T_{\min}\Bigr{)}=1-\epsilon.$
(22)

We notice in the limit case of unbounded pulse amplitude a Dirac $\delta$
pulse, i.e. of infinite amplitude and zero duration with a finite area
$\Omega_{0}T_{\text{min}}$ given by (21). The dynamics from the ground state
is shown in Fig. 1 for $\epsilon=0.002$.

The minimum time $T_{\text{min}}$ can be used as the definition of the so-
called quantum speed limit in this system as suggested in Frey . In linear
systems, the mimum time is given by $T_{\min,\text{lin.}}$:
$\cos\bigl{(}\Omega_{0}T_{\min,\text{lin.}}/2\bigr{)}=\frac{1}{2}\sqrt{1-\eta_{3\text{i}}}\sqrt{1-\eta_{3\text{f}}}+\frac{1}{2}\sqrt{1+\eta_{3\text{i}}}\sqrt{1+\eta_{3\text{f}}}$
Boscain ; Frey  with $\eta_{3\text{i}}=-1$ for the ground state and
$\eta_{3\text{f}}=1$ for the excited state. In the non-linear case we obtain
$\tanh(\Omega_{0}T_{\text{min}}/2)=\sqrt{\frac{1}{2}+\eta_{\text{3f}}}$ for
$\eta_{3\text{i}}=-1/2$ from Eq. (21).

We consider $\eta_{\text{3i}}=-1/2$ (south pole). The pulse area
$\mathcal{A}_{\min}=\Omega_{0}T_{\min}$ is given from Eq. (21) with its
leading order for $\epsilon\to 0$:

$\mathcal{A}_{\min}=2\,\text{atanh}\sqrt{1-\epsilon}\sim-\ln\Bigl{(}\frac{\epsilon}{4}\Bigr{)},\text{
i.e. }\epsilon\sim 4e^{-\mathcal{A}_{\min}}.$ (23)

We conclude that the minimum time (or pulse area) for the (almost complete)
inversion increases in a (slow) logarithmic way with respect to small
deviations $\epsilon$, or reciprocally that the deviation $\epsilon$ decreases
exponentially with respect to the minimum time (or pulse area). Figure 2
compares the pulse area cost between the linear,
$p_{\text{lin}}=\sin^{2}\left(\Omega_{0}T_{\min}/2\right)$, and the non-linear
probability (22). One can observe that the non-linearity affects weakly the
transfer for small transfers ($\epsilon\to 1$).

Figure 2: Population transfer probability $1-\epsilon$ as a function of the
optimal pulse area for the linear
$p_{\text{lin}}=\sin^{2}\left(\Omega_{0}T_{\min}/2\right)$ (green dotted-
dashed line) and nonlinear (red full line) two-level models, respectively.
Asymptotics (23) of the probability (blue dashed line) shows its accuracy when
$\epsilon\to 0$.

###  II.3 Energy optimal control

We consider the cost functional for the energy optimal control

$J\equiv E=\hbar\int^{t_{f}}_{t_{i}}\Omega^{2}(t)dt.$ (24)

In this case, we can rewrite the control Hamiltonian (in units such that
$\hbar=1$)

$\displaystyle h_{c}$
$\displaystyle=\tilde{\Delta}(\lambda_{1}\eta_{2}-\lambda_{2}\eta_{1})$
$\displaystyle+\Omega\Bigl{[}\frac{\lambda_{2}}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}+\lambda_{3}\eta_{2}\Bigr{]}-p_{0}\Omega^{2}$
(25)

with the standard choice $p_{0}=1/2$. In this case, the costate has the
angular frequency unit. After applying the PMP,

$\displaystyle\frac{\partial
h_{c}}{\partial\tilde{\Delta}}=0,\quad\frac{\partial
h_{c}}{\partial\Omega}=0,$ (26)

we obtain $\lambda_{1}\eta_{2}-\lambda_{2}\eta_{1}=0$ and

$\displaystyle\Omega=\frac{\lambda_{2}}{2}\Bigl{(}3\eta^{2}_{3}-\eta_{3}-\frac{1}{4}\Bigr{)}+\lambda_{3}\eta_{2}.$
(27)

The control Hamiltonian can be rewritten as

$\displaystyle h_{c}=\frac{1}{2}\Omega^{2},$ (28)

leading to a constant coupling $\Omega=\Omega_{0}=\sqrt{2h_{c}}$. Similarly as
for the time-optimum case, we obtain $\eta_{1}=0$, $\lambda_{1}=0$,
$\tilde{\Delta}=0$, and thus the same dynamics as in Fig. 1. We also derive
Eq. (21) but interpreted differently, i.e. for a given time of interaction
$T=t_{f}-t_{i}$, we determine the minimum $\Omega_{0}$:

$\displaystyle\Omega_{0,\min}=\frac{2}{T}\left|\text{atanh}\sqrt{\frac{1}{2}+\eta_{\text{3f}}}-\text{atanh}\sqrt{\frac{1}{2}+\eta_{\text{3i}}}\right|$
(29)

leading to the minimum energy

$\displaystyle
E_{\min}=\hbar\Omega_{0,\min}^{2}T=\hbar\mathcal{A}_{\min}^{2}/T$ (30)

corresponding to the minimum area given by (29):
$\mathcal{A}_{\min}=\Omega_{0,\min}T=2\left|\text{atanh}\sqrt{\frac{1}{2}+\eta_{\text{3f}}}-\text{atanh}\sqrt{\frac{1}{2}+\eta_{\text{3i}}}\right|$.

##  III 1:2 nonlinear three-level $\Lambda$ model

###  III.1 The model

The equations of motion for the three-level Raman model (forming a $\Lambda$
system), including the second- and third-order nonlinearities, read

$\displaystyle i\dot{\psi}_{1}=K_{1}\psi_{1}+\Omega_{p}\psi^{*}_{1}\psi_{2},$
(31a) $\displaystyle
i\dot{\psi}_{2}=K_{2}\psi_{2}+\Delta_{P}\psi_{2}+\frac{\Omega_{p}}{2}\psi_{1}^{2}+\frac{\Omega_{s}}{2}\psi_{3},$
(31b) $\displaystyle
i\dot{\psi}_{3}=K_{3}\psi_{3}+\frac{\Omega_{s}}{2}\psi_{2}+(\Delta_{P}-\Delta_{S})\psi_{3},$
(31c)

where $\Omega_{p}$ and $\Omega_{s}$ are the time-dependent pump and Stokes
fields, respectively, $\Delta_{P}$ is the one-photon detuning associated to
the pump coupling, and $\Delta_{P}-\Delta_{S}$ is the two-photon detuning
associated to the Raman process. The second-order 1:2 nonlinearity appears
here through the pump coupling. This is typically the situation for the two-
color photoassociation process pra65_Drummond . The third-order nonlinearities
$K_{j}$, $j=1,2,3$, can be absorbed in the definition of the detunings and the
change of phases $\psi_{1}\to\psi_{1}e^{-i\gamma}$,
$\psi_{2,3}\to\psi_{2,3}e^{-2i\gamma}$: $\dot{\gamma}=K_{1}$,
$\Delta_{P}=2K_{1}-K_{2}$, $\Delta_{S}=K_{3}-K_{2}$, in order to lock the
(one- and two-photon) resonances prl119_Stephane :

$\displaystyle i\dot{\psi}_{1}=\Omega_{p}\psi^{*}_{1}\psi_{2},$ (32a)
$\displaystyle
i\dot{\psi}_{2}=\frac{\Omega_{p}}{2}\psi_{1}^{2}+\frac{\Omega_{s}}{2}\psi_{3},$
(32b) $\displaystyle i\dot{\psi}_{3}=\frac{\Omega_{s}}{2}\psi_{2}.$ (32c)

The amplitude probabilities satisfy
$|\psi_{1}|^{2}+2(|\psi_{2}|^{2}+|\psi_{3}|^{2})=1$. We decompose the
components into real and imaginary parts, $\psi_{i}=x_{i}+iy_{i}$ $(i=1,2,3)$,
and assume real Rabi frequencies, which allows one to separate the original
problem into two disjoint dynamics that emerge according to the real and
imaginary parts of the initial state, respectively. We consider the initial
state $\psi_{1}$ being real, $x_{1}(t_{i})=1$, and we get the equations of
motion from (32):

$\displaystyle\dot{x}_{1}=\Omega_{p}x_{1}y_{2},$ (33a)
$\displaystyle\dot{y}_{2}=-\frac{\Omega_{s}}{2}x_{3}-\frac{\Omega_{p}}{2}x^{2}_{1},$
(33b) $\displaystyle\dot{x}_{3}=\frac{\Omega_{s}}{2}y_{2}.$ (33c)

This system of equations can be analyzed using an isomorphism with the non-
linear two-level problem, similarly to the linear problem ShoreVitanov  (see
Appendix A). One can also show the incomplete transfer between the two ground
states for finite pulse areas. We however prefer to keep the original
coordinates for solving the problem.

Without loss of generality, we can parameterize the dynamics with the two
dynamical angles $\theta(t)\in[0,\pi]$ and $\phi(t)\in[0,2\pi[$:

$\displaystyle x_{1}=\cos\phi\cos\theta,$ (34a) $\displaystyle
y_{2}=-\frac{1}{\sqrt{2}}\sin\phi,$ (34b) $\displaystyle
x_{3}=-\frac{\cos\phi\sin\theta}{\sqrt{2}},$ (34c)

which satisfies the normalization condition
$x_{1}^{2}+2(y_{2}^{2}+x_{3}^{2})=1$. Inserting the definition (34) into Eqs.
(33) leads to

$\displaystyle\dot{\phi}=\frac{\Omega_{p}\cos\phi\cos^{2}\theta}{\sqrt{2}}-\frac{\Omega_{s}\sin\theta}{2},$
(35a)
$\displaystyle\dot{\theta}=\frac{\sin\phi\sin\theta}{\cos\phi\cos\theta}\left(\frac{\Omega_{s}}{2\sin\theta}+\frac{\Omega_{p}\cos\phi\cos^{2}\theta}{\sqrt{2}}-\frac{\Omega_{s}\sin\theta}{2}\right),$
(35b)

which by inversion provide the shape of the fields as functions of the angles:

$\displaystyle\Omega_{s}=2\left(\dot{\theta}\cot\phi\cos\theta-\dot{\phi}\sin\theta\right),$
(36a)
$\displaystyle\Omega_{p}=\sqrt{2}\left(\frac{\dot{\phi}}{\cos\phi}+\frac{\dot{\theta}\tan\theta}{\sin\phi}\right).$
(36b)

We define the target state in the vicinity of the state 3:

$\displaystyle|x_{3}(t_{f})|^{2}=\frac{1}{2}(1-\epsilon),$ (37)

where $\epsilon$ is, as in the two-level case, a small deviation due to the
incomplete population transfer in such nonlinear system prl119_Stephane . To
determine the expression of the optimal trajectories from the ground state to
the target state (for a given $\epsilon$), we apply the PMP taking for the
cost the time or the energy.

###  III.2 Time optimal control

####  III.2.1 Definition

We consider the time minimizing cost functional (9) and impose similarly to
the two-level case the constraint of bounded pulses amplitudes:

$\Omega^{2}_{p}+\Omega^{2}_{s}\leq\Omega^{2}_{0}.$ (38)

The control Hamiltonian reads (where we have added a constant $p_{0}$)

$\displaystyle H_{c}$
$\displaystyle=\lambda_{\phi}\dot{\phi}+\lambda_{\theta}\dot{\theta},$ (39)
$\displaystyle=\lambda_{\phi}\left(\frac{\Omega_{p}\cos\phi\cos^{2}\theta}{\sqrt{2}}-\frac{\Omega_{s}\sin\theta}{2}\right)$
$\displaystyle+\lambda_{\theta}\left(\frac{\Omega_{s}\cos\theta\tan\phi}{2}+\frac{\Omega_{p}\cos\theta\sin\theta\sin\phi}{\sqrt{2}}\right),$
(40)

where $\lambda_{\phi,\theta}$ are the component of the (dimensionless) costate
$\Lambda=[\lambda_{\phi},\lambda_{\theta}]^{\mathsf{T}}$ of dynamics

$\displaystyle\dot{\lambda}_{\phi}=-\frac{\partial H_{c}}{\partial\phi}$
$\displaystyle=\frac{\lambda_{\phi}{\Omega}_{p}\sin\phi\cos^{2}\theta}{\sqrt{2}}-\lambda_{\theta}\left(\frac{{\Omega}_{s}\cos\theta}{2\cos^{2}\phi}+\frac{{\Omega}_{p}\sin
2\theta\cos\phi}{2\sqrt{2}}\right),$ (41a)
$\displaystyle\dot{\lambda_{\theta}}=-\frac{\partial H_{c}}{\partial\theta}$
$\displaystyle=\lambda_{\phi}\left(\frac{{\Omega}_{p}\cos\phi\sin
2\theta}{\sqrt{2}}+\frac{\Omega_{s}\cos\theta}{2}\right)$
$\displaystyle+\lambda_{\theta}\left(\frac{{\Omega}_{s}\sin\theta\tan\phi}{2}-\frac{{\Omega}_{p}\cos
2\theta\sin\phi}{\sqrt{2}}\right).$ (41b)

Without loss of the generality, we can consider the control functions
satisfying the contraint (38) as follows:

$\displaystyle\Omega_{p}=\Omega_{m}(t)\cos\beta(t),\quad\Omega_{s}=\Omega_{m}(t)\sin\beta(t)$
(42)

with the condition

$\Omega^{2}_{s}+\Omega^{2}_{p}=\Omega_{m}^{2}\leq\Omega_{0}^{2}.$ (43)

The constraint (38) is thus transferred to the condition (43) (which is
independent on $\beta$). The PMP maximization of $H_{c}$ is thus reduced to
the necessary condition

$\frac{\partial H_{c}}{\partial\beta}=0,$ (44)

which gives

$\frac{\partial
H_{c}}{\partial\beta}=\frac{\partial\Omega_{p}}{\partial\beta}H_{1}+\frac{\partial\Omega_{s}}{\partial\beta}H_{2}=-H_{1}\sin\beta+H_{2}\cos\beta=0,$
(45)

where

$\displaystyle
H_{1}=\frac{\lambda_{\phi}\cos\phi\cos^{2}\theta}{\sqrt{2}}+\frac{\lambda_{\theta}\sin
2\theta\sin\phi}{2\sqrt{2}},$ (46a) $\displaystyle
H_{2}=\frac{\lambda_{\theta}\cos\theta\tan\phi}{2}-\frac{\lambda_{\phi}\sin\theta}{2}.$
(46b)

We then deduce

$\displaystyle\cos\beta=\frac{H_{1}}{\sqrt{H^{2}_{1}+H^{2}_{2}}},\quad\sin\beta=\frac{H_{2}}{\sqrt{H^{2}_{1}+H^{2}_{2}}}.$
(47)

Substituting Eq. (47) in Eq. (39), we obtain

$H_{c}=\Omega_{m}\sqrt{H^{2}_{1}+H^{2}_{2}}.$ (48)

from which we conclude that the control Hamiltonian is maximum for
$\Omega_{m}=\Omega_{0}$, i.e. when the maximum of (38) is reached at all
times.

Figure 3: Contour plot for $\log_{10}(T-T_{\min})$ as a function of
dimensionless $\lambda_{\phi,\text{i}}$ and $\lambda_{\theta,\text{i}}$ when
the dynamics reaches the targeted transfer, for a given $\epsilon$. The
transfer time is determined in a given range $T\in[0,\hat{T}]$ (taken here as
$\hat{T}=15/\Omega_{0}$). The white areas mean that there is no transfer time
found in the prescribed interval. We obtain $T_{\min}\approx 7.40/\Omega_{0}$
for $\epsilon=0.002$ (upper frame), $T_{\min}\approx 6.78/\Omega_{0}$ for
$\epsilon=0.005$ (lower frame) along the dark straight lines, respectively.

####  III.2.2 Numerics

We first determine the systematic landscape of the time $T$ to reach the
target (for a given $\epsilon$) by solving the set of equations (35,41) with
the controls (42) and $\Omega_{m}=\Omega_{0}$ rewritten as functions of the
angles and the costate components via (47,46), as a function of the two
parameters $\lambda_{\phi,\text{i}}$ and $\lambda_{\theta,\text{i}}$, which
with $\theta_{\text{i}}=0$, $\phi_{\text{i}}=0$ form the set of initial
conditions. The landscape is shown in Fig. 3, where the white areas
corresponds to the absence of solution reaching the given target in the
prescribed interval. Figure 3 shows an infinite set of initial
$\lambda_{\phi\text{i}}$, $\lambda_{\theta\text{i}}$ forming two straight
(symmetric) lines that lead to the same minimum $T_{\min}\approx
7.4/\Omega_{0}$.

The four quadrants give all the possible respective signs of the controls. The
controls are both positive when the initial values $\lambda_{\phi\text{i}}$,
$\lambda_{\theta\text{i}}$ are taken positive. In order to determine a more
accurate value of the optimum, we choose a certain value
$\lambda_{\phi\text{i}}$ (e.g. $\lambda_{\phi\text{i}}=1.85$) and run an
optimal procedure leading to the (positive) value
$\lambda_{\theta\text{i}}\approx 0.45266$ corresponding to the minimum time
(using a Nelder-Mead simplex algorithm as described in NMSM_matlab ).

Figure 4: Optimal pulse area (blue full line) as a function of $\epsilon$ (in
a logarithmic scale) for $\lambda_{\phi\text{i}}=1.85$ and the resulting
optimum $\lambda_{\theta\text{i}}$ and the logarithmic fit of the asymptotics
for small $\epsilon$: ${\cal
A}_{\min}=\Omega_{0}T_{\min}\sim-(\ln\epsilon)/\sqrt{2}+3$ (red dashed line).
Figure 5: Optimal pulse shapes (36) in units of $\Omega_{0}$ with
$\epsilon=0.002$, $\lambda_{\phi\text{i}}=1.85$ and
$\lambda_{\theta\text{i}}\approx 0.45266$ operating in minimum time
$T_{\min}\approx 7.40/\Omega_{0}$ (upper frame), and the resulting population
history and the ansatz (49) (lower frame).

Figure. 4 shows as a function of $\epsilon$ the minimum time, via the minimum
generalized pulse area defined as ${\cal
A}_{\min}=\int_{0}^{T_{\min}}dt\sqrt{\Omega_{p}^{2}+\Omega_{s}^{2}}=T_{\min}\Omega_{0}$.
It exhibits a logarithmic decreasing behavior similar to the two-level case.
This suggests the following ansatz for the population inspired by the two-
level problem and the behavior of the optimal pulse area in Fig. 4 (taking
$t_{i}=0$):

$\displaystyle
y_{2}^{2}+x_{3}^{2}\approx\frac{1}{2}\tanh^{2}\Big{(}\frac{\Omega_{0}}{\sqrt{2}}t\Bigr{)},$
(49)

which fits well the numerics of the dynamics shown in Fig. 5, but slightly
overestimates the accuracy (of the order of $\epsilon$). The boundary at
$t_{f}=T_{\min}$ gives the asymptotic expansion for small $\epsilon$ (with
$y_{2\text{f}}^{2}\ll x_{3\text{f}}^{2}$):

$\displaystyle
y_{2\text{f}}^{2}+x_{3\text{f}}^{2}\approx\frac{1}{2}(1-\epsilon)\sim\frac{1}{2}\bigl{(}1-4e^{-\sqrt{2}{\cal
A}_{\min}}\bigr{)},$ (50)

i.e.

$\displaystyle{\cal A}_{\min}\sim-\frac{1}{\sqrt{2}}(\ln\epsilon-\ln 4).$ (51)

We notice that the resulting logarithmic scaling with respect to $\epsilon$ is
the same as the scaling of the fit determined from Fig. 4. The absolute value
is different due to the systematic error of (49) mentioned above.

The dynamics and the controls of Fig. 5 show an intuitive sequence of pulses
with a large transient population in the upper state, similarly to the linear
case which features pulses of explicit form $\cos-\sin$ Boscain ; TNXu .
Figure 6 shows the trajectory in the angles $\phi,\theta$ space corresponding
to the optimal non-linear solution for $\epsilon=0.002$, it is compared to the
optimal solution of the linear problem.

Figure 6: Trajectories in the $\phi,\theta$ space of the optimal non-linear
solution corresponding to the dynamics shown in Fig. 5 (full line) and the
optimal solution of the linear model (dotted line).

###  III.3 Energy optimal control

For the energy optimization, we use the cost functional

$J\equiv E=\hbar\int^{t_{f}}_{t_{i}}(\Omega^{2}_{p}+\Omega^{2}_{s})dt.$ (52)

We derive the control Hamiltonian (using the standard value $p_{0}=1/2$):

$\displaystyle H_{c}$
$\displaystyle=\lambda_{\phi}\dot{\phi}+\lambda_{\theta}\dot{\theta}-p_{0}(\Omega_{p}^{2}+\Omega_{s}^{2}),$
$\displaystyle=\lambda_{\phi}\left(\frac{\Omega_{p}\cos\phi\cos^{2}\theta}{\sqrt{2}}-\frac{\Omega_{s}\sin\theta}{2}\right)$
$\displaystyle+\lambda_{\theta}\left(\frac{\Omega_{s}\cos\theta\tan\phi}{2}+\frac{\Omega_{p}\cos\theta\sin\theta\sin\phi}{\sqrt{2}}\right)$
$\displaystyle-\frac{1}{2}(\Omega_{p}^{2}+\Omega_{s}^{2}).$ (53)

The PMP $\frac{\partial H_{c}}{\partial\Omega_{p}}=0$, $\frac{\partial
H_{c}}{\partial\Omega_{s}}=0$ leads to the pulse shape:

$\displaystyle\Omega_{p}=\frac{\lambda_{\phi}\cos\phi\cos^{2}\theta}{\sqrt{2}}+\frac{\lambda_{\theta}\sin
2\theta\sin\phi}{2\sqrt{2}},$ (54a)
$\displaystyle\Omega_{s}=\frac{\lambda_{\theta}\cos\theta\tan\phi}{2}-\frac{\lambda_{\phi}\sin\theta}{2},$
(54b)

and

$H_{c}=\frac{1}{2}({\Omega^{2}_{p}+\Omega^{2}_{s}}),$ (55)

which is constant as before, i.e.
$\Omega^{2}_{p}+\Omega^{2}_{s}=\Omega_{0}^{2}$. The dynamics of the components
of the (angular frequency unit) costate
$\Lambda=[\lambda_{\phi},\lambda_{\theta}]^{\mathsf{T}}$ is still given by
(41). We thus obtain the same dynamics shown in Fig. 5 as in the time-optimal
control but for a given time of interaction $T=t_{f}-t_{i}$ (instead of
$T_{\min}$). The corresponding optimal generalized pulse area ${\cal
A}_{\min}$ is determined from Fig. 4 for a given $\epsilon$: ${\cal
A}_{\min}\approx-\ln(\epsilon)/\sqrt{2}+3$, and we deduce the corresponding
$\Omega_{0,\min}={\cal A}_{\min}/T$. The minimum energy is then

$\displaystyle
E_{\min}=\hbar\Omega_{0,\min}^{2}T=\hbar\mathcal{A}_{\min}^{2}/T.$ (56)

##  IV Conclusion

In this paper, we have determined the ultimate bounds in terms of optimal time
and optimal energy for the two- and Raman three-level problems featuring a 1:2
non-linear resonance when an accurate (but not strictly complete) population
transfer is targeted. In both cases, we have incorporated the third-order Kerr
terms in the detuning locking the dynamics to the resonance at all times. In
the two-level system, we have shown the equivalence of the dynamics for the
optimal time or energy, given by a resonant and constant pulse, as it the case
for the linear problem. The behavior of the resonant nonlinear dynamics is
qualitatively different with the linear one: the complete inversion in only
(exponentially) asymptotic instead of the Rabi oscillations of the linear
problem (see Fig. 2). The optimal time features an asymptotic logarithmic
increasing as a function of the accuracy. For the three-level problem, the
optimal solution can be obtained only numerically. However, we have fitted it
using the resembling results of the two-level problem. In this case, the
generalized pulse area is constant. We have determined the shape of the
individual pulses featuring an intuitive pump-Stokes sequence, as in the
linear case but with different shapes. We have also obtained an asymptotic
logarithmic increasing of the optimal time as a function of the accuracy.

The finding of the ultimate bounds (time or energy) for nonlinear systems
provide an important benchmark. The issue of robustness of the process will
have to be considered in future analysis involving optimal inverse engineering
Dridi , since the (almost) complete transfer is very unstable when the
resonance is not perfectly satisfied (see Fig. 1 of Ref. pra_SG ). This is
also the case for two-level adiabatic transfer pra102_Jingjun . We notice that
this instability does not exist when only third-order nonlinearities apply,
neither for the usual nonlinear Raman three-level processes, where the target
state is linked with the Stokes coupling which is linear.

###### Acknowledgements.

We acknowledge support from the EUR-EIPHI Graduate School (17-EURE-0002) and
from the European Union’s Horizon 2020 research and innovation program under
the Marie Sklodowska-Curie Grant No. 765075 (LIMQUET). X.C. acknowledges EU
FET Open Grant EPIQUS (Grant No. 899368), the Basque Government through Grant
No. IT1470-22, the project grant PID2021-126273NB-I00 funded by
MCIN/AEI/10.13039/501100011033 and by “ERDF A way of making Europe” and “ERDF
Invest in your Future” and ayudas para contratos Ramon y Cajal–2015-2020
(RYC-2017-22482).

##  Appendix A Isomorphism between non-linear three-level and two-level models

The three-level problem (33) can be rewritten as

$\displaystyle\dot{x}_{1}=-Px_{1}z_{2},$ (57a)
$\displaystyle\dot{z}_{2}=Sz_{3}+Px^{2}_{1},$ (57b)
$\displaystyle\dot{z}_{3}=-Sz_{2},$ (57c)

with $z_{2}=-y_{2}\sqrt{2}$, $z_{3}=x_{3}\sqrt{2}$, $P=\Omega_{p}/\sqrt{2}$,
$S=\Omega_{s}/\sqrt{2}$, and the normalization condition
$x_{1}^{2}+z_{2}^{2}+z_{3}^{2}=1$. We assume that $P$ and $S$ are real and
that $x_{1}(t_{i})=1$. It can be reinterpreted as a density “matrix”
formulation on the Bloch sphere:

$\frac{d}{dt}\left[\begin{array}[]{ccc}\rho_{z}\\\ \rho_{y}\\\
\rho_{x}\end{array}\right]=\left[\begin{array}[]{ccc}0&-P\rho_{z}&0\\\
P\rho_{z}&0&S\\\ 0&-S&0\end{array}\right]\left[\begin{array}[]{ccc}\rho_{z}\\\
\rho_{y}\\\ \rho_{x}\end{array}\right]$ (58)

with $\rho_{z}=x_{1}=\rho_{11}-\rho_{22}=|a_{1}|^{2}-|a_{2}|^{2}$,
$\rho_{y}=z_{2}=i(\rho_{21}-\rho_{12})=2\Im(a_{1}\bar{a}_{2})$,
$\rho_{x}=z_{3}=\rho_{21}+\rho_{12}=2\Re(a_{1}\bar{a}_{2})$,
$\rho_{ij}=a_{i}\bar{a}_{j}$ of a non-linear two-level problem

$H_{2,n\ell}=\frac{1}{2}\left[\begin{array}[]{cc}-S&P(|a_{1}|^{2}-|a_{2}|^{2})\\\
P(|a_{1}|^{2}-|a_{2}|^{2})&S\end{array}\right]$ (59)

with $i\frac{d}{dt}\left[a_{1}\ a_{2}\right]^{t}=H_{2,n\ell}\left[a_{1}\
a_{2}\right]^{t}$ and the normalization $|a_{1}|^{2}+|a_{2}|^{2}=1$. As a
consequence, the non-linear three-level problem (57) is isomorphic to the
above non-linear two-state problem (59). The non-linearity which appears here
is not the one usually encountered (1). We obtain a similar isomorphic
relation for linear problems ShoreVitanov .

As a consequence, the transfer is, at the final time $t_{f}$, complete,
$|z_{3}(t_{f})|=1$, i.e. $\rho_{x}(t_{f})=\pm 1$, when, in the counterpart
two-state problem, from the initial state $\rho_{11}(t_{i})=1$, the
superposition of state of maximal coherence is produced:
$\rho_{12}(t_{f})=\rho_{21}(t_{f})=\pm 1/2$. This shows a similar qualitative
behavior than for its linear analog.

The general solution of the two-state problem (59) can be parametrized by
three angles in general:

$\left[\begin{array}[]{cc}a_{1}\\\
a_{2}\end{array}\right]=\left[\begin{array}[]{cc}\cos(\theta/2)\\\
\sin(\theta/2)e^{-i\varphi}\end{array}\right]e^{-i\gamma},$ (60)

and the Schrödinger equation leads to the set of equations

$\displaystyle\dot{\theta}=P\cos\theta\sin\varphi,$ (61a)
$\displaystyle\dot{\varphi}=S+P\frac{\cos^{2}\theta}{\sin\theta}\cos\varphi,$
(61b)
$\displaystyle\dot{\gamma}=-\frac{S}{2}+\frac{P\cos\theta\tan(\theta/2)\cos\varphi}{2}.$
(61c)

One can solve Eq. (61a) exactly [for any $P(t)$ and $S(t)$]:

$\tan(\theta/2)=\tanh\left[\frac{1}{2}\int_{t_{i}}^{t}P(s)\sin\varphi(s)ds\right].$
(62)

This shows that, in order to have a complete transfer from state 1 to state 3
in the original model, i.e. $\theta(t_{i})=0$ and $\theta(t_{f})=\pi/4$, one
needs an infinite pulse area of $P(t)$. This contrasts with the linear model
where a complete population transfer is possible for finite pulse areas
Boscain .

## References

  * (1) M.H. Levitt, Composite pulses, Prog. Nucl. Magn. Reson. Spectrosc. 18, 61 (1986).
  * (2) S. Wimperis, Broadband, Narrowband, and Passband Composite Pulses for Use in Advanced NMR Experiments, J. Magn. Reson. 109, 221 (1994).
  * (3) B.T. Torosov, S. Guérin, and N.V. Vitanov, High-Fidelity Adiabatic Passage by Composite Sequences of Chirped Pulses, Phys. Rev. Lett. 106, 233001 (2011).
  * (4) J.A. Jones, Designing short robust NOT gates for quantum computation, Phys. Rev. A 87, 052317 (2013).
  * (5) S.S. Ivanov, B.T. Torosov, and N.V. Vitanov, High-Fidelity Quantum Control by Polychromatic Pulse Trains, Phys. Rev. Lett. 129, 240505 (2022).
  * (6) Z. Shi, C. Zhang, D. Ran, Y. Xia, R. Ianconescu, A. Friedman, X.X. Yi, and S. Zheng, Composite pulses for high fidelity population transfer in three-level systems, New J. Phys. 24, 023014 (2022).
  * (7) N.V. Vitanov, A.A. Rangelov, B.W. Shore, and K. Bergmann, Stimulated Raman adiabatic passage in physics, chemistry, and beyond, Rev. Mod. Phys. 89, 015006 (2017).
  * (8) K. Bergmann et al., Roadmap on STIRAP applications, J. Phys. B 52, 202001 (2019).
  * (9) G. Dridi, S. Guérin, V. Hakobyan, H.R. Jauslin, and H. Eleuch, Ultrafast stimulated Raman parallel adiabatic passage by shaped pulses, Phys. Rev. A 80, 043408 (2009).
  * (10) X. Chen, I. Lizuain, A. Ruschhaupt, D. Guéry-Odelin, and J.G. Muga, Shortcut to Adiabatic Passage in Two- and Three-Level Atoms, Phys. Rev. Lett. 105, 123003 (2010).
  * (11) A Ruschhaupt, X. Chen, D. Alonso and J. G Muga, N. J. Phys. 14, 093040 (2012).
  * (12) D. Guéry-Odelin, A. Ruschhaupt, A. Kiely, E. Torrontegui, S. Martínez-Garaot, J.G. Muga, Shortcuts to adiabaticity: Concepts, methods, and applications, Rev. Mod. Phys. 91, 045001 (2019).
  * (13) J. Zhu and X. Chen, Fast-forward scaling of atom-molecule conversion in Bose-Einstein condensates, Phys. Rev. A. 103, 023307 (2021).
  * (14) D. Daems, A. Ruschhaupt, D. Sugny, and S. Guérin, Robust Quantum Control by a Single-Shot Shaped Pulse, Phys. Rev. Lett. 111, 050404 (2013).
  * (15) L. Van-Damme, D. Schraft, G.T. Genov, D. Sugny, T. Halfmann, and S. Guérin, Robust not gate by single-shot-shaped pulses, Phys. Rev. A 96, 022309 (2017).
  * (16) X. Laforgue, Xi Chen, and S. Guérin, Robust stimulated Raman exact passage using shaped pulses, Phys. Rev. A. 100, 023415 (2019).
  * (17) C.P. Koch, U. Boscain, T. Calarco, G. Dirr, S. Filipp, S.J. Glaser, R. Kosloff, S. Montangero, T. Schulte-Herbrüggen, D. Sugny, F.K. Wilhelm, Quantum optimal control in quantum technologies. Strategic report on current status, visions and goals for research in Europe, EPJ Quantum Technol. 9, 19 (2022).
  * (18) T. Caneva, M. Murphy, T. Calarco, R. Fazio, S. Montangero, V. Giovannetti, and G. E. Santoro, Optimal Control at the Quantum Speed Limit, Phys. Rev. Lett. 103, 240501 (2009).
  * (19) G.C. Hegerfeldt, Driving at the Quantum Speed Limit: Optimal Control of a Two-Level System, Phys. Rev. Lett. 111, 260501 (2013).
  * (20) M.R. Frey, Quantum speed limits - primer, perspectives, and potential future directions, Quantum Inf Process 15, 3919 (2016).
  * (21) F-Q. Dou, J. Liu, and L-B. Fu, Fast quantum driving in two-level systems with interaction and nonlinear sweep, Phys. Rev. A 98, 022102 (2018).
  * (22) W. Zhu, J. Botina, and H. Rabitz, Rapidly convergent iteration methods for quantum optimal control of population, J. Chem. Phys. 108, 1953 (1998).
  * (23) W. Zhu and H. Rabitz, Noniterative algorithms for finding quantum optimal controls, J. Chem. Phys. 110, 7142 (1999).
  * (24) J. Ruths and J. Lia, A multidimensional pseudospectral method for optimal control of quantum ensembles, J. Chem. Phys. 134, 044128 (2011).
  * (25) M. Lapert, R. Tehini, G. Turinici, and D. Sugny, Monotonically convergent optimal control theory of quantum systems under a nonlinear interaction with the control field, Phys. Rev. A 78, 02340 (2008).
  * (26) A.I. Konnov and V.A. Krotov, Automation and Remote Control 60, 1427 (1999).
  * (27) N. Khaneja, T. Reiss, C. Kehlet, T. Schulte-Herbrueggen, and S.J. Glaser, Optimal control of coupled spin dynamics: Design of NMR pulse sequences by gradient ascent algorithms, Journal of Magnetic Resonance 172, 296 (2005).
  * (28) L.S. Pontryagin, V.G. Boltyanskii, R.V. Gamkrelidze, E.F. Mischenko, The Mathematical Theory of Optimal Processes (Wiley, New York, 1962).
  * (29) U. Boscain, G. Charlot, J.-P. Gauthier, S. Guérin, and H. R. Jauslin, Optimal control in laser-induced population transfer for two- and three-level quantum systems, J. Math. Phys. 43, 2107 (2002).
  * (30) L. Van Damme, Q. Ansel, S.J. Glaser, and D. Sugny, Robust optimal control of two-level quantum systems, Phys. Rev. A 95, 063403 (2017).
  * (31) U. Boscain, M. Sigalotti , and D. Sugny, Introduction to the Pontryagin Maximum Principle for Quantum Optimal Control, PRX Quantum 2, 030203 (2021).
  * (32) J. Zeng, C.H. Yang, A. S. Dzurak, and E. Barnes, Geometric formalism for constructing arbitrary single-qubit dynamically corrected gates, Phys. Rev. A 99, 052321 (2019).
  * (33) G. Dridi, K. Liu, and S. Guérin, Optimal robust quantum control by inverse geometric optimization, Phys. Rev. Lett. 125, 250403 (2020).
  * (34) X. Laforgue, G. Dridi, and S. Guérin, Optimal robust quantum control against pulse inhonomogeneities: Analytic solutions, Phys. Rev. A 106, 052608 (2022).
  * (35) P.D. Drummond and K.V. Kheruntsyan, Stimulated Raman adiabatic passage from an atomic to a molecular Bose-Einstein condensate, Phys. Rev. A. 65, 063619 (2002).
  * (36) V.I. Arnold, V.V. Kozlov, and A.I. Neishtadt, Mathematical Aspects of Classical and Celestial Mechanics (Springer, Berlin, 2006), 3rd ed.
  * (37) J. Henrard, The Adiabatic Invariant Theory and Applications, in Hamiltonian Dynamics - Theory and Applications, Lecture Notes in Mathematics, edited by G.Benettin, J. Henrard, S.B. Kuksin, A. Giorgilli (Springer-Verlag, Berlin, Heidelberg, 2005), Vol. 1861.
  * (38) A.P. Itin and S. Watanabe, Integrability, Stability, and Adiabaticity in Nonlinear Stimulated Raman Adiabatic Passage, Phys. Rev. Lett 99, 223903 (2007).
  * (39) S. Guérin, M. Gevorgyan, C. Leroy, H. R. Jauslin, and A. Ishkhanyan, Efficient adiabatic tracking of driven quantum nonlinear systems, Phys. Rev. A 88, 063622 (2013).
  * (40) M. Gevorgyan, S. Guérin, C. Leroy, A. Ishkhanyan, and H.R. Jauslin, Adiabatic tracking for photo- and magneto-association of Bose-Einstein condensates with Kerr nonlinearities, Eur. Phys. J. D 70, 253 (2016).
  * (41) J.J. Zhu, X. Chen, H.R. Jauslin, and S. Guérin, Robust control of unstable nonlinear quantum systems, Phys. Rev. A 102, 052203 (2020).
  * (42) V. Dorier, M. Gevorgyan, A. Ishkhanyan, C. Leroy, H.R. Jauslin, and S. Guérin, Nonlinear Stimulated Raman Exact Passage by Resonance-Locked Inverse Engineering, Phys. Rev. Lett. 119, 243902 (2017).
  * (43) T. Huang, B.A. Malomed, and X. Chen, Shortcuts to adiabaticity for an interacting Bose-Einstein condensate via exact solutions of the generalized Ermakov equation, Chaos 30, 053131 (2020).
  * (44) X. Chen, Y. Ban, and G.C. Hegerfeldt, Time-optimal quantum control of nonlinear two-level systems, Phys. Rev. A 94, 023624 (2016).
  * (45) K. Efstathiou, Metamorphoses of Hamiltonian Systems with Symmetries, Lecture Notes in Mathematics (Springer-Verlag, Berlin, Heidelberg, 2005), Vol. 1864.
  * (46) N.V. Vitanov and B.W. Shore, Stimulated Raman adiabatic passage in a two-state system, Phys. Rev. A 73, 053402 (2006).
  * (47) J.C. Lagarias, J.A. Reeds, M.H. Wright, and P.E. Wright, Convergence Properties of the Nelder-Mead Simplex Method in Low Dimensions, SIAM Journal of Optimization, 19, 112 (1998).
  * (48) T.-N. Xu, K. Liu, X. Chen, and S. Guérin, Invariant-based optimal composite stimulated Raman exact passage, J. Phys. B 52, 235501 (2019).