# Optimal Control of Quantum State Transfer under Control Constraints

Nahid Binandeh Dehaghani  A. Pedro Aguiar  Research Center for Systems and
Technologies (SYSTEC), Electrical and Computer Engineering Department, FEUP -
Faculty of Engineering, University of Porto, Rua Dr. Roberto Frias sn, i219,
4200-465 Porto, Portugal (e-mail<EMAIL_ADDRESS>pedro.aguiar@fe.up.pt).

###### Abstract

We consider a control constrained optimal control problem of quantum state
transformation from an initial given state to a desired target state
satisfying Liouville-von Neumann equation. The cost functional to be optimized
is viewed as a trade off between maximizing fidelity and minimizing energy
consumption. As a new approach, we drive the optimality conditions in the form
of Pontryagin Maximum Principle (PMP) for the related matrix-valued dynamics,
and next we present a time-discretized computational scheme to solve the
proposed optimal control problem by using an indirect method based on the PMP.
The algorithm is applied for a spin $\frac{1}{2}$ particle interacting with a
magnetic field.

††thanks: The authors acknowledge the support of FCT for the grant
2021.07608.BD, the ARISE Associated Laboratory, Ref. LA/P/0112/2020, and the
R$\&$D Unit SYSTEC-Base, Ref. UIDB/00147/2020, and Programmatic, Ref.
UIDP/00147/2020 funds, and also the support of projects SNAP, Ref.
NORTE-01-0145-FEDER-000085, and RELIABLE (PTDC/EEI-AUT/3522/2020) funded by
national funds through FCT/MCTES. The work has been done in the honor and
memory of Professor Fernando Lobo Pereira.

Quantum technology aspires to develop practical applications based on
properties of the systems obeying the laws of quantum mechanics. This
objective requires efficient manipulation of quantum objects in order to
obtain desired behaviors. Quantum control holds within a number of techniques
to obtain the time evolution of control parameters and enable useful
performance in applications ranging from quantum computing, (Palao and
Kosloff, 2002; Xu et al., 2020) to sensing, (Rembold et al., 2020; Poggiali et
al., 2018), simulation, (Holland et al., 2020), and metrology, (Lin et al.,
2021). Several quantum control methods, (Stefanatos and Paspalakis, 2021),
such as brute-force optimization of a few pulse parameters, (Cheng et al.,
2020), Brumer-Shapiro coherent control, (Gruebele, 2001), pulse-timing
control, (McDermott and Vavilov, 2014), stimulated-Raman-Adiabatic-Passage,
(Vitanov et al., 2017), genetic algorithms, (Lahoz-Beltra, 2016), and optimal
control theory (OCT), (Glaser et al., 2015; James, 2021; Werschnik, 2006),
have been exploited in order to discover the optimal pulse sequences. The
studies on Quantum Optimal Control (QOC) began in the late 1980’s, (Peirce et
al., 1988), and has undergone continuous developments up to now.

A wide range of problems arising in quantum technology, as in quantum
computing or nuclear magnetic resonance spectroscopy, are compatible to be
formulated in the framework of OCT. Amongst the most dominant advances of QOC,
we can point out the introduction of rapidly converging iterative algorithms,
(Maday and Turinici, 2003; Zhu and Rabitz, 1998), and its generalization for
dissipative systems, (Ohtsuki et al., 1999), while taking several control
criteria into account, (Ohtsuki et al., 2001). In fact, QOC aims at developing
an organized and rigorous design methodology in order to control the behavior
of quantum systems so that a desired set of objectives can be obtained in an
optimal way. A control field $u(t)$ is able to produce the global maximum or
minimum value of a performance index $J(u(t))$, e.g. maximum fidelity,
(Dehaghani and Pereira, 2022a) or minimum time, (Sugny et al., 2007), while
being able to overcome decoherence and dissipation.

Despite the considered research works that have been done in quantum optimal
control theory, Pontryagin maximum principle of optimal control is still far
from being fully exploited in quantum context. In the current state of the
art, the quantum systems to be controlled are usually simple closed systems,
in which the quantum state is expressed in the form of unit vectors, and
evolves according to the Schrödinger equation, (Boscain et al., 2021;
Werschnik and Gross, 2007; D’alessandro and Dahleh, 2001). However, in
practical applications, the quantum systems to be controlled are usually not
simple closed systems. They may be quantum ensembles, and their states cannot
be expressed in the form of unit vectors. In this paper, we address precisely
this problem by considering the evolution of the matrix-valued probability
density function, and describing the system dynamics in terms of a density
operator, by means of, e.g., master equations. As far as the authors’
knowledge, such a study has not been done, with exception of our previous
recent work in (Dehaghani and Pereira, 2022b). In this paper, we extend and
generalize our previous results by applying the Maximum Principle of
Pontryagin for the matrix valued quantum dynamic control system, where the
state of the quantum system is described through the density operator evolving
according to Liouville-von Neumann equation. We aim a trade off between the
goal of attaining the maximum fidelity in order to signify a security level
for quantum state transformation, which is of high importance notably in
quantum information theory, and keep the energy of the field small. The first-
order necessary optimality conditions obtained from the application of PMP
results in a two-point boundary value problem, which we managed to solve by
proposing a shooting algorithm, where we also consider the problem of control
constraints. Simulation results illustrate the effectiveness of the proposed
methodology.

The paper is organized as follows: We first review the set of PMP optimality
conditions for a general system. Then, we turn to the context of quantum
systems, and we study the basis of a quantum optimal control problem,
including the current state of the art. Next, we address the questions of
existence of an optimal control and also controllability of a quantum system.
Then, we present the physical description of a simple two-level quantum-
mechanical system, which we have used in our work. In the next section, we
describe the system under study, and, then, we formulate the optimal control
problem under control constraints and obtain the necessary conditions of
optimality. We later present the application of PMP by means of an indirect
method through an algorithm, including simulation results. The paper ends with
brief conclusions and an overview on prospective research challenges.

###  0.1 Notation.

