==== Front Sci Rep Sci Rep Scientific Reports 2045-2322 Nature Publishing Group UK London 33335278 79309 10.1038/s41598-020-79309-8 Article Quantum control operations with fuzzy evolution trajectories based on polyharmonic magnetic fields http://orcid.org/0000-0003-2787-8123 Fuentes Jesús j.fuentesaguilar@ugto.mx grid.412891.7 0000 0001 0561 8457 División de Ciencias e Ingenierías, Departamento de Física, Universidad de Guanajuato, Loma del Bosque 103, 37150 León, Guanajuato Mexico 17 12 2020 17 12 2020 2020 10 222566 11 2020 7 12 2020 © The Author(s) 2020, corrected publication 2021 https://creativecommons.org/licenses/by/4.0/ Open AccessThis article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. We explore a class of quantum control operations based on a wide family of harmonic magnetic fields that vary softly in time. Depending on the magnetic field amplitudes taking part, these control operations can produce either squeezing or loop (orbit) effects, and even parametric resonances, on the canonical variables. For these purposes we focus our attention on the evolution of observables whose dynamical picture is ascribed to a quadratic Hamiltonian that depends explicitly on time. In the first part of this work we survey such operations in terms of biharmonic magnetic fields. The dynamical analysis is simplified using a stability diagram in the amplitude space, where the points of each region will characterise a specific control operation. We discuss how the evolution loop effects are formed by fuzzy (non-commutative) trajectories that can be closed or open, in the latter case, even hiding some features that can be used to manipulate the operational time. In the second part, we generalise the case of biharmonic fields and translate the discussion to the case of polyharmonic fields. Using elementary properties of the Toeplitz matrices, we can derive exact solutions of the problem in a symmetric evolution interval, leading to the temporal profile of those magnetic fields suitable to achieve specific control operations. Some of the resulting fuzzy orbits can be destroyed by the influence of external forces, while others simply remain stable. Subject terms Physics Quantum physics Quantum mechanics Single photons and quantum effects Consejo nacional de ciencia y tecnología, México558745 Fuentes Jesús issue-copyright-statement© The Author(s) 2020 ==== Body Introduction A typical quantum control programme consists of a setup of oscillating external fields irradiating a micro-object (qubit) confined in a very small cavity either to induce a class of driven motion or to manoeuvre its degrees of freedom toward a special configuration. One shall be particularly cautious with some other details that cannot be disregarded, in that the settled presence of environmental noise, the radiative pollution or the emergence of electric shocks, amongst other issues, could result in an inoperative control protocol. Then, how to achieve the desired effects, at least averagely, in a simple but realisable way? Customarily, some of these drawbacks can be successfully circumvented if a reasonable approximation is considered. There are control techniques based on magnetic nuclear resonance1,2, that accomplish acceptable results even if one ignores the intrinsic inhomogeneities as well as any electrical component conveyed by the sequence of variable magnetic fields, although beforehand the approximation relies on the quality of the generated magnetic pulses. Other schemes addressed to systems with discrete spectrum3 can suppress unwanted effects if the pulses of electrical forces depend only on the temporal domain. It is also the case for the effective description of an atom that moves across an optical device and interacts with a photon, which can be represented as an oscillating electric field4,5 as long as the dimensions of the cavity are small enough such that the approximation is not drastically far from the accurate description. Nonetheless, whenever the particle is not strictly localised, another class of control techniques shall be considered. Think about the motion of an electron, in an unbounded state, that is controlled with a protocol of simple Rabi rotations. Hence, the implementation of ample cavities would be necessary due to the electron’s free propagation. These operations are usually implemented via ion traps6–9 or optical devices10,11, and receive special attention in a variety of applications where the unbounded motion shall be guided in a particular way, e.g. particle accelerators. The interest in this kind of problems has motivated a number of works not limited to ion traps, but also in the realm of squeezed states12—useful to produce efficient non-demolition measurements13,14. Although the method is, once again, based on an approximation. It basically relies on the application of harmonic oscillator potentials with time-dependent frequencies, producing an effective oscillator field inside the walls of a hyperbolic structure, which is in turn connected to a pulsating electric potential6,15. Our discussion below is addressed to these types of control problems, focusing our attention on time-dependent magnetic fields that are modulated in intensity by way of control currents. We shall design magnetic field pulses whose temporal profile is such that a charged particle under their influence evolves in time with either a stable (ion traps) or unstable (squeezing) motion, even neither of them, which yields a particular operation that produces parametric resonance. Whenever the magnetic pulse is smooth enough, we will be allowed to look for exact solutions of the evolution problem. However, the application of soft operations will not be absent from very tiny discontinuities in the evolution process16 due to the concatenation of pulses. Indeed the presence of control currents produces a circular electric field enclosing the magnetic field, which could yield a sudden electrical jump. In our proposal, such discontinuities will be partially avoided by observing that the unitary matrices of evolution satisfy a simple Toeplitz algebra. Even more, these unitary matrices of evolution truly reproduces the evolution loop phenomenon within the structure of fuzzy orbits17,18, which makes the particle to return to its initial values after a certain period of time (i.e. the interval of operation), albeit according to the type of motion taking part the fuzzy orbits can be open or broken. From our viewpoint, these orbits have a semiclassical character in that we are to consider only periodic, quadratic Hamiltonians, yet the loop phenomenon in itself also arises in problems with non-quadratic Hamiltonians19,20. This interesting aspect entices a variety of quantum control applications16,21,22 as well as some theoretical aspects in quantum tomography23–25. In agreement with causality, the change in the field intensity due to time-varying currents will not be registered instantaneously in every part of the laboratory, instead the effects will be sensed depending on the dimensions of the laboratory and on the velocity of the signals. In that regard, an accurate treatment of the problem should include the retarded field effects, nevertheless, for simplicity we have restricted our analysis to the non-relativistic regime neglecting the amount of time in which the fields propagate. Our simplification deserves special attention for a practical implementation, and we shall survey the consequences as well. For a self-contained discussion, let us start summarising the known facts26. Consider a spinless particle of charge e and fixed mass m moving in a uniform magnetic field B(x)∈R3, generated by a class of vector potential A(x) through the relation B(x)=∇×A(x). In particular, if the magnetic field is homogeneous a possible choice of A(x) is1 A(x)=12B(x)×x, where B(x) will be aligned in the direction Oz, for the sake of convenience. As well, it is worthy of note how this gauge condition does not define uniquely A(x) when B(x) is given. Let us start with a time-independent, quadratic Hamiltonian of the form:2 H=m2V2, where V=1mP-ecA(X) is the operator that accounts for the velocity of the particle and P is the operator of momentum. At this level A(X) is also an operator and depends on the observable of position X. Besides, we see that the observables X and P obey the usual commutation relations [Xj,Pk]=iħδjk, [Xj,Xk]=0 and [Pj,Pk]=0. Whereas the components of the operator V satisfy [Vj,Vk]=ieħm2cϵjklBl(X), where ϵjkl is 1 (-1) if j, k, l is and even (odd) permutation of 1,2,3, otherwise it is zero. Likewise, the commutation relations between the components of X and V are effortlessly obtained given that A(X) commutes with Xj, it follows that [Xj,Vk]=imħδjk. There are special features that shall not escape attention. Unlike the classical problem where the trajectory described by the particle is fully localisable, in the quantum analogue a type of abstract space localisations will enter into the picture. The quantum particle will describe a trajectory whose centre (X¯,Y¯) is a fuzzy point, i.e. it is not fully localisable given that the operators associated with its position do not commute [X¯,Y¯]≠0 —even either commuting or not there is nothing at that point, despite its geometric origin27,28. For instance, take the gauge (1), then the Hamiltonian will separate into two components H=H⊥+H‖, the transversal motion (the one of our interest) will occur in the xOy plane whereas the parallel motion will be aligned to Oz. The components of the operator V are:3 Vx=Pxm+ωc2YVy=Pym-ωc2XVz=Pzm, where ωc=emcB is the cyclotron frequency and B=|B(X)|. The respective rotation centre results in a fuzzy point with coordinates:4 X¯=X+1ωcVy,Y¯=Y-1ωcVx, that satisfies the commutation relation [X¯,X¯]=-iħmωc. Note that if we define ρ2=(X-X¯)2+(Y-Y¯)2, then πρ2=2πmωc2H can be interpreted as a fuzzy surface spanned by the transversal orbits. One can conclude that the surface πρ2 is a constant of the motion, yet it cannot change continuously29 since it is quantised. This fact has a number of applications in the framework of loop quantum gravity30,31, for instance. Indeed the fuzzy centres are not limited to appear in quantum Hall-effect models32 under the influence of a fixed magnetic field, but still occur in time-dependent problems. In our case, we will discuss how the time-varying magnetic operations become responsible for a number of noncommutative aspects supported on closed or open semiclassical trajectories. In turn, we are to define the time-dependent version of the Hamiltonian (2), actually it can be written in a similar fashion but considering a magnitude-controlled current I(t) that modulates A(X) and B(X), in this way their corresponding time-dependent versions are written as A(t,X)=I(t)A(X) and B(t,X)=I(t)B(X). Such definitions imply that the time-dependent, quadratic Hamiltonian describing the motion of the non-relativistic, spinless, charged particle reads as33:5 H(t)=m2V2(t)=12mP2-ecA(t,X)2=12mPx2+Py2+eB(t)2c2(X2+Y2)⏟Hosc(t)-eB(t)2mcLz⏟Hrot(t)+12mPz2⏟H‖. Given that H‖ leads to the well known free particle eigenfunction 12πħeipzz/ħ, we are to separate this part from our analysis, and only pay attention to H⊥(t)=Hosc(t)+Hrot(t), where the oscillatory term Hosc(t) is truly equivalent to a two dimensional harmonic oscillator, whereas the easily integrable term Hrot(t) represents the rotations caused by Lz—the component of the angular momentum projected onto Oz—which is a conserved quantity. The former picture resembles, for example, a cylinder whose symmetry axis is parallel to Oz and carries a homogeneous current density on its surface. To proceed it is convenient to express the perpendicular terms in (Eq. 5) in dimensionless variables, in order, we symbolise with T=2πω the time scale or period of operation. Our new variables become:6 t→tT,P→TħmP,X→mħTX, in consequence the intensity of the oscillatory field is written as:7 β(t)=eTB(t)2mc. The function β(t) is, in some extent, the cornerstone of our study: We shall implement a continuous, real function that solves the evolution problem. One is much in need to be specially careful with this task, for if β(t) is bounded and piecewise continuous, the evolution of X and P will be continuous as well, although the same is not true for the components of the kinetic momentum provided there are sudden jumps that can be interpreted as electric shocks34,35 E=-ec∂∂tA(t,X). The substitution of Eqs. (6) and (7) into the perpendicular terms of the Hamiltonian (Eq. 5) yields the simplified, dimensionless expression:8 H⊥(t)=12P2+β2(t)X2⏟Hosc(t)-β(t)Lz⏟Hrot(t), in the special case that B(t) oscillates periodically (a regular Floquet problem) with frequency ω, the substitution H⊥(tω)→H⊥(t)ħω corresponds also to a dimensionless representation. One must be aware that in such case the stability regions of the motion described by the older Hamiltonian will no be valid for the dimensionless Hamiltonian any more. Even though our Report is wholly concerned with magnetic operations, we would like to remark that despite the elementary and minimalistic structure of the Hamiltonian (Eq. 8), it can also describe the evolution of charged particles in hyperbolically shaped Paul’s traps6 of radius r0. Under such circumstances, the time-dependent electric potentials of the type Φ(t,X)=eϕ(t)2r02X2+Y2-2Z2 or Φ(t,X)=eϕ(t)2r02X2-Y2 define a modulated intensity β(t)=eT2ϕ(t)mr02. One immediately notes that the evolution problem can be handled in a similar fashion to the magnetic problem—in both cases yielding evolution matrices u(t,t0) identical for classical and quantum dynamics. Below we are to solve the evolution problem described by Eq. (8). The already studied biharmonic amplitudes9,34 β(t) occupy the first part of our survey, but we shall reformulate the stability map where the magnetic control operations are comprised. In the second part, we shall generalise the biharmonic approach in terms of polyharmonic amplitudes β(t), that at the same time play the role of exact solutions. For this case we cannot generate a stability map such as the Strutt diagram that governs the harmonic and biharmonic cases. The systematic classification of fuzzy orbits in the polyharmonic case remains open. Results We start with the operator equations satisfied by the unitary evolution operators U(t,t0):9 ddtU(t,t0)=-iH(t)U(t,t0),ddt0U(t,t0)=iU(t,t0)H(t0),U(t0,t0)=1, nonetheless, for convenience and without loss of generality, we chose deliberately t0=0 in order to introduce the abbreviations u(t)=u(t,0) and U(t)=U(t,0). As well let us represent with the column-vector Q of dimension 2N the set of ordered pairs of observables X and P as Q=(Q1,…,Q2N)T=((X,P)1,…,(X,P)2N)T (the notation AT indicates the transpose of the matrix A). In our particular case Q=(X,Px,Y,Py)T, recall that we are exclusively interested in the motion projected into the plane xOy, described by the Hamiltonian H⊥(t). Then a set of time-dependent observables X(t) and P(t) in the Heisenberg picture will evolve as a linear combination of their initial values X and P according to the rule:10 Q(t)=U†(t)QU(t)=u(t)Q, where u(t) is a 2N×2N time evolution matrix that determines uniquely the evolution operator U(t), and vice versa, iff Q spans a complete set of observables in a Hilbert space H. We can write down the j-th coordinate of the fuzzy centre (Eq. 4):11 X¯j(t)=u¯j1X1+u¯j2P1+⋯+u¯j,2N-1XN+u¯j,2NPN,(X1,X2,X3)=(X,Y,Z),u¯jk=1T∫0Tdtujk(t). The direct application of the matrix u(t) to the initial vector Q will generate the whole trajectory of motion. Furthermore, since H⊥(t) is quadratic in its variables, the evolution operator U(t) will be the same in the classical and quantum regimes. Thus, the classical motion trajectories U(t) will be analogously interpreted in the quantum case, where, the mean position ⟨X⟩ will evolve according to iħddt⟨X⟩=⟨[X,H⊥(t)]⟩ (Ehrenfest’s theorem). Note that if two unitary operators U1(t) and U2(t) lead to the transformationU1†(t)QU1(t)=U2†(t)QU2(t)thenU2(t)U1†(t)Q=QU2(t)U1†(t) and U2(t)U1†(t) commutes simultaneously with any function depending on X and P. Moreover, the observables X and P generate an irreducible algebra in L2(R), therefore U1(t)U2†(t) must be a phase factor provided these are unitary operators. As a result we have U1(t)U2†(t)=eiφ, or equivalently U1(t)=eiφU2(t) where φ is a real number. Such detail becomes relevant in quantum control problems, for if any two unitary operators differ by a phase factor, they generate the same transformation of quantum states, so both are equivalent in this sense. The condition of evolution loop will be satisfied inasmuch as the whole set of observables Q return to their initial conditions after a time interval T. In consequence the analysis is considerably simplified if we assume that Eq. (8) possesses a periodic temporal dependence, such that H⊥(T+t)=H⊥(t), that is a Floquet Hamiltonian. As well, the periodicity of Eq. (8) implies that U(t)=U(T+t), concluding that any observable Q in the Heisenberg picture will evolve periodically even if it does not depend explicitly on time. It follows from Eq. (10) that Q(T+t)=U†(T+t)QU(T+t)=Q(t), or for a fixed period of time evolution Q(T)=U†(T)QU(T)=Q, which means that the loop condition indicates U(T)=eiφ1. Biharmonic fields We shall now analyse the evolution loop trajectories (flowing inside of a solenoid with a time-varying, homogeneous current density on its surface) as well as the squeezing effects achieved by biharmonic oscillations of the type34:12 β(t)=β0+β1sin(ω1t)+β2sin(ω2t),ω1,ω2∈R, where the selection of either β2=0 or ω2=0 enables the harmonic case. The temporal profile of β(t) is not restricted to a particular form, in fact it can be quite arbitrary. Nevertheless any reasonable physical implementation of β(t) shall preferably avoid any sudden jump such as the ones furnished by kicked protocols16,21,34,35, since the resulting electric delta shocks could compromise an even evolution. Accordingly, an alternative is to choose a sufficiently smooth pulse in terms of harmonic functions. An interesting consequence of our cylindrical model is that the commutation relation [Hosc(t),Hrot(t′)]=0 is invariably satisfied, thus the evolution operator U(t) can be split into two distinct components, U(t)=Urot(t)Uosc(t), describing the evolution around Oz and the evolution of the magnetic oscillator Uosc(t), respectively, each one satisfying Eq. (9). The resulting evolution matrix u(t), therefore, will consist of two steps of evolution starting from the initial condition Q to Qosc(t) to Q(t) unfolded by the consecutive application of the oscillatory and rotational evolutions, in that order. The matrix of rotations urot(t) will be generated by the operator Urot(t), which is straightforwardly integrated as13 Urot(t)=e-iLz∫0tdt′β(t′), whereas the matrix uosc(t) of dimension 4×4 will be constructed from Uosc(t). Fortunately, the integration of the oscillatory part can be further simplified due to uosc(t) can be reduced to a single matrix h(t) that will make each pair of observables Qx=(X,Px)T and Qy=(Y,Py)T to evolve simultaneously, i.e.14 uosc(t)=h(t)h(t) will act at once in each of the subspaces spanned by Qx and Qy. In that regard, for Qj with j=x,y, we have:15 dh(t)dtQj=iUosc†(t)[Hosc(t),Qj]Uosc(t)=Λ(t)Uosc†(t)QjUosc(t)=Λ(t)h(t)Qj, where16 Λ(t)=01β2(t)0, concluding that h(t) simply obeys the differential equation:17 dh(t)dt=Λ(t)h(t),h(0)=12×2, while the integration of this equation becomes standard for stationary fields β(t)=β0, yielding solutions of the form eΛt, in general we want a computer to do the job. Observe that the determinant of the symplectic matrix h(t) is an integral of the evolution in that det(h(t))=1. This attribute turns essential to classify the particle’s motion as well as to find exact solutions, as we shall show later. The matrix h(t) is at the core of any evolution process generated by the class of periodic fields β(T+t)=β(t), and we must discuss the ascribed dynamical features to such matrix. First, h(t) is symplectic, it has two eigenvalues λ+,λ- such that λ+=1λ-, hence the algebraic structure of h(t) is wholly determined by the scalar Σ=tr(h(t)), as suggested by the characteristic polynomial D(λ)=λ2-Σλ+1 and its roots λ±=12(Σ±Δ), where Δ=Σ2-4. We are granted, therefore, to classify the outperformed motion through |Σ| into three categories34: |Σ|<2. The matrix h(t) has the eigenvalues λ+=e+iσ and λ-=e-iσ, σ∈[0,π2]. The motion is completely stable and restricted to an oscillatory evolution. The annihilation and creation operators, a- and a+, can be obtained from the eigenvectors of h(t). Control operations of this type are useful in ion traps. |Σ|=2. The eigenvalues of h(t) are merely ±1, defining a separatrix between the stable and unstable regions. The motion described by field amplitudes in the threshold region might originate an effective parametric resonance, which supports a variety of interesting applications, e.g. a charged particle could be attracted by repulsive forces. |Σ|>2. The matrix h(t) has a pair of real eigenvalues, λ+=e+σ and λ-=e-σ, σ∈[0,π2]. This class of motion is unstable, producing a squeezing effect over the annihilation and creation operators a- and a+, i.e. expanding a- at the cost of a+, or conversely as well. Despite the classification given above is entirely based on the Heisenberg’s evolution of canonical observables, regarding material particles, the parametric resonance region, |Σ|=2, can be understood as an analogue of the parametric amplification studied by Mollow and Glauber36 in the context of coherent photon states. Even though the trajectory picture commonly attracts the concern of many researches as for the design of ion traps.Figure 1 Strutt diagram of a biharmonic field (Eq. 12). The stability regions, |Σ|<2, correspond to the clear areas, whereas the instability regions, |Σ|>2, appear in colour. The threshold belts, |Σ|=2, can be tracked with the aid of the matrix entries h11 or h21. If both entries are simultaneously equal to zero we have a pure squeezing effect, otherwise, the squeezing operations can be generated by taking the points inside the coloured areas. That is the case of the tracked subregions of amplitudes that would achieve the effects of h11=λ=2 or h11=λ=4. A point (β1,β2) lying in the clear areas could be used to produce a class of stable motion, showing a loop effect, such as the required for ion traps. As an example, a well known ion trap is the radio-frequency, quadrupole device engineered by Paul6 in which the oscillating elastic forces β(t) are written in terms of Mathieu functions. Indeed, the stability regions for such trap are simply identified from a Strutt map, useful in a variety of time-dependent problems with elastic potentials, e.g. tomography23. We shall now illustrate how this diagram works, not only for ion traps, but for another quantum control operations. To this aim, we have devised a computer routine to integrate numerically Eq. (17) in terms of a biharmonic field (Eq. 12) with fixed β0=0. A scanning process shall be meticulously undertaken, to track the different values of |Σ| in the amplitude domain. We have chosen β1,β2∈[-15,15] with angular frequencies ω1=2π, ω2=4π and t∈[0,1]. As a result, we have drawn the biharmonic map in Fig. 1, where the three regions according to |Σ|, can be located in the amplitude space spanned by {β1,β2}, in other words, depending on the values of such amplitudes, a type of quantum control operation will be furnished. First we are to survey the stable motion assured by the set of amplitudes {β1,β2} lying in the clear areas of the map in Fig. 1. To explain how the evolution loops are originated, we have integrated Eq. (17) up to T=6 field periods, finding the closed trajectory shown in Fig. 2-A, whose corresponding dimensionless amplitudes are β1=π4 and β2=-10. In this example, the particle starts its evolution with velocities (px,py)=(-5,20) at the point (x,y)=(10,-20) and finishes at the same point with velocities (px,py)=(6.25,-2.41). Any evolution loop with symmetry under parity reflection, such as this one, will have a vanishing fuzzy point (X¯,Y¯) indicating immunity to the effect of external driven forces. In fact, whenever the elastic field β(t) is biharmonic, the potential A(t,X) will evanesce at the beginning t=0 and at the end t=T of the process, then the canonical momenta and the kinetic momenta coincide, simplifying the interpretation of the wave packet at these points. The oscillatory motion is not restricted to circuits of evolution, it truly can be broken into open trajectories under some circumstances. These cases are particularly interesting since they involve a kind of free evolution embedded in the particle’s history of motion. We recall that in the Schrödinger picture the free evolution along a time interval τ corresponds to the operator e-iτħP22, which transforms the canonical observables into Xj→Xj+τPj and Pj→Pj via a matrix of evolution of the form (Eq. 14), although in this case h(t)→h(τ) has the simple structure18 h(τ)=1τ01, therefore, should h(τ) represents a segment of evolution loop the rest of the trajectory merely corresponds to the inverse operation e+iτħP22, that is h(-τ), such that h(τ)h(-τ)=1. In other words, the action of h(-τ) reverts the effects produced by h(τ) over the set of observables, concentrating the wave packet instead of contributing to its spread across the interval of operation.Figure 2 (A) Semiclassical closed trajectory in the plane xOy generated from the initial conditions (x,y)=(10,-20) and (px,py)=(-5,20). The corresponding dimensionless amplitudes of the elastic field are β1=π4 and β2=-10 which lie in the region of stability as for Fig. 1. The loop closes after a period T=6. The orange dot represents its initial and final positions. (B) Broken loop of inverted free evolution. It takes the initial conditions (x,y)=(0,0) and (px,py)=(-5,20) with the amplitudes β1=-11.86 and β2=-0.4 which lie exactly in the separatrix between the stability and instability regions. For any operation that is realised with amplitudes of this kind, it will take a period of T=2 to complete the control operation. The orange (red) dot represents its initial (final) position. Such effects are produced by allowing the field β(t) to take amplitudes {β1,β2} living in the threshold region, |Σ|=2, allowing three classes of temporal manipulations in accordance with the effective time of operation τ: (1) Accelerated free evolution, τ>T, the system will get older faster than the actual time of operation; (2) Retarded free evolution, 0<τ0. As an example, let us consider the separatrix point (β1,β2)=(-11.86,-0.4) on the map of Fig. 1. As for initial conditions we set (x,y)=(0,0) and (px,py)=(-5,20), thereon after two periods the particle will arrive at the point (x,y)=(-4.47,17.85) with the same velocities, i.e. it will end at a distance 2τ|v| behind its original location, see Fig. 2B. Finally, we are to review the squeezing effects granted by those points pertaining to the coloured areas of the Strutt map, Fig. 1. Our comments ahead, however, are not addressed to the squeezed-states philosophy, rather we focus our attention on the source of the squeezing operations regarding the periodic-evolution model (Eq. 8). Once again the symplectic structure of the matrix h(T), permits to simplify the analysis, and formulate the squeezing condition in terms of its eigenvalues as:19 λ+λ-=1, fulfilled inside the region |Σ|>2. However, the time evolution outlined by Eq. (17) has to be calculated within a time interval where β2(t) is antisymmetric around the interval centre, or the squeezing effects would not occur37. Moreover, depending on the actual value of Σ, there will be two possible squeezing transformations: the purely positive transformation for Σ>2 and the parity transformation for Σ<2. One can freely choose between these two options via the points (β1,β2), based on the control purposes. In the light of the condition (Eq. 19), it follows that the matrix elements h12 and h21 are equal to zero, and therefore any canonical observable Q will be transformed as:20 Q′→λQ′,Q′′→1λQ′′, which read as the amplification of Q′ at the cost of compressing Q′′, or vice versa. Such transformations are, in fact, embedded in the eigenvectors of h(T). It might also be noted that if a fuzzy point (Eq. 4) does not vanish, the squeezing operations will sense the influence of any external force but in an orthogonal direction to its exertion, otherwise the control operations become resistant to external perturbations, such as radiation pollution38. Nonetheless, in general the kind of transformations (Eq. 20) regarding eigenvalues xj→λxj, pj→1λpj in intervals [nT,(n+1)T], can only be produced at those separatrix points where h12=h21=0 is satisfied, namely at the intersection points of the matrix trajectories in the separatrix belt. The squeezing transformations (Eq. 20) are realisable if the amplitudes β1,β2 take values inside the unstable regions specified by the Strutt map. We have traced some of these points at which an intense squeezing λ=h11 is achieved, particularly for λ=2 and λ=4. For instance, a biharmonic field with the amplitudes (β1,β2)=(-10.3,-6.9) would attain a squeezing/amplification factor λ=4. Polyharmonic fields: exact operations There exists a miscellany of methods to generate quantum control operations of the form (Eq. 20), amongst them the programmes of kicked (or discrete) pulses that disrupt a continuous evolution process16,21,34. Unfortunately, as we have already mentioned, this kind of operations result technically impractical leading to the necessity of a more method based on soft evolution operations. If the discrete pulses are replaced by smooth pulses, e.g., biharmonic fields β(t), the imperfections are partially circumvented. Indeed, the design of such operations could be devised considering two different smooth pulses35, both in the stable region of the Strutt diagram, but demanding that the product between them belongs to the instability region, yet the successive application of both pulses makes the particle to absorb an amount of energy39, that shall be less than the difference between the neighbouring energy levels of the particle. How can we really avoid any sudden jump in the composition of quantum control operations? We shall not answer that question with utter certainty, rather we would like to sketch an attempt in the following discussion. From our viewpoint, the composition approach consists of taking portions of the time evolution course described by the Hamiltonian (Eq. 8), with the exception that this time the elastic fields (Eq. 7) shall be transformed as β2(t)→β(t) to facilitate the quest of exact solutions. Hence, in the subspace of oscillations, the matrix h(t) becomes:21 h(t)=cos(ωt)1ωsin(ωt)-ωsin(ωt)cos(ωt), noting that it adopts the form of a symplectic matrix of rotations. In particular, for any ωt=nπ2 (n is an integer) we obtain the squeezing transformations:22 hnπ2ω=0±1ω∓ω0, now, the successive application of two matrices of this type leads to:23 hsqueezing=0±1ω1∓ω100±1ω2∓ω20=λ001λ,λ=-ω2ω1,ω1≠ω2, leading again to the set of operations in Eq. (20) with the requirement of two distinct frequencies ω1 and ω2 at different times, say t1 and t2, such that ω1t1=ω2t2=π2. The sudden jumps have not been removed already, for example, should the protocol be applied in the void background there would be at least three jumps involved: 0→ω1→ω2→0. What if we translate this procedure to the continuum? We have shown that the set of operations (Eq. 20) are actually generated by symplectic matrices h(t) with the property h11=h22=12Σ, i.e. h(t) has the form of a Toeplitz matrix. This agreeable property grants that if η and ξ are Toeplitz matrices, their anti-commutator algebra ηξ+ξη, as well as their symmetric products ηξη and ξηξ, furnish matrices of the same class. It follows that the squeezing transformations (Eq. 20) can be built from the contribution of a big number of symmetric products between symplectic matrices η(t) of the form (Eq. 21), each one at different time subintervals tj with a definite elastic field amplitude β(tj), where j=0,1,2,…, namely:24 h(t)=η(tk)⋯η(t1)η(t0)η(t1)⋯η(tk), once again, conveying the property h11=h22=12Σ. To obtain the continuous analogue we shall assume that the infinitesimal jumps dh(t) are assembled by the infinitesimal contributions dη(t)=Λ(t)dt, from the right to the left sides as in Eq. (24), each one depending on a symmetric field β(t)=β(-t) around t=0, arriving at the matrix differential equation in the expanded interval [-t,t]:25 dh(t)dt=Λ(t)h(t)+h(t)Λ(t),h(0)=1,Λ(t)=01-β(t)0, which is explicitly written as26 dh(t)dt=h21-h12β(t)Σ-Σβ(t)h21-h12β(t)=(h21-h12β(t))1+Σ01-β(t)0, nonetheless, the actual trajectory determined by the whole evolution process, requires the integration of Eq. (17) along an asymmetric interval, since the character of Eq. (25) is rather auxiliary as will be clear below. On top of that the anti-commutative structure of Eq. (25) defines a β(t) in terms of a smooth enough, real function θ(t). What we are to discuss now is precisely how to obtain an exact solution of the inverse evolution problem. To this aim, first note that the diagonal elements of h(t) satisfy dh11dt=dh22dt=h21-h12β(t), but according to the initial condition h11=h22=1 at t=0, allowing to write27 h11=h22=12dθ(t)dt, consistently, it follows that det(h(t))=(12dθ(t)dt)2-h21θ(t)=1, therefore:28 h21=12dθ(t)dt2-1θ(t), substituting into (27) and rearranging terms, yields h12β(t)=h21-dh11dt. Even more, since θ(t)=h12 it directly means dh11dt=12d2θ(t)dt2, in this fashion we get:29 β(t)=-d2θ(t)dt22θ(t)+12dθ(t)dt2-1θ2(t), one can immediately notice that the special case θ(t)=1ωsin(2ωt) conducts to the basic harmonic oscillator β(t)=ω2=const. Even without any specialised consideration on adiabatic invariants40,41, β(t) as expressed in Eq. (29) is an exact solution of the inverse evolution problem for h(t) given a function θ(t), which is, in principle, arbitrary though constricted to satisfy non-trivial conditions at singular points. Nevertheless, for any asymmetric interval [t0,t] the dependence of h(t) on β(t) must be determined by the integration of Eq. (17), as we have already shown for the case of biharmonic fields. Some simple relations between β(t) and θ(t) deserve observation. Any field amplitude given by Eq. (29) in a symmetric interval [-t,t], is accomplished by a sufficiently smooth function θ(t) to assure continuity and differentiability. Specifically, at any point t where θ(t)=0, there must be dθ(t)dt=±2. As well, if θ(t)≠0 but dθ(t)dt=0 then Eq. (25) becomes a matrix of squeezing transformations such as the one in Eq. (23). Moreover, if d3θ(t)dt3=0 it follows necessarily that dβ(t)dt=0.Figure 3 Temporal profile of three polyharmonic pulses β(t) in the symmetric interval -π2,π2, with b=2,c=-3 (solid line), b=95,c=-72 (dashed), and b=2,c=-5 (dotted). The three pulses vanish at the endpoints of the interval, even though only the amplitudes in solid and dashed curves can be utilised to achieve squeezing transformations of the form (Eq. 20), since β(t)>0 in both cases, whereas the field amplitude depicted with a dotted curve could generate an evolution loop effect over the canonical variables—useful for ion traps, for instance—in that the matrix h(t) will fulfil the condition |Σ|<2 for full stability. The only missing piece is the construction of θ(t). As we said, this function is truly arbitrary, although from an empirical viewpoint a natural choice would be a θ(t) in terms of harmonic functions to induce soft control operations and prevent the evolution process from any leap. Amongst the most elementary cases, let us consider the polyharmonic function30 θ(t)=a1sin(ω1t)+a3sin(ω3t)+a5sin(ω5t)+a7sin(ω7t), since this function is antisymmetric around t=0 the corresponding β(t) as for Eq. (29) will be symmetric around that point. Particularly, if we fix the frequencies as ω1=1,ω3=3,ω5=5 and ω7=7, at the endpoints of the symmetric interval -π2,π2 a matrix of squeezing (Eq. 23) will emerge with h12=±ω=b. Accordingly, the following initial conditions must be fulfilled:31 θπ2=a1-a3+a5-a7=b,θ′(0)=a1+3a3+5a5+7a7=2,θ′′π2=-a1+9a3-25a5=-2b,θ′′′(0)=-a1-27a3-125a5-343a7=c, where c≠0 is a real parameter that we shall carefully fix to achieve the desired control operation. This set of rules gives additional information regarding the nontrivial relation between β(t) and θ(t). The condition θπ2=b determines the magnitude of the squeezing transformation depending on the whole trajectory, whereas θ′(0)=2 endows β(t) with non-singularity at t=0. In turn, the coefficients of θ(t) are:32 a1=(105b+c+58)b-10128b,a3=-(35b-c-74)b+2128b,a5=18-(21b+c-22)b384b,a7=(15b-c-26)b-6384b, see the numerical examples of the generated pulses β(t) in Fig. 3, in which we have selected β-π2=βπ2=0. Unlike the method based on biharmonic fields, where we could generate a map such as the one presented in Fig. 1, in the current case we are unable to pursue an analogue survey to scan the type of control operations produced by a set of amplitudes {β1,β2}. Even though, in the light of polyharmonic pulses the functions (Eq. 29) and (Eq. 30) truly offer additional information regarding the effects achieved at the end of the symmetric interval of operation: If β(t)>0 the magnetic field would induce a squeezing effect over the canonical observables in the concatenated intervals -π2,π2 and π2,3π2, otherwise, the magnetic field would reproduce an evolution loop effect, either closed or broken. We shall remark that the trajectory smeared in the interval -π2,π2∪π2,3π2 is not yet determined at this stage, which calls for a separate computer routine to integrate (Eq. 17) in the asymmetric interval where the initial time will be t0=-π2. That is, the integration of Eq. (17) will be split into two subintervals, each one with a definite pulse β(t). The evolution process will start at the time -π2 with the first pulse β(t), then continuing the integration in the following subinterval at π2 with the second pulse β(t) and finishing at 32π, in order to draw a congruent trajectory.Figure 4 Semiclassical trajectories generated by the consecutive application of two polyharmonic pulses of the type (Eq. 29) over a particle with initial conditions (x,y)=12,10 and (px,py)=(-5,20). (A) The evolution process starts at t0=-π2 in terms of a polyharmonic amplitude (Eq. 29) defined by the parameters b=2,c=-3 finalising its application at t=π2—where the amplitude vanishes. Immediately, at the very same time, the evolution process continues its development now in terms of the pulse defined by b=95,c=-72 that will end its application at t=3π2. The resulting operation attains a squeezing transformation with λx=-0.403 and λy=-1.125. (B) To invert the operation and recover the initial observables’s values, the second pulse has to be applied once again from t=3π2 to t=5π2, where the evolution process is now conducted by the first pulse during from t=5π2 to t=7π2. The dotted line represents the distorted evolution generated by the presence of an external harmonic force sin(t) in the direction Ox, even though note the original configuration is entirely recovered by reversing the operation. Please note, as for the pulses in Fig. 3, the resulting control operation would be absent of any abrupt discontinuity even when it was regulated by the concatenation of two different fields, for example, the amplitudes β(t) with solid and dashed curves can achieve a squeezing transformation in the fashion of (Eq. 20), whose actual trajectory in the xOy plane is plotted in Fig. 4(A) for a charged particle with initial conditions (x,y)=12,10 and (px,py)=(-5,20). The canonical observables of position are transformed as X→λxX and Y→λyY in agreement with (Eq. 20), our numerical computation indicates that such particle would be subject to a squeezing transformation with λx=-0.403 and λy=-1.125, meaning that the momenta have been amplified and compressed, respectively. These effects, furthermore, can become significantly boosted by external driven forces either time-dependent or not. In such scenario the evolution process must account for the unperturbed Hamiltonian (Eq. 8) plus the external field. For instance, we want to consider a force aligned in the Ox direction, F=(F,0,0), leading to the perturbed Hamiltonian H~(t)=H⊥(t)+FX, such that the corresponding evolution operator is constituted as U~(t)=U(t)W(t), where U(t) is the unperturbed evolution operator, whereas W(t) represents the contribution of the perturbative force that obeys the differential equation33 dW(t)dt=iFX(t)W(t),W(0)=1, where X(t) is the time-dependent observable fostered by the Heisenberg evolution (Eq. 10). Even more, since the evolved observable X(t) is linear in the canonical variables, the commutators [X(t′),X(t)] are reduced to numbers. The continuous Baker–Campbell–Hausdorff formula42 then implies that Eq. (33) has a solution of the form eiF∫0tdt′X(t′) up to a phase factor eiϕ, where ϕ is real. It follows that if t=T is the loop period in the unperturbed evolution operator, then any initial coordinate X¯j of the fuzzy centre (Eq. 4) will be reconstructed from the integral ∫0tdt′X(t′), therefore U~(T)=U(T)W(T)=eiϕeiTFX¯j. Suggesting that X¯j will be unaffected by the transit of W(t), namely U~†(T)X¯jU~(T)=e-iTFX¯jXjeiTFX¯j=X¯j. However any other coordinate of the fuzzy centre could suffer a drift U~†(T)X¯kU~(T)=X¯k-iTF[X¯j,X¯k], of course j≠k. We conclude that any fuzzy point (X¯j,…,X¯k) with noncommutative coordinates will result in a broken loop every time it is manoeuvred by an external force and the coordinates of its centre will drift in a transversal direction with respect to the force. On the contrary, if the fuzzy centre is such that its centre [X¯j,X¯k]=0, then at the end of the period T the point (X¯j,…,X¯k) will return to its initial value. To illustrate our discussion, we have examined the inverted squeezing operation referred in Fig. 4(A) in the presence of a time-dependent perturbation sin(t) aligned in the Ox direction to distort the operational trajectory, see Fig. 4(B). As we mentioned before, the evolution process (17) first sweeps in the interval t∈-π2,π2 with a given polyharmonic pulse β1(t), and continues with a different β2(t) in the interval t∈π2,3π2. To invert the operation it is required, therefore, a second application of β2(t) with t∈3π2,5π2 and then a second application of β1(t) with t∈5π2,7π2. In our example, note that the perturbed evolution returns to its initial position after completing the period of application. Also the effect of the driven external force sin(t) in this case produces a set of transformations (Eq. 20) at t=3π2 with λx=3.92 and λy=-0.96—unlike the unperturbed case, this time we have amplified the observable X and compressed Y at the cost of compressing Px and amplifying Py.Figure 5 (A) Evolution loops performed via a polyharmonic field amplitude β(t) characterised with the parameters b=2,c=-5, (see dotted pulse in Fig. 3). The charged particle has the initial conditions (x,y)=(0,0) and (px,py)=(10,5). The operation completes the loop after a period T=52. (B) Distorted evolution loop in terms of the polyharmonic field β(t) with parameters b=32,c=-7 and initial conditions (x,y)=(1,-1), (px,py)=(-30,30) subject to an external force of magnitude F=12 in the direction Ox. In solid (dotted) line, the unperturbed (perturbed) loop, closing after a period T=39. The perturbed loop exhibits a drifting effect in the opposite direction of the externally applied force. The dot at the centre of both figures represents the initial and final positions. We have commented before that those polyharmonic pulses (Eq. 29) greater than zero convey the adequate oscillatory motion to produce squeezing operations, given that the corresponding matrix h(t) fulfils |Σ|>2. However, it has not been surveyed yet a more general type of pulses that can be either greater or less than zero, as the amplitude β(t) with dotted line in Fig. 3, and whose associated matrix h(t) furnishes a stable motion in the region |Σ|<2. Accordingly, this class of pulses β(t) does bring the alternative of evolution loop operations by merely implementing the programme outlined for biharmonic fields. In this regard, please refer to the loops shown in Fig. 5. In panel A, we present a symmetric polyharmonic loop that is resistant to the application of an external potential. By resistant we mean that the loop is absolutely stable—as the examples shown in Fig. 2—and cannot be broken under the influence of forces coming from the outside, in fact, under these circumstances any solution of Eq. (33) will account as a phase factor. On the other hand, the loop in panel B has been deformed by an external field of constant magnitude (see the figure caption for details). It truly has a fuzzy point (Eq. 4) with noncommutative coordinates, then its fuzzy centre will suffer a drift in the transversal direction and eventually the loop could be broken. Discussion Our discussion above is fundamentally focused on the time dependent family of matrices h(t) responsible for a number of quantum control operations, although our scheme does merely constitute an approximation. Usually this type of imperfections are not entirely disappointing: In Paul’s ion traps6, for instance, the propagation of the electromagnetic signals through the device’s interiors is usually disregarded due to the pretty small size of the trap. However, adopting a rigorous attitude, even those softly changing potentials sensed on the trap surfaces must produce field corrections starting from 1c (post Newtonian) terms in the Einstein–Infeld–Hoffmann (EIH) approximation43. Below we are to briefly discuss a similar problem for softly changing magnetic fields β(t) such as Eq. (12) or (29). Throughout this report, we have been considering a time dependent, homogeneous magnetic field B(t,x) in a cylindric solenoid. It truly does not fulfil the Maxwell’s equations, yet it obeys a sequence of EIH approximations. To evaluate the conveyed errors, let us look for the exact time dependent vector potentials of a cylindrical solenoid in the form (Eq. 1):34 A(t,x)=12B(t,r)n×x=12B(t,r)-yx, where the magnetic field B instead of depending only on t, could also depend on the radius r of the solenoid’s cross section, and n is a unit vector in the field direction. To assure the relativistic sense of Eq. (34) we must assume □A(t,x)=4πcJ(t,x)=0, with the notation □=1c2∂2∂t2-∇2 (the d’Alembert operator). As easily seen, the application of the Laplacian to the right hand side of Eq. (34) is equivalent to the acts of the operator D=∂2∂r2+3r∂∂r to B(t, r) alone. Hence, □A(t,x) means:35 1c2∂2∂t2-DB(t,r)=0. The magnetic field B(t, r) must be analytic around Oz, then we preferably look for a solution in the form:36 B(t,r)=B0(t)+B2(t)r2+B4(t)r4+⋯ where B0(t)=B(t) is the homogeneous quasi-static approximation. Since Dr2n=4n(n+1)r2(n-1), the substitution of Eq. (36) into Eq. (35) yields:37 B2n(t)=14nn!(n+1)!1c2n∂2n∂t2nB(t). We now introduce a dimensionless time τ=tT, where T is some conventional time unit corresponding to the laboratory observation, and by writing Eq. (37) in terms of the time derivatives ∂∂τ one can reduce it to:38 B(t,r)=B(t)+18rcT2∂2B(t)∂2τ+124rcT4∂4B(t)∂4τ+⋯. We kept here the B(t) depending on the actual time t, in order to assure that all derivatives ∂n∂τnB(t) will be expressed in magnetic field units. The curious property of this formula is the absence of terms proportional to 1c —in fact, the field propagation law (Eq. 35) is solved exclusively in terms of extremely small contributions proportional to 1c2. It suggests that the superpositions of delicate wave fronts running towards the solenoid centre create a good approximation of the quasi-static theory. Even though, this does not explain how to create (or at least approximate) the first magnetic step B(t), to wake up the whole iterative series (Eq. 38). (That is, how to induce the homogeneous surface currents which do not depend on z, but depend on time in any desired way.) In the static case, the magnetic field inside the solenoid is generated by a stationary current I circulating around the surface, i.e., B=4πcΔIΔz. The question is, how to produce the circulating currents depending on time, but homogeneous (independent from z) on every section of the surface? Should the solenoid was constructed as a single spiral wire around the cylindrical surface, connected at its both extremes to the potential difference Φ(t), then even a subtle change δΦ(t) would be propagated along the solenoid as a current pulse, creating a softly changing but z-dependent field, instead of the desired quasi-static B(t). There is an alternative. One rather could consider a cylindrical surface of non-conducting material (e.g., glass) of radius a, charged uniformly with a surface density σ, so that each circular belt of 1cm of height contains a charge aσ. The experimental challenge is not extraordinary. Suppose the cylinder has a radius a=20cm, and it is rotating with an angular velocity ω=1s-1 around its symmetrical axis aligned to Oz and has 1C of charge at every horizontal belt of 1cm, then inside the cylinder will be generated a homogeneous magnetic field nB of intensity:39 B=4πcωRσ≃1.25gauss, at least in the post-post-Newtonian approximation. Hence, by employing the softly changing angular velocity ω→ω(t) one truly can generate the practically homogeneous magnetic field nB(t) of the quasi-static environment described by Eq. (38). Will such technique work? As a remark, the actual operation time tT→t given in Eq. (6) can be arbitrarily large (or short) and the fields (Eq. 7) arbitrarily strong (or weak). To form an idea of the real orders of magnitude required by our control operations, we shall present the Table 1 below. To this aim, we have also decided to check the magnitude of the Abraham-Lorentz radiative force. While the ordinary force of the variable oscillator trajectory is simply Fosc(t)=mx¨(t), the hypothetical radiative force is expressed as Frad(t)=mγx⃛(t), where γ=23e2mc3 is the particle dependent characteristic time. Using now the definitions of dimensional quantities (Eq. 6) we can compare the magnitudes of the conventional and radiative forces for the squeezing operations (see the last row in Table 1) finding the radiative ones extremely small albeit slowly increasing as T decreases (higher frequencies).Table 1 Physical conditions (in cgs) to achieve the squeezing transformations λx=-0.43 and λy=-1.125 on a proton moving inside a cylindrical solenoid. The reported quantities have been calculated regarding the example in Fig. 4. The physical magnitudes of q, p and v, with q=x=y and p=px=py correspond to the dimensionless values x′=y′=px′=py′=1 according to Eq. (6). Note how the incredibly modest strength of the magnetic field grows as the operation time T becomes shorter, to form an idea, a control operation of T=10-2s. would require half the strength of Earth’s magnetic field on the equator. In the last row it is reported the average ratio of the Abraham-Lorentz radiative force to the time-dependent oscillator forces, at different operation intervals T. T (s) 10-2 1 102 q (cm) 2.5×10-3 2.5×10-2 2.5×10-1 p (g cm s-1) 4.2×10-25 4.2×10-26 4.2×10-27 v (cm s-1) 2.5×10-1 2.5×10-2 2.5×10-3 Bmax (gauss) 1.5×10-2 1.5×10-4 1.5×10-6 Frad/Fosc 5×10-25 5×10-27 5×10-29 Our preliminary evaluations become valid whenever the physical size of the solenoid is realisable. In many laboratories the techniques to cool and trap ions are adequate to hardly study the atomic structure but insufficient for wider purposes, since the time dependent oscillator potential is created only in a strictly local scale, such is the case of a quadrupole trap whose immediate vicinity from the central axis is conformed, in some cases, just by four metal bars15. In our calculations we have assumed pure states, evolving under the influence of slowly changing external fields, without taking corrections for traces of matter across the ion traps or solenoids. Moreover, we have neglected the possible packet reflection or absorption by the laboratory walls. Suggesting an apparatus whose dimensions make any of these considerations to become insignificant. Actually, in case that particle would be absorbed by the (surface) walls, the problem leads to the fundamental question of the time operator which persists without a truly convincing solution even in case of flat surfaces. We thus hope that for wide traps, our proposal brings something of interest to the quantum control theory. As well, we have noticed that in works41,44 using the Ermakov-Milne invariants the operations are supposed to be faster that the times given here. In that regard, our contribution is slightly different: Its main goal is the use of a simple Toeplitz algebra to obtain the exact solutions (Eq. 29) whereas other approaches use the suggestive idea of frictionless driving, which transports the states without changing the eigenvalues of certain invariants. But could the more general variational methods be applied at the level of a θ(t) function? The question remains open. Acknowledgements This work is dedicated to the memory of Prof. Dr. Bogdan Mielnik, whom has deeply inspired this study. The results reported here are nothing but the materialisation of our early discussions. The financial support granted by CONACYT (National Council of Science and Technology, Mexico) is well appreciated. Author contributions J.F. has prepared the entire manuscript. Competing interests The author declares no competing interests. The original online version of this Article was revised: The original version of this Article contained a repeated error, where the Blackboard Bold 1 symbol did not display correctly in Equations 9, 17, 25, 26, and 33, and in the Results section. In addition, two corrections were made in the Introduction and the Results section. Full information regarding the corrections made can be found in the correction for this Article. Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Change history 8/30/2021 A Correction to this paper has been published: 10.1038/s41598-021-96665-1 ==== Refs References 1. Ramakrishna V Flores KL Rabitz H Ober RJ Quantum control by decompositions of su(2) Phys. Rev. A 2000 62 053409 10.1103/PhysRevA.62.053409 2. Ramakrishna V Ober RJ Flores KL Rabitz H Control of a coupled two-spin system without hard pulses Phys. Rev. A 2002 65 063405 10.1103/PhysRevA.65.063405 3. Schirmer SG Greentree AD Ramakrishna V Rabitz H Constructive control of quantum systems using factorization of unitary operators J. Phys. A Math. Gen. 2002 35 8315 8339 10.1088/0305-4470/35/39/313 4. Haroche S Entanglement, decoherence and the quantum classical boundary Phys. Today 1998 51 36 10.1063/1.882326 5. Nogues G Seeing a single photon without destroying it Nature 1999 400 239 242 10.1038/22275 6. Paul W Electromagnetic traps for charged and neutral particles Rev. Mod. Phys. 1990 62 531 540 10.1103/RevModPhys.62.531 7. Hu, C.-K. et al. Quantum thermodynamics in adiabatic open systems and its trapped-ion experimental realization. npj Quantum Inf. 6, 73, 10.1038/s41534-020-00300-2 (2020). 8. Flühmann C Encoding a qubit in a trapped-ion mechanical oscillator Nature 2019 566 513 517 10.1038/s41586-019-0960-6 30814715 9. Mielnik B Ramírez A Ion traps: Some semiclassical observations Phys. Scr. 2010 82 055002 10.1088/0031-8949/82/05/055002 10. Blais A Girvin SM Oliver WD Quantum information processing and quantum optics with circuit quantum electrodynamics Nat. Phys. 2020 16 247 256 10.1038/s41567-020-0806-z 11. Gessner M Smerzi A Pezzé L Multiparameter squeezing for optimal quantum enhancements in sensor networks Nat.Commun. 2020 11 3817 10.1038/s41467-020-17471-3 32733031 12. Burd, S. C. et al. Quantum amplification of mechanical oscillator motion. Science 364, 1163–1165, 10.1126/science.aaw2884 (2019). https://science.sciencemag.org/content/364/6446/1163.full.pdf. 13. Thorne KS Drever RWP Caves CM Zimmermann M Sandberg VD Quantum nondemolition measurements of harmonic oscillators Phys. Rev. Lett. 1978 40 667 671 10.1103/PhysRevLett.40.667 14. Braginsky, V. B., Vorontsov, Y. I. & Thorne, K. S. Quantum nondemolition measurements. Science 209, 547–557, 10.1126/science.209.4456.547 (1980). https://science.sciencemag.org/content/209/4456/547.full.pdf. 15. Thompson RI Harmon TJ Ball MG The rotating-saddle trap: A mechanical analogy to rf-electric-quadrupole ion trapping? Can. J. Phys. 2002 80 1433 1448 10.1139/p02-110 16. Emmanouilidou A Zhao X-G Ao P Niu Q Steering an eigenstate to a destination Phys. Rev. Lett. 2000 85 1626 1629 10.1103/PhysRevLett.85.1626 10970574 17. Mielnik B Global mobility of Schrödinger’s particle Rep. Math. Phys. 1977 12 331 339 10.1016/0034-4877(77)90031-3 18. Mielnik B Evolution loops J. Math. Phys. 1986 27 2290 2306 10.1063/1.527001 19. Fernández C, D. J. Geometric phases and mielnik’s evolution loops. Int. J.Theor. Phys. 33, 2037–2047, 10.1007/BF00675169 (1994). 20. Chen T Higher-Order Supersymmetry, in Quantum Mechanics, 187–188 2004 Netherlands, Dordrecht Springer 21. Harel G Akulin VM Complete control of hamiltonian quantum systems: Engineering of floquet evolution Phys. Rev. Lett. 1999 82 1 5 10.1103/PhysRevLett.82.1 22. Viola L Lloyd S Knill E Universal control of decoupled quantum systems Phys. Rev. Lett. 1999 83 4888 4891 10.1103/PhysRevLett.83.4888 23. Mancini S Manko V Tombesi P Symplectic tomography as classical approach to quantum systems Phys. Lett. A 1996 213 1 6 10.1016/0375-9601(96)00107-7 24. Castaños O López-Peña R Manko MA Manko VI Squeeze tomography of quantum states J. Phys. A Math. Gen. 2004 37 8529 8544 10.1088/0305-4470/37/35/009 25. Asorey M Generalized tomographic maps Phys. Rev. A 2008 77 042115 10.1103/PhysRevA.77.042115 26. Johnson MH Lippmann BA Motion in a constant magnetic field Phys. Rev. 1949 76 828 832 10.1103/PhysRev.76.828 27. Connes, A. Noncommutative Geometry, 1 edn (Springer, New York, 1994). 28. Bellissard, J., van Elst, A. & Schulz- Baldes, H. The noncommutative geometry of the quantum Hall effect. J. Math. Phys. 35, 5373–5451, 10.1063/1.530758 (1994). 29. Vagner ID Gvozdikov VM Wyder P Quantum mechanics of electrons in strong magnetic field HIT J. Sci. Eng. 2006 3 5 55 30. Ashtekar A Gravity and the quantum N. J. Phys. 2005 7 198 198 10.1088/1367-2630/7/1/198 31. Rovelli C A dialog on quantum gravity Int. J. Mod. Phys. D 2003 12 1509 1528 10.1142/S0218271803004304 32. Chalopin T Probing chiral edge dynamics and bulk topology of a synthetic hall system Nat. Phys. 2020 16 1017 1021 10.1038/s41567-020-0942-5 33. Landovitz LF Levine AM Schreiber WM Time-dependent harmonic oscillators Phys. Rev. A 1979 20 1162 1168 10.1103/PhysRevA.20.1162 34. Mielnik B Ramírez A Magnetic operations: A little fuzzy mechanics? Phys. Scr. 2011 84 045008 10.1088/0031-8949/84/04/045008 35. Hong-Yi F Zaidi HR Squeezing and frequency jump of a harmonic oscillator Phys. Rev. A 1988 37 2985 2988 10.1103/PhysRevA.37.2985 36. Mollow, B. R. & Glauber, R. J. Quantum theory of parametric amplification. I. Phys. Rev. 160, 1076–1096, 10.1103/PhysRev.160.1076 (1967). 37. Wolf, K. B. Geometric Optics on Phase Space, 1 edn (Springer, Berlin, 2004). 38. Vepsäläinen AP Impact of ionizing radiation on superconducting qubit coherence Nature 2020 584 551 556 10.1038/s41586-020-2619-8 32848227 39. Grübl G Dynamical squeezing in quantum mechanics J. Phys. A Math. Gen. 1989 22 3243 3252 10.1088/0305-4470/22/16/015 40. Suslov SK Dynamical invariants for variable quadratic hamiltonians Phys. Scr. 2010 81 055006 10.1088/0031-8949/81/05/055006 41. Berry, M. V. Quantal phase factors accompanying adiabatic changes. in Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, Vol. 392, 45–57 (1984). 42. Mielnik B Plebański J Combinatorial approach to baker-campbell-hausdorff exponents Ann. l’IHP Phys. Théor. 1970 12 215 254 43. Infeld L Motion and relativity 1960 New York Pergamon 44. Chen X Fast optimal frictionless atom cooling in harmonic traps: Shortcut to adiabaticity Phys. Rev. Lett. 2010 104 063002 10.1103/PhysRevLett.104.063002 20366818