Second order phase dispersion by optimised rotation pulses

TL;DR

This paper introduces a method to design shorter broadband control pulses for two-level systems by optimizing phase dispersion functions, improving control efficiency in quantum systems.

Contribution

It proposes a novel quadratic phase dispersion approach and an optimal control method to create shorter, broadband universal rotation pulses avoiding local extrema traps.

Findings

Pulse duration is halved with quadratic phase dispersion control.

The method achieves broadband universal rotations efficiently.

Optimal control avoids convergence to local extrema.

Abstract

We show that the duration of broadband universal control pulses can be halved by choosing control targets with a quadratic function of phase dispersion. This class of control pulses perform a broadband universal rotation around an axis, in the Bloch sphere representation of two-level systems, given by this phase dispersion function. We present an effective optimal control method to avoid the problem of convergence to local extrema traps.

Figures

Equations

$$\hat{H}(t)=\sum _{k=1}^K{\omega }_k{\hat{\sigma }}_{\text{z}}^{(k)}+A\cos (\phi (t)){\hat{\sigma }}_{\text{x}}^{(k)}+A\sin (\phi (t)){\hat{\sigma }}_{\text{y}}^{(k)}$$

$${\hat{\sigma }}_{\text{x,y,z}}^{(k)}=E^{(k)}\otimes {\hat{\sigma }}_{\text{x,y,z}}$$

$$\displaystyle\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}=\exp{\big{[}-i\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{L}}}}_{n}\Delta t\big{]}},$$

$$\max_{\varphi}\big{\{}\mathcal{F}\big{\}}=\frac{1}{d}\max_{\varphi}\Big{\{}\text{Re}\big{\langle}\mathcal{R}\big{|}\mathcal{U}_{N}\big{\rangle}\Big{\}}$$

$$\frac{\partial\mathcal{F}}{\partial\varphi_{n}}=\text{Re}\big{\langle}\mathcal{V}_{n}\big{|}\Bigg{[}A\cos{(\varphi_{n})}\frac{\partial\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}}{\partial\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}_{\text{y}}}-A\sin{(\varphi_{n})}\frac{\partial\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}}{\partial\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}_{\text{x}}}\Bigg{]}\mathcal{U}_{n}\big{\rangle}$$

$$\exp{\Bigg{[}-i\begin{bmatrix}\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{L}}}}_{n}&\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}\\ \bm{0}&\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{L}}}}_{n}\end{bmatrix}\Delta t\Bigg{]}}=\begin{bmatrix}\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}&\dfrac{\partial\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}}{\partial\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}}\\ \bm{0}&\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}\end{bmatrix}$$

$${\Omega }_1{T}_1={\Omega }_2{T}_2=b$$

$$\mathcal{R}_{k}(\beta)=\exp{\big{[}-i\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{n}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{n}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{n}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{n}}}}_{k}\beta\big{]}}$$

$${\hat{\hat{n}}}_k=\cos ({\alpha }_k){\hat{\hat{\sigma }}}_{\text{x}}^{(k)}+\sin ({\alpha }_k){\hat{\hat{\sigma }}}_{\text{y}}^{(k)}$$

$$\alpha_{k}\triangleq\pi bQ\bigg{(}1-\frac{\omega^{2}_{k}}{\Omega^{2}}\bigg{)}$$

$$\displaystyle\varphi(t)=\pi\Omega_{\text{c}}T\bigg{(}\frac{t}{T}-\frac{1}{2}\bigg{)}^{2},$$

Full text

Optimal Control, Unitary Propagators, Universal Rotations, Phase Dispersion, Broadband Pulses

I Introduction

The critically important milestones on the road-map of quantum technology development are expected to include quantum optimal control [1, 2]. The main result of this communication is to reduce the time and energy required to control a quantum system. As an application of optimal control theory, quantum optimal control can be tasked with finding control pulses, usually in the form of electromagnetic radiation, to perform a desired operation on an arbitrary quantum state. Already, this method has been applied successfully to a wide range of experiments including magnetic resonance spectroscopy [3, 4], error-correction for quantum computing [5], quantum information registers [6], robust atom interferometry [7, 8], and high-resolution medical imaging satisfying legal irradiation constraints [9].

The growing number of reported applications of quantum optimal control are expanded by the concurrent development of optimal control algorithms. There are a number of approaches to optimal control, such as Lagrangian methods [10], time optimal control [11], annealing quantum systems to a desired effective Hamiltonian [4], sophisticated gradient-free searches [12], rapidly converging Krylov-Newton methods [13], and optimal control using analytic controls [14]. This development is essential as desired solutions push to the limit of what is physically possible and control problems become computationally, and numerically challenging [15].