For a general continuous-time trajectory $x$, the term $x(t)$ indicates the
trajectory assessed at a specific time $t$. For writing partial differential
equations (PDEs), we denote partial derivatives using subscripts. In the
general situation where $f$ denotes a function of $n$ variables including $x$,
then $f_{x}$ denotes the partial derivative relative to the $x$ input.
Throughout the paper, we have used $i$ as the imaginary unit. For a matrix
$A$, $A^{T}$, $A^{*}$, and $A^{\dagger}$ represent the transpose, complex
conjugate, and conjugate transpose of matrix A, respectively. To denote the
wavefunctions as vectors, we use the Dirac notation such that
$\left|\psi\right\rangle=\sum\limits_{j=1}^{n}{{{\alpha}_{j}}\left|{{{\hat{\psi}}}_{j}}\right\rangle}$,
where $|\psi\rangle$ indicates a state vector, ${\alpha}_{j}$ are the complex-
valued expansion coefficients, and $|{\hat{\psi}}_{j}\rangle$ are basis
vectors that are fixed.
We denote a finite-dimensional separable Hilbert space by $\mathbb{H}$, and
define it over the complex field $\mathbb{C}$, so
$\mathbb{H}\simeq\mathbb{C}^{N}$, where $N$ indicates the dimension of the
space. We consider the set $\mathcal{B}(\mathbb{H})$ as the set of linear
operators on the Hilbert space and define it as
$\mathcal{B}(\mathbb{H})\simeq\mathbb{C}^{N\times N}$. The inner product in
the set $\mathcal{B}(\mathbb{H})$ is the Hilbert-Schmidt inner product,
defined as $\left\langle A,B\right\rangle=tr\left({{A}^{\dagger}}B\right)$,
where $tr(A)$ indicates the trace of an square matrix $A$. The commutator of
two elements $A$ and $B$, as linear operators on the Hilbert space, is
indicated by $[A,B]:=AB-BA$.

##  1 Pontryagin’s maximum principle

In general terms, a fundamental optimal control problem is formulated as the
following (Werschnik, 2006):

Given the dynamical system
$\dot{x}\left(t\right)=f\left(x\left(t\right),u\left(t\right)\right)$ and a
set of admissible controls $u(t)\in\mathcal{U}$, we have to determine an
admissible control signal such that the objective functional

$J=\Phi(x(T))+\int\limits_{0}^{T}{L}\left(x\left(t\right),u\left(t\right)\right){dt}$
(1)

is minimized. Let $u^{*}(t)\in\mathcal{U}$ and $x^{*}(t)$ represent the
optimal control and state trajectory for the defined optimal control problem.
Then, there is a time-varying adjoint trajectory $\lambda(t)$ that together
with $u^{*}(t)$ and $x^{*}(t)$ satisfy

system equation and initial state condition (2)
$\displaystyle\dot{x}^{*}\left(t\right)=f\left(x^{*}\left(t\right),u^{*}\left(t\right)\right)$
$\displaystyle x^{*}\left(t=0\right)={{x}_{0}}$ adjoint equation and
transversality condition (3)
$\displaystyle-{{{\dot{\lambda}}}^{\dagger}}\left(t\right)={{\lambda}^{\dagger}}\left(t\right){{f}_{x}}\left(x^{*}\left(t\right),u^{*}\left(t\right)\right)-{L_{x}}\left(x\left(t\right),u\left(t\right)\right)$
$\displaystyle{{\lambda}^{\dagger}}\left(T\right)=-{\Phi_{x}}\left(x^{*}\left(T\right)\right)$

For all $t\in[0,T]$ and $u(t)\in\mathcal{U}$

Maximum condition (4)
$\displaystyle\mathcal{H}\left(\lambda\left(t\right),x^{*}\left(t\right),u\left(t\right)\right)\leq\mathcal{H}\left(\lambda\left(t\right),x^{*}\left(t\right),u^{*}\left(t\right)\right)$

where $\mathcal{H}$ is the Pontryagin Hamiltonian defined as

$\mathcal{H}\left(\lambda\left(t\right),x\left(t\right),u\left(t\right)\right):={{\lambda}^{\dagger}}\left(t\right)f\left(x\left(t\right),u\left(t\right)\right)-L\left(x\left(t\right),u\left(t\right)\right)$
(5)

##  2 An overview on quantum optimal control

Quantum Optimal Control, (Glaser et al., 2015), can intuitively be formulated
in the above-mentioned setting. For such problems, the state of the system may
be described by the pure state vector, density operator (for both mixed and
pure quantum states), or we can consider the dynamics of the evolution
operator. The evolution of a pure state, which is not entangled with the
environment, can be described by a wave function $\left|\psi(t)\right\rangle$,
which evolves in time according to a control-dependent Schrödinger equation,
(Boscain et al., 2021; Cong, 2014; Werschnik, 2006),

$i\hbar{\left|\dot{\psi}\left(t\right)\right\rangle}=H\left(u\left(t\right)\right)\left|\psi\left(t\right)\right\rangle,\quad\left|\psi(t=0)\right\rangle=\left|{{\psi}_{0}}\right\rangle$
(6)

where ${H(u(t))}$ is the quantum-mechanical Hamiltonian of the system, and
$\hbar$ is the reduced Planck constant, usually set as $\hbar=1$ for
convenience. The system control can be realized by an admissible set of
external control signals ${{u}_{k}}(t)\in\mathbb{R}$, which are coupled to the
quantum system via time independent interaction Hamiltonians ${{H}_{k}}$.
Therefore, the total quantum Hamiltonian defined as

$H\left(u\left(t\right)\right)={{H}_{0}}+\sum\limits_{k=1}^{m}{{{u}_{k}}(t)}{{H}_{k}}$
(7)

determines the controlled evolution, in which ${{H}_{0}}$ indicates the time
independent internal (free) Hamiltonian. In (6), both the system state and
Hamiltonian are complex quantities. However, we can express the problem in
terms of only real quantities by introducing
$x:={{\left[\vec{\psi}_{R}^{T},\vec{\psi}_{I}^{T}\right]}^{T}}$, where
$\left|\psi_{R}(t)\right\rangle,\left|\psi_{I}(t)\right\rangle\in\mathbb{R}^{n}$,
and separating the real and imaginary parts of
$-iH\left(u(t)\right)=R\left(u(t)\right)+iI\left(u(t)\right)$, where
$R\left(u(t)\right),I\left(u(t)\right)\in{{\mathbb{R}}^{n\times n}}$ are skew-
symmetric and symmetric matrices for all values of $u(t)$, respectively. Now,
we can rewrite the differential equation describing the dynamics of the system
implying only real values by

