Quantum simulations of time-dependent Hamiltonians beyond the quasi-static approximation

Existing approaches to analogue quantum simulations of time-dependent quantum systems rely on perturbative corrections to quantum simulations of time-independent quantum systems. We overcome this restriction to perturbative treatments with an approach based on flow equations and a multi-mode Fourier expansion. The potential of the quantum simulations that can be achieved with our approach is demonstrated with the pedagogical example of a Lambda-system and the quench in finite time through a quantum phase transition of a Chern insulator in a driven non-interacting Hubbard system. The example of the Lambda-system demonstrates the ability of our approach to describe situations beyond the validity of adiabatic approximations.


I. INTRODUCTION
There is an abundance of open questions in quantum physics, that we will most likely not be able to solve with classical computational means.Only the use of quantum simulators seems to allow us to overcome the computational complexity of many quantum mechanical manybody problems [1,2].
The hardware that is required to accurately mimic the dynamics induced by a given Hamiltonian is sufficiently advanced to use quantum simulators for problems that are outside the reach of classical computational hardware.Notable platforms include atomic gases in optical lattices [3], crystals of trapped ions [4], arrays of Rydberg atoms [5][6][7], and superconducting qubits [8] that can be used to quantum simulate strongly interacting Hubbard models [3], topologically non-trivial phases of matter [9], interacting quantum spin models [4,10] and quantum chemistry [11].While many such problems are defined in terms of a time-independent Hamiltonian, there is also a broad range of problems resultant from time-dependent Hamiltonians, such as laser-driven dynamics of electrons in molecules [12], quenches across boundaries between quantum phases [13], time crystals [14], diabatic switching between different Hamiltonians [15,16] or cycles of quantum thermodynamical machines [17].
The theory and experiments on quantum simulation, so far, have mostly focused on time-independent Hamiltonians [18,19].A crucial reason for this restriction in theoretical work is the rigorous footing that Floquet theory provides for the definition of an effective, timeindependent Hamiltonians, whereas the definition of an effective time-dependent Hamiltonian in a driven quantum system is much more problematic.
While generalizations of the Floquet theorem to aperiodically driven systems have proven difficult to find, they are also not necessary for purposes of time-dependent quantum simulations.An effective Hamiltonian can be defined in terms of any finite interval of system dynamics, and flow equations [20][21][22] provide a solid basis for this [23].In particular, they allow to ensure that the system dynamics is covered exactly by the effective Hamiltonian at periodic instances, even though the timedependence of the actual system has no such periodicity.
Despite the solid foundations that the flow equations provide for the definition of a time-dependent effective Hamiltonian [23,24], any practical construction requires a separation of time-scales, with fast time-dependencies resulting in desired effective processes, and a slow timescale for the time-dependence of the effective Hamiltonian.This, in turn, implies either very fast and strong driving or long duration of an experiment.The former unavoidably induces undesired processes, such as heating in the case of atomic gases [25][26][27], or leakage beyond the levels that define individual qubits [28,29], and the latter typically results in conflicts with coherence time [30].
The goal of this paper is thus to develop a framework that allows for quantum simulations of time-dependent quantum systems without the requirement of such a separation of time scales.

