OPENMMF: a library for multimode driven quantum systems

OPENMMF is a numerical library designed to evaluate the dynamics of quantum systems with a discrete spectrum and driven by an arbitrary combination of harmonic couplings. The time-evolution operator of such systems is calculated as a multifrequency Fourier expansion, which results from expressing the time-dependent Schrödinger equation in the frequency domain [1]. The library provides a generic tool to study systems with arbitrary spectral composition, limited only by the available computational resources. OPENMMF includes functionalities to build dense and sparse matrix representations of the Hamiltonian and functions to calculate the micromotion operator and time/phase average of state populations. The library introduces a generalised notion of dressed state for systems with polychromatic driving. In this paper, we describe the design and functionality of OPENMMF, provide examples of its use and discuss its range of applicability in problems of current interest in quantum dynamics. The library is written in object-oriented style Fortran90 and includes a set of wrappers for C++ and Python..


Motivation and significance
The dynamics of quantum systems driven by polychromatic (or multimode) electromagnetic radiation is of interest to many domains of physical sciences. Over the last few decades, there has been a rapid development of experimental platforms where coherent dynamics of driven quantum systems has been explored with unprecedented precision and sensitivity. Examples of this include experiments with ultracold atomic ensembles [2,3], trapped ions [4], Nitrogen-Vacancy centres [5,6,7] and superconducting circuits [8,9,10]. In all these platforms, the state of an ensemble of quantum systems is manipulated by applying a combination of monochromatic driving fields, produced by lasers or electric/electronic devices. Technical advances of such controlling techniques have led to practical applications (e.g. highprecision spectroscopy [11] or Nuclear Magnetic Resonance [12]), and are instrumental to developing quantum technologies for sensing and computing [13].
A growing number of applications exploit the possibility to modify the properties of a quantum system through its interaction with electromagnetic radiation. The phenomena of electromagnetic induced transparency (EIT) is a textbook example of this situation [14]: The absorption spectrum of a three-level system is strongly modulated by coupling them through a monochromatic field. Expanding this idea to more complex quantum systems has lead to applications such as robust atomic clocks [15,16], ultra-sensitive field sensors with NV-centres [11], and to discovering new physical phenomena, e.g. topological features in driven many-body systems [17] and time crystals [18], which can be implemented with solid state platforms [19] and ultracold-atomic ensembles loaded into optical lattices [20].
As schematically shown in Fig. 1, the physical systems mentioned above have a similar structure: they are described by a discrete set of states {|i } and polychromatic (or multimode) couplings. The OPEN-MMF library, introduced in this work, provides tools to evaluate the time-evolution operator of such driven systems. The library implements the multimode Floquet method described originally by Ho, Chu and Tiet in [1], where the semiclassical time-dependent Schrödinger equation is substituted by a time-independent matrix representation in the frequency-domain. The library also provides functions to evaluate physical quantities that cannot be easily accessed from direct integration of the Schrödinger equation in the time-domain, such as the micromotion operator and dressed energies and states [21].
The OPENMMF library includes algorithms that enable its use for a wide variety of problems currently relevant in quantum physics. For example, in the context of quantum technologies, the library can be used to design sequences of electromagnetic pulses for qubit gates [22,7] as well as elements for manipulating matter-waves [23]. The library is also ready to assist emerging areas of research such time-domain quantum simulations [24,25,26] and (multimode) Floquet engineering [27], where harmonic driving of the system plays an essential role.