$\dot{x}=\tilde{H}\left(u(t)\right)x$ (8)

in which

$\tilde{H}\left(u(t)\right)=\left(\begin{matrix}R\left(u(t)\right)&-I\left(u(t)\right)\\\
I\left(u(t)\right)&R\left(u(t)\right)\\\ \end{matrix}\right)$ (9)

is both symplectic and skew-symmetric for all values of $u(t)$. The cost in
(1) can also be rewritten by introducing appropriate functions $\tilde{\Phi}$
and $\tilde{L}$. Hence, the Pontryagin Hamiltonian takes the form

$\mathcal{H}\left(\lambda\left(t\right),x\left(t\right),u\left(t\right)\right)={{\lambda}^{T}}\left(t\right)\tilde{H}\left(u\left(t\right)\right)x\\!-\\!\tilde{L}\left(x\left(t\right),u\left(t\right)\right)$
(10)

from which the optimal control $u(t)$ has to satisfy the maximum condition
(4). Since $\tilde{H}\left(u(t)\right)$ is skew-symmetric, (3) can be
rewritten as (Dehaghani and Pereira, 2022c),

$\displaystyle\dot{\lambda}\left(t\right)=\tilde{H}\left(u\right)\lambda\left(t\right)+\tilde{L}_{x}^{T}$
(11)
$\displaystyle\lambda\left(T\right)=-\tilde{\Phi}_{x}^{T}\left(x(T)\right)$

##  3 Existence of the Optimal Control

Since we use the necessary optimality conditions, we also need to address the
problem of existence of the optimal control. There are already some standard
results, well explained in (Fleming and Rishel, 2012) Chapter III, that are
applicable to the case considered in our work. According to (Fleming and
Rishel, 2012), the existence theorem guarantee the existence of an optimal
control in the Lebesgue Integrable ($LI$) set if a control steering a given
initial state to a desired target state exists. The control functions
$u_{k}(t)$ in (8) are supposed to be $LI$ as well. We also will suppose in
this paper that the optimal control are in the set $LI$, according to the
above-mentioned theorem. In order to guarantee the existence of a control, one
needs to address the question of controllability. The controllability concerns
the possibility of steering the system from one state to another for every
pair of states.

The complex sphere $\mathbb{S}^{2N-1}\subset\mathbb{H}$, representing the pure
quantum states, is a homogeneous space of the Lie group
$U\left(N\right)=\left\\{U\in
GL\left(N,\mathbb{C}\right)\left|U{{U}^{\dagger}}\right.={{U}^{\dagger}}U=I\right\\}$,
and its proper subgroup $SU\left(N\right)=U\left(N\right)/U\left(1\right)$.
The Lie algebras of $U(N)$ and $SU(N)$ are

$u\left(N\right)=\left\\{A\in{{C}^{N\times
N}}\left|{{A}^{\dagger}}=-A\right.\right\\}$ (12)

and

$su\left(N\right)=\left\\{A\in
u\left(N\right)\left|tr\left(A\right)=0\right.\right\\},$ (13)

respectively. The Schrödinger equation stated in (6) can be lifted to the Lie
group $SU(N)$ to obtain the Schrödinger equation for the unitary propagator as

$i\frac{d}{dt}U\left(t\right)=\left({{H}_{0}}+\sum\limits_{k=1}^{m}{{{u}_{k}}\left(t\right){{H}_{k}}}\right)U\left(t\right),\quad
U(0)=I$ (14)

where $U\in SU(N)$. Equation (14) is a right invariant control system on the
compact Lie group $su(n)$, so

$Lie\left\\{i{{H}_{0}},\ldots,i{{H}_{m}}\right\\}=su\left(n\right)$ (15)

which represents a necessary and sufficient condition for the controllability
of the system indicated in (14), (Boothby and Wilson, 1979). Overall, if a
right invariant system is controllable, then the bilinear system is also
controllable, (Sachkov, 1997). Therefore, if the evolution operator of the
system satisfies a controllable equation, then the bilinear system described
by the certain quantum-mechanical Hamiltonian $H\left(u\left(t\right)\right)$
is controllable as well, so it is possible to design controls to steer the
bilinear system from one state to another state in the state-space.
Consequently, the system indicated in (6) is controllable provided that (15)
holds. The controllability of two-level quantum systems, as the case
considered in this paper, is analysed in more details in (D’Alessandro, 2000).

##  4 Physical Description of a Control System

In Nuclear Magnetic Resonance (NMR) experiments, a single spin $\frac{1}{2}$
particle is controlled by means of an electromagnetic field
$\vec{B}\left(t\right)$, in which one component is kept constant in $z$
direction while $x$ ans $y$ components vary in time in order to change the
direction of the spin, (D’Alessandro, 2021). The total quantum mechanical
Hamiltonian is then described by the interaction of the external magnetic
field
$\vec{B}\left(t\right)=\left({{B}_{x}}\left(t\right),{{B}_{y}}\left(t\right),{B}_{z}\right)^{T}$
with the spin angular momentum
$\hat{S}:=\left({{{\hat{S}}}_{x}},{{{\hat{S}}}_{y}},{{{\hat{S}}}_{z}}\right)^{T}$
as

$\displaystyle H\left(t\right)$
$\displaystyle=\gamma\hat{S}^{T}\vec{B}\left(t\right)$ (16)
$\displaystyle=\gamma\left({{{\hat{S}}}_{x}}{{B}_{x}}\left(t\right)+{{{\hat{S}}}_{y}}{{B}_{y}}\left(t\right)+{{{\hat{S}}}_{z}}{{B}_{z}}\right)$