Quantum optimal control problems can be solved numerically by using a piecewise-constant control pulse approximation [16, 17, 18, 19], and the trajectory gradient can be calculated and followed using the gradient ascent pulse engineering (grape) method [3]. The initial grape method showed linear convergence to a desired optimum destination, and was later extended to give superlinear convergence [20], then quadratic convergence [21, 22]. This communication introduces an ordered optimal control method, based on the grape method, to avoid local convergence traps encountered in the optimal control problem set out below. This method is termed morphic-grape by the authors.

In addition, and vital to effective experimental applications, optimal control has advanced theoretical studies calculating universal gates for quantum computing [23, 24, 25, 26, 27]. A subset of universal quantum gates include control targets designed as an effective rotation in the Bloch sphere representation of a two-level system, with a defined rotation axis and rotation angle, for any initial basis-state vector [28, 6]. The developments reported here will be concerned with this universal rotation class of control pulses.

One of the drawbacks associated with universal rotation solutions is the long duration of control pulses when compared with the easier control problem of optimising state-to-state problems [28, 1]. The optimal control method presented here will show that the previously required pulse duration can be halved by defining target propagators as a function of phase dispersion.

II Broadband optimal control to effective propagators

An ensemble of noninteracting two-level systems should include frequency dispersion terms from local environmental conditions. Uniform manipulation of the ensemble is a difficult practical task but by describing the ensemble as a bilinear control problem [16, 18], numerical optimisation can be used to find solutions that are difficult to find analytically. In formulating this bilinear control problem, the controllable parts of the Hamiltonian are the pulses, and the uncontrollable part is the local frequency dispersion term. Restricting control pulses to phase modulation, $\varphi(t)$, with constant amplitude $A$, the time-dependent Hamiltonian for this ensemble of $K$ noninteracting two-level systems can be written as

[TABLE]

where the angular frequency $\omega_{k}$ is variously termed the resonant frequency offset, chemical shift, or detuning, and describes the local frequency dispersion within the ensemble. Practically, this ensemble of resonant frequency offsets form a discrete grid of noninteracting two-level systems spread over a relevant bandwidth [29, 30]. $\hat{\sigma}_{\text{x,y,z}}^{(k)}$ are operators of the $k^{\text{th}}$ two-level system, operating in a Hilbert space and related to the Pauli matrices $\hat{\sigma}_{\text{x,y,z}}$ by

[TABLE]

where $E^{(k)}$ is a $K\times K$ single-entry matrix with the diagonal element $E_{kk}^{(k)}=1$.

The grape (gradient ascent pulse engineering) method of optimal control [3] proceeds to describe the time-dependent control pulses as piecewise constant over a small time interval $\Delta t$ [18, 19]. This approximation allows for a numerical solution of Eq. (1) through time-ordered propagation,

[TABLE]

with the Liouvillian $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{L}}}}_{n}$ being an adjoint representation of a Hamiltonian: $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{L}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{L}}}}=\openone\otimes\hat{H}-\hat{H}^{\dagger}\otimes\openone$ [20]. At a time increment $n$, the effective propagator $\mathcal{U}_{n}$ evolves the system forward in time over the interval $[0,t_{n}]$, and the effective propagator of the adjoint control problem, $\mathcal{V}_{n}$, evolves backwards from a desired target $\mathcal{R}$ over the interval $[T,t_{n}]$. The effective propagator $\mathcal{U}_{N}$ is the solution to the ensemble Hamiltonian Eq. (1) over the total pulse duration $T$. The effective propagator of the adjoint problem, $\mathcal{V}_{n}$, is needed for a fidelity gradient calculation, outlined below.

The Hilbert-Schmidt inner product of the desired effective propagator $\mathcal{R}$ and the effective propagator over the control duration, $\mathcal{U}_{N}$, gives a measure to numerically maximise and is termed the fidelity [3, 4, 28, 6]:

[TABLE]

This fidelity measure $\mathcal{F}$ is normalised by the dimension of the Hilbert space $d$ to give sensible and predictable bounds $\mathcal{F}\in[-1,+1]$. The desired target $\mathcal{R}$ is interpreted as a rotation in the Bloch sphere representation in all that follows.