II. THE FORMAL FRAMEWORK
Hardly any Hamiltonian that one would want to quantum simulate can be exactly realized experimentally.It is rather necessary to realize a time-dependent Hamiltonian H(t) that induces dynamics which approximates the dynamics of interest.The central, underlying mechanism is the fact that a Hamiltonian defined in terms of a time-dependent unitary transformation U (t) can describe different physics than the underlying Hamiltonian H(t).It can thus be possible to realize a Hamiltonian in the frame defined by U (t), even if this Hamiltonian is practically out of reach in the laboratory frame [9,31,32].The propagator V U (t) induced by H U (t) reads V U (t) = U (t)V H (t) in terms of the propagator V H (t) induced by H(t).Since the propagators V U (t) and V H (t) differ by a factor U (t), it is essential for U (t) to reduce to the identity 1 when observations are being taken.The identification of a transformation U (t) that achieves this reduction arXiv:2305.17097v3[quant-ph] 12 Feb 2024 periodically is formalized in terms of the Floquet theory for time-independent effective Hamiltonians [18,33].In the case of aperiodically driven systems, the unitary U (t) will typically not be periodic, but it is essential that it reduces to the identity at well-defined points in time.This can be achieved with the framework of flow equations [22,23] that considers a family of unitaries U s (t) parametrised by a parameter s in terms of the differential equation ∂Us(t) ∂s = iη s (t)U s (t) with an hermitian generator η s (t).Associated with each such unitary U s (t) is a Hamiltonian H s (t) (following Eq. ( 1)) that satisfies the differential equation [21,22] ∂H s (t) The generator η s (t) is typically an explicit function of the flowing Hamiltonian H s (t), such that Eq. ( 2) is actually non-linear in H s (t).It typically has a stationary solution in the limit s → ∞, and the flowing Hamiltonian H s (t) with the initial condition H s=0 = H(t) approaches the effective Hamiltonian in this limit.
Through the explicit choice of generator η s (t) one can specify general properties that the effective Hamiltonian shall have.While, typically the effective Hamiltonian is expected to be time-independent, we will require that the effective Hamiltonian does not have any time-dependence associated with a fundamental driving frequency, but that it can still have time-dependence associated with some other frequencies.
The definition of the generator η s (t) is facilitated in term of a multi-mode Fourier series [34][35][36][37] H s (t) = m h m s e im•ωt for the flowing Hamiltonian H s (t), where m is a vector of integers, ω is a vector of mutually incommensurate frequencies (i.e.frequencies whose ratios are not rational bumbers), and the operators h m s are generalized Fourier coefficients.The number of frequencies, i.e. the dimension of m and ω can be chosen in accordance with the problem to-be-quantum-simulated; subsequent examples in this paper are based on driving with two fundamental frequencies.
The generator with m 0 = [0, m 2 , m 3 , . ..], ω 0 = [0, ω 2 , ω 3 , . ..] and real scalars f m contains the operators h m s and is thus a function of the flowing Hamiltonian H s (t).The factor (e im1ω1t − 1) ensures that the generator η(s, t) vanishes at any time t that is an integer multiple of the period 2π/ω 1 .This, in turn, guarantees that the unitary U (t) in Eq. (1) periodically reduces to the identity.
Besides the condition f m = −f −m that ensures that the generator η s (t) is hermitian, and the condition f m0 = 0 required for convergence, there is substantial freedom in the choice of the scalars f m , and this freedom of choice can be used to specify which of the frequencies in ω are contained in the effective Hamiltonian, and which frequencies are meant to play the role of driving in order to realize effective processes.In the following, we will use ω 1 as driving frequency such that all the other elements of ω correspond to time-dependencies in the effective Hamiltonian.
In order to construct the effective Hamiltonian explicitly, it is helpful to express the flow equation Eq. ( 2) in terms of the generalized Fourier amplitudes h m s .The explicit equation of motion for the terms h m0 s and for h m s with m 1 ̸ = 0 read They can be solved in the wellestablished high-frequency expansion [23,38] with the expansion coefficient 1/ω 1 .
Crucially, however, the factors m 0 •ω/ω 1 and m•ω/ω 1 do not need to be taken into account perturbatively.While the frequency ω 1 needs to be large as compared to the amplitudes in the system Hamiltonian for the highfrequency expansion to be valid, it does not need to be large as compared to the frequencies ω j (j > 1).Provided that the inequality m • ωf m > 0 is satisfied for m 1 ̸ = 0, the components h m s with m ̸ = m 0 will suffer from an exponential attenuation in the dynamics described by Eq. ( 4), such that they vanish in the limit s → ∞.The resulting effective Hamiltonian H e , thus only has components h m0 e = h m0 s→∞ , i.e. time-dependence in terms of ω 1 is no longer present.Similar to time-independent effective Hamiltonians, the high-frequency expansion of solutions of Eq. ( 4) can be specified in terms of nested commutators of the generalized Fourier amplitudes h m = h m s=0 of the system Hamiltonian.The lowest order h m0 e,0 of the effective Hamiltonian H e reads and it is independent of the choice of the constants f m .The next order, that also shares this independence is specified in Eq. (A1) in the Appendix.