where $\gamma$ is the gyromagnetic ratio. Let express the state
$\left|\psi\left(t\right)\right\rangle={{c}_{1}}\left(t\right)\left|\frac{1}{2}\right\rangle+{{c}_{2}}\left(t\right)\left|-\frac{1}{2}\right\rangle$.
From (6) and (16), the differential equation for
$\vec{c}\left(t\right)={{\left({{c}_{1}}\left(t\right),{{c}_{2}}\left(t\right)\right)}^{T}}$
is expressed as

$i\frac{d}{dt}\vec{c}=\frac{\gamma}{2}\left({{\sigma}_{x}}{{B}_{x}}\left(t\right)+{{\sigma}_{y}}{{B}_{y}}\left(t\right)+{{\sigma}_{z}}{{B}_{z}}\right)\vec{c}$
(17)

in which the matrix representation of the operators $\sigma_{x}$,
$\sigma_{y}$, and $\sigma_{z}$ is done by the so-called Pauli matrices. By
appropriate scaling of time and setting the magnetic field arguments as
controls we have

$-iH\left(u\left(t\right)\right)=\left({{{\bar{\sigma}}}_{z}}{{u}_{z}}+{{{\bar{\sigma}}}_{x}}{{u}_{x}}\left(t\right)+{{{\bar{\sigma}}}_{y}}{{u}_{y}}\left(t\right)\right)$
(18)

where

$\\!\\!{{\bar{\sigma}}_{x}}=\frac{1}{2}\left(\begin{matrix}0&i\\\ i&0\\\
\end{matrix}\right),\\!\\!\quad{{\bar{\sigma}}_{y}}=\frac{1}{2}\left(\begin{matrix}0&-1\\\
1&0\\\
\end{matrix}\right),\\!\\!\quad{{\bar{\sigma}}_{z}}=\frac{1}{2}\left(\begin{matrix}i&0\\\
0&-i\\\ \end{matrix}\right)\\!$

span $SU(2)$ of skew-Hermitian matrices with trace equal to zero, and satisfy
the following commutation relations

$\left[{{{\bar{\sigma}}}_{x}},{{{\bar{\sigma}}}_{y}}\right]={{\bar{\sigma}}_{z}},\quad\left[{{{\bar{\sigma}}}_{y}},{{{\bar{\sigma}}}_{z}}\right]={{\bar{\sigma}}_{x}},\quad\left[{{{\bar{\sigma}}}_{z}},{{{\bar{\sigma}}}_{x}}\right]={{\bar{\sigma}}_{y}}\cdot$
(19)

Equation (20) can be implemented according to (8), to obtain

$\tilde{H}\left(u\left(t\right)\right)={{T}_{z}}{{u}_{z}}+{{T}_{y}}{{u}_{y}}\left(t\right)+{{T}_{x}}{{u}_{x}}\left(t\right)$
(20)

where

$\displaystyle{{T}_{x}}=\frac{1}{2}\left(\begin{matrix}0&0&0&-1\\\ 0&0&-1&0\\\
0&1&0&0\\\ 1&0&0&0\\\
\end{matrix}\right),{{T}_{y}}=\frac{1}{2}\left(\begin{matrix}0&-1&0&0\\\
1&0&0&0\\\ 0&0&0&-1\\\ 0&0&1&0\\\ \end{matrix}\right),$ (21)
$\displaystyle{{T}_{z}}=\frac{1}{2}\left(\begin{matrix}0&0&-1&0\\\ 0&0&0&1\\\
1&0&0&0\\\ 0&-1&0&0\\\ \end{matrix}\right)\cdot$

From now on, we show the two components of the time varying control by the
vector
$\mathbf{u}(t)=\left[\begin{matrix}{{u}_{x}}\left(t\right)&{{u}_{y}}\left(t\right)\\\
\end{matrix}\right]^{T}$.

##  5 System Description

Every state vector $\left|{{\psi}_{j}}\left(t\right)\right\rangle$ at time $t$
can be obtained by

$\left|{{\psi}_{j}}\left(t\right)\right\rangle=U(t)\left|{{\psi}_{j}}\left(0\right)\right\rangle$
(22)

in which $U(t)$ is the solution of (14) expressed as

$U\left(t\right)=\exp\left(-i\int\limits_{0}^{t}{H\left(s\right)ds}\right)$
(23)

called Dyson’s series in the physics context and is similar to a Volterra
series in control theory, (Altafini and Ticozzi, 2012). Hence, the outer
product
$\left|{{\psi}_{j}}\left(t\right)\right\rangle\left\langle{{\psi}_{j}}\left(t\right)\right|$
results in

$\left|{{\psi}_{j}}\left(t\right)\right\rangle\left\langle{{\psi}_{j}}\left(t\right)\right|=U\left(t\right)\left|{{\psi}_{j}}\left(0\right)\right\rangle\left\langle{{\psi}_{j}}\left(0\right)\right|{{U}^{\dagger}}\left(t\right)$
(24)

and the same for any convex sum of
$\left|{{\psi}_{j}}\left(t\right)\right\rangle\left\langle{{\psi}_{j}}\left(t\right)\right|$.
More precisely, let $p_{j}$ be the fraction of population of an ensemble
$\left\\{\left.{{p}_{j}},\left|{{\psi}_{j}}\left(t\right)\right\rangle\right\\}\right.$,
so the corresponding quantum density operator is expressed as

$\rho=\sum\limits_{j}{{{p}_{j}}}\left|{{\psi}_{j}}\right\rangle\left\langle{{\psi}_{j}}\right|,\quad{{p}_{j}}\geq
0,\quad\sum\limits_{j}{{{p}_{j}}}=1$ (25)

which belongs to the set of Hermitian, semi definite, and positive matrices
with trace equal to one on the system’s Hilbert space $\mathbb{H}$. Hence, we
can rewrite (25) as

$\rho\left(t\right)=U\left(t\right)\rho\left(0\right){{U}^{\dagger}}\left(t\right)$
(26)

The infinitesimal version of (26) is the quantum Liouville-von Neumann
equation, (Berman and Kosloff, 1991), expressed by

$\dot{\rho}\left(t\right)=-i\left[H\left(u\left(t\right)\right),\rho\left(t\right)\right],\quad\rho(0)=\rho_{0}$
(27)

The controllability of (27) can consequently be obtained from the necessary
and sufficient condition indicated in (15). The main characteristic of the
Liouville-von Neumann equation is that it generates isospectral evolutions,
meaning that