Time-evolution of multimode driven quantum systems
The main task of the OPENMMF library is to calculate the timeevolution operator, U (t , t), t > t, of systems whose Hamiltonian has Figure 1: (a) Schematic energy level structure of a generic quantum system. The basis of states consists of a discrete set of energy states, which define several bands according to the level energy spacing. Inter-and intra-band coupling is induced by electromagnetic radiation tuned at the corresponding frequencies, as indicated by the coupling terms. The wide variety of physical systems described by this model includes (b) trapped ions [4], (c) superconducting qubits [10] and (d) diamond NV-centres [5]. the form: where D is the dimension of the Hilbert space, E i,j defines the static component of the Hamiltonian H, V ,n i,j is the coupling between the states i and j oscillating at frequency nω (i.e. the n-th harmonic of the -th fundamental frequency ω ) and N is the number of non-commensurable frequencies.
To calculate the time-evolution operator, we generalise the Rotating (or Resonant) Wave Approximation (RWA) taking into account the complex time dependence of Eq. (1). For this, we rephrase the problem in terms of building a time-dependent unitary transformation, U F (t) to a new basis {|ī }, that leads to a time-independent and diagonal Hamiltonian,H. The standard quantum-mechanical transformation rule to the Schrödinger equation [28,29] leads to the condition: In the basis of states defined by this transformation, {|ī }. the time evolution operator is diagonal and has the form: which let us to calculate the time evolution operator in the original basis, {|i }, just by inverting the transformation U F (t) according to [29]: To formulate a fully defined computational problem, we express the micromotion operator U F (t) as the multifrequency Fourier series [1]: where ω = (ω 1 , ω 2 , . . . , ω N ) and n is a N -dimensional vector of integers.
Using this expansion in Eq. (2) and performing an integral over time, we obtain a fully defined eigenproblem for the eigenvaluesĒī and Fourier components of the unitary transformation U n i,Ē : where the couplings V ,n i,j define the values of V n i,j with the vector m = (0, . . . , m, . . . , 0), where the value m locates at the −th position. To obtain a finite matrix representation of this problem we truncate the sum over the number of modes of the Fourier expansion Eq. (5), which corresponds to the number of Floquet manifolds for each fundamental frequency. In Appendix A, we show an specific example of the shape of the matrix for a bichromatic driven problem.
This formulation to calculate the time-evolution operator is equivalent to the multimode Floquet representation of the Hamiltonian that introduces the extended Hilbert space |E i , n [1,21,30]. However, the semiclassical description presented here makes emphasis in the experimentally accessible states, which are commonly used to express the static part of the Hamiltonian Eq. (1).

Software description
The OPENMMF library is written in object-oriented style FORTRAN 90 and includes wrappers to enable its use with C++. The software architecture resembles a typical sequence of operations enabled by this library, which is shown in Fig. 2.
The starting point to use the OPENMMF library is defining parameters to build a matrix representation of the Hamiltonian components of Eq. (1). Then, the user can execute a series of tasks that lead to the time-evolution operator.
The components of the library can be classified into four groups: 1. Fortran modules that define as global variables relevant physical constants (physical _ constants) and properties of atomic systems (ATOMIC _ PROPERTIES).
2. Two Fortran derived data types (and corresponding C++ structures) that define the parameters required to build a matrix representation of the system's Hamiltonian (TYPE::ATOM and TYPE::MODE).

3.
A set of computational routines that perform distinct computational tasks, as in Fig. 2.

4.
A set of driver routines, each one to call a sequence of computational routines for solving standard problems with a single subroutine call.
To give a brief overview of the library structure and functionalities, we describe the last three groups. We refer the reader to the user

Setting the system configuration
The parameters required to build a matrix representation of the Hamiltonian Eq. (1) must be initialised using the following three data structures: 1. The Fortran derived type TYPE::ATOM, which contains one component to set the dimension of the Hilbert space (D _ Bare) and another one (declared as an array of type double), which specifies the energy spectrum of the system.
2. An array of integer type, with as many components as the number of fundamental driving frequencies plus one (i.e N +1). The entries of this array indicate the number of harmonics of each driving frequency, with the first component corresponding to the static component of the field.
3. An array of the Fortran derived type TYPE::MODE, with as many elements as the total number of driving frequencies plus one (i.e. 1 + N =1 n ). The components of this derived data type include the driving angular frequency (omega), complex matrices to store the explicit representation of the coupling matrix (V ), and the number of modes of the Fourier expansion N _ Floquet.