III. EXAMPLES
A. The driven Λ-system A pedagogical example for the present framework is given by the realization of a time-dependent effective coupling between two low-lying states in a Λ-system with driving of two fundamental frequencies ω 1 and ω 2 .
The system Hamiltonian reads with the balanced superposition |+⟩ = (|1⟩ + |2⟩) / √ 2 of two degenerate ground states, and each of the driving functions Γ(t) and Ω(t) is a Fourier sum with fundamental frequency ω 2 .
At the lowest order in the high-frequency expansion, the effective Hamiltonian reads with the Fourier coefficients Ω p and Γ p of Ω(t) and Γ(t) and the ratio η = ω 2 /ω 1 of the two driving frequencies.
The dependence of H (0) e on η reflects the fact that ω 1 is not assumed to be large as compared to the second driving frequency ω 2 .If this assumption was made, the leading order of the effective Hamiltonian would be independent of η (in the quasi-static approximation), or it would contain perturbative corrections in η (beyond the quasi-static approximation), in contrast to the actual dependence in Eq. (7).
Resulting from the non-negligible variations in the Rabi-frequency, there is a direct coupling between the low-lying state |+⟩ and the excited state |3⟩, in contrast to the to regular case of the monochromatically driven Λ-system.This regular Floquet result is naturally contained in Eq. ( 7) in the limit η → 0, and in the absence of any resonant driving, i.e.Γ(t) = 0.
The explicit dependence of H (0) e on the driving parameters can also be used to identify driving profiles that ensure that no undesired excitations to the excited state |3⟩ occur.For any component Ω m of the off-resonant drive, the corresponding component Γ m of the resonant drive can be chosen such that H (0) e vanishes.Given this choice, the first order contribution to the effective Hamiltonian reads with the effective Rabi-frequency (9) In the limit η → 0, the effective Rabi frequency reduces to e (t) and underlying driving profile Ω(t).Insets a) and b)correspond to parameter values η = 1/(5 √ 3) and η = 1/ √ 7 respectively.The thin dashed green line depicts a normalized Gaussian profile as targeted time-dependence.The solid red line depicts the actual effective Rabi-frequency (normalized) following Eq.( 9) realized with 7 Fourier components (between −3 and 3); the solid blue line depicts the quasi-static approximation (η → 0) of the effective Rabi-frequency.The dotted dark red and dashed dotted dark green lines depict the real and imaginary part of the driving profile Ω(t).
i.e. to the quasi-static solution.
Eq. ( 9) permits to identify driving profiles Ω(t) that realize a desired time-dependent effective Rabi frequency Ω e (t).Fig. 1 depicts the exemplary case of a Gaussian time-dependence for the effective coupling between the two low-lying states.Insets (a) and (b) correspond to the parameter values η = 1/(5 √ 3) and η = 1/ √ 7 i.e. in one case the adiabatic approximation is expected to hold approximately while in the other the adiabatic approximation to be violated.The thin dashed green line depicts the desired time-dependence; the actual effective Rabi-frequency Ω e (t) can be made to approximate the desired time-dependence arbitrarily well, but a restriction to a finite number of 7 Fourier components Ω m results in a small deviation from the desired behavior.
The solid orange-red line depicts the quasi-static approximation Ω qs (t); in inset (a), this approximation is indeed good, but inset (b) shows that a clear separation of time-scales (i.e.η ≪ 1) is required for the quasi-static approximation to hold.
The dotted dark red and dashed dotted dark green lines depict the real and imaginary parts of the actual driving function Ω(t) in Eq. (6).Given the validity of the quasi-static approximation (Fig. 1(a)), the driving function Ω(t) is approximately mirror symmetric/antisymmetric around the mid-point of the depicted time-window, and in the quasi-static limit, this symmetry of the targeted Gaussian dependence is given exactly.Outside the regime of validity of the quasi-static approximation (Fig. 1(b)) this symmetry is clearly violated by Ω(t), highlighting that a diabatic increase of Ω(t) requires different driving than a diabatic decrease.