$sp\left(\rho\left(t\right)\right)=sp\left(\rho\left(0\right)\right)=\Phi\left(\rho\right)=\left\\{{{\mu}_{1}},\ldots,{{\mu}_{N}}\right\\}$
(28)

where $\mu_{1},\dots,\mu_{N}$ are the eigenvalues of $\rho(t)$. As a
consequence of the isospectrality of (27), the set $\Phi\left(\rho\right)$
forms a complete set of constants of motion of (27). Let consider the set
$\mathcal{D}\left(\mathbb{H}\right)=\left\\{\rho\in\mathcal{B}\left(\mathbb{H})\right)\left|\rho={{\rho}^{\dagger}}\geq
0,tr\left(\rho\right)=1\right.\right\\}$, which is foliated into leaves
uniquely determined through the set $\Phi(\rho)$, and
$\mathcal{C}\subset\mathcal{D}\left(\mathbb{H}\right)$ as such as leaf, where
$\rho_{0}\in\mathcal{C}$. Hence,
$\mathcal{C}=\left\\{U{{\rho}_{0}}{{U}^{\dagger}},U\in
SU\left(N\right)\right\\}$ corresponds to the orbit of $SU(N)$ under the
action of conjugation passing through the initial density operator. Assume
$j_{1},\dots,j_{l}$ as the geometric multiplicities of the eigenvalues
$\Phi(\rho_{0})$, and $j_{1}+\dots+j_{l}=N$, $2\leq l\leq N$, then
$\mathcal{C}$ is the homogeneous space as
$C=U\left(N\right)/\left(U\left({{j}_{1}}\right)\times\ldots\times
U\left({{j}_{l}}\right)\right)$. As the eigenvalues of the density operator
vary, the geometric multiplicities $j_{1},\dots,j_{l}$ form a flag, so the
$\mathcal{C}$ are called complex flag manifolds, (Bengtsson and Życzkowski,
2017). The flag also determines the dimension of flag manifolds $\mathcal{C}$,
varying from $2N-2$, for pure states, to $N^{2}-N$, for all various
eigenvalues, (Altafini and Ticozzi, 2012). Equation (27) will be adopted as
the investigated model in the following subsections.

##  6 Formulation of the Optimal Control Problem under Control Constraints

In optimal control of quantum state transfer problems, a possible cost in the
form of (1) to be minimized is to consider a trade off between the goal of
attaining the maximum fidelity in order to signify a security level for
quantum state transformation, while simultaneously keeping the energy of the
field small. Motivated by this consideration, we propose the following optimal
control problem $(P_{1})$

Minimize $\displaystyle
J=-\mathcal{F}(\rho(T),\sigma)+\eta\int\limits_{0}^{T}{\mathbf{u}^{T}(t)\mathbf{u}(t)}dt$
$\displaystyle t\in[0,T],\quad\eta\in[0,1]$ subject to
$\displaystyle\textit{dynamics}:\dot{\rho}(t)=F\left(\rho(t),u(t)\right),$
$\displaystyle\qquad\qquad\quad\rho\left(t\right)\in{{\mathbb{C}}^{n\times
n}}$ $\displaystyle\textit{initial condition}:\rho\left(0\right)={{\rho}_{0}}$
$\displaystyle\textit{control constraint}:u(\cdot)\in\mathcal{U},$
$\displaystyle\qquad\qquad\textit{i.e.,
}u:[0,T]\to\Omega\subset\mathbb{R}^{2}$

The coefficient $\eta\geq 0$ is considered to signify the importance of energy
minimization. The density operator $\rho(t)$ is the quantum state variable
supposed to satisfy the differential constraints according to the Liouville-
von Neumann equation (27), and $\rho_{0}$ is the so-called initial quantum
state. The set of admissible controls $\Omega$ is defined as

$\Omega=\left\\{\left[\begin{matrix}{{u}_{x}}(t)&{{u}_{y}}(t)\\\
\end{matrix}\right]{{\left|\left\|u\right\|\right.}^{2}}\leq\sqrt{2}{{u}_{\max}}\right\\}$

where ${{u}_{\max}}$ is a given positive parameter. In P1, we have used the
well-known Uhlmann-Jozsa definition for fidelity, (Jozsa, 1994), representing
the maximal transition probability between the purification of a pair of
density matrices, $\rho(T)$ and the desired target state $\sigma$, (Liang et
al., 2019), defined as

$\mathcal{F}\left(\rho(T),\sigma\right):={{\left(tr\sqrt{\sqrt{\rho(T)}\sigma\sqrt{\rho(T)}}\right)}^{2}}\cdot$
(29)

##  7 Necessary Conditions of Optimality in the form of a Maximum Principle

The Pontryagin-Hamilton function ${\mathcal{H}}$ is defined for almost all
$t\in[0,T]$ by introducing the matrix-valued time-varying multiplier
$\Lambda$, designated by costate or adjoint variable of the system, (Dehaghani
and Pereira, 2022b). Thus,

${\mathcal{H}}\left(\rho,u,\Lambda\right)=tr\left(\Lambda^{\dagger}F\left(\rho,u\right)\right)-L(u(t))$
(30)

According to the Pontryagin’s Maximum Principle (4), for the optimal state
trajectory $\rho^{*}$ and the corresponding adjoint variable $\Lambda$, the
optimal control $u^{*}(t)$ maximizes the Pontryagin-Hamiltonian function
${\mathcal{H}}$ for almost all $t\in[0,T]$ and all admissible control values
$u\in\Omega$ such that

${\mathcal{H}}\left(\rho^{*}(t),u(t),\Lambda(t)\right)\leq{\mathcal{H}}\left(\rho^{*}(t),u^{*}(t),\Lambda(t)\right)$
(31)

Consequently, the adjoint equation implies that

$\displaystyle-\dot{\Lambda}^{\dagger}(t)$
$\displaystyle={\mathcal{H}}_{\rho}(\rho^{*}(t),u^{*}(t),\Lambda(t))$ (32)
$\displaystyle=i\left[H(u^{*}(t)),\Lambda^{\dagger}(t)\right]$