Typical computing sequences
Using the previous three data structures and an explicit representation of the components of the Hamiltonian, the user can proceed to build a finite matrix representation of the eigenproblem Eq. (6). This task can be done using the subroutine MULTIMODEFLOQUETMATRIX, which creates a square-matrix representation of the right-hand-side of Eq. (5).
The eigenvalues and eigenstates of the extended representation of the Hamiltonian can then be obtained using standard diagonalisation routines. The library includes wrapers to call the relevant subroutines of the LAPACK and the MKL-Intel libraries. However, the user is free to use other tools. The resulting eigenvalues and eigenvectors are then employed to calculate the time-evolution (U (t , t), Eq. (3)) and micromotion (U F (t), Eq. (5)) operators, using the routines MULTIMODETIMEEVOLUTIONOPERATOR and MICROMOTIONOPERATOR, respectively.
The openMMF library includes the driver subroutine TIMEEVOLUTIONOPERATOR, which allows the user to perform these three tasks (building the Hamiltonian representation, diagonalisation and evaluation of U (t , t)) with a single subroutine call after defining the components of the Hamiltonian.
To take advantage of the sparse nature of Eq. (6), the library also allows us to create a sparse representation of this equation and diagonalisation with the MKL-intel library. The names of the corresponding subroutines have appended the letters _ SP, and the user must ensure the correct allocation of memory.

Time-evolution in a dressed basis
The openMMF library provides tools to evaluate the dynamics of driven quantum systems in the so-called dressed basis. Figure 3 presents a situation where it is convenient to express the time evolution operator on a time-dependent basis that differs from the one used to represent the Hamiltonian. In this case, the driving modes applied to the quantum system can be classified into two groups: dressing and probing modes. Initially, the dressing modes are adiabatically switched on over a period of time T adiabatic , after which their amplitudes and frequencies remain constant.
This procedure is shown schematically in Fig. 3, with the red line and the colour gradient indicating the slow change of the amplitudes of the dressing modes. The adiabatic change of the dressing modes ensures that the basis of the static system (on the left-hand-side of the Fig. 3) deforms into the time-dependent basis defined by the transformation U F (t), Eq. (6) defined by only the dressing modes. Thus, for example, when at T = 0 the system is prepared in the state with the i-th energy, |i , after T adiabatic , this state becomes the corresponding linear superposition associated with theī-th eigenvalueĒī.
After a period of evolution driven by the dressing modes, the system is probed by pulses of the probing subset, while the dressing modes are kept constant. During the duration of the probing pulses, t ∈ [T Probe , T Probe + ∆T ], the time-evolution of the system is determined by a Hamiltonian of the form Eq. (1), which must include both the dressing and probing modes. The OPENMMF evaluates the time evolution operator initially in the bare basis and includes functions to define the micromotion operator associated with the subset of dressing fields and the corresponding transformation of the time-evolution operator. An example of this situation is shown in the following section.

Functionality for spin systems
Spin systems are of particular interest due to the simplicity of their energy spectrum and their experimental feasibility with a diverse set of experimental platforms. For these systems, the couplings to harmonic drivings are of the form: where S i is the i-th Cartesian component of angular momentum operator, whose matrix representation follows the standard definition in Sec. 3.5 of [32]. Each one of these terms is defined by six real parameters: three phases (φ ,n ) and three amplitudes (V ,n i ), which are included as components of the Fortran data structure TYPE::MODE. The OPENMMF library includes the subroutine SETHAMILTONIANCOMPONETS that uses these parameters to define Hamiltonian components Eq. (7).
The library can also handle the energy configuration of systems composed by two spins S 1 and S 2 , coupled via the dipole interaction H dipole = AS 1 · S 2 and each spin interactions as above in Eq. (7). The library evaluates the Hamiltonian and time-evolution operator in the eigenbasis of the dipole interaction, with the energy scale defined by the dipole coupling constant A. The libary defines the matrix representation of the coupled system using the uncoupled basis |S 1 , S 2 , m 1 , m 2 and the definitions in the Sec. 3.5 of [32]. For details of the implementation, werefeer the reader to the User guide [31]. This type of system of coupled system is present in the electronic ground states of alkali atoms [33] and single-charged atomic ions [34], which are currently the focus of intense experimental research in quantum technology. The interaction between the nuclear and electronic magnetic dipole moments defines two sub-sets of states separated by an energy gap of order A, which has enabled the use of radio-frequency and microwaves for coherent control of quantum states [34,33]. The library let us to readily evaluate the time-evolution of the atomic states in such situations [35].