Given an initial guess for control pulses, optimal control pulses can be found by stepping in the direction of the fidelity gradient vector over control manifold, until a maximum fidelity is found. The distance to step is set by an appropriate Newton-type optimisation method [31]. The gradient-following grape method [3] requires a fidelity gradient vector $\nabla\mathcal{F}$ and, for phase modulated control, this is constructed from the elements

[TABLE]

The two partial derivatives on the right-hand side are directional propagator derivatives, and can be interpreted at the derivative of the propagator, $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{P}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{P}}}}_{\!n}$, in the direction of the operators $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}_{\text{x}}$ and $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}_{\text{y}}$. The fidelity gradient is calculated by using the same fidelity definition of Eq. (6), with the directional propagator derivatives sandwiched by the effective propagator and its adjoint in Eq. (5)

For each time slice, the time propagator in Eq. (5) and each of the directional derivatives in Eq. (7) can be calculated analytically with one exponentiation of a block-triangular matrix [32, 33],

[TABLE]

in each direction $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}\in\{\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}_{\text{x}},\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{\sigma}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{\sigma}}}}_{\text{y}}\}$. Considering the large and sparse Hamiltonian formulation of Eq. (1), the main computational cost is $2N$ exponential operations, which can be performed in parallel and can include a efficient propagator recycling scheme [22].

Pulses are linearly scaleable with the desired bandwidth [1], with the transform from one pulse $\text{p}_{1}$, with duration $T_{1}$, amplitude $A_{1}$, and covering a bandwidth $\Omega_{1}$, to a second pulse $\text{p}_{2}$, with duration $T_{2}$, amplitude $A_{2}$, and covering a bandwidth $\Omega_{2}$, having the same effect if

[TABLE]

where the subscripts refer to the according property of pulse $\text{p}_{1}$ and $\text{p}_{2}$. The dimensionless bandwidth factor $b$ and scaling factor $s$ are defined here to characterise the pulse, rather than the system-specific quantities of bandwidth, pulse amplitude, and pulse duration. This is common in the discipline of magnetic resonance, and the rotation angle $\beta$ is included in the scaling factor to normalise with this use in magnetic resonance. In this definition, bandwidths are measured in Hertz, pulse duration in seconds, pulse amplitudes in radians per second, and rotation angles in radians. In the investigation that follows the scaling factor is chosen to be $s={}^{2\beta}\!/_{\pi b}$ which, in effect, gives the same pulse amplitude over all bandwidth factors for a defined system.

The choice of constant amplitude, phase modulated, control pulses is one of inference: both chirped pulses [34, 35, 36, 37, 38] and universal rotation pulses [28] are close to constant amplitude, and it is expected that pulses produced in this work will have a similar form.

The piecewise-constant time slices used in calculations are a very good approximation in nuclear magnetic resonance (nmr) spectroscopy. In cases where this piecewise-constant approximation is not valid, such as applications to electron paramagnetic resonance (epr) spectroscopy, a feedback control method [39, 40] or a transfer-matrix method [41, 42, 43] can be used to calibrate pulses due to these hardware-specific cavity effects. Further to this, optimal control can be designed to include a robustness to pulse amplitude miscalibration [29].

III Rotation axes with quadratic phase dispersion

Derived from Landau-Zener-Stückelberg-Majorana theory [37], the ensemble of Eq. (1) can be controlled with adiabatic evolution [44, 45, 46] and is described as a linear frequency sweep over the two-level systems. The phase dispersion resulting from a linear frequency sweep is quadratic, and can be used to define the local targets of an optimal control method. Previous work on this theme has optimised state-to-state problems with linear phase dispersion [47, 48, 49] and quadratic phase dispersion [50, 44, 51, 52].

Although state-to-state solutions are useful, they are only functionally defined for a specified state preparation: the desired control is only effective for a given initial state of the system. A universal control method, independent of the initial state, attempts to find desired unitary propagators which represent a rotation of the Bloch sphere; rotating all components of the coordinate system about an axis through the origin [28, 6].

The effective unitary propagator representing a universal rotation can be formulated for each member of the ensemble, $k$,

[TABLE]

for a rotation angle $\beta$. The desired axes of rotation, $\mathchoice{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{n}}}}{\hat{\vphantom{\rule{1.0pt}{6.36115pt}}\smash{\hat{n}}}}{\hat{\vphantom{\rule{1.0pt}{6.77017pt}}\smash{\hat{n}}}}{\hat{\vphantom{\rule{1.0pt}{7.06079pt}}\smash{\hat{n}}}}_{k}$, can be defined as