in which $H^{*}(t)$ denotes the quantum mechanical Hamiltonian evaluated at
each time along the optimal control $u^{*}$. Equation (32) has the formal
solution

$\Lambda^{\dagger}(t)=e^{i\int_{t}^{T}H^{*}(s)ds}\Lambda^{\dagger}(T)e^{-i\int_{t}^{T}H^{*}(s)ds}\cdot$
(33)

The boundary condition at the final time for adjoint variable implies that

$\displaystyle\Lambda^{\dagger}(T)$
$\displaystyle=\nabla_{\rho}\left(tr\sqrt{\sqrt{\rho(T)}\sigma\sqrt{\rho(T)}}\right)^{2}$
(34)
$\displaystyle=2tr\sqrt{\rho(T)\sigma}\sum\limits_{k=0}^{n-1}{\alpha}_{k}\sum_{i=0}^{k-1}\bar{\rho}(T)^{i}\sqrt{\sigma}\bar{\rho}(T)^{k-i-1}$

in which $n$ indicates the dimension of the density matrix and
$\bar{\rho}(T)=\rho(T)-I$. The derivation of (34) can be find in (Dehaghani
and Pereira, 2022b), where the coefficients $\alpha_{k}$, $k=0,\ldots,n-1$ are
obtained from the application of the Cayley–Hamilton theorem.

##  8 Application of the Pontryagin Maximum Principle

In this section, we solve the optimal control problem defined in $(P_{1})$ by
means of an indirect method based on the PMP. In the presented algorithm, we
have discretized the time interval $[0,T]$ into $N$ sub intervals such that
$\displaystyle t_{k}=\frac{k}{N}$ for $k=0,\ldots,N-1$, and setting
$j=0,\ldots$ as the iterations counter. Hence, $j^{th}$ iteration of the
function $f$ at time $t_{k}$ is represented by $f_{k}^{j}$. Here, both the
dynamics of the system and adjoint equation are considered by a first order
Euler approximation. The proposed algorithm is as follows:

  * Step 1 - Initialization.

Initialize the values of ${u_{x}}_{k}^{j}$ and ${u_{y}}_{k}^{j}$ for
$k=0,\ldots,N-1$, and go to the first iteration (j=1).

  * Step 2 - Computation of the state trajectory.

For $k=0,\dots,N-1$, compute
$U_{k}^{j}={{e}^{\sum\limits_{l=0}^{k-1}{\frac{T}{N}\tilde{H}_{l}^{j}}}}$,
where $\tilde{H}_{k}^{j}$ is the implemented Hamiltonian of the system
dynamics associated to the control $u$ at time $t_{k}$ and $(j-1)th$
iteration. Obtain

$\rho_{k}^{j}=U_{k}^{j}\rho_{0}^{j}U_{k}^{j\dagger}$

  * Step 3 - Computation of the adjoint trajectory

Compute $\Lambda_{N}^{j}$ by using $\rho_{N}^{j}$ computed in Step 2.

Compute $\Lambda_{k}^{j}$ by using the discretized version of (33), that is to
compute $V_{k}^{j}={e}^{\sum\limits_{k}^{N-1}{\frac{T}{N}\tilde{H}_{k}^{j}}}$,
and obtain

$\Lambda_{k}^{j\dagger}=V_{k}^{j\dagger}\Lambda_{N}^{j\dagger}V_{k}^{j}$

  * Step 4 - Computation of the Pontryagin Hamilton function.

For $k=0,\dots,N-1$, let

$\displaystyle\mathcal{H}_{k}^{j}(u_{k}^{j})=tr\left({\Lambda^{j}_{k}}^{\dagger}F_{k}^{j}(\rho_{k}^{j},u_{k}^{j})\right)-L^{j}_{k}$

  * Step 5: Computation of the control function.

For $k=0,\ldots,N-1$, compute the temporary control values
${u_{temp,x}}_{k}^{j}$, and ${u_{temp,y}}_{k}^{j}$ that maximize the map

$({u_{x}},{u_{y}})\to\mathcal{H}_{k}^{j}({u_{x}},{u_{y}}).$

  * Step 6: Apply the control constraints.

${u_{x,}}_{k}^{j}=\min(|{u_{temp,x,}}_{k}^{j}|,u_{max})\text{sign}({u_{temp,x}})$

${u_{y,}}_{k}^{j}=\min(|{u_{temp,y,}}_{k}^{j}|,u_{max})\text{sign}({u_{temp,y}})$

  * Step 7: Stopping test.

For a determined tolerance error $\varepsilon>0$, check the algorithm
convergence by verifying if

$\displaystyle\max_{k=0,\ldots,N-1}\\{|{u_{x}}_{k}^{j}-{u_{x}}_{k}^{j-1}|\\}<\varepsilon$
$\displaystyle\max_{k=0,\ldots,N-1}\\{|{u_{y}}_{k}^{j}-{u_{y}}_{k}^{j-1}|\\}<\varepsilon$

holds true. If yes, let ${u_{x}}^{*}(t_{k})={u_{x}}_{k}^{j}$ and
${u_{y}}^{*}(t_{k})={u_{y}}_{k}^{j}$ for $k=0,\ldots,N-1$, and exit the
algorithm.

Otherwise, update the temporary values for the control function according to

${u_{temp,x}}_{k}^{j}={u_{x}}_{k}^{j-1}+\delta({u_{x}}_{k}^{j}-{u_{x}}_{k}^{j-1})$
${u_{temp,y}}_{k}^{j}={u_{y}}_{k}^{j-1}+\delta({u_{y}}_{k}^{j}-{u_{y}}_{k}^{j-1})$

where $\delta>0$ is the learning rate coefficient.

Check the control constraints according to Step 6.
Then, let $j=j+1$, go to Step 2.

Here, we check the algorithm convergence by verifying whether the control
functions obtained in the current iteration approach the ones of the previous
iteration with an acceptable tolerance or not. If not, we repeat the above
steps until we achieve the desired convergence.

##  9 Simulation Results