Illustrative examples
The main features of the library are illustrated with two examples, which are relevant across many experimental platforms. We present a selection of examples that illustrate the potential of OPENMMF for a wide range of applications. The example in Sec. 4.1 focuses on a two-level systems, which can be extended to study quantum control of quantum systems with a small number of consituents. In Section we study a modified 4.1. Rabi oscillations of a dressed qubit As our first example, we consider the time-evolution of a two-level system (or qubit) driven by two harmonic drives with the Hamiltonian: where the subindice D (P ) indicates that the mode is a Dressing (Probing) the qubit, in the sense defined in section 3.3.
In this case, the data structures required to setup the matrix representation of the Hamiltonian components are declared with the instructions: where the first and second lines produced the time-dependent micromotion operators of the dressed basis, and the third line executes the basis transformation. As a concrete example of the convinience of the dressed basis, in Fig. 4 we plot the time evolution of the probability of remaining in an initially pure state as a function of time and the frequency of the second driving field. Fast oscillations of this probability are displayed in the original qubit basis, while a slower dynamics is observed in the dressed basis. Note the resonant feature at the frequency ω P = Ω x /2.0 [21].

Harmonically driven 1-D lattice
Technical advances in atomic and optical physics have enabled experiments with ultracold atomic ensembles loaded in spatially periodic potentials landscapes defined by optical lattices [27]. Controlled temporal modulation of the lattice parameters can induce non-equilibrium phenomena of interest for quantum simulations, and the experimental demonstration of such phenomena is leading to rapid advances in this area. Many important dynamical features of such driven systems are diagnosed via the dependence of eigenvalues and eigenvectors of the Hamiltonian in the frequency-domain with respect to the variation of the driving parameters [36,37].
To illustrate the use of the OPENMMF library for these problems, here we consider a spatially periodic one-dimensional lattice with a time-periodic modulation of the next-neighbour hopping. The lattice has D/2 unit cells, each containing two sites with energy offsets ±V 0 .
The corresponding Hamiltonian reads: where the first terms defines the energy offset at each lattice site, the second line defines a site-dependent tunneling rate between neighbouring sites and the last term defines a uniform (i.e. site-independent) modulation of the near-neighbour tunneling.
The three Fortran data structures containing the parameters of the Hamiltonian are initialised as before, but here only two modes are needed (the static component plus one driving mode). The static and driving components of the Hamiltonian are declared with the following set of instructions: The second line calls a wrapper to a diagonalisation subroutine from the LAPACK library. The resulting eigenvectors are stored in the matrix U _ F and the eigenvalues in the array E _ FLOQUET. A typical example of this calculation is shown in Fig. 5, where the set of eigenenergies are plotted as a function of the driving frequency.
In addition to evaluation of the Floquet spectrum, the numerical calculation of the time-evolution operator is of interest for Hamiltonian engineering with atomic ensembles in atomic lattices [38]. This capability of the OPENMMF can be used to evaluate corrections to analytical approximations of the Effective Hamiltonian [38], and their effects on the on the dynamics of the atomc cloud [39].