[TABLE]

where the angle $\alpha_{k}$ is a phase on the transverse plane of the Bloch sphere and is described by a phase dispersion function. This study uses a quadratic phase dispersion similar to that of the linear frequency swept chirped pulses [34, 35, 36, 37, 38]:

[TABLE]

This is the phase that the transverse components of states acquire during the pulsing time [47]. The introduction of $Q$ is defined here as the quadratic coefficient of a phase dispersion function, which is closely related to previous work using a linear phase dispersion function [47]. The case when $Q=0$ has been previously published as Broadband Universal Rotation By Optimised Pulses (burbop) [28].

It should be emphasised that the desired pulses in the communication are an improvement towards universal rotation pulses with reduced control duration, using a phase dispersion function similar to chirped pulses, rather than a generalisation of chirped pulses to universal rotation pulses. Furthermore, the function in Eq. (12) is a quadratic function of $Q$ but could be formulated, in principle, by any function of $Q$. This new class of pulses are named Second ORder phase Dispersion by Optimised Rotation (sordor) pulses by the authors.

IV Morphic optimal control

In the authors’ experience, the standard procedure of starting grape from a random guess or an informed guess, e.g. an adiabatic passage pulse, produces a fidelity profile periodically reducing to zero over the bandwidth, particularly when setting $Q\gtrapprox 0.1$ in Eq. (12). The reason for this is that, as the optimisation progresses, fidelity sections either side of the zeroed fidelity increase, also increasing the average fidelity, trapping these local fidelity minima into the control manifold [15, 22]. This can make the process of finding the desired optimal solutions set out in Eq. (10) a difficult task. Furthermore, a high average fidelity over the bandwidth is not sufficient for acceptable pulse performance – the pulse should uniformly manipulate the entire bandwidth.

An ordered optimisation procedure of morphic-grape is set out below to avoid these traps. This method is effective, yet easy to use, and proceeds by morphing one optimisation problem into an incrementally modified optimisation problem; using the solution of one as the new starting point of the next.

A grid of coefficients, $Q\in[0,1]$, and bandwidth factors, $b\in(0,18]$, sets a morph area for the results that follow. Sordor pulses are created with four directional morphs: forward morph, incrementing $Q$ with $+\Delta Q$; backward morph, decrementing $Q$ with $-\Delta Q$; compressed morph, interpolating a control pulse from $b$ to $b-\Delta b$; and expanded morph, interpolating a control pulse from $b$ to $b+\Delta b$. An additional smoothing stage can be used to remove discontinuities in the fidelity profile, using sequences of forward or backward morphs from the largest extrema in the $\frac{\mathrm{d}\mathcal{F}}{\mathrm{d}Q}$ profile. The ordered recipe of finding sordor pulses is outlined as follows:

Start with a universal rotation pulse ($Q=0$) and use morphic-grape:

1. (a)

Forward, increasing from $Q=0$. 2. (b)

Backward, decreasing from $Q=1$. 2. 2.

Smoothing morphic-grape for $b=18$, from extrema in the $\frac{\mathrm{d}\mathcal{F}}{\mathrm{d}Q}$ profile. 3. 3.

Start with a sordor pulse of high bandwidth factor ($b=18$) and use morphic-grape:

1. (a)

Compressed, decreasing from $b=18$. 2. (b)

Expanded, increasing from $b=\Delta b$.

Empirical investigation finds an adequate grid spacing of $\Delta b=0.2$ and $\Delta Q=0.01$. Convergence is accepted when $\|\nabla\mathcal{F}\|\leqslant\min{\left(10^{-4},10^{-5\left(\frac{210}{583}+\frac{b}{36}\right)}\right)}$. The smoothing morph starts from the eleven largest extrema, with the first five of these shown in Fig. 1. Figure 1(a) shows the places on the $Q$-coefficient grid where the smoothing starts, after the forward and backward stages, and Fig. 1(b) shows how the fidelity increases and decreases over these stages. Both are plotted as a function of the number of gradient calculations which shows the computational cost of morphic-grape at $b=18.0$, which can be regarded as the seed for the compressed and expanded stages. The forward and backward stages at different values of $b$ show a similar profile to those stages in Fig. 1.

A direct comparison of the computational cost of this morphic-grape method to the previously published universal rotation pulses [28] is not appropriate, because morphic-grape essentially includes that previous work.

V Results