Let consider the quantum state transfer problem given the initial state
$\rho_{0}=\left(\begin{matrix}1&0\\\ 0&0\\\ \end{matrix}\right)$ and the
desired target state $\sigma=\left(\begin{matrix}0&0\\\ 0&1\\\
\end{matrix}\right)$. By means of implementing the algorithm explained in the
previous section, we solve the quantum state transfer problem $P_{1}$ to drive
the given initial quantum state to a desired set while simultaneously
minimizing the required control power with the factor $\eta=0.01$. In this
simulation, the control is constrained to $u_{x},u_{y}\in[-1.5,1.5]$, and the
algorithm has been considered to evolve in $t\in[0,1.6]$ to accomplish the
desired fidelity. The control signal in $z$ direction is considered constant
as $u_{z}=0.001$, and in $x$ and $y$ directions is initialized by
${u_{x}}_{k}^{j}=0$ and ${u_{y}}_{k}^{j}=0$ for $k=0,\ldots,N-1$. By setting
the time-slicing $N=20$, the learning coefficient $\delta=0.1$, and the
stopping threshold $\varepsilon=10^{-3}$, the algorithm converges after 114
iterations. The residual graph showing the convergence of the algorithm is
displayed in Fig. 1.

Figure 1: Residuals: The algorithm converges after 114 iterations.

In Fig. 2 and Fig. 3, the evolution of the 4 elements of the density matrix
during the time in the last iteration is presented. The figures show the
imaginary and real parts of the state trajectory elements separately. As seen,
a smooth evolution targeting the desired state is obtained at the final
iteration. The density matrix has unit trace and is Hermitian during its
evolution over the time. The resulted final state is
$\rho\left(T\right)=\left(\begin{matrix}0.0005&\text{0}\text{.0157 +
0}\text{.0145i}\\\ \text{0}\text{.0157 +
0}\text{.0145i}&\text{0}\text{.9995}\\\ \end{matrix}\right)$.

Figure 2: Evolution of Diagonal Elements of the State Trajectory at Final
Iteration Figure 3: Evolution of Off Diagonal Elements of the State Trajectory
at Final Iteration

The $x$ and $y$ components of the control signal are depicted in Fig. 4 and
Fig. 5. The control signals hold the time evolution of the optimal trajectory
to reach the maximum fidelity. Without these signals and by only considering
the drift part of quantum mechanical Hamiltonian, fidelity would stay zero.

Figure 4: Control Evolution in $x$ Direction Figure 5: Control Evolution in
$y$ Direction

The evolution of fidelity along the iterations is demonstrated in Fig. 6. The
graph shows the fidelity reaches 99.8 percent after 44 iterations, and it
continues to reach at more than 99.9 percent, which proves that the algorithm
has performed well in optimizing a cost functional combining fidelity and
control power required to drive the state of the system from a given point to
the given final target state.

Figure 6: Fidelity

##  10 Conclusion

In this paper, we have shown the application of Maximum Principle of
Pontryagin in order to compute the constrained optimal control at the cost of
compromising between the goal of maximizing fidelity and keeping the energy of
the field relatively small. For the studied optimal control problem, we
described the state of the quantum system by the density operator, evolving
through Liouville-von Neumann equation. As a new approach in quantum optimal
control, we obtained the first-order necessary optimality conditions for the
matrix-valued dynamics. The application of the Maximum Principle of Pontryagin
resulted in the proposed shooting algorithm, which has been used to solve the
two point boundary value problem. The suggested procedure can be easily
applied for pure state vectors, or one can also consider the dynamics of the
evolution operator. Future challenges consists in exploiting the versatility
of the optimal control paradigm further by considering open quantum system,
and also additional constraints, e.g., state constraints.