Impact
Considering that the Hamiltonian in Eq. (1) describes a wide diversity of physical systems, the OPENMMF library constitutes a generalpurpose tool for studying quantum dynamics of driven systems. The library is open-source and provides a ready-to-use interface, which can be exploited by users with a broad range of expertise ranging from graduate students to established researchers. An important feature of the library is the set of functions that allow the user to evaluate the time-evolution on a dressed basis, which is of current interest across different fields. To illustrate the degree of generality of the library, and therefore its potential impact, here we mention three specific research areas where the OPENMMF library can find the direct application.

Coherent control of quantum systems
Precise control over the state of a quantum system is key for quantum technologies. An important operational step in such applications is the precise creation of quantum superposition, which is generally achieved by interleaving periods of free and driven evolution described by a Hamiltonian of the form Eq. 1 [40]. Decoherence induced by the interaction of the system with the environment is a major obstacle to achieve the levels of control required for the realisation of practical applications. However, superposition of dresses states can be more robust to environmental noise than bare states [41,42]. This mechanism to reduce the sensitivity of a quantum system to noise is known as Continuous Dynamical Decoupling and, to be implemented, it requires an experimental sequence that follows the path shown in Fig. 3. This technique has found application improving the coherence of atomic clocks [15,43] and qubit gates [34,44], as well as for the sensitivity and spatial resolution of low-frequency magnetometers operating [9,11,45].

Time-domain quantum simulations
The frequency-domain representation of the Hamiltonian Eq. (1) with N mutually irrational frequencies is given by Eq. (6). This representation can be interpreted as the Hamiltonian of an N -dimensional lattice with sites label by the index n, a DC force equal to ω and couplings between sites controlled by the driving modes. The typical signatures of dynamics of the simulated lattice (e.g. band structure, conductivity, edge-states, etc.) arise in the time-evolution of the system, which is directly accessible with the OPENMMF library. Topologically non-trivial states [25,46] and Anderson localisation in quasi-crystals [26] are two phenomena accessible with this scheme for quantum simulations, which could be realised with a significantly simplified experimental setup.
Out-of-equilibrium phases of driven interacting systems The behaviour of driven systems can be qualitatively very different from their static counterparts. Various driving schemes have been recently used to engineer effective Hamiltonian displaying out-of-equilibrium phases that are impossible to realise with static systems [47]. Such applications are subject of intense experimental and theoretical research with ultracold atomic gases loaded in optical lattices [38], where the Hamiltonian and driving schemes are similar to those in the example in Sec. 4.2. In the single-particle case, the dynamics in the limits of high frequency or strong and resonant driving are relatively well understood [39], and novel phases of matter in such regimes have been proposed and realised in recent years [36,2]. However, there are several fundamental questions still to be addressed in the context of many-body driven systems, which emerge from the interplay between interactions, driving and symmetries/disorder [26,48]. The OPENMMF library provides a tool to investigate such phenomena, as we show in our last illustrative example. Though current computational limitations allows only to investigate systems with relatively small Hilbert space dimension, such studies provide a solid understanding of the non-equilibrium dynamics of driven systems [49,50,51] and the library provides a platform of simple usability.

Conclusions
In this paper, we introduce the OPENMMF library as a tool to calculate the time-evolution operator of quantum systems with discrete spectrum driven by harmonic couplings. The library implements the multimode Floquet expansion of a time-dependent Hamiltonian and provides several functionalities relevant for its analysis. In particular, the library allows the user to evaluate the micromotion operator associated with arbitrary drivings built as a superposition of harmonic terms. In passing, this procedure generalises the notion of dressed states to such situations.
The general formulation implemented in the library enables its use for applications where determining the time-evolution operator is rel-evant, such as designing optimal control sequences, tailoring the response of quantum systems and studying non-equilibrium quantum dynamics and Hamiltonian engineering. The code is written in Fortran 90, and we provide wrappers for C++ and python, making it widely accessible and easy usage.
and the full matrix results in:

Python wrapper
As an experimental feature, the library includes wrappers to use with Python through calls to the C++ wrapper functions using ctypes.