The main results of this work are shown in Fig. 2. The Spinach optimal control toolbox [20, 33, 22] is used to simulate the ensemble of two-level systems, defined by Eq. (1), with minimal modification to enable optimisation to the desired unitary propagators in Eq. (10). A linearly spaced resonance offset profile in Eq. (1), with $K=1+\lceil 10b\rceil$, and a piecewise constant approximation of Eq. (5), with $N=50b$, is optimised with the $\ell$-bfgs method [20], forming a Hessian approximation from the previous twenty gradients.

It is not only interesting to see how the fidelity increases over the grid of $Q$-coefficients but it is also practically useful because experiments requiring combinations of pulses with different rotation angles $\beta$ should be matched, i.e. a sequence consisting of pulse $\text{p}_{1}$ followed by pulse $\text{p}_{2}$, should have $b_{\text{p}_{1}}=mb_{\text{p}_{2}}$ with $Q_{\text{p}_{2}}=mQ_{\text{p}_{1}}$, where $m$ is an arbitrary multiplier. The maximum fidelities found over the grid of $Q$-coefficients and bandwidth factors, from the forward, backward, compressed, and expanded morphic optimisations, are shown in Figs. 2(a) and 2(b).

It should be noted that the morphic optimisations for $\beta={}^{\pi}\!/_{2}$ are difficult, numerically and computationally; the backward, compressed, and expanded stages of the morphic optimisations give little improvement when $\beta=\pi$, and were designed primarily for the case when $\beta={}^{\pi}\!/_{2}$.

Figure 2(c) shows the maximum fidelity at each bandwidth factor compared with the fidelity of an equivalent universal rotation pulse, broadband universal rotations by optimised pulses (burbop): the shortest duration universal rotation pulses known to the authors [28, 1]. The burbop in this comparison have the same constant amplitude $A$ and scaling factor $s$ as the sordor pulses. Burbop pulses start from the same ones as published previously [28], then are further optimised with the same exact gradients [Eq. (8)] and optimisation tolerances. It is clear that sordor pulses show significant improvement over burbop.

VI Example of an application in magnetic resonance spectroscopy

This section will proceed to apply the phase modulated pulses above, the subject of this article, to the experimentally relevant system of nuclear magnetic resonance spectroscopy.

As an example, a spin-${}^{1}\!/_{2}$ 13C nucleus in a static magnetic field of $B_{0}=14.0954~{}$\mathrm{T} forms the ensemble of two-level systems in Eq. ([1](#S2.E1)). For the analysis of small molecules, pulses should cover a bandwidth of $\Omega=37.5~{}$\mathrm{kHz} [53] at this magnetic field due to chemical shielding effects from varying molecular environments. This chemical shift range is linearly spaced with $\omega_{k}\in[-\pi\Omega,+\pi\Omega]$ covering a bandwidth of $\Omega=40~{}$\mathrm{kHz}, slightly more than the required $37.5~{}$\mathrm{kHz}, and represents a 13C spin ensemble with $K=451$ members.

Pulses produced with the method outlined are morphed from one grid point to another in Fig. 2, and this morphing is key to obtaining high fidelities. This indicates that pulses should be related to each other and asks the question of whether pulses close on the grid look like each other. Twelve pulses produced by the expanded stage of the morphic optimal control method are shown in Fig. 3, six of each have $\beta={}^{\pi}\!/_{2}$ in Fig. 3(a) and $\beta={\pi}$ in Fig. 3(c), at $(200,250,300,350,400,450)$\mathrm{\SIUnitSymbolMicro s}$$.

Interestingly, sordor ${}^{\pi}\!/_{2}$ and sordor ${\pi}$ both are of very similar constitution. As shown in Fig. 3 a) and (c), pulse shapes in all cases are symmetric in time, while no symmetry constraints are imposed during the optimisation procedure. In essence, a parabolic time course of the pulse phase is observed, as is reminiscent of linearly frequency swept chirped pulses [34, 54, 55, 35, 36], but with regular spike-like features that also frequently come with approximate $\pi$, $2\pi$ ,or $4\pi$ jumps in phase. The features are even better visualised if the quadratic phase contribution is subtracted [Figs. 3(b) and 3(d)]. Similar phase behaviour has previously been reported for bip [56] and bibop [57, 58] $\beta=\pi$ inversion pulse shapes, but not for $\beta={}^{\pi}\!/_{2}$ rotation or excitation pulses.

The comparison of the sordor pulses to a chirped pulse show an additional phase modulation, where the chirped phase is from the formula

[TABLE]

The bandwidth of this chirped phase, $\Omega_{\text{c}}$, is chosen to give the same pulse amplitude as the sordor pulses, $A=(2\pi)10~{}$\mathrm{kHz}. Clearly, this does not match the quadratic component of all sordor pulses in Figs. [3](#S6.F3)(b) and [3](#S6.F3)(c), particularly for the shorter duration pulses. However, this is a good match for the $450~{}$\mathrm{\SIUnitSymbolMicro s} pulses, and the comparison chosen in this way has no difference in pulse amplitude, only the extra phase modulation on top of the chirped pulse.

Sordor pulses are very close to symmetric, which should not be a surprise because a rotation should be a reversible operation and a time-symmetric pulse is a common form for this [59, 28]. Although not implemented in this work, including this symmetry as a constraint should halve the computation time because only half the time propagators and their directional derivatives need to be calculated with Eq. (8).

Two sordor pulses are further investigated in combination: a ${}^{\pi}\!/_{2}$ pulse from Fig. 2(a), and a $\pi$ pulse from Fig. 2(b), each with $Q=0.78$ and $b=18$. This corresponds to a pulse duration of $T=450~{}$\mathrm{\SIUnitSymbolMicro s}, a time increment of $\Delta t=0.5~{}$\mathrm{\SIUnitSymbolMicro s} in Eq. (5), and the constant pulse amplitude of $A=(2\pi)10~{}$\mathrm{kHz}$$ in Eq. (1). Both of the sordor pulses were obtained from the expanding stage of morphic optimal control method outlined above, these being the best performing combination of pulses produced in Fig. 2.

Illustrative examples of sequences using these two sordor pulse are shown in Fig. 4. Figures 4(a) and 4(b) depict the performance of sordor ${}^{\pi}\!/_{2}|_{\text{x}}$ and ${\pi}|_{\text{x}}$ pulses (subscripts denote the axis about which an $\omega=0$ rotation is defined) projected onto a Bloch sphere for three different initial states. Figure 4(c) depicts these two pulses in the sequence ${}^{\pi}\!/_{2}|_{\text{x}}\to{\pi}|_{\text{x}}$, named the Hahn echo [60], Fig. 4(d) is the sequence ${}^{\pi}\!/_{2}|_{\text{x}}\to{\pi}|_{\text{x}}\to{}^{\pi}\!/_{2}|_{\text{x}}$, a building block within the inept sequence [61], and Fig. 4(e) is the perfect echo sequence ${}^{\pi}\!/_{2}|_{\text{x}}\to{\pi}|_{\text{y}}\to{}^{\pi}\!/_{2}|_{\text{y}}\to{\pi}|_{\text{y}}\to{}^{\pi}\!/_{2}|_{-\text{x}}$ [62]. The fidelities presented for each of these last three examples clearly show that sordor pulses can perform effectively in sequence.

VII Conclusions

Broadband universal rotation pulses with second order phase dispersion define a new class of pulses, named sordor pulses by the authors, with the universality of rotation pulses and an evolution character of adiabatic pulses. The difficulty in producing universal rotation pulses with second order phase dispersion, similar to adiabatic passage pulses, poses an interesting question regarding methods of optimal control.

A simple but effective morphic-grape optimal control method has been designed to avoid local extrema traps and find shorter duration control pulses when compared with the current best universal rotation (burbop) pulses. The trend of this comparison shows that sordor pulses are approximately $50\%$ of burbop durations: a significant improvement.

The advance in pulse performance, compared with burbop, derives from a phase dispersion function used to define target unitary propagators. The particular phase dispersion function used in this work consists of a constant term, a quadratic term, and a variable multiplier. The exact form of this function is not derived within this work and is used in this form simply because the function produces high fidelity pulses. Studies using other phase dispersion functions, e.g. of higher or arbitrary order, are left to future investigations.

Acknowledgements.

The authors are grateful to Andrin Doll, Steffen Glaser, Maria Grazia Concilio, Jens Haller, Malcolm Levitt, and Herbert Ullrich for their for insightful comments and discussion. Thanks are given to Mohammadali Foroozandeh, Martin Goodwin, Ilya Kuprov, and Jason Ralph for help finalising this paper. B.L. acknowledges funding by the Deutsche Forschungsgemeinschaft (DFG LU 835/ 13-1), the HGF programme BIFTM (47.02.04), and the EU Horizon 2020 programme (667192 HYPERDIAMOND).