B. Quench across phase-transitions of a Chern insulator
A more involved example is given by the problem of crossing of a phase transition in a Chern insulator [13,39].It is defined by a driven non-interacting Hubbard Hamiltonian H = H B + H S with (11) c † j and c i are creation and annihilation operators on a hexagonal lattice and ⟨...⟩ denotes the nearest neighbours.This lattice is given by a triangular Bravais lattice and two-site basis, or, equivalently by two triangular sublattices, depicted in Fig. 2 by red empty (sub-lattice A) and black full (sub-lattice B) circles respectively.The onsite energies ∆ i are chosen such that all sites on sublattice A have the same onsite-energy ∆ and all sites on sub-lattice B have the opposite onsite-energy −∆.
There are three inequivalent directions of nearest neighbor tunnelling processes depicted by a 1 , a 2 and a 3 .The tunnelling rate J(t) for those processes does not depend on the direction in real space, but it is time-dependent.All tunnelling processes beyond nearest neighbors are neglected.The Hamiltonian H S in Eq. ( 11) describes shaking [9,[40][41][42] in terms of time-dependent onsite energies with driving amplitudes q k , driving frequencies ω k and phases δ k and δ ′ k [43], and the positions [x i , y i ] of the lattice sites.In the frame defined by the shaking term H S , the resulting effective Hamiltonian has the same type of tunnelling processes as H B in Eq. ( 11), but the tunnelling rates are renormalized, and they are generally complex.For suitable phase relations between the different tunnelling processes, the effective Hamiltonian captures the Haldane Model [9,44] with two topological non-trivial phases (with Chern number +1 and −1) and one topologically trivial phase (with Chern number 0).Deviations from these phase relations result in a deformation of the system's phase diagram [43] (typically depicted in terms on onsite energy and phase of next-nearest-neighbor tunnelling rate), but, they do not affect the existence of several phases with different topological properties.
An effective Hamiltonian with time-dependent onsite energy can be used to investigate the creation of topological defect generations [13] as the system is quenched 2: Geometry of a honeycomb lattice (left) and quenches through a phase diagram (right).The honeycomb lattice is comprised of two triangular sub-lattices, consist of black and hollow red dots.Directions of nearest neighbour and next nearest neighbour tunnelling processes are depicted with black and dashed red vectors.An exponential decrease of the tunnelling rate J(t) satisfying the conditions Eq. ( 14) allows exploration the two-dimensional phase diagram spanned by J 2 κ + /ω 1 ∆ and J 2 κ − /ω 1 ∆.Grey regions corresponds to Chern number being 0, while purple and blue regions correspond to C = ±1.Three quenches profiles in dashed lines are plotted based on three solutions of Eq. ( 14) with the sequence {q k , ω k , δ k − δ ′ k } inherited from [43].
through the phase boundary between topological trivial and non-trivial phases.In addition to the effective onsite energies ∆(t) and − ∆(t) for sub-lattices A and B, the effective Hamiltonian has nearest-neighbor tunnelling rates Jk (t) (from sub-lattices A to B along a k (k = 1, 2, 3) in Fig. 2), and the rates Gk (t) (− Gk (t)) of next-nearest neighbor tunnelling processes (along b k (k = 1, 2, 3) in Fig. 2) within sublattice A (B).All these quantities can be obtained without assuming the separation of timescales discussed above, but for the sake of clarity, the following discussion is focused on the dominant corrections to the quasi-static approximation.
The effective Hamiltonian in the frame induced by the shaking Hamiltonian H S (t) is characterized by the system parameters where C ∆ , C∆ , D k , Dk , E k and Ẽk (given in Eq. (B2) and Eq.(B4)) are time-independent scalars that are specific to the hexagonal lattice geometry and that depend on the amplitudes q k , frequencies ω k and phases δ k , δ ′ k of the shaking profile in Eq. (12).With a suitable time-dependence of the tunnelling rate J and the shaking profile, one can realize a broad range of time-dependencies in the effective Hamiltonian, as exemplified in the following.The lattice of the underlying Hamiltonian (Eq.( 11)), and many lattice models of interest such as the Haldane model or the Kitaev model [45], are invariant under a rotation of 2π/3.The shaking necessarily breaks this invariance so that the resulting effective Hamiltonian does not satisfy this symmetry for general driving patterns.The symmetry can, however, be partially recovered with suitably chosen shaking profiles.For the given timedependent tunnelling rate J(t) = J 0 e −γt for example, the condition ensures the desired symmetry for the nearest neighbor tunnelling, i.e.J1 = J2 = J3 (though next-to-nearest tunnelings are still anisotropic), and this condition can indeed be satisfied for tri-chromatic shaking profile, i.e.N = 3 in Eq. ( 12).The Chern number of the resultant system is given by ), where α k (k = 1, 2, 3) are the phases of complex NNN hopping rates [47].The parameters h ± that the Chern number depends on can also be expressed as h ± = ∆ − J 2 (t)κ ± /ω 1 , where κ ± = τ 1 + 2 k τk cos (α k ± 2π/3) with τ 1 and τk given in Eq. (B6) in the Supp.Mat. are timeindependent.The system can thus be characterized in terms of a phase diagram spanned by J 2 κ + /ω 1 ∆ and J 2 κ − /ω 1 ∆, as depicted in the right panel of Fig. 2.
Also shown are three exemplary solutions for quenches with exponential time-dependence that can be realised in terms of suitably modulated lattice shaking.The quenches depicted by dashed brown and dotted magenta lines correspond to initial conditions in topologically nontrivial phases, and a final point in a topologically trivial phase, whereas the start and the end point of the quench depicted with a solid green line lies in domains of topologically trivial phase, but the quench takes the system through an ordered phase.In all these cases, the rate γ in the tunnelling rate J(t) can be varied from a regime of adiabatic to diabatic quenches.

