Optimal control with a multidimensional quantum invariant

Optimal quantum control of continuous variable systems poses a formidable computational challenge because of the high-dimensional character of the system dynamics. The framework of quantum invariants can significantly reduce the complexity of such problems, but it requires the knowledge of an invariant compatible with the Hamiltonian of the system in question. We explore the potential of a Gaussian invariant that is suitable for quadratic Hamiltonians with any given number of motional degrees of freedom for quantum optimal control problems that are inspired by current challenges in ground-state-to-ground-state shuttling of trapped-ions.


I. INTRODUCTION
Quantum optimal control is widely accepted as one of the central tools in the development of quantum technological applications [1,2].While control of discrete degrees of freedom, such as qubits, is well established both in theory and experiment [3], the control over motional degrees of freedom is still in an early stage of development [3,4].In particular, the identification of optimal solutions for the control of motional quantum states is rather computationally challenging due to the strictly infinite-dimensional, or practically high-dimensional, Hilbert space of the underlying control problem.
While most proof-of-principle demonstrations of elementary building blocks of quantum information processing did not require control over motional quantum states beyond a mechanism that holds qubits in place, the current challenges in developing scalable technologies lead to a growing need to control motional quantum states [5].In trapped-ion hardware, for example, scalability is expected to require modular architectures [6][7][8] that allow the ions, which act as qubits, to be separated into many small groups.Interconnecting spatially separated clusters with quantum logical operations then requires the ability to move individual ions from one cluster to another [9].Since trapped-ion quantum logic requires the motional states of the ions to be of sufficiently low energy, it is essential that such shuttling operations do not result in motional excitations at the end of the process.This can be achieved with adiabatically slow shuttling, but, in practice, diabatic protocols [6,8,10] are sought-after because of the requirement to realize all operations of a quantum algorithm within the system's coherence time.
Transport along a line can be performed diabatically, and diabatic ground-state to ground-state transport is within experimental capabilities [11,12].Diabatic transport beyond such one-dimensional problems, has been experimentally demonstrated in favorable geometries [13].In typical trapping geometries, however, actual realiza-tions of diabatic transport protocols result in excess motional energy in the final state.
Quantum invariants have been applied to a wide variety of theoretical control problems, from atom cooling [19], fast separation of two trapped-ions [25], unidimensional atomic transport in harmonic-traps [21,29] or expansions and compressions of trapped-ion chains with minimal motional excitation [24].Combined with variational methods or included in the master equation formalism, they have also been shown to be useful in addressing problems that do not admit a standard treatment, such as anharmonic-potential transport [30], experimental noise minimization [31] or atomic state population inversion [32].
While a substantial part of the conceptual development of quantum invariants is devoted to the construction of invariants for systems with only one translational degree of freedom, many practical applications of this framework require invariants for systems with more translational degrees of freedom.Here, we focus on a recently developed Gaussian invariant for multidimensional systems [27], and explore its suitability for applications in quantum optimal control with theoretical examples inspired by the current ion shuttling experiments.
Sec. IV.In particular, we sketch how a time-dependent Hamiltonian can be constructed such that it induces dynamics with desired properties.

A. Hamiltonian model
Real-world potentials are never strictly harmonic, but if an object's quantum state is localized within a potential well, such as a sufficiently cooled ion, even an ostensibly anharmonic potential landscape can be well approximated as harmonic around the expected position of the ion.During the shuttling protocol, this quadratic approximation is valid as long as the localization in realspace remains sufficiently narrow.
The Hamiltonian for an ion of mass m, with a quantum state localized around a d-dimensional trajectory z(t), can be expressed as where the upper-script T represents the transpose operation, x and p ≡ −i ∇ x are the canonical displacement and momentum operators with respect to the trajectory z(t), and the vector F and matrix M correspond to the first and second derivatives of the trapping potential taken along z(t).
Up to an irrelevant scalar term, this quadratic Hamiltonian can also be cast as with the trap center C(t) = 1 m M −1 F (t).If the initial quantum state of the ion is Gaussian, such as the ground-state of the Hamiltonian with a harmonic trapping potential, or a thermal state, the Gaussian character of the quantum state is preserved during the shuttling protocol within this quadratic approximation.The dynamics are then characterized completely in terms of the classical trajectory z(t) and momentum p(t) (i.e., the expected value of the displacement and momentum operators x and p) and the covariance matrix Σ(t) of all variances and covariances of x and p.

B. Multi-dimensional Gaussian Quantum Invariant
The goal of shuttling an ion such that it ends up in the ground-state of its final Hamiltonian can be understood as an optimization problem, with motional excitation as figure of merit.Instead of finding a suitable protocol as result of an optimization, one can also design suitable shuttling protocols more directly with the framework of quantum invariants.
A quantum invariant Î(t) is an operator that satisfies the equation of motion The property of such invariants that is crucial for quantum control is that the eigenstates of Î(t) -which are generally time-dependent -are solutions of the timedependent Schrödinger equation with the Hamiltonian Ĥ(t).
A non-degenerate eigenstate of an invariant Î(t) that commutes with the Hamiltonian Ĥ(t) at the beginning and at the end of the dynamics thus corresponds to a shuttling protocol in which an ion starts in an eigenstate of the initial Hamiltonian Ĥ(t = 0) and ends up in an eigenstate of the final Hamiltonian Ĥ(t = T ) even though it does not need to follow the dynamics of an eigenstate as it would be the case in an adiabatic protocol.
A shuttling protocol can thus be defined in terms of a time-dependent invariant, and Eq. ( 3) is the defining relation determining a time-dependent Hamiltonian Ĥ(t) that achieves the desired shuttling protocol.
The following analysis is based on an invariant [27,28] that is defined in terms of the actual classical trajectory z(t) that the ion is meant to take, and a time-dependent, positive (i.e.positive semi-definite), d-dimensional matrix R(t).
The quadratic component M(t) of the potential is obtained from R(t) via the relation [27 where [X, Y ] Z = XZY − Y ZX is the generalized commutator, {X, Y } = XY + Y X the anticommutator and In practice, for a given a choice of R(t), the above equations must be solved in reverse order, i.e., first find J with Eq. ( 6), then compute A using Eq. ( 5) and finally M with Eq. ( 4).Once M(t) is obtained in this fashion, the linear force term F is given by in terms of the classical trajectory z(t) [27].
For any choice of z(t) and R(t), one can thus find a time-dependent Hamiltonian of the form given in Eq. ( 2), or equivalently Eq. ( 1), such that the eigenstates of the invariant are solutions of the Schrödinger equation.Since a general choice of z(t) and R(t), however, does not ensure that the invariant commutes with Ĥ(t) at the start or end of the shuttling protocol, it is not yet ensured that any eigenstate of the invariant is also an eigenstate of Ĥ(t) at the start or at the end of the shuttling protocol.
Commutativity of the invariant and Hamiltonian at any instant t is achived, if the set of relations The solid lines correspond to the optimal solutions based on a diagonal form (D) of R, while the dashed lines result from the optimisation of all entries of R (F).The color code distinguishes the number of control functions (i.e., free parameters) considered in the optimisation.
are satisfied [27].If the boundary conditions Eqs. ( 8) are satisfied at t = 0 and at t = T , with C(0), C(T ), M(0) and M(T ) determined by the initial and final Hamiltonian, then the time-dependent ground-state of the invariant indeed defines a shuttling protocol in which an ion evolves from the ground-state of an initial Hamiltonian towards the ground-state of a final Hamiltonian.

III. OPTIMAL CONTROL
The framework sketched above in Sec.II enables the construction of time-dependent Hamiltonian resulting in ground-state to ground-state transfer.Indeed, any choice of the time-dependent functions z(t) and R(t) results in a suitable Hamiltonian.Since the shuttling protocol is thus anything but unique, one can aim at finding the protocol that is optimal in a sense to be specified.To this end, one would need to perform a variational analysis over z(t) and R(t) respecting the boundary conditions Eqs. (8).
For that purpose, it is convenient to use a set of functions { z i } such that z 0 (t) satisfies the boundary conditions Eqs. (8a) and such that any function z i (t) with i > 0 satisfies the homogeneous boundary conditions, i.e., Eqs.(8a) with vanishing right-hand-sides.
The parametisation then satisfies the boundary conditions for any set of N a expansion coefficients a i .A parametrisation of the matrix R(t) can be chosen in a similar way, with one matrix R 0 (t) satisfying the inhomogeneous boundary conditions, and a set of matrices R i (t) satisfying the homogeneous version of Eqs.(8b), such that the parametrized matrix satisfies the boundary conditions for all values of the N b free parameters b i .The condition that R(t) be positive, however, requires some extra care.Demanding that all the matrices R i (t) (including i = 0) are positive, and that all the expansion coefficients b i are non-negative, does not give access to all positive matrices within the spanning set.In all the following optimizations, there are thus no direct restrictions on the expansion coefficients b i , but positivity of R(t) is assessed numerically, either in terms of sub-determinants (Sylvester's criterion, [33]), moments or numerically obtained eigenvalues.

IV. RESULTS
The formalism introduced in the previous section is applicable to any control problem that requires the identification of ground-state to ground-state transport protocols.In particular, it can be used to derive the time-dependent trapping potentials for ground-state-toground-state transport that is optimal according to a figure of merit to be specified.The following section shows three examples of such optimized transport protocols.

A. 3D shuttling with reduced displacement from the potential center
Realizing a fast, diabatic shuttling requires strong acceleration and deceleration at the beginning and the end of the protocol.Such processes are most easily realised in terms of a large displacement of the ion from the center of the trapping potential.Since in practice such large displacements imply that the ion traverses domains of substantial anharmonicity in the trapping potential, it is desirable to minimize the displacement.The following discussion exemplifies the control problem of minimizing the maximum displacement of the ion while it moves around a corner in a potential landscape.
The initial and final position are denoted by the vectors C 0 = (0, r, h) and C T = (r, 0, h).An additional boundary condition r J (t = T /2) = (r, r, h ) reflects the fact that realistic trapping potential do not always ensure dynamics with a constant value of the z-component of the ion trajectory (see Fig. 1a).
The trapping potential along the ion trajectory is characterized in terms of a matrix M(t) following Eq.( 1), that satisfies the boundary conditions with and axial frequency ω t that is substantially smaller than the radial frequency ω r , so that the preferred direction of motion is initially along the x-axis and finally along the y-axis.Since the displacement required for an acceleration or deceleration can be reduced to any desired value by increasing the confining potential along the axis of acceleration or deceleration, the following optimization includes the constraints that With these boundary conditions and constraints, adiabatic shuttling yields close-to-perfect transport (with fidelities F > 99%, [34]) in regimes where T > 200/ω t .Any transport with T < 200/ω t is thus considered diabatic.The following results correspond to T = 10/ω t , i.e., protocols an order of magnitude faster than in the adiabatic regime.
As described above in Sec.III, both the classical trajectory of the ion z and the matrix R can be parameterized such that the boundary conditions are satisfied for any value of the expansion coefficients.In order to highlight the ability of the present framework to find good solutions of a control problem in the presence of restrictions on achievable potentials, the following discussion will compare two different parametrizations of R. In case 1, the matrix R is taken to be diagonal, and only the time-dependent diagonal elements are subject to optimization.The physical implication is that the principal axes of the trapping potential coincide with the x-, y-and z-axis.The potential can thus not be rotated, but its strength along the principal axes can be modulated.In case 2, each matrix element of R is subject to optimization, such that rotations of the principal axes are also possible.lutions with R in diagonal (D) or full (F) form.The red curve depicts the solution without any optimization.For τ 0.15 the ion is being separated from the trap center, until a maximum separation of about 0.14r is reached.Subsequently, the ion moves towards the trap center before it goes through another interval of separation.Due to the symmetry in trap geometry and time-dependence of the invariant, the trajectory is symmetric around the instant τ = 1/2.The large separation in the early and late part of the shuttling protocol is not unexpected, because the ion needs to be strongly accelerated and decelerated, which is easily achieved far away from the trap center where the trapping forces are strong.Achieving a similarly strong acceleration and deceleration without substantial displacement from the trap center requires a suitably designed time-dependence of the shuttling protocol.The green curves depict a numerically optimised protocol with 3 control functions for each independent component of z and R.Both forms of R lead to a reduction of the maximum displacement compared to the first solution.Nevertheless, rotating the trapping potentials (i.e., using the full R matrix, dashed line), slightly improves the protocol.Similarly, the optimization based on 5 control functions (light blue curves) shows even greater improvement, as the displacement does not exceed 0.04r.Likewise, releasing the ma-jor axes of the potential (dashed line) leads to a better solution.Further minimizing the displacement, however, requires a much more significant increase in the number of control functions.Thus, the optimisation of up to 25 control functions (i.e., 225 free parameters when the full R matrix is considered and 150 otherwise) leads to a maximum displacement smaller than 0.02r (dashed purple line).The impact of employing rotating potentials is minor in this case.
Unsurprisingly, at any given number of control functions (per matrix element), the maximum displacement obtained in an optimization with a general matrix R is slightly smaller than in the corresponding optimization with a diagonal matrix R. Crucially, however, optimizations with a diagonal matrix R can outperform optimizations with a general matrix R if the number of control functions in the former problem is only moderately larger than in the latter problem.This means that restrictions in experimentally realizable rotating trapping potentials can be overcome in terms of additional temporal degree of freedom in the tuneable trapping parameters.
The discussion so far was focused on the trajectory on the ion, but a motional ground-state-to-ground-state protocol also requires the covariances of the Gaussian state to evolve towards their appropriate value.Fig. 2 depicts the time-dependent components M xx , M yy and M zz of the trapping potential in (a) and the linear force in (b).In between the initial weak confinement along the x-axis and the final weak confinement along the y-axis, there are two periods of time of enhanced confinement, reflecting the fact that strong confinement is necessary to rapidly accelerate or decelerate the ion.Fig. 2b shows that such an acceleration a force (see Eq.( 7)) that is particularly strong in the two directions subjected to change (y, z in the first half and x, z in the second).Merely enhancing the confinement and the force in order to ease control over the trajectory is, however not possible, since this would result in substantial undesired dynamics of the covariances.Therefore, control of the covariances requires that strong confinement be applied only during some time windows, while weak confinement is necessary for the covariances to evolve in such a way that the ion ends up in the ground-state of the final Hamiltonian.

B. Narrow wave-packets in weak confinement
The spatial uncertainty, i.e.
the width of the wave-packet in diabatic dynamics is determined by the strength of the confining potential.A weak potential implies a broad wavepacket that is susceptible to anharmonicties in the potential landscape.In diabatic shuttling, however, it is possible that the wave-packet propagates through a domain of weak confinement without broadening to the natural width of its trapping potential.
Beyond the fact that a fast shuttling protocol leaves little time for the wave-packet to broaden, it is possible to initialise the wave-packet in a shape such that the x as a function of normalized time.The red curve depicts the dynamics experienced by σ 2 x when the evolution follows the optimal protocol, whereas the green curve corresponds to an adiabatic transport.The dashed black lines indicate upper and lower frequency limits and all quantities are expressed in units of h mω 0 , with ω0T = 1 in this case.
natural dynamics makes it more narrow.If this happens shortly before the wave-packet enters the domain of weak trapping potential, i.e. when it can still be controlled well by the stronger potential, the time-window in which the wave-packet remains narrow can be prolonged.
In the following, this is exemplified with the problem of reducing the spatial width σ x of a wave-packet in a time-dependent potential trapping frequency ω(t).The trapping frequency has to satisfy the boundary conditions ω(0) = ω(T ) = ω 0 , and the wave-packet has to evolve from the ground-state of the initial Hamiltonian to the ground-state of the final Hamiltonian.During the time intervals [0, T /4] and [3T /4, T ], the trapping frequency can be modulated with a maximally allowed frequency √ 2ω 0 , but in the time-window [T /4, 3T /4], the maximally allowed trapping frequency is given by ω 0 .The goal of the optimisation involves minimizing a weighted average of the maxima of σ x (t) and ω(t), i.e., the cost function at the central window T /4 < t < 3T /4 of the dynamics.
The protocol shown in Fig. 3 is based on the optimization of 25 control functions.The left panel depicts the squared trapping frequency ω 2 (in blue) as a function of normalized time.Its optimal evolution consists of increasing the strength of the potential to approximately twice the initial value and then gradually release the confinement so that in the central timewindow (T /4 < t < 3T /4) the confinement is weak.The squared trapping frequency reaches a minimum near zero at t = T /2 and then exhibits a symmetric time evolution with respect to that time instance.The red curve in the figure represents the temporal evolution of the wave-packet spread in a potential evolving with the optimal ω 2 (t).Not only does the uncertainty not increase throughout the protocol, but it even decreases slightly 0.0 0.2 0.4 0.6 0. in the central region, when the potential is weak (see inset).Interestingly, the minimum uncertainty is obtained at t = T /2, coinciding with the point of minimum confinement.
In weak confining potentials, wave-packets with narrow widths can only be obtained in the non-adiabatic regime.In fact, following the optimal time-evolution with a total protocol time 10 times longer, T = 10T , radically different results are obtained.In this case, uncertainty evolves inversely proportional to confinement, as shown by the green curve in Fig. 3.In the first quarter of the protocol, it slightly decreases as a consequence of a greater confinement.As the potential strength is reduced, the wave packet spreads out more widely, reaching a maximum in the central window.Again, increasing the confining strength in the second half of the protocol results in a relocalization of the wave packet before it evolves to the final motional ground-state.

C. Shuttling based on null control over the curvature of the potential
The problem of shuttling around a corner discussed above in Sec.IV A highlights that restrictions in realizable trapping potential can be compensated by temporal degrees of freedom.The ability to compensate lack of control over some degrees of freedom by controlling other degrees of freedom can be demonstrated in the following example, in which only linear force terms, but not confinement, are subject to controllable time-dependence.
Crucially, the dynamics of the controlled particle depends on the curvature along the classical trajectory taken by the particle.If the actual trapping potential is anharmonic -as it is the case in any realistic scenario -then the particle can experience a position-dependent curvature, which can be controlled in a time-dependant way via the trajectory of the particle.The choice of tra-jectory can therefore be used to compensate for restrictions on the tuneability of the trap confinement strength.
This feature can be exemplified with the potential of a tapered trap, in which the confinement in the xdirection varies along the y-axis.The trapping parameters ω x , ω y and y c are time-independent, and only the trap-center y 0 can be temporally modified.The actual confinement in the x-direction is characterized by the frequency ω x y(t) yc , and since the trajectory y(t) can be controlled in terms of the trap center y 0 (t), is is possible to effectively modulate the confinement in the x-direction.
Because there is no direct control over the quadratic components of the trapping potential, the control framework with the present invariant does not provide a direct way to construct a protocol that guarantees groundstate-to-ground-state shuttling.Applying the framework straight-forwardly results in an over-determination of the quadratic component of the potential, i.e. the matrix M in Eq. ( 2).On the one hand, M is determined via Eq.( 4) following the regular framework; but, on the other hand, it is also determined via curvatures of the trapping potential (Eq. ( 12)) along the ion trajectory y(t).The prescription for M(t) following Eq.( 4) ensures ground-stateto-ground-state shuttling, but the prescription following Eq.( 12) is what is experienced by the ion.Without a means to ensure that those two prescriptions coincide, the requirement that the ion end up in the ground-state of the final potential is not automatically met by the basic framework.Nevertheless, the deviation between the two prescriptions can be numerically minimized over the degrees of freedom in selecting the time dependent trajectory y(t), and a successful minimization will result in a ground-state-to-ground-state shuttling protocol.
In the following, this will be exemplified with the con-trol problem of moving the trap center from its initial position y 0 (0) = y i to its final position y 0 (T ) = y f , such that the particle evolves from ground-state to groundstate of initial and final Hamiltonian.Specific parameter values used in the example are y i = 10y c , y f = 1000y c and ω y = 10ω x .With those parameter values, adiabatic transport is considered to occur when T > 70/ω y (for which F > 99%, or, equivalently, a final motional heating n < 10 −2 ).The results shown here correspond to T = 3/ω y (regime in which extending the adiabatic protocol leads to n = O(10 3 )).
Fig. 4(a) depicts the solution obtained for the trap center y 0 (t) in black.It exhibits oscillations far outside the interval through which the ion is meant to be shuttled.The ion trajectory (depicted in red), however, remains within this interval, and the motion of the ion is comparatively slow in the first 60% of the time-window.Most of the transport occurs within a short time-window (about 20% of the full duration), and in the last 10% of the time-window the ion is, again, moving very slowly.
As a result of this dynamics, the confinement along the x-direction remains rather weak in the first half of the shuttling protocol, and it increases towards its final value only later, as depicted in Fig. 4(b) in red.The time-dependent confinement that ensures perfect groundstate-to-ground-state transfer (i.e., derived from the invariant framework) is depicted in blue.The deviation between the two is fairly minor (in fact, the infidelity of the protocol is below 0.02%, i.e., the final motional excitation n < 10 −4 ), meaning that the dynamics experienced by the ion will be very similar to that which guarantees ground-state-to-ground-state transfer in the case of direct control of the curvature of the potential.
The close-to-perfect fidelity of the shuttling protocol also implies that the covariances of the ion evolve towards the values corresponding to the ground-state.This can also be seen in Fig. 4(c).The green curve in Fig. 4(c) depicts the instantaneous spatial covariance, i.e. the covariance corresponding to adiabatic transport along the actual trajectory of the ion (as depicted in red in Fig. 4(a)).The red curve depicts the dynamics of Σ xx following the same trajectory, but with the actual velocity.Despite the strong acceleration of the ion (corresponding to the high amplitudes of the trap center in Fig. 4(a)), and the periods of slow dynamics of the ion, that results in the plateau of the instantaneous covariance (green in Fig. 4c)), the actual dynamics of Σ xx is rather unspectacular, i.e. the width of wave-packet simply shrinks monotonically and adopts its desired final value at the end of the shuttling protocol.The covariance thus behaves exactly the way that one would have designed it in an adiabatic protocol, even though the protocol is far outside the regime of adiabaticity.

V. CONCLUSIONS
Ground-state-to-ground-state shuttling remains one of the crucial challenges in the realization of scalable quantum information processing of trapped ions.Invariant-based quantum control is a promising route towards finding shuttling protocols that facilitate the suppression on undesired motional excitations.Given the multiple experimental imperfections that can contribute to reduced quality of shuttling operations, the ability to optimize shuttling protocols without compromising the goal to have an ion end up in its ground state can become a central step towards a practical technology.
The exemplary tasks of reducing the displacement from the center of the trapping potential and the spatial width of a wave-packet without strong confining potential can give a flavor of what can be achieved with suitable temporal shaping of trapping potentials.
While the framework of invariants can ensure that the resultant Hamiltonians satisfy physically motivated constraints, such as the decomposition into a pre-defined kinetic energy term and a time-dependent potential term, any practical implementation typically requires compliance with additional constraints resultant e.g. from the geometry of the trap electrodes.Given the notorious difficulty to construct an invariant that is consistent with a Hamiltonian of a specified set of properties, it is hopeless to expect that suitable invariants for given trap geometries can be found.The ability to minimize the deviations between solutions of the general invariant framework and solutions that can be obtained in the presence of experimental constraints, however, enables the incorporation of such constraints in practice.As explained in Section III, a quantum invariant compatible with a quadratic Hamiltonian is determined by the choice of z(t) and R(t), subject to satisfy a set of boundary conditions.Any solution for z(t) can be parametrised in terms of a set function z i (t) (with i ≥ 0): where z 0 (t) satisfies the inhomogeneous boundary conditions (Eqs.(8a)) and z i>0 (t) fulfil the corresponding homogeneous equations.Thus, z 0 (t) can be chosen to be the lowest-degree polynomial satisfying the associated boundary equations in each component, i .e. parametrizing each component of z 0 (t) with the following function: where C I , C F refer to the initial and final values respectively and τ ≡ t T is the normalized time.If extra constraints are imposed beyond the boundary conditions, polynomials of higher degree must be considered.For example, in the case of constraining f (τ = 1/2) = C J , the function that fulfills all the constraints is The rest of functions z i>0 satisfy the homogeneous boundary conditions.A similar polynomial expansion may be used for each of the components of z i>0 , but the following parametrisation is more suitable for problems with fast dynamics: where the coefficients ω 1 , ω 2 are free provided that ω 1 , ω 2 ∈ Z.Therefore, g(τ ) defines a set of vectors z i>0 that result in a valid solution for z(t) independently of the expansion coefficients a i .As in the case of z 0 , it may also necessary to impose constraints on the functions g(τ ) for some particular control problems.Of particular interest is to impose that g(τ = 1/2) = 0.In that case, a set of functions that satisfy simultaneously the homogeneous boundary conditions and the added constraint is where we introduced the notation S i ≡ sin(πω i /2) and the frequencies ω 1 , ω 2 , ω 3 ∈ Z subject to ω 1 = ω 2 .
Similarly, the positive semi-definite matrix R(t) can be parametrised as where R 0 (t) is a d × d diagonal matrix that fulfils the non-homogeneous boundary conditions and can be conveniently expressed as where f i (t) are polynomials of p-degree satisfying that f i (t) > 0, ∀i, t.The boundary conditions are imposed as F,ii , ḟi (0) = fi (0) = ḟi (T ) = fi (T ) = 0.These conditions are fulfilled by the control function defined in Eq. (A2).Note that the nonnegativity is also guaranteed because C I , C F > 0 and the control function in Eq. (A2) has no minimum for τ ∈ (0, 1).
On the other hand, the set of matrices R i (t) (i > 0) can be chosen following two different alternatives.The first option is to assume each R i (t) as a matrix with functions in each entry that satisfy the homogeneous boundary conditions.Then, R(t) > 0 is numerically assured at every instant of time either in terms of the eigenvalues, the Sylvester's criterion or the calculation of moments.
Alternatively, each matrix R i (t) (i > 0) can be decomposed in terms of lower-triangular matrices L i such that R i (t) = L i L † i satisfy the homogeneous boundary conditions and ensure that R(t) > 0. As it is sufficient to require all the functions g ij (t) to satisfy the homogeneous boundary conditions to guarantee that L i L † i does so.In addition, in order for L i L † i to be positive semi-definite, the diagonal must be nonnegative throughout the time domain.Both constraints are fulfilled if the off-diagonal functions are chosen from Eq. (A4) and the diagonal is parametrised in term of the non-negative functions h(τ ) ≡ g(τ ) 2 , where g(τ ) is defined in Eq. (A4).Therefore, restricting the expansion coefficients b i to be positive guarantees that R(t) is positive semi-definite in the whole time domain.Note that, although this method guarantees the positivity of R(t) at any time, also limits the spectrum of valid solutions.

FIG. 1 :
FIG. 1: (a) Example of a corner trajectory (red solid line) from C0 to CT restricted to pass through rJ .Dashed black lines are drawn to facilitate the three-dimensional view of the problem.(b) Displacement as a function of normalized time (τ = t/T ).The solid lines correspond to the optimal solutions based on a diagonal form (D) of R, while the dashed lines result from the optimisation of all entries of R (F).The color code distinguishes the number of control functions (i.e., free parameters) considered in the optimisation.

FIG. 3 :
FIG.3: Squared trapping frequency ω 2 and spatial uncertainty of the quantum state σ 2x as a function of normalized time.The red curve depicts the dynamics experienced by σ 2x when the evolution follows the optimal protocol, whereas the green curve corresponds to an adiabatic transport.The dashed black lines indicate upper and lower frequency limits and all quantities are expressed in units of

FIG. 4 :
FIG.4:(a) Trajectory of the ion along the y axis as a function of time, following the optimal temporal evolution (red) and the corresponding trajectory of the trap center y0 (black).(b) In red, curvature along the x-axis experienced by the ion as a function of time, achieved purely indirectly via the position-dependent curvature of the tapered trap.In blue, the curvature that guarantees perfect ground-state to ground-state transfer, obtained with direct dynamic control of the curvature independently of position.(c) Spread of the wave-packet in x (in normalized units) as a function of time.The red line shows the dynamics resulting from following the trajectory y(t) (with the appropriate velocity), while the green curve represents the instantaneous spread corresponding to the trajectory's curvature Mxx, ignoring the dynamics of the full trajectory.
Appendix A: Parametrisation of z(t) and R(t)