## References

  * Altafini and Ticozzi (2012) Altafini, C. and Ticozzi, F. (2012).  Modeling and control of quantum systems: An introduction.  _IEEE Transactions on Automatic Control_ , 57(8), 1898–1917.
  * Bengtsson and Życzkowski (2017) Bengtsson, I. and Życzkowski, K. (2017).  _Geometry of quantum states: an introduction to quantum entanglement_.  Cambridge university press.
  * Berman and Kosloff (1991) Berman, M. and Kosloff, R. (1991).  Time-dependent solution of the liouville-von neumann equation: Non-dissipative evolution.  _Computer physics communications_ , 63(1-3), 1–20.
  * Boothby and Wilson (1979) Boothby, W.M. and Wilson, E.N. (1979).  Determination of the transitivity of bilinear systems.  _SIAM journal on control and optimization_ , 17(2), 212–221.
  * Boscain et al. (2021) Boscain, U., Sigalotti, M., and Sugny, D. (2021).  Introduction to the pontryagin maximum principle for quantum optimal control.  _PRX Quantum_ , 2(3), 030203.
  * Cheng et al. (2020) Cheng, J., Deng, H., and Qia, X. (2020).  Accqoc: Accelerating quantum optimal control based pulse generation.  In _2020 ACM/IEEE 47th Annual International Symposium on Computer Architecture (ISCA)_ , 543–555. IEEE.
  * Cong (2014) Cong, S. (2014).  _Control of quantum systems: theory and methods_.  John Wiley & Sons.
  * D’Alessandro (2000) D’Alessandro, D. (2000).  Topological properties of reachable sets and the control of quantum bits.  _Systems & control letters_, 41(3), 213–221.
  * D’Alessandro (2021) D’Alessandro, D. (2021).  _Introduction to quantum control and dynamics_.  Chapman and hall/CRC.
  * D’alessandro and Dahleh (2001) D’alessandro, D. and Dahleh, M. (2001).  Optimal control of two-level quantum systems.  _IEEE Transactions on Automatic Control_ , 46(6), 866–876.
  * Dehaghani and Pereira (2022a) Dehaghani, N.B. and Pereira, F.L. (2022a).  High fidelity quantum state transfer by pontryagin maximum principle.  _arXiv preprint arXiv:2203.04361_.
  * Dehaghani and Pereira (2022b) Dehaghani, N.B. and Pereira, F.L. (2022b).  High fidelity quantum state transfer by pontryagin maximum principle.  _IFAC-PapersOnLine_ , 55(16), 214–219.  https://doi.org/10.1016/j.ifacol.2022.09.026.  18th IFAC Workshop on Control Applications of Optimization CAO 2022.
  * Dehaghani and Pereira (2022c) Dehaghani, N.B. and Pereira, F.L. (2022c).  Optimal control of quantum systems by pontryagin maximum principle.  _U. Porto Journal of Engineering_ , 8(2), 194–201.
  * Fleming and Rishel (2012) Fleming, W.H. and Rishel, R.W. (2012).  _Deterministic and stochastic optimal control_ , volume 1.  Springer Science & Business Media.
  * Glaser et al. (2015) Glaser, S.J., Boscain, U., Calarco, T., Koch, C.P., Köckenberger, W., Kosloff, R., Kuprov, I., Luy, B., Schirmer, S., Schulte-Herbrüggen, T., et al. (2015).  Training schrödinger’s cat: quantum optimal control.  _The European Physical Journal D_ , 69(12), 1–24.
  * Gruebele (2001) Gruebele, M. (2001).  Fully quantum coherent control.  _Chemical Physics_ , 267(1-3), 33–46.
  * Holland et al. (2020) Holland, E.T., Wendt, K.A., Kravvaris, K., Wu, X., Ormand, W.E., DuBois, J.L., Quaglioni, S., and Pederiva, F. (2020).  Optimal control for the quantum simulation of nuclear dynamics.  _Physical Review A_ , 101(6), 062307.
  * James (2021) James, M. (2021).  Optimal quantum control theory.  _Annual Review of Control, Robotics, and Autonomous Systems_ , 4, 343–367.
  * Jozsa (1994) Jozsa, R. (1994).  Fidelity for mixed quantum states.  _Journal of modern optics_ , 41(12), 2315–2323.
  * Lahoz-Beltra (2016) Lahoz-Beltra, R. (2016).  Quantum genetic algorithms for computer scientists.  _Computers_ , 5(4), 24.
  * Liang et al. (2019) Liang, Y.C., Yeh, Y.H., Mendonça, P.E., Teh, R.Y., Reid, M.D., and Drummond, P.D. (2019).  Quantum fidelity measures for mixed states.  _Reports on Progress in Physics_ , 82(7), 076001.
  * Lin et al. (2021) Lin, C., Ma, Y., and Sels, D. (2021).  Optimal control for quantum metrology via pontryagin’s principle.  _Physical Review A_ , 103(5), 052607.
  * Maday and Turinici (2003) Maday, Y. and Turinici, G. (2003).  New formulations of monotonically convergent quantum control algorithms.  _The Journal of Chemical Physics_ , 118(18), 8191–8196.
  * McDermott and Vavilov (2014) McDermott, R. and Vavilov, M. (2014).  Accurate qubit control with single flux quantum pulses.  _Physical Review Applied_ , 2(1), 014007.
  * Ohtsuki et al. (2001) Ohtsuki, Y., Nakagami, K., Fujimura, Y., Zhu, W., and Rabitz, H. (2001).  Quantum optimal control of multiple targets: Development of a monotonically convergent algorithm and application to intramolecular vibrational energy redistribution control.  _The Journal of Chemical Physics_ , 114(20), 8867–8876.
  * Ohtsuki et al. (1999) Ohtsuki, Y., Zhu, W., and Rabitz, H. (1999).  Monotonically convergent algorithm for quantum optimal control with dissipation.  _The Journal of chemical physics_ , 110(20), 9825–9832.
  * Palao and Kosloff (2002) Palao, J.P. and Kosloff, R. (2002).  Quantum computing by an optimal control algorithm for unitary transformations.  _Physical review letters_ , 89(18), 188301.
  * Peirce et al. (1988) Peirce, A.P., Dahleh, M.A., and Rabitz, H. (1988).  Optimal control of quantum-mechanical systems: Existence, numerical approximation, and applications.  _Physical Review A_ , 37(12), 4950.
  * Poggiali et al. (2018) Poggiali, F., Cappellaro, P., and Fabbri, N. (2018).  Optimal control for one-qubit quantum sensing.  _Physical Review X_ , 8(2), 021059.
  * Rembold et al. (2020) Rembold, P., Oshnik, N., Müller, M.M., Montangero, S., Calarco, T., and Neu, E. (2020).  Introduction to quantum optimal control for quantum sensing with nitrogen-vacancy centers in diamond.  _AVS Quantum Science_ , 2(2), 024701.
  * Sachkov (1997) Sachkov, Y.L. (1997).  Controllability of right-invariant systems on solvable lie groups.  _Journal of dynamical and control systems_ , 3(4), 531–564.
  * Stefanatos and Paspalakis (2021) Stefanatos, D. and Paspalakis, E. (2021).  A shortcut tour of quantum control methods for modern quantum technologies.  _Europhysics Letters_ , 132(6), 60001.
  * Sugny et al. (2007) Sugny, D., Kontz, C., and Jauslin, H.R. (2007).  Time-optimal control of a two-level dissipative quantum system.  _Physical Review A_ , 76(2), 023419.
  * Vitanov et al. (2017) Vitanov, N.V., Rangelov, A.A., Shore, B.W., and Bergmann, K. (2017).  Stimulated raman adiabatic passage in physics, chemistry, and beyond.  _Reviews of Modern Physics_ , 89(1), 015006.
  * Werschnik and Gross (2007) Werschnik, J. and Gross, E. (2007).  Quantum optimal control theory.  _Journal of Physics B: Atomic, Molecular and Optical Physics_ , 40(18), R175.
  * Werschnik (2006) Werschnik, J. (2006).  _Quantum Optimal Control Theory: Filter Techniques, Time-Dependent Targets, and Time-Dependent Density-Functional Theory_.  Cuvillier Verlag.
  * Xu et al. (2020) Xu, J., Li, S., Chen, T., and Xue, Z.Y. (2020).  Nonadiabatic geometric quantum computation with optimal control on superconducting circuits.  _Frontiers of Physics_ , 15(4), 1–8.
  * Zhu and Rabitz (1998) Zhu, W. and Rabitz, H. (1998).  A rapid monotonically convergent iteration algorithm for quantum optimal control over the expectation value of a positive definite operator.  _The Journal of Chemical Physics_ , 109(2), 385–391.