IV. RANGE OF APPLICABILITY
In order to gauge the range of applicability of the involved approximations, this section provides a comparison with numerically exact simulations utilising the QuSpin package [48,49].

A. Frequency regime
The present expansion is derived under the premise of the high-frequency expansion (i.e.ω 1 exceeds all relevant rates in the Hamiltonian and the inequalities ω 1 ≳ ω j . Since the time-dependence in the underlying Hamiltonian H(t) can include higher harmonics of the frequencies ω j , this implies that there can be frequency components in the time-dependent effective Hamiltonian that exceed the fundamental driving frequency ω 1 .This is corroborated in Fig. 3 that depicts the infidelity of the dynamics induced by the time-dependent effective Hamiltonian for a spin chain comprised of 16 interacting spins and for a one-dimensional Fermi-Hubbard system with 16 sites.
The underlying Hamiltonian for the spin chain reads with the time-dependent functions d x (t), d zz (t) and d yy (t).The underlying Hamiltonian for the Fermi-Hubbard chain with L sites reads with annihilation operator c i,σ of a Fermion of spin σ at site i and corresponding creation operator c † i,σ , and with time-dependent tunnelling amplitude J(t) and interaction rate U I (t).Both models are understood with periodic boundary conditions.The resulting translational invariance implies conservation of quasi-momentum, and the following discussion is focused on the subspace with zero quasi-momentum.
With the propagator U 0 (t) induced by the underlying Hamiltonian and the propagator U e (t) induced by the time-dependent effective Hamiltonian in first order in 1/ω 1 , the fidelity of the effective dynamics at t = 2π/ω 1 is defined as [50] where the factor dim(H) in terms of the dimension of the Hilbert space ensures a maximal value of 1 of the fidelity.Fig. 3 depicts the infidelity 1 − F (in log-scale) as a function of η = ω 2 /ω 1 .The data shown in Fig. 3, is based on random choices of the driving functions d x (t), d zz (t), d yy (t), J(t) and U I (t) in terms of bi-chromatic Fourier sums where the real parameter γ is the ratio between driving strength and driving frequency ω 1 ; i.e. the highfrequency approximation is valid for γ ≪ 1.The integers M 1 and M 2 specify the spectral width of the driving functions, and the complex numbers A pq are chosen at random from within the interval [−0.5, 0.5].All the driving parameters are available on [51].
The left-most data points (for η = 0) in each of the three panels correspond to the regular Floquet case with a time-independent effective Hamiltonian, and the infidelities obtained for η = 0 give an idea of what infidelity one can reasonably expect for a given value of γ.In large parts of the parameter space, the infidelities obtained for the time-dependent effective Hamiltonian are comparable or even lower than in the regular Floquet case.Only when η is approximated well by a rational number p/q with small integers p and q, there is a sizeable increase of the infidelity, i.e. a decrease in fidelity, as expected from the breakdown of multi-mode Floquet theory that requires incommensurate frequencies.The increase in infidelity is most pronounced for η ≃ 1/2 and for η ≃ 1, but also visible for η ≃ 1/4, η ≃ 1/3, η ≃ 2/3 and η ≃ 3/4 in insets (a) and (b).Inset (c) shows fewer instances of increased infidelity, since it is based on driving patterns with fewer ω 2 frequency components (M 1 = 4 and M 2 = 1 as opposed to M 1 = 3 and M 2 = 2 in insets (a) and (b)), highlighting that the accuracy of effective Hamiltonians at given values of η can be controlled through the spectral properties of the driving functions and a careful design of driving functions is essential for a challenging quantum simulation such as the quench dynamics discussed in Sec.III B.

B. Long-time dynamics and heating
A crucial issue with driven quantum systems is heating, and as shown in the following the heating resultant from the present driving schemes is comparable to the heating obtained in regular Floquet engineering of timeindependent effective Hamiltonians.
Heating is most suitable characterized in term of a comparison between the energy expectations of the driven system and the quantum-simulated system over several driving periods.The expectation value of the driven Hamiltonian for both of these dynamics is shown in Fig. 4 for the spin chain and the Fermi-Hubbard chain (with the same system parameters as above in Sec.IV A, and with the system initialized in ground state of the bare Hamiltonian for η = 0 at t = 0 and γ = 0.06).The spectral width of the driving functions are characterized by M 1 = 3 and M 2 = 1 for both the spin chain and the Fermi-Hubbard model.
Energy expectation values are depicted by dots (with different colors for different values of η.The dots are connected with straight lines to guide the eye and to distinguish between energy expectations of the driven system (solid) and the quantum-simulated system (dashed).
In both insets one can see that the energy difference between the driven and quantum-simulated systems are small as compared to the energy expectations, and that these differences to not grow noticeably with time or with η.This suggests that the heating caused by polychromatic driving is sufficiently low over several driving periods that is does not jeopardize an accurate quantum simulation.

V. CONCLUSIONS
While the field of quantum simulations of timeindependent quantum systems has demonstrated the readiness of quantum mechanical hardware for problems that are far outside the range of classically achievable simulations, ideas for quantum simulations of timedependent Hamiltonians [23] are only in its infancy.The ability to realize time-dependencies that do not need to be slow as compared to the driving time-scale and that can be designed with several fundamental frequencies enables the experimental realization of quantum simulations of a broad range of physical problems with explicit time-dependence in numerous state of the art platforms, such as sweeping between different quantum phases in the lattice Z 2 gauge theories [52] and the Hofstadter model [33] or realizing time-dependent Hamiltonians that exhibits Floquet symmetry-protected topological phases [53,54].The first-order component of the effective Hamiltonian is which along with the zeroth-order expression in the same limit, would be identical to the results using the framework in [23] in the Floquet gauge by further decomposing h (n1) (t) = m0 h (n1,m0) e im0•ω0t .This shows that the flow equation in the main text is consistent with the flow equation constructed in [23] under the same limit.Taylor expanding Eq. ( 5) and Eq.(A1) to the first order in ω/ω 1 ≡ η yields [h n , h (m1−n1,m0−n0) ]. (B1) Eq. (B1) then yields the explicit form of the coefficients (B6)