Optimal control of the transport of Bose-Einstein condensates with atom chips

Using Optimal Control Theory (OCT), we design fast ramps for the controlled transport of Bose-Einstein condensates with atom chips’ magnetic traps. These ramps are engineered in the context of precision atom interferometry experiments and support transport over large distances, typically of the order of 1 mm, i.e. about 1,000 times the size of the atomic clouds, yet with durations not exceeding 200 ms. We show that with such transport durations of the order of the trap period, one can recover the ground state of the final trap at the end of the transport. The performance of the OCT procedure is compared to that of a Shortcut-To-Adiabaticity (STA) protocol and the respective advantages/disadvantages of the OCT treatment over the STA one are discussed.

We show that with such transport durations of the order of the trap period, one can recover the ground state of the final trap at the end of the transport. The performance of the OCT procedure is compared to that of a shortcut-to-Adiabaticity (stA) protocol and the respective advantages/disadvantages of the oCt treatment over the stA one are discussed.
The measurement's outcome of a phase-sensitive sensor probing forces exerted on neutral atoms by inertial, material or electromagnetic sources depends dramatically on the initial conditions, i.e. on the position, velocity and size of the input matter-wave. A lack of knowledge or scattering in these initial properties inevitably leads to systematic effects or statistical errors harming the sensor's performance. An example of the degree of control needed can be grasped if one considers making a test of the Universality of Free Fall (UFF) with two different atomic species to put bounds on a possible violation of the UFF at the femto-level in the Eötvös ratio 1 , level at which state-of-the-art experiments perform with material test masses 2 . Such a precise experiment requires that the initial positions, center-of-mass velocities and expansion rates are defined at a level better than 1 μm, 1 μm/s and 100 μm/s (35 pK in 3D), respectively 3 .
To meet these stringent requirements, the temperatures of the atomic ensembles have to be drastically reduced (down to a sub-nK level) and their size must remain compact (not exceeding a few mm after several seconds of free expansion) clearly indicating the necessity of using Bose-Einstein Condensates (BEC). Such a direction is taken by several metrology groups worldwide 4-10 , including the QUANTUS 11 and MAIUS 12 consortia, which reached important milestones in controlling quantum gases dynamics in microgravity conditions using atom chips 13,14 .
In a recent work 15 , we considered an approach based on Shortcut-To-Adiabaticity (STA) protocols to obtain analytic solutions for the transport of the BEC in an atom chip setup with realistic anharmonic and rotating trapping potentials. This approach based on the reverse engineering technique allows for a full control of the translational degrees of freedom of the BEC. It is, however, exciting several collective modes of the quantum gas, an effect which could eventually compromise the expected metrological gain if such a source is used without any precaution as an input of an atom interferometer. It is in this context that the use of optimal control theory (OCT) can reveal an unchallenged potential of targeting a given final state in timescales shorter than the trivial adiabatic manipulation, which is of no practical use in the metrology context since it is associated with poor cycling rates.
The aim of optimal control theory is to bring a dynamical system from one state to another, while minimizing a cost functional, such as the control time or the energy of the pulse used. The modern version of OCT is born with the Pontryagin's Maximum Principle (PMP) in the late 1950s 16,17 . Originally applied to problems of space dynamics, OCT is nowadays a key tool to study a large spectrum of applications both in classical 18,19 and quantum physics [20][21][22] . In the Pontryagin formulation, solving an optimal control problem is equivalent to finding extremal trajectories which are solutions of a generalized Hamiltonian system. These trajectories satisfy the maximization condition of the PMP as well as specific boundary conditions [18][19][20] . The implementation of the PMP is far from being trivial and numerical control algorithms have been developed to approximate the optimal solution 23 .
Among others, we can mention the gradient 19,24 and the Krotov 23,25 algorithms, which are nowadays standard tools in physics.
OCT has been applied with success to quantum systems since the 1980s in domains extending from molecular physics and nuclear magnetic resonance to quantum information science (see refs 22,26 , for recent reviews, and references therein). The application of OCT to BEC dynamics has also been explored in different contexts. Using the Gross-Pitaevskii equation, the optimal coherent manipulation of an atomic BEC has been investigated in a series of studies (see refs [27][28][29][30][31][32][33] to cite a few, and references therein). The transport of cold atoms has also been optimized for simple models in combination with invariant-based inverse methods 15,[34][35][36][37] . It should be mentioned here that OCT and STA are usually compatible in the sense that an OCT methodology can be built on top of a basic STA frame of solutions 34,[38][39][40][41][42] . One can also note that recently, new methods have been tested successfully to bridge the gap between an ideal STA and a realistic experimental implementation for the optical transfer of a degenerate gas, demonstrating fast highly non-adiabatic transfer with almost no residual sloshing using corrected STA trajectories 43 .
In this paper, we discuss the application of optimal control theory for the fast transport of Bose-Einstein condensates with atom chips while simultaneously controlling the quantum degrees of freedom of the problem to target the ground state of the final trap as the optimization result. The article is organized as follows: We describe in the next section the chip model used to transport the BEC. The chosen cost functional and the associated transport ramp are presented in the next section, which is followed by a comparison of our findings to the results of the STA technique applied in a similar context 15 . We finally illustrate the impact of the OCT ramp duration on the internal degrees of freedom of the final BEC state. We conclude by discussing the limits of the methodology we have developed, and by mentioning potential experimental implementations.

Chip model
We consider the case of a Z-shaped chip configuration used to trap and manipulate cold Rb atoms in micro-gravity (See reference 15 for a detailed description of the numerical model and reference 44 for the description of an experimental implementation). The three spatial dimensions are denoted by the three coordinates X, Y and z. z is the direction perpendicular to the chip. X and Y are two orthogonal directions in the plane of the chip. The trap is naturally rotating in the (XY) plane when the physical parameters governing the trap potential change. The diagonalization of the associated Hessian matrix allows to define two new eigen-coordinates x and y of the trap, rotated compared to the fixed X and Y coordinates 15 . The physical parameters which govern the trap potential are the chip intensity I w and the bias magnetic field B bias . For the present study, I w is fixed at 5 A and the control parameter for the implementation of the transport ramp is the time-dependent bias magnetic field B bias (t), which varies between B bias (0) = B i = 21.5 G at the initial time t = 0 and B bias (t f ) = B f = 4.5 G at the end of the transport corresponding to t = t f .
In such a configuration already described in our previous study 15 , the minimum of the trap is at the origin in x and y, and it is located at a distance z 0 (t) from the chip surface. At t = 0 we have .  z (0) 0 45 0 mm and at the end of the transport .  z t ( ) 1 65 f 0 mm. If we limit ourselves, in a first approximation, to the simplest case of a time-dependent harmonic trap, the center-of-mass of the condensate z A (t) in the direction normal to the surface follows Newton's equations of motion where ω z (t) denotes the frequency of the trapping potential along z at time t. Moreover, in the Thomas-Fermi approximation 45 , the evolution of the size of the BEC is accurately described by a scaling approach 46,47 . The size of the BEC is defined by the three time-dependent radii r x (t), r y (t) and r z (t) of the paraboloid associated with the BEC wave function, with It was shown 46,47 that the time-dependent scaling factors λ x (t), λ y (t) and λ z (t) obey the three coupled second order differential equations where ω x (t) and ω y (t) denote the frequencies of the trapping potential along x and y at time t. The full behavior of the trapping frequencies as a function of the control parameter B bias can be found in 15 . Initially the trapping frequencies are ω π ⋅  (0) 2 15 Hz . At the end of the transport ω π ⋅  t ( ) 2 10 Hz . The largest time scale associated with the trap is therefore of the order of 100 ms. An adiabatic transport would thus require transport durations larger than 1 s. Here we want to design a simple, fast and efficient transport ramp for the BEC. The OCT technique being very powerful, we have decided to optimize a single control parameter, B bias (t), in order to control the final position of the BEC z A (t f ), its final speed v A (t f ), and its final size defined by the three final scaling factors λ x (t f ), λ y (t f ) and λ z (t f ). We also wish to control the final expansion rates given www.nature.com/scientificreports www.nature.com/scientificreports/ ( ) z f . Finally, since we want the harmonic approximation to hold during the entire transport, we also wish to limit the time-dependent offset between the position of the center of mass of the BEC and the center of the trap |z A (t) − z 0 (t)| as well as the the time-dependent offset between their respective speeds . To be compatible with metrology applications with an integration over thousands of experimental cycles, we want this transport to be realized quickly, i.e. in a duration of the order of the largest time scale associated with the trap, that is of the order of 100 ms with the present chip configuration.

Cost functional
To implement such an optimal control scheme, we first introduce the "classical" point-wise translational energy of the condensate in the reference frame of the trap as well as the "quantum" energy associated with the 3D Thomas-Fermi wave function is the scattering amplitude, a s is the s-wave scattering length of Rb-87 and N = 10 5 denotes the number of condensed atoms. The first term in Eq. (5) describes the potential energy associated with the finite size of the condensate, the second term is the kinetic energy associated with the size dynamics, and the third and last term is the average mean-field interaction energy between the atoms of the condensate. The numerical factors (1/14) and (15/28) seen in Eq. (5) come from the specific definition given in Eq. (2) of the size of the condensate using a Thomas-Fermi expression for the probability density.
The goal we want to achieve is the minimization of a total cost functional C tot , defined by the sum and a running cost The terminal cost was designed to insure the formation of the ground state of the trap at time t f . It imposes the minimization of the total energy of the condensate at the end of the transport. The running cost is introduced in order to limit the transient excitation of the condensate in the moving harmonic trap. Here we fix λ 1 = 1 and the two other dimensionless parameters λ 2 and λ 3 are chosen to express the relative weights between the three terms of the cost functional. Changing the values of λ 2 and λ 3 affects the progress of the optimization procedure by changing the path it takes during optimization. This can lead in practice to different final transport ramps, which will take into account the relative weight assigned to each of the terms of the cost functional.

transport ramp
The initial and final traps are defined by the initial and final values B i = 21.5 G and B f = 4.5 G of the bias magnetic field B bias (t). Since in experiments one can be limited by the switch on/off speed of the magnetic field we circumvent this problem by insuring a smooth variation of B bias (t) at t = 0 and at t = t f . For this reason we have chosen to impose i  B Note that a consequence of these boundary conditions imposed on B bias (t) is that similar relations hold for all trap parameters such has the trap position z 0 (t) and the trap frequencies in all directions ω x (t), ω y (t) and ω z (t). The optimization procedure we have adopted is therefore using the dimensionless control function u(t), from which we can calculate the optimal bias magnetic field using Eq. (9).
www.nature.com/scientificreports www.nature.com/scientificreports/ optimal control We now reformulate our optimization problem in the framework of optimal control theory. We refer the interested reader to standard textbooks for details [18][19][20][21] . The state of the system is described by a state vector x, with λ λ λ λ λ λ As suggested by Eqs (1) and (3), the time evolution of all components of the state vector x is governed by a set of coupled first order differential equations controlled by u(t) through the time dependence of the trap position and frequencies. Once u(t) is chosen and for well defined initial conditions at t = 0, these equations are easily solved using a Runge-Kutta algorithm 48,49 or the Verlet method 50 , for instance.
According to the Pontryagin maximum principle 16,17 , the extremal solutions of the problem, candidates to be optimal, satisfy the equations of Hamiltonian p cl f 3 From Eq. (12), it can be easily shown that the dynamics of the adjoint state is governed by the following set of coupled first order differential equations   = − .  p p (14h) 8 7 In addition, the transversality conditions for the adjoint state read We use a standard first-order gradient algorithm which is adapted to the control problem under study. The optimization procedure proceeds as follows: (i) First we fix an initial control ramp u(t) arbitrarily, such as the linear ramp u(t) = t/t f for instance, or the STA ramp obtained from ref. 15 . (ii) We then compute the magnetic field B bias (t) using Eq. (9) and we deduce the trap dynamics by calculating the trap motion z 0 (t) and the trap frequencies ω x (t), ω y (t) and ω z (t); (iii) Using the Verlet method 50 , we then solve Eqs (1) and (3)  This procedure is repeated until convergence is reached. Figure 1 shows a typical example of convergence of this algorithm. The condensate is assumed to be initially at rest in the ground state of the initial trap. The initial control ramp is the shortcut-to-adiabaticity solution (see ref. 15 for details). The weight parameters are λ 1 = 1, λ 2 = 5.10 5 and λ 3 = 0.001. We have chosen in this example a large value for λ 2 in order to impose a fast convergence for the control of the final size of the condensate. In practice the correction parameter ε has to be chosen small enough to insure the convergence of the optimization algorithm. Since the correction to the control ramp is introduced at first order only, decreasing the value of ε beyond a reasonable limit does not improve the accuracy of the optimization procedure but it slows down the convergence. In the present example we have chosen ε = 10 −11 . In Fig. 1, panel ( Since the total cost functional given in Eq. (6) is characterized by a very large weight λ 2 associated with the final quantum energy, we see that E qu (t f ) is very quickly minimized, in about 1,000 iterations. This limit of 1,000 iterations is emphasized in Fig. 1 with a vertical dashed red line. Once this convergence is reached, the final 3D size of the condensate adopts the size of the ground state of the final trap and the size dynamics of the BEC is frozen. This convergence was obtained at the cost of a transient degradation of the final classical energy, which reaches a maximum of about 20 nK after about 60 iterations, but the final classical energy is then minimized very quickly to reach a near-zero value in about 1000 iterations. It is only when this first stage of convergence is reached (iteration number >1000) that the last cost functional, associated with a smaller weight λ 3 , starts to decrease. One can note that the convergence of the average classical energy during the transport [in panel (c)] is rather slow since it requires more than 10 7 iterations before it starts to stabilize at values close to 30 nK. This value can be compared with the energy of the condensate in the initial trap, which is close to 120 nK, and with the energy of the condensate in the final trap, close to 10 nK. The transient excitation during the transport is therefore relatively limited. www.nature.com/scientificreports www.nature.com/scientificreports/

Comparison of different optimization procedures
In Fig. 2 the Shortcut-To-Adiabaticity (STA) transport ramp obtained in ref. 15 (dotted blue line) is compared with two results obtained with the present optimal control technique (OCT). Note that strictly speaking, these two methods assume slightly different constraints and that, in principle, STA can be combined with the minimization of a cost functional. However, a quantitative study on this specific point is beyond the scope of this article where we concentrate mainly on developing the transport method with the OCT protocol. The correction parameter is ε = 10 −10 . The dashed green line labeled as "cl-OCT" shows the result obtained for the weight factors λ 1 = 1, λ 2 = 0 and λ 3 = 5.510 −4 . The solid red line labeled as "qu-OCT" is for λ 1 = 1, λ 2 = 3.3 and λ 3 = 5.510 −4 . The difference between these two OCT results lies in the fact that qu-OCT takes into account the influence of the finite size of the BEC in the cost functional, while cl-OCT considers the BEC as a classical point-wise particle. The BEC model used for cl-OCT is therefore similar to the model used in STA and these two approaches can be compared directly. We see in panels (b) and (c) that the three methods are very efficient for the control of the final average position and velocity of the BEC since the condensate is fully at rest in the center of the trap at the end of the transport and for all times t > t f = 150 ms. In addition, the transient position and velocity offsets during the transport reach similar values using these three different optimization methods. One can note in panels (b) and (c) that in terms of maximum transient offset in position and speed, from the two methods that we can compare directly, cl-OCT is a little better than STA (maximum offsets of 4.5 μm vs. 5.3 μm in position and 14 μm/ms vs. 22 μm/ms in speed) but this difference is not very significant in practice. The transient offsets of the qu-OCT approach are slightly larger than those of the cl-OCT method (with maximum offsets of 6.2 μm in position and 15 μm/ms in speed). Again this increase would be very benign in a practical implementation. Note finally that the three control fields B bias (t) shown in panel (a) are relatively similar, with a fast initial decrease during the first half of the ramp, before www.nature.com/scientificreports www.nature.com/scientificreports/ 75 ms, followed by a much slower decrease afterward. A first conclusion of this study is therefore that, if one is mainly interested in the control of the average translational degree of freedom of the BEC, the STA approach, whose numerical implementation is much simpler than OCT, is sufficient.
It is in the size dynamics shown in panels (d), (e) and (f) that there is a striking difference between qu-OCT and the two other optimization methods. In terms of size dynamics, cl-OCT and STA give very similar results which consist in a persistent size excitation of the condensate after the transport. This result was already seen in ref. 15 where it was shown that it was mainly the first quadrupole mode Q 1 which was excited, thus explaining that the size oscillation along x, y and z is almost periodic after the transport. The qu-OCT approach is able to suppress efficiently this quadrupole-mode excitation and, at the end of the transport, the sizes Δx, Δy and Δz remain constant. We can therefore conclude that the introduction of a minimization goal for the quantum energy associated with the finite size dynamics of the condensate allows the qu-OCT transport ramp to prepare the true ground state of the final trap at t = t f . When the size dynamics is not accounted for, as in the STA and cl-OCT approaches, it is impossible to insure the preparation of the lowest energy state in the final trap using short transport ramps.
The optimized OCT transport ramps were obtained using a Thomas-Fermi approximation in a 3D harmonic trap. We have therefore verified, by solving the 3D mean-field time-dependent Gross-Pitaevskii equation for the evolution of the time-dependent macroscopic condensate wave function ψ(x, y, z, t), that this control is robust when taking into account the anharmonicities and the rotation of the trap. The numerical method used for this calculation is described in ref. 15 . This result is illustrated in Fig. 3, showing the time evolution of the average atomic densities

Influence of the transport duration
What remains to be seen is the efficiency of these various optimization procedures for different transport durations. Figure 4 shows in panel (a), for the three optimized ramps, the variation of the average translational energy cl f t cl 0 f as a function of the ramp duration t f . Panel (b) shows, in the same conditions, the maximum position offset Max|z A (t) − z 0 (t)| during the transport. We see here that whatever the transport duration STA and cl-OCT are characterized by a very similar performance in terms of transient excitations. This confirms the advantage of the STA protocol in practical implementations, due to its overall simplicity when compared to cl-OCT. We also see that, on one hand, when the transport duration is larger than 140 ms (i.e. about 1.4 times the largest time scale associated with the trap), the transient excitations realized by the improved qu-OCT procedure are very close to the ones of cl-OCT and STA. On the other hand, for transport durations smaller than 140 ms larger transient excitations are obtained when using qu-OCT.
We could however verify that for all transport durations in the range 100 ms ≤ t f ≤ 200 ms, the qu-OCT method is able to minimize very efficiently the residual size excitations after the transport, a goal which is not achievable with the STA or cl-OCT procedures. This can be seen in Fig. 5, which shows the residual oscillation amplitudes strate that for t f ≥ 140 ms the residual size excitations of the condensate can be limited efficiently by optimal control and that this limitation does not introduce any detrimental effect on the transient excitation of the BEC. The same result can also be obtained by optimal control with shorter transport ramps, but at the cost of an increased transient excitation of the condensate.

Conclusion
In conclusion, we engineered optimal control theory protocols allowing for the fast, excitation-less transport of BECs over large distances compatible with a precision atom interferometric use. The ramps presented in this work relied on a single-parameter (bias magnetic field) optimization to shift the trap minimum position of the atom chip, promising a straightforward experimental implementation. The results of the OCT procedure relied on a scaling approach assuming a harmonic trapping. Real-life implementations on atom chips come with anharmonic corrections, mainly cubic in the direction of the transport, that scale with the position offset between the atoms and the trap minimum during the transport and with an inherent rotation of the trap. We demonstrated in this study, by solving 3D Gross-Pitaevskii equations for typical anharmonic and rotating chip traps, that the proposed OCT protocol does not compromise the target state solution even for very competitive ramp times of 150 ms. This also suggests a successful transfer to experiments. Moreover, we indicated by studying the impact of different transport durations, the methodology to follow in order to device the shortest ramps possible. Indeed, by quantifying the maximum offset induced by each ramp duration, every experimental implementation would be characterized by an anharmonicity range explored according to the specific trap configuration considered. This range determines, ultimately, the success of the ramp in reaching the ground state of the final trap. The positive outcome of this study suggests a natural generalization to the dual-species transport case. For this latter, no analytic neither intuitive solutions do exist. The STA approach generally fails since the two species experience different potential frequencies due their mass difference. A comparable OCT approach to the one adopted in this study, based on