Solvable Model of a Generic Driven Mixture of Trapped Bose–Einstein Condensates and Properties of a Many-Boson Floquet State at the Limit of an Infinite Number of Particles

A solvable model of a periodically driven trapped mixture of Bose–Einstein condensates, consisting of N1 interacting bosons of mass m1 driven by a force of amplitude fL,1 and N2 interacting bosons of mass m2 driven by a force of amplitude fL,2, is presented. The model generalizes the harmonic-interaction model for mixtures to the time-dependent domain. The resulting many-particle ground Floquet wavefunction and quasienergy, as well as the time-dependent densities and reduced density matrices, are prescribed explicitly and analyzed at the many-body and mean-field levels of theory for finite systems and at the limit of an infinite number of particles. We prove that the time-dependent densities per particle are given at the limit of an infinite number of particles by their respective mean-field quantities, and that the time-dependent reduced one-particle and two-particle density matrices per particle of the driven mixture are 100% condensed. Interestingly, the quasienergy per particle does not coincide with the mean-field value at this limit, unless the relative center-of-mass coordinate of the two Bose–Einstein condensates is not activated by the driving forces fL,1 and fL,2. As an application, we investigate the imprinting of angular momentum and its fluctuations when steering a Bose–Einstein condensate by an interacting bosonic impurity and the resulting modes of rotations. Whereas the expectation values per particle of the angular-momentum operator for the many-body and mean-field solutions coincide at the limit of an infinite number of particles, the respective fluctuations can differ substantially. The results are analyzed in terms of the transformation properties of the angular-momentum operator under translations and boosts, and as a function of the interactions between the particles. Implications are briefly discussed.


Introduction
Quantum dynamics is central to nuclear, atomic, molecular, and optical physics [1][2][3][4][5][6][7][8]. In order to theoretically describe the dynamics of a quantum many-particle system, it is generally required to solve and analyze the time-dependent Schrödinger equation. Usually, a numerical treatment of a faithful representation of the time-dependent system under investigation is the only available tool, whereas the existence of solvable time-dependent many-particle models is scarcer.
A particularly exciting arena for quantum dynamics across many disciplines is that in which the system's Hamiltonian is periodic in time . The time-periodic Hamiltonian defines the Floquet Hamiltonian and governs the Floquet eigenvalue equation in an extended Hilbert space [11,12], in analogy to the static Hamiltonian, which governs the Schrödinger eigenvalue equation in the standard Hilbert space of a time-independent system.
In the field of Bose-Einstein condensates [4,6], the relation between the many-body and mean-field treatments of bosonic systems at the limit of an infinite number of particles has drawn much attention [59][60][61][62][63][64][65][66][67][68][69][70][71][72][73][74]. At the infinite-particle-number limit, the interaction parameter (i.e., the product of the number of bosons times the interaction strength) is held fixed while the number of particles is increased to infinity. Under well-defined conditions, it has been proved for the ground state of trapped bosons that the many-body reduced one-particle and two-particle density matrices [75,76] per particle are 100% condensed and coincide with the respective Gross-Pitaevskii results [60,61]. Similarly, the many-body energy per particle as well as the one-particle and two-particle densities per particle coincide in this limit with the Gross-Pitaevskii results. Analogous results exist for the dynamics of interacting bosons with time-independent Hamiltonians at a definite propagation time [62,63]. Generalizations for trapped mixtures of two kinds of identical bosons cover the many-body energy per particle, intraspecies and interspecies densities per particle, and reduced density matrices per particle, and show coincidence in the limit of an infinite number of particles with the Gross-Pitaevskii results, 100% condensation of each of the species, and separability of interspecies quantities in the mixture [67][68][69].
However, differences between many-body and mean-field results of fully condensed bosonic systems can emerge in other quantities, in particular, variances of operators, such as the position, momentum, and angular-momentum operators [64,65,72,74], and stem from the differences between the many-body and mean-field wavefunctions themselves [66,70]. The differences in the variances of observables are instrumental in defining and quantifying correlations at the infinite-particle-number limit; for instance, distinguishing between 100% condensed trapped bosons, either in the ground state or undergoing dynamics in two spatial dimensions, in which the anisotropy along the x and y directions of their position variance per particle and density per particle are alike, and 100% condensed bosons, in which this anisotropy is opposite [71,73].
In this work, we present a solvable model of a periodically driven trapped mixture of Bose-Einstein condensates. Mixtures of Bose-Einstein condensates have been a vivid topic . Here, we introduce the driven harmonic-interaction model for mixtures, that is, N 1 bosons of mass m 1 interacting by harmonic interparticle interaction of strength λ 1 and driven by a force of amplitude f L,1 , and N 2 bosons of mass m 2 interacting by harmonic interparticle interaction of strength λ 2 and driven by a force of amplitude f L,1 . All bosons are trapped in a harmonic potential of frequency ω and, furthermore, bosons of species 1 and species 2 interact with each other by harmonic interparticle interaction of strength λ 12 . The model generalizes the harmonic-interaction model for mixtures to the time-dependent domain, and builds on a long-standing interest in harmonic-interaction models [69,[105][106][107][108][109][110][111][112][113][114][115][116][117][118][119][120][121].
Besides solving the model and discussing properties, our main purpose is to investigate the infinite particle limit of a Floquet many-boson state, namely, the connection between the many-body and mean-field time-dependent reduced density matrices, time-dependent densities, quasienergies, and, eventually, variances. In this way, we expand on studies from the literature on the ground state and dynamics with time-independent Hamiltonians of bosonic systems to the realm of Floquet bosonic systems. Both the many-body and mean-field solutions of the trapped driven mixture are required and derived in closed form, making the proofs of properties explicit. As an application, the expectation value and variance of the many-boson angular momentum operator when a Bose-Einstein condensate is steered by interacting bosons consisting of a different species are explored fully analytically within our model. The modes of rotations are investigated as a function of the steering frequency and interspecies interaction, and the differences between the many-body and mean-field variances are explored.
The structure of the paper is as follows. In Section 2, we consider the driven single-species harmonic-interaction model, solve it at the many-body and mean-field levels, and derive working tools. In Section 3, we present the driven harmonic-interaction model for mixtures, solve it at the many-body and mean-field levels of theory, discuss its properties, and analyze the limit of an infinite number of particles. In Section 4, we analyze the imprinting of angular momentum when one of the species is steered and embedded in a second species that is not steered, and discuss the many-body and mean-field angular-momentum variances at the infinite-particle-number limit. In Section 5, a summary and outlook are given. Finally, Appendices A, B and C collect supplemental results, details of derivations, and further tools.

The Driven Harmonic-Interaction Model
Consider the driven classical harmonic oscillator in one spatial dimension. The classical equation of motion and its particular solution are given by Here and hereafter, the amplitude of the driving force satisfies f L = 0, the driving frequency ω L > 0, and the trapping frequency w > 0; for convenience, the off-resonance condition ω L = ω holds.
For ω L < ω, the amplitude of oscillations a is in phase with f L , whereas for ω L > ω and the high-frequency limit, a has the opposite phase to f L . The particular solution x(t) appears in the driven quantum harmonic oscillator problem and in its many-boson generalization, the driven harmonic-interaction model. The driven quantum harmonic oscillator [122,123] obeys the time-dependent Schrödinger equation where, throughout this work,h = 1. It has a set of Floquet solutions corresponding to the eigenfunctions of the non-driven system. The solution corresponding to the ground state of the harmonic oscillator, hereafter termed, for brevity, the ground Floquet solution, is given by where is the time-dependent phase and [123] is an eigenvalue of the Floquet Hamiltonian, ĥ (x, t) − i ∂ ∂t ψ (x, t) = ε Fψ (x, t), and called quasienergy.
The Floquet quasienergy stateψ(x, t) is the periodic part of ψ(x, t) = e −iε F tψ (x, t). There are three contributions to the time-dependent phase (4) and, hence, to the quasienergy (5): The first is the energy of the stationary problem, the second originates, upon acting with the kinetic-energy operator, from the spatial-dependent phase of the time-dependent wavefunction, and the third results from completing the square in the time-dependent Hamiltonian. The latter's contribution and its renormalization by interaction in particular would play an interesting role in the driven many-boson problem. In the present work, we study interacting many-boson generalizations of the driven harmonic oscillator-explicitly, time-dependent harmonic-interaction models-and analyze their ground Floquet solutions at the many-body and mean-field levels of theory.
Consider now the single-species harmonic-interaction model [108] plus a driving term, describing N harmonically interacting trapped bosons driven by a time-periodic force of amplitude f L . Using the Jacobi coordinates the time-dependent Hamiltonian is readily diagonalized, and describes a single driven center-of-mass oscillator and N − 1 time-independent relative coordinates' degrees of freedom. The frequency Ω = ω 2 + 2λN m of the relative coordinates is dressed by the interparticle interaction. The (driven) system is bound provided λ > − mω 2 2N , i.e., the interaction is either attractive (λ > 0) or not too repulsive (λ < 0).
The normalized ground Floquet solution of the time-dependent many-boson Schrödinger equation is given in terms of the Jacobi coordinates by where is the time-dependent phase, and is the many-boson quasienergy. Comparison of the many-boson (11) and single-particle (5) quasienergies shows that they have a similar structure. The later emanates from the fact that for both systems the driving force is coupled to the center of mass only, the frequency of which being the trapping frequency ω. The many-boson quasienergy stateΨ(Q 1 , . . . , Q N , t) is an eigenfunction of the many-body Floquet Hamiltonian, To proceed, we transform from the Jacobi coordinates, for which the wavefunction is separable, to the laboratory frame, and get In obtaining the wavefunction in the laboratory frame, we made use of the relations that and that the relative Jacobi coordinates are invariant under the translation of all coordinates by x(t).
The determination of the time-dependent reduced density matrices follows in essentially the same way as for the undriven system [108]. Importantly, since the driving is coupled to the center-of-mass coordinate only, all coordinates are translated by the same amplitude, the classical amplitude x(t), and the coherence properties of the driven system, like the eigenvalues of the time-dependent reduced one-particle density matrix, are exactly those of the undriven system; in particular, these eigenvalues are time independent. Explicitly, using the normalized time-dependent wavefunction and from the N-particle density matrix, the reduced one-particle density matrix ρ (1) is defined. Performing the integrations consecutively over x j , x j = x j , j = N, . . . , 2, for which the spatial-dependent phase factor falls, and while changing the variable x j − x(t) → x j , the final result reads The diagonal of the time-dependent reduced one-particle density matrix is the time-dependent density and is given by For the respective two-particle quantities, See Appendix A.
So far, we have presented the solution of the driven harmonic-oscillator model at the many-body level of theory. The similar structure of the time-dependent reduced one-particle density matrix per particle of the ground Floquet state (15) and the static reduced one-particle density matrix per particle of the ground state [108] already tells one that in the limit of an infinite number of particles, the former, like the latter, is 100% condensed. To discuss properties of the driven system at the infinite-particle-number limit, we now derive the solution at the mean-field level of theory. This would allow us to explicitly compare the time-dependent reduced density matrices per particle and quasienergies per particle, and prescribe the connection between the many-body and mean-field levels of theory for driven Bose-Einstein condensates, at least within the driven harmonic-interaction model.
Starting from the Gross-Pitaevskii mean-field ansatz for the time-dependent many-boson wavefunction, and using the time-dependent variational principle, in either of its forms, the Gross-Pitaevskii mean-field equation for the driven system takes on the form where Λ = λ(N − 1) is the mean-field interaction parameter. Solving the Gross-Pitaevskii Equation (18) combines the approach used in the treatment of the non-interacting driven oscillator and that of the undriven static mean-field system, and is detailed in Appendix B. The final result can be cast as where Ω GP = ω 2 + 2Λ m is the mean-field interaction-dressed frequency and is a time-dependent phase. Here, µ F can, in analogy to the static case, be called the quasichemical potential, and is an eigenvalue of the mean-field Floquet Hamiltonian, Analogously,φ(x, t), which may be called the mean-field quasichemical-potential state, is the periodic Since |φ(x, t)| 2 = |φ(x, t)| 2 , the mean-field Floquet Hamiltonian (21) can be defined solely in terms ofφ(x, t). Of course, the Gross-Pitaevskii solution φ(x, t) boils down to ψ(x, t) for non-interacting particles, Λ = 0. The determination of the mean-field wavefunction (17) allows one to prescribe the reduced density matrices. Thus, the time-dependent reduced one-particle density matrix at the mean-field level reads and its diagonal, the density, is ρ (1) Next, letΦ(x 1 , . . . , x N , t) = ∏ N j=1φ j (x j , t) be the periodic part of the mean-field many-boson wavefunction (17), i.e., Φ(x 1 , . . . , x N , t) = e −iNµ F tΦ (x 1 , . . . , x N , t). Then, the mean-field Floquet energy per particle is defined by sandwiching the many-boson Floquet Hamiltonian, in analogy to a Bose-Einstein condensate in its ground state for which the mean-field energy per particle is obtained by sandwiching the static many-boson Hamiltonian with the static mean-field wavefunction. We are now in the position to examine the relations between the many-body and mean-field solutions of the driven interacting bosons at the limit of an infinite number of particles. Thus, the following relations readily hold: for the time-dependent reduced one-particle density matrix per particle, lim N→∞ GP (x, t) for the density per particle, and for the quasienergy per particle, thus generalizing the literature results for the ground state of trapped interacting bosons at the limit of an infinite number of particles.

The Driven Harmonic-Interaction Model for Mixtures
Consider a driven mixture of two distinguishable kinds of identical bosons, which we denote as species 1 and species 2. The bosons are trapped in a one-dimensional harmonic potential of frequency ω and interact via harmonic particle-particle interactions. In the present work, we only deal with the ground Floquet solution of the driven mixture. We solve for the case of a generic driven mixture, i.e., a mixture comprising N 1 bosons of species 1 and mass m 1 and N 2 bosons of species 2 and mass m 2 . The total number of bosons is denoted by N = N 1 + N 2 . The species are driven independently by forces of amplitude f L,1 and f L,2 , respectively. Furthermore, the two intraspecies interaction strengths are denoted by λ 1 and λ 2 , and the interspecies interaction strength by λ 12 . Positive values of λ 1 , λ 2 , and λ 12 mean attractive particle-particle interactions, whereas negative values imply repulsive interactions. In what follows, we solve the driven harmonic-interaction model for mixtures at the many-body and mean-field levels, obtain the respective quasienergies, construct the time-dependent intraspecies and interspecies reduced density matrices and densities, discuss properties of these quantities as a function of the interactions between the bosons and strengths of the driving forces, and establish connections between the many-body and mean-field Floquet quantities at the infinite-particle-number limit.

Many-Body Solution and Time-Dependent Reduced Density Matrices
The Hamiltonian of the driven mixture in the laboratory frame is then that of the harmonic-interaction model for mixtures [69] plus two driving terms, and readŝ Here, the coordinates x 1,j with subscript 1 stand for species 1 bosons and x 2,k with subscript 2 for species 2 bosons. In the present section, we solve the one-dimensional driven mixture, where the amplitudes of the two driving forces f L,1 and f L,2 are general. In Section 4, we examine an application in two spatial dimensions with f L,2 = 0, i.e., when only species 1 in the interacting mixture is driven.
Using the Jacobi coordinates [109] x 2,j , the Hamiltonian (26) transforms to the diagonal form where the masses M 12 = m 1 m 2 M , M = N 1 m 1 + N 2 m 2 , and inverse relations ∑ N 1 The transformed Hamiltonian of the mixture (28) is that of N uncoupled harmonic oscillators, two driven oscillators built on the center-of-mass Q N and relative center-of-mass Q N−1 coordinates , and N − 2 time-independent oscillators for the relative coordinates of the individual species. The frequencies of the relative coordinates are dressed by the interactions, , and have to be positive for the driven mixture to be bound, implying the following conditions for the interspecies λ 12 > − M 12 ω 2 2 and the intraspecies λ 1 > − N 2 We can now proceed and prescribe the normalized ground Floquet wavefunction of (28), where are the time-dependent oscillating amplitudes of the relative center-of-mass and center-of-mass coordinates, respectively, is the time-dependent phase, and is the many-particle quasienergy. The quasienergy for the driven mixture (33) is, of course, more intricate than the above quasienergy of the driven Bose-Einstein condensate alone (see (11)). In particular, the appearance of the interspecies interaction-dependent resonance is appealing. For simplicity, in what follows, we assume the off-resonance condition ω L = Ω 12 with the relative center-of-mass frequency, in addition to the above-used off-resonance condition ω L = ω with the center-of-mass frequency. It is instructive to analyze the solution of the driven mixture. We begin with the poles. There are two poles in the quasienergy; the first is associated with driving the center-of-mass coordinate Q N of the mixture, and the second with driving the relative center-of-mass coordinate Q N−1 . In the absence of interspecies interaction, the two poles coincide, Ω 12 = ω, and the solution boils down to two independent driven Bose-Einstein condensates. Let us move to the zeros. For general forces f L,1 and f L,2 , both the Q N and Q N−1 oscillations are activated. There are, however, particular relations between the forces and parameters of the mixture for which either of the oscillations is not activated, i.e., there is zero contribution to the quasienergy. Explicitly, for the relation f L,1 m 1 = f L,2 m 2 between the forces and masses, the relative center-of-mass coordinate is not activated, whereas for the relation f L,1 N 1 + f L,2 N 2 = 0 between the forces and numbers of particles, the center-of-mass coordinate is not activated. Importantly, both modes of oscillations are (non-trivially) activated if only species 1 is driven, i.e., when f L,2 = 0, due to the interspecies interaction. We shall be investigating this particular case in two spatial dimensions in Section 4, focusing on angular-momentum matters.
To move from the Jacobi coordinate to the laboratory-frame representation for the driven mixture, one must express the time-dependent oscillation amplitudes Q N−1 (t) and Q N (t) in the former representation in terms of corresponding amplitudes for each of the species in the latter. Since each species is made of identical particles, only two such amplitudes are required, which we shall denote by x 1 (t) and x 2 (t) for species 1 and 2, respectively. Recall that the Jacobi coordinates of the relative motions, Q j , j = 1, . . . , N − 2, are invariant under translations of each species. Using the definition of the relative center-of-mass and center-of-mass Jacobi coordinates (27), we require Q N−1 − The solutions to these conditions are the translations which are the relations between the individual species' center-of-mass coordinates and the respective Jacobi coordinates [121]. We see that, in the laboratory frame, both contributing poles (31) are intermixed, and the time-dependent translation of each species is more intricate.
We can now translate the Floquet wavefunction (30) from the Jacobi to Cartesian coordinates, where the parameters are and x 2,j is used. From the time-dependent wavefunction (35,36), the (N 1 + N 2 )-particle density matrix of the driven mixture is defined and given by This allows us to compute the reduced density matrices of the driven mixture. It turns out that the determination of the time-dependent reduced density matrices of the driven mixture follows in practically the same way as for the undriven mixtures [69], although this is not a straightforward result. The driving forces are coupled to the relative center-of-mass Q N−1 and center-of-mass Q N coordinates only. Consequently, as shown above, all coordinates of species 1 are translated by the same amplitude x 1 (t), and likewise, all coordinates of species 2 are translated by x 2 (t). Performing the integrations consecutively over x 1,j , x 1,j = x 1,j , j = N 1 , . . . , 2, 1 and x 2,k , x 2,k = x 2,k , k = N 2 , . . . , 2, 1, for which the spatial-dependent phase factors are canceled, and while changing the variables x 1,j − x 1 (t) → x 1,j and x 2,k − x 2 (t) → x 2,k , one finds the same recurrence relations [69] determining the coefficients of the now time-dependent reduced density matrices (see below). As a result of this, the coherence properties of the driven mixture described by the ground Floquet wavefunction (35,36), like the eigenvalues of the time-dependent intraspecies reduced one-particle density matrices and interspecies reduced two-particle density matrix, are exactly those of the undriven mixture; in particular, these eigenvalues are time independent. Thus, the time-dependent intraspecies reduced one-particle density matrices are explicitly given as x 2 , t) are the intraspecies time-dependent densities and are given explicitly by We see that the two species generally oscillate with different amplitudes due to the intermixing of the two poles, with the specific cases of a single-pole contribution when only the center-of-mass (same-amplitude motion; for f L,1 Finally, the time-dependent interspecies reduced two-particle density matrix of the mixture, which is the lowest-order interspecies quantity, is given by where , Then, the diagonal part of ρ (2) 12 (x 1 , x 1 , x 2 , x 2 , t) is the time-dependent two-particle interspecies density and reads The coupling term in ρ (2) 12 (x 1 , x 2 , t) stems from the interspecies interaction and renders the time-dependent two-particle density not separable at the many-body level of theory. Furthermore, it mixes the time-dependent oscillation amplitudes x 1 (t) and x 2 (t) of both species. This could raise interesting questions pertaining to the intensities of ω L -harmonics in many-body quantities describing driven mixtures, which are left for investigations elsewhere.

Mean-Field Solution and the Limit of an Infinite Number of Particles
We are now in the position to solve the driven mixture at the mean-field level, and thereafter to compare the many-body and mean-field outcomes at the limit of an infinite number of particles. The related structures of the time-dependent reduced one-particle density matrices per particle of the mixture's ground Floquet state (38) and of the static reduced one-particle density matrices per particle of the mixture's ground state [69] already tell us that at the infinite-particle-number limit, each of the driven species is 100% condensed. Furthermore and on similar grounds, the time-dependent interspecies reduced two-particle density matrix per particle (40,41) is separable and given by the product of the time-dependent intraspecies reduced one-particle density matrices per particle at the limit of an infinite number of particles. To obtain further properties of the driven mixture at this limit, we now derive the solution of the driven harmonic-interaction model for mixtures at the mean-field level of theory. This would allow us to explicitly compare the respective intraspecies and interspecies time-dependent reduced density matrices per particle and quasienergies per particle and identify the connection between many-body and mean-field levels of theory for driven interacting mixtures, at least within the driven harmonic-interaction model for mixtures. Beyond the above, the Floquet solution of the coupled non-linear Schrödinger equations is interesting and instrumental for itself.
The time-dependent many-body wavefunction is given at the mean-field level by All species 1 bosons are described by one and the same time-dependent orbital φ 1 (x 1 , t), and all species 2 bosons are similarly described by another time-dependent orbital φ 2 (x 2 , t), which are to be determined self-consistently. Employing the time-dependent variational principle in any of its forms leads to the coupled time-dependent Gross-Pitaevskii equations where the mean-field interaction parameters are given by Only the interaction parameters appear within the mean-field treatment of the (driven) mixture. To compare the many-body and mean-field treatments, the mean-field interaction parameters are held fixed, while the numbers of particles N 1 and N 2 are increased towards the infinite-particle-number limit. In this way, the total number of particles N = N 1 + N 2 is increased to infinity while the ratios N 1 Λ 12 +Λ 21 , defined using the (non-zero) intraspecies interaction strength λ 12 , are held fixed.
The solution of the time-dependent coupled non-linear Schrödinger Equation (44) is somewhat involved, and for brevity and completeness it is detailed in Appendix B. The final result for the time-dependent mean-field orbitals can be written as where and Ω 2,GP = ω 2 + 2(Λ 2 +Λ 12 ) m 2 are interaction-dressed frequencies. The main finding at this point (see Appendix B for details) is that the time-dependent mean-field orbitals (45) oscillate with the same amplitudes, x 1 (t) and x 2 (t), found in the many-body solution (see (34)), and are otherwise similar to the ground-state mean-field orbitals of the static problem being dressed by both the intraspecies and intersepcies interactions [69]. These results suggest that the respective many-body and mean-field time-dependent reduced density matrices per particle of the driven mixture are to coincide at the limit of an infinite number of particles (see explicitly below). More delicate and intricate are the respective relations of the quasienergies.
As in the driven single-species case, the time-dependent phases µ 1 (t) and µ 2 (t) in (45) consist of periodic and linear in t parts. The former do not enter the quasienergy and can be read off (A20) in Appendix B, and the latter are given explicitly as The quantities µ 1,F and µ 2,F may be called quasichemical potentials of the driven mixture, and are simultaneous eigenvalues of the coupled non-linear Floquet equations are the corresponding periodic parts of the time-dependent orbitals, Thus, when solved at the mean-field level of theory, a driven mixture gives rise to coupled non-linear Floquet equations whose eigenvalues are analogs of the chemical potentials of a static mixture.
As can be seen from (46), the quasichemical potentials µ 1,F and µ 2,F are built from two parts. The first is the contribution of the static problem to the solution and is dressed by the intraspecies and interspecies interactions, and the second originates from the driving forces and depends on these forces and the interspecies interaction only. In particular, the poles' structure of µ 1,F and µ 2,F would play a central role in the mean-field quasienergy and its relation to the many-body quasienergy (33).
With the mean-field wavefunction (43) computed, we proceed to the time-dependent reduced density matrices and thereafter to the limit of an infinite number of particles. Thus, the lowest-order intraspecies and interspecies reduced density matrices are given at the mean-field level by Recall that, at the mean-field level, ρ (2) 12,MF (x 1 , x 1 , x 2 , x 2 , t) is naturally built as a product of intraspecies reduced one-particle density matrices.
We are now in the position to examine the relations between the many-body and mean-field solutions of the driven bosonic mixture at the limit of an infinite number of particles. Thus, the following relations readily hold: for the time-dependent reduced one-particle density matrices per particle, with 2,GP (x 2 , t) for the time-dependent densities per particle, and 12 (x 1 ,x 2 ,t) 2,GP (x 2 , t) for the time-dependent interspecies two-particle density.
The above relations between the time-dependent many-body and mean-field reduced density matrices generalize for time-periodic (Floquet) mixtures results existing in the literature for the ground state of trapped mixtures and, when available, results for the dynamics of interacting bosons with time-independent Hamiltonians. On top of and contrary to that, an appealing and different result emerges for the quasienergy per particle. Here, we find implying that the quasienergy per particle computed at the many-body level of theory (33) is different from that at the mean-field level (49), even at the infinite-particle-number limit. Comparatively, the result (52) for the quasienergy per particle of the driven mixture differs from the literature results for the ground-state energy per particle of trapped interacting bosons and mixtures at the limit of an infinite number of particles. Here, the difference can be seen as the additional factor , the ratio of the two resonance frequencies' poles, appearing in the first term inside the square brackets of the mean-field quasienergy (49). This additional factor emerges because the non-linear interaction terms in the mean-field treatment cannot compensate both poles in the driven amplitudes x 1 (t) and x 2 (t). More precisely, the mean-field treatment renormalizes the corresponding phases only to the center-of-mass frequency ω 2 and not also to the relative center-of-mass frequency Ω 2 12 (see Appendix B for the derivation).
When are the quasienergies equal? Trivially, if there is no interaction between the two species, implying that Ω 12 = ω and, therefore, ω 2 −ω 2 L Ω 2 12 −ω 2 L = 1, i.e., we have then two harmonically trapped non-coupled driven Bose-Einstein condensates; otherwise, as long as the relative center-of-mass coordinate is not activated by the forces, the respective quasienergies per particle are equal, namely, The latter occurs independently of the number of particles and strengths of interactions, when the forces and masses are interrelated by f L,1 This brings the present section on many-body and mean-field quasienergies and Floquet reduced density matrices to an end.

Angular Momentum and Its Fluctuations in a Bose-Einstein Condensate Steered by an Interacting Bosonic Impurity
So far, we considered periodically driving a mixture in one spatial dimension. A natural extension and application is to consider a two-dimensional mixture and drive it correspondingly along the x and y directions such that it is steered, namely, circularly driven. Then, angular momentum is imprinted on the bosons. We use the harmonic-interaction model for the driven mixture to investigate the roles of intraspecies and interspecies interactions and resulting angular-momentum fluctuations [72,74] at the many-body and mean-field levels. The availability of a solvable model is instructive.
Our strategy in this section is to begin with steering a single particle in two spatial dimensions and computing the angular momentum and fluctuations, namely, the expectation value of the angular-momentum operator and its variance. Then, we steer the N-boson problem and investigate the role of interaction at the many-body and mean-field levels of theory. Finally, we steer species 1 bosons, which serve as an interacting bosonic impurity, embedded in species 2 bosons, where the latter are not steered, and then investigate the modes of rotation, distribution of angular momentum between the species, and fluctuations.
Consider the circularly driven two-dimensional quantum harmonic oscillator [r = (x, y)], h(r, t) = − 1 2m The ground Floquet solution of the time-dependent Schrödinger equation is given by where ε F is given in (5).
The time-dependent wavefunction (55) can be seen, apart from the time-dependent phase e −i2ε F t , which is irrelevant for the properties discussed in this section, as the ground state of the two-dimensional isotropic harmonic oscillator being translated by x(t) and y(t) and boosted by mẋ(t) and mẏ(t) along the x and y directions, respectively. A useful way to evaluate the expectation value of the angular-momentum operator and its variance carried by a translated and boosted wavefunction is to use the transformation properties of the angular-momentum operator under translations and boosts. In Appendix C, we derive and discuss the corresponding expressions for the most general wavefunction treated in this work, that of driven mixtures, where each species can be translated and boosted differently.
Using (A25) for N = 1 particles, we can readily evaluate the expectation value of the angular-momentum operatorl z =xp y −ŷp x , which reads The value (56) is that of a 'classical' particle encircling the origin with angular velocity ω L and radius f L m(ω 2 −ω 2 L ) . The radius is determined by the driving force and diverges as one approaches the resonance. On the other end, the radius is always non-zero (for f L = 0), meaning that the expectation value of the angular momentum for the driven particle is non-zero (and positive).
The Floquet solution (55) of the steered harmonic oscillator is not an eigenfunction of the angular-momentum operator. Using the transformation properties of the latter under translations and boosts, one can readily compute the angular-momentum variance, which is found to depend on the position and momentum variances along the x and y directions (see Appendix C). Thus, employing (A27) for N = 1 particles, we find The dependence of the variance (57) on the radius f L m(ω 2 −ω 2 L ) is in conjunction with the expectation value (56). For instance, for tight confinements (ω 1), the radius diminishes, and so does the variance. Note, however, that for ω → 0, the variance diverges, whereas the expectation value of the angular momentum does not. This signifies that the angular-momentum variance depends on the size of the wavepacket, and for weaker confinements, this size increases. We will return to these points when interactions set in.
Let us move to the circularly driven isotropic harmonic-interaction model in two spatial dimensions, The ground Floquet solution is given by where E F is given in (11). Now, the expectation value of the many-particle angular-momentum operator j=1 (x jpy,j −ŷ jpx,j ) and the respective variance are given by where (A25) and (A27) are employed. We see that the interaction between particles enters neither the expectation value nor the variance of the angular-momentum operator at the many-body level, and the values per particle are given by those of the steered single-particle system (see (56) and (57)). Let us examine the respective quantities at the mean-field level, where the Floquet solution was analytically derived in Section 2 and reads in two spatial dimensions: where µ F is given in (20). Computing now the expectation value and variance ofL Z , we find where the variance exhibits dependence on the interaction strength. This behavior is because the position and momentum variances (of the non-boosted/translated wavefunction Φ) themselves depend, at the mean-field level, on the interaction strength , 1 In particular, for repulsive [Λ → (−mω 2 /2) + ] as well as for attractive (Λ 1) interactions, the angular-momentum variance (62) can diverge, unlike the situation at the many-body level of theory (60). On the other end, the variance does not vanish for any interacting strength. The situation becomes more interesting when we steer species 1 bosons embedded in species 2 bosons (see below). Finally, we note that the expectation values per particle of the angular momentum at the many-body and mean-field levels of theory coincide (for finite systems and) at the limit of an infinite number of particles, as is expected for one-body operators [64,72].
Let us proceed to the scenario of a steered mixture. Consider the harmonic-interaction model in two spatial dimensions for circularly-driven species 1 bosons embedded in un-steered species 2 bosons: H(r 1,1 , . . . , r 1,N 1 , r 2,1 , . . . , r 2,N 2 , t) = We expect the interspecies interaction λ 12 to mediate the rotation from species 1 to species 2, and the question to answer is in what way. The ground Floquet solution to (63) is given by Ψ(r 1,1 , . . . , r 1,N 1 , r 2,1 , . . . , r 2,N 2 , t) = m 1 Ω 1 where the quasienergy E F is given in (33), x 1 (t) = a 1 cos(ω L t), y 1 (t) = a 1 sin(ω L t), x 2 (t) = a 2 cos(ω L t), and y 2 (t) = a 2 sin(ω L t), and the amplitudes a 1 and a 2 , determining the respective radii of rotations, are: The resulting amplitudes a 1 and a 2 deserve a discussion. Despite that only species 1 is being steered, there are contributions from two distinct poles, at ω 2 and at Ω 2 12 = ω 2 + 2 Λ 12 m 2 + Λ 21 m 1 , signifying the activation of the center-of-mass Q N and relative center-of-mass Q N−1 coordinates in the presence of coupling between the two species. The radius of rotation of species 2, a 2 , is linear in the intraspecies interacting parameter, Λ 12 = N 1 λ 12 , and vanishes, as expected for diminishing coupling between the steered species 1 and un-steered species 2 Bose-Einstein condensates.
It is appealing to see that there is also a zero for a 1 when which is 'born' only due to the interaction between species 1 and species 2 bosons. Physically, a vanishing value for a 1 means that species 1 does not rotate or move, despite being driven by force f L,1 at the angular frequency ω 0 L . Alternatively, we may say that the effect of the driving force f L,1 exerted on species 1 is then completely transferred to species 2. The value of a 2 at ω 0 L (see (66)) is inversely proportional to the intraspecies interaction parameter Λ 21 = N 2 λ 12 . Further analysis of a 1 and a 2 is given below when analyzing the distribution of angular momentum between the two species.
The mean-field Floquet wavefunction in two spatial dimensions is needed to compute the respective properties and is given by Φ (r 1,1 , . . . , r 1,N 1 , r 2,1 , . . . , r 2,N 2 where the quasichemical potentials µ 1,F and µ 2,F are given in (46). We are now in the position to compute the expectation value of the angular-momentum operator and its variance and to compare the many-body and mean-field results. For a mixtureL Z =L Z 1 +L Z 2 , we shall also investigate the connection between the mixture's and individual species' angular momenta. We start with the expectation values at the many-body level of theory, which take on the form: The amount of angular momentum carried by each species is proportional to the square of the radius of its rotation. The poles and zeros of the radii govern the sizes of angular momenta, yet only the zeros dictate the relative distribution of angular momentum between the two species; recall that only species 1 is steered. Explicitly, the ratio of the angular momenta per particle of the two species is and can vary from arbitrarily large for small interspecies interaction parameters |Λ 12 | to arbitrarily small due to the zero of a 1 at ω 0 L = ω 2 + 2Λ 12 m 2 (see (66)). Another interesting regime occurs in the vicinity of the center-of-mass pole 1 (ω 2 −ω 2 L ) 2 (see (68)), where the ratio of the angular momenta per particle approaches the mass ratio m 1 m 2 and is independent of other parameters. Of course, in the vicinity of any of the poles, the rotation radii increase ad infinitum.
To compute the variances ofL Z ,L Z 1 , andL Z 2 of the time-dependent mixture, we need, according to Appendix C, the respective position, momentum, and angular-momentum variances of the individual species, computed at the many-body and mean-field levels with respect to the wavefunctions Ψ and Φ without translations and boosts (see (A28) and (A30)). Using the definition of the center-of-mass and relative center-of-mass Jacobi coordinates in terms of the individual species' center-of-mass coordinates, the final result at the many-body level reads [121]: The respective quantities at the mean-field level are simpler to evaluate, since Φ is a product state, and the final result reads [121]: It is worthwhile to list the main properties of and differences between the variances (69) and (70), starting with the position and momentum variances, as far as they are needed for our investigations: (i) The individual species' variances depend only on the interspecies interaction parameters Λ 12 and Λ 21 at the many-body level, whereas at the mean-field level, they additionally depend on the intraspecies interaction parameters Λ 1 and Λ 2 . (ii) For attractive interactions, the momentum variances can diverge, and for repulsive interactions, the position variances can diverge. Combining (i) and (ii), one can have situations where a position or momentum variance at the many-body level does not diverge, whereas it does diverge at the mean-field level, and vice versa, i.e., mean-field variances tend to zero, whereas many-body variances remain finite. Furthermore, for the individual species' angular-momentum variances: (iii) There is a non-vanishing value (69) at the many-body level of theory, which solely originates from the interaction between the species, and which can diverge for both repulsive and attractive interspecies interaction. Interestingly, it is marginal in the number of particles, that is, the angular-momentum variance per particle, unlike the position and momentum variances per particle, diminishes for large particle numbers and at the limit of an infinite number. Finally, at the mean-field level of theory, the angular-momentum variance (70) is zero due to separability of Φ. We are now in the position to put the pieces together.
We find for the angular-momentum variance at the many-body level 1 N ∆ 2 L Z Ψ(r 1,1 ,...,r 2,N 2 ,t) where the individual species' values are (r 1,1 ,...,r 2,N 2 ,t) Examining the expressions (71) for the angular-momentum variances , the first point to note is that the individual species' variances (71b) do not sum up to the mixture's variance (71a). This is due to the position Ψ|X 1X2 |Ψ = Ψ|Ŷ 1Ŷ2 |Ψ and momentum Ψ|P X 1P X 2 |Ψ = Ψ|P Y 1P Y 2 |Ψ cross-terms appearing in (A30). For the individual species (see (A28)), one of course does not have contributions from any cross-term. On the other hand, for the mixture's angular-momentum variance, the angular-momentum cross-term Ψ|L Z 1L Z 2 |Ψ cancels the individual species' contributions Ψ|L 2 Z 1 |Ψ and Ψ|L 2 Z 2 |Ψ due to conservation of angular momentum, L Z 1 +L Z 2 |Ψ = 0, and hence, this cross-term does not appear in (A30). All in all, the angular-momentum variance of the driven mixture (71a) does not vanish, depends on the interspecies interaction parameters Λ 12 and Λ 21 only, can diverge either due to repulsive or attractive interspecies interaction, and is governed by the poles of the center-of-mass and relative center-of-mass coordinates.
To proceed, the angular-momentum variances at the mean-field level of theory are evaluated and read: (r 1,1 ,...,r 2,N 2 ,t) (r 1,1 ,...,r 2,N 2 ,t) (r 1,1 ,...,r 2,N 2 ,t) where (A28) and (A30) are employed. At the mean-field level, in contrast with the outcome (71) at the many-body level, the individual species' angular-momentum variances add up to the mixture's angular-momentum variance, since all cross- On top of that, the mean-field quantities depend on all interaction parameters, intraspecies and interspecies, not only on the latter ones. Altogether, the angular-momentum variances per particle at the mean-field level (72) are different from those at the many-body level (71) for finite systems as well as at the limit of an infinite number of particles where each of the species is 100% condensed, i.e., lim N 1 →∞ ..,r 2,N 2 ,t) for the mixture's angular-momentum variance and lim N 1 →∞ 1 ,...,r 2,N 2 ,t) for the individual species' angular-momentum variances. A final point: What happens to the angular-momentum expectation values and variances at the frequency ω L = ω 0 L when the effect of the steering force exerted on species 1 is completely transferred to species 2? First, since a 1 ω L =ω 0 L = 0, the expectation value of the mixture's angular momentum and its fluctuations at the many-body and mean-field levels are now given by the respective expressions of the un-steered species 2 only. Second, one gets for the angular-momentum variance of species 1 that the only contribution at the many-body level is the marginal part, which diminishes to zero, i.e., to the mean-field value, as the number of particles increases. Third, the respective many-body and mean-field expressions for the mixture's and species 2 variances remain different, even at the limit of an infinite number of particles when both species are 100% condensed. This wraps up the present work.

Summary and Outlook
We have presented in this work a solvable time-dependent model of a many-particle system and used it to amalgamate together several themes: (i) Mixtures of trapped Bose-Einstein condensates; (ii) Periodically driven many-particle systems; (iii) Many-body versus mean-field theory of Floquet bosonic systems at the limit of an infinite number of particles; (iv) The imprinting and fluctuations of angular momentum when a Bose-Einstein condensate is steered by interacting bosons consisting of a different species.
We introduced the driven harmonic-interaction model for mixtures, namely, N 1 bosons of mass m 1 interacting by harmonic interparticle interaction of strength λ 1 and driven by a time-periodic force of amplitude f L,1 and frequency ω L , and N 2 bosons of mass m 2 interacting by harmonic interparticle interaction of strength λ 2 and driven by a time-periodic force of amplitude f L,1 and the same frequency ω L . All bosons are trapped in a harmonic potential of frequency ω, and, furthermore, bosons of species 1 and species 2 interact with each other by harmonic interparticle interaction of strength λ 12 . The model generalizes the harmonic-interaction model for mixtures to the time-dependent domain.
Using Jacobi coordinates, the time-dependent model is diagonalized, and governed by two time-dependent many-particle oscillators, which are associated with the center of mass and relative center of mass of the driven mixture, as well as two poles of respective frequencies. The frequency of the center-of-mass motion does not depend on any interaction strength, whereas the frequency of the relative center-of-mass motion depends on the interspecies interaction strength. We concentrated in this work on the ground Floquet solution of the driven mixture, and began by computing its quasienergy, the eigenvalue of the time-periodic many-particle Floquet Hamiltonian, as a function of the above parameters.
Back in the laboratory frame, all Cartesian coordinates in the time-dependent many-particle quasienergy state are coupled, and the coordinates of each species oscillate with a different amplitude that combines the two poles associated with the center-of-mass and relative center-of-mass motions. This structure of the time-dependent many-boson wavefunction allows one to explicitly obtain the intraspecies and interspecies time-dependent reduced one-particle and two-particle density matrices of the driven mixture, which possess the same coherence properties as the corresponding static quantities computed for the ground state of the trapped mixture.
One of the main goals of the present work is to ask and answer questions on the relation between many-body and mean-field theories of (trapped) interacting bosons. These have been dealt with in the literature for the ground state or for the dynamics with time-independent Hamiltonians, and we would like to extend the topic to the realm of Floquet many-boson systems. To discuss properties of the driven mixture at the limit of an infinite number of particles-specifically, the connection between the respective many-body and mean-field results for the time-dependent reduced density matrices, densities, the quasienergies that are unique for Floquet systems, as well as the variances-the solution of the coupled Floquet non-linear Schrödinger equations of the driven mixture has been explicitly derived.
With closed-form expressions at the many-body and mean-field levels of theory, we have proven explicitly that, at the limit of an infinite number of particles, the intraspecies time-dependent reduced one-particle and two-particle density matrices per particle are 100% condensed and given by the respective mean-field quantities. Furthermore, in this limit, we have proven that the interspecies time-dependent reduced two-particle density matrix per particle is separable, i.e., given by the product of the time-dependent reduced one-particle density matrices per particle of species 1 and species 2 bosons, and is given by the respective mean-field quantity. These results generalize the literature results for the ground state and the dynamics with time-independent Hamiltonians to a Floquet bosonic system.
We have found that the relation between the many-body and mean-field quasienergies per particle is intriguing, namely, that the many-body and mean-field results need not coincide. This is unlike the situation for the many-body and mean-field energies per particle of the ground state of trapped bosons and mixtures. Specifically, the coupled non-linear Schrödinger equations cannot appropriately renormalize the intensity of both poles of the oscillating mixture in the time-dependent mean-field wavefunction, and hence, in the mean-field quasienergy, in comparison with the many-body treatment. Thus, the many-body quasienergy per particle coincides with the mean-field value at the infinite-particle-number limit in the specific case where only the center-of-mass motion of the driven mixture is activated by the driving forces. In the general case when the relative center-of-mass motion is also activated, the many-body quasienergy per particle differs from the mean-field one even at the limit of an infinite number of particles.
We then investigated an application in two spatial dimensions, in which only species 1 bosons are circularly driven by a time-periodic force, and through their interspecies interaction with the un-steered species 2 bosons, induce rotation and imprint angular momentum on the latter. The dynamics are solved analytically and analyzed within many-body and mean-field theories, as presented above. We first concentrated on the average of the angular-momentum operators of each of the species. Despite steering species 1 bosons only, we show that the ratio of the angular momenta per particle of species 1 and species 2 bosons can be marginally small or arbitrarily large because of the interspecies interaction. Explicitly, the latter induces a zero in the radius of rotation of species 1 bosons at a frequency just between the center-of-mass and relative center-of-mass resonances. We have also shown that the expectation values per particle of the angular-momentum operators for the many-body and mean-field solutions of the time-periodic mixture coincide (for finite systems and) at the limit of an infinite number of particles.
Finally, we investigated the variances of the angular-momentum operators of each of the species and of the mixture. To this end, it is useful to express the Floquet quasienergy state of the steered mixture using translation and boost operators along the x and y directions. Consequently, we can express the angular-momentum variance computed from a translated and boosted wavefunction using the position, momentum, and angular-momentum variances of the un-translated/boosted wavefunction. The angular-momentum variances of species 1 and species 2 bosons as well as of the mixture were computed explicitly and compared at the many-body and mean-field levels of theory. The mean-field variances depend on the resonances and all interactions, and can diverge due to either the attractive or repulsive intraspecies and interspecies interactions. Interestingly, at the mean-field level of theory, the individual species' angular-momentum variances are additive and sum up to the mixture's angular-momentum variance. On the other hand, the many-body quantities depend on the resonances and interspecies interaction only, can diverge due to either the attractive or repulsive interspecies interaction, and are non-additive because of the correlation term in the time-dependent many-boson wavefunction. These differences persist at the limit of an infinite number of particles, where each of the species comprising the Floquet time-dependent mixture is 100% condensed.
Last but not least, a brief outlook: There are several directions worth pursuing of which we list four. We have discussed, in the present work, the driving of a bosonic mixture by a monochromatic force, and can envision that using a polychromatic force would be instructive. We have studied the driving of structureless bosons, and foresee that solvable models for driving structured targets embedded in a Bose-Einstein condensate would be interesting. We have investigated the imprinting of angular-momentum in a steered two-species mixture in two spatial dimensions, and suggest that exploring three-species mixtures would be valuable. Finally, the present model can serve to benchmark numerical algorithms developed to compute real-space many-boson Floquet states, which would exhibit intriguing and both fundamental and applicative facets. where C 1 can be read from (15). The diagonal of (A1) is the time-dependent two-particle density and reads From the mean-field wavefunction (17), the reduced two-particle density matrix at the mean-field level reads where ρ (1) (22)].
The two-particle mean-field density reads ρ (2) . Thus, at the limit of an infinite number of particles, one finds that for the time-dependent reduced two-particle density matrix per particle and lim N→∞ GP (x , t) for the time-dependent two-particle density per particle. For the driven mixture, we have for the ground Floquet solution (35,36) the following closed-form expressions: where the time-dependent two-body densities are given explicitly by From the mean-field wavefunction (43), the time-dependent intraspecies reduced two-particle density matrices are given at the mean-field level of theory by [Equation (48)]. The two-particle mean-field intraspecies densities read ρ (2) . Therefore, with the mean-field quantities explicitly prescribed, one obtains at the limit of an infinite number of particles that for the reduced two-particle density matrices per particle and lim N 1 →∞ 2,GP (x 2 , t) for the time-dependent two-particle intraspecies densities per particle themselves.

Appendix B. Further Details of Solving the Time-Dependent Coupled Gross-Pitaevskii Mean-Field Equations
In this appendix, we derive the solution of the Floquet non-linear Schrödinger equations of the driven mixture. To this end, we first obtain the respective solution of the driven single-species harmonic-interaction model and thereafter generalize the techniques used.
The Gross-Pitaevskii equation for the driven harmonic-interaction model (see (18)) is ∂t . We suspect that and examine whether the solution can be cast as where the interaction-dressed frequency, Ω GP , time-dependent phase, µ(t), and amplitude of oscillations, a, are to be determined self-consistently. Substituting the function (A11) into the square brackets on the left-hand side of the time-dependent Gross-Pitaevskii equation, expanding the interaction term using the driven oscillations, and collecting terms, we find, for the square brackets, − 1 2m ∂ 2 ∂x 2 + 1 2 m ω 2 + 2Λ m x 2 − x( f L + 2aΛ) cos(ω L t) + Λa 2 cos 2 (ω L t) + Λ 2Ω GP . This implies Ω GP = ω 2 + 2Λ m to hold, just like in the static problem. Next, as the solution of a driven oscillator with the interaction-dressed frequency Ω GP , the amplitude a must satisfy i.e., resulting in the same amplitude of oscillations as for the bare driven harmonic oscillator and for the many-body solution in the laboratory frame. To verify consistency, we re-expand 1 2 mΩ 2 GP x 2 = 1 2 mΩ 2 GP [x − a cos(ω L t)] 2 + mΩ 2 GP xa cos(ω L t) − 1 2 mΩ 2 GP a 2 cos 2 (ω L t) and indeed find mΩ 2 for the linear-in-x term and − 1 2 mΩ 2 GP a 2 cos 2 (ω L t) + Λa 2 cos 2 (ω L t) = − 1 2 mω 2 a 2 cos 2 (ω L t) = − f 2 L ω 2 2m(ω 2 −ω 2 L ) 2 cos 2 (ω L t) for the contribution to the phase, as is the case with the bare system. Finally, when operating from the left-hand side of the Gross-Pitaevskii equation with the kinetic-energy operator and from the right-hand side with the derivative in time on (A11), the remaining pieces of the time-dependent phase are collected, , determining the time-dependent mean-field solution (19). Summing up, despite the interaction dressing of the frequency of the (driven) oscillator from ω to Ω GP , the second and third contributions (see main text) to the time-dependent phase of φ(x, t) are renormalized 'back' to those of the bare system.
We employ the same strategy to derive the more involved solution of the coupled time-dependent Gross-Pitaevskii equations for the driven-mixture system.
The coupled Gross-Pitaevskii equations for the harmonic-interaction model for a driven mixture are (see (44)). We guess that and check whether the solution can be cast as where the interaction-dressed frequencies, Ω 1,GP and Ω 2,GP , time-dependent phases, µ 1 (t) and µ 2 (t), and amplitudes of oscillations, a 1 and a 2 , are to be determined self-consistently. Substituting the functions (A14) into the curly brackets on the left-hand sides of the two coupled time-dependent Gross-Pitaevskii equations, expanding the interaction terms with the respective oscillation of each orbital, and collecting the terms in each set of curly brackets, we find These imply that Ω 1,GP = ω 2 + 2(Λ 1 +Λ 21 ) m 1 and Ω 2,GP = ω 2 + 2(Λ 2 +Λ 12 ) m 2 hold, just as in the static mixture's problem.
Next, as the solutions of driven coupled (many-particle) oscillators with the interaction-dressed frequencies Ω 1,GP and Ω 2,GP , the amplitudes a 1 and a 2 have to satisfy the coupled linear relations the solution of which is readily found: Interestingly, these are the same amplitudes of oscillations, with the same poles, as for the many-body solution (see (34)). We may say that the non-linear mean-field terms renormalize the poles from the values of the interaction-dressed oscillators Ω 2 1,GP and Ω 2 2,GP (see (A16)) to the oscillators' frequencies of the many-body problem, ω 2 and Ω 2 12 (see (A17)).
To proceed for the terms linear in x 1 and x 2 , after re-expanding 1 2 m 1 Ω 2 1,GP x 2 1 = = 1 2 m 1 Ω 2 1,GP [x 1 − a 1 cos(ω L t)] 2 + m 1 Ω 2 1,GP x 1 a 1 cos(ω L t) − 1 2 m 1 Ω 2 1,GP a 2 1 cos 2 (ω L t) and likewise for 1 2 m 2 Ω 2 2,GP x 2 2 , we find and for the contributions to the phases The meaning of the above is that the linear terms (A18) are appropriately renormalized in comparison with the many-body treatment by the non-linear interaction terms, whereas the phases (A19) are not (also see below). In summary, the mean-field treatment correctly captures the amplitudes of oscillations and contributions of the linear terms (associated with operation of the kinetic-energy operator) to the wavefunction, but not the contribution of the time-dependent phases. These contributions do not match both poles' frequencies, ω 2 and Ω 2 12 , appearing in the oscillations' terms at the mean-field level. Adding up all terms, we have, for the time-dependent phases in (A14), where are the quasichemical potentials, determining the time-dependent mean-field solution of the driven mixture in the main text (see (45)). After substitution of the amplitudes a 1 and a 2 , the final expressions for µ 1,F and µ 2,F are obtained (see (46)).
The final expression for the quasienergy per particle at the mean-field level of theory simplifies, since the interspecies interaction terms in (A21) cancel each other upon the addition N 1 µ 1 (t) + N 2 µ 2 (t) within the time-dependent phase of the mean-field Floquet wavefunction, and reads: Upon substitution of the oscillations' amplitudes a 1 and a 2 , the expression for the mixture's quasienergy given in the main text is obtained (see (49)).

Appendix C. Angular-Momentum Variances in Presence of Translations and Boosts
In this appendix, we derive expressions for the angular-momentum variance of translated and boosted many-particle wavefunctions in terms of quantities computed from the un-translated and un-boosted wavefunctions. We start with the analysis for single-species bosons (for the case of translations only, see [74]), then generalize it to the two-species mixture or interacting-impurity system.
Consider the translated-boosted many-particle wavefunction in two spatial dimensions, e −i(P X a+P Y b) e +i(Xα+Ŷβ) Ψ(X, Y) ≡ Ψ(a, b; α, β), whereP X = ∑ N j=1px,j ,P Y = ∑ N j=1py,j andX = ∑ N j=1xj ,Ŷ = ∑ N j=1ŷ j , and the shorthand collective notation for the coordinates Ψ(X, Y) ≡ Ψ(x 1 , y 1 , . . . , x N , y N ) is implied. The order of the translation and boost operators does not matter for the computation of matrix elements. The wavefunction (before translations and boosts) as well as the translation and boost operators themselves can be time dependent. What are the implications for the variances of observables when computed with respect to the translated-boosted wavefunction Ψ(a, b; α, β)?
For the position operatorX (and equivalently forŶ), we have e −i(Xα+Ŷβ) e +i(P X a+P Y b)X × ×e −i(P X a+P Y b) e +i(Xα+Ŷβ) =X + Na, and using its square, one has 1 N ∆ 2X Ψ(a,b;α,β) Similarly, for the momentum operatorP X (and equivalently forP Y ), we have e −i(Xα+Ŷβ) × ×e +i(P X a+P Y b)P X e −i(P X a+P Y b) e +i(Xα+Ŷβ) =P X + Nα, and using its square leads to 1 N ∆ 2P X Ψ(a,b;α,β) i.e., both the position and momentum variances are invariant under translations, boosts, and their combined operation. For the angular-momentum operatorL Z = ∑ N j=1 x jpy,j −ŷ jpx,j , the situation is more intricate because it consists of the position and momentum operators, which do not commute with each other. From e −i(Xα+Ŷβ) e +i(P X a+P Y b)L Z e −i(P X a+P Y b) e +i(Xα+Ŷβ) = =L Z + aP Y − bP X + βX − αŶ + N(aβ − bα) (A25) and its square, we have 1 N ∆ 2 L Z Ψ(a,b;α,β) Equation (A26) warrants a discussion. We see that the angular-momentum variance of a translated-boosted wavefunction can be decomposed into various terms of physical meaning computed from the un-translated/boosted wavefunction. These include, for a start, the angular-momentum variance of the wavefunction and the momentum and position variances in each direction (see the first line in (A26)). Then, there are terms that originate from translations only, as well as 'interference' between them, followed by terms that originate from boosts only, as well as 'interference' between them, and, finally, 'interference' terms between translations and boosts. For the most general wavefunction, all involved expectation values contribute. When symmetries set in (see below), or when the wavefunction is real (excluding a trivial phase factor), corresponding terms fall.
For our needs, let us examine translated and boosted rotationally symmetric systems. Consider the ground state Ψ of an interacting many-boson system in a rotationally symmetric trap (to be more precise, we do not consider the ground state of a rotationally symmetric system in the rotating frame). Clearly, 1 N ∆ 2 L Z Ψ = 0 holds. Now, when the ground state is translated and boosted, the angular-momentum variance can be written as and seen to solely be described by the momentum and position variances along the y and x directions multiplied by the respective values of the translations and boosts along the x and y directions. Interestingly, the angular-momentum variance of a translated-boosted many-boson system 1 N ∆ 2 L Z Ψ(a,b;α,β) differs at the many-body and mean-field levels of theory, i.e., when a, b = 0, α, β = 0, and λ = 0. This, as can be seen in (A27), is because of the respective many-body and mean-field momentum, 1 N ∆ 2P , variances. For instance, in harmonic traps, these variances do not dependent on the interaction strength at the many-body level of theory, but are dressed by the interaction at the mean-field level of theory (see, in this context, [64]). The analytical result (A27) is employed to compute and analyze the findings for the angular-momentum variance of a many-boson Floquet state in the main text. Generally, as mentioned above, in the absence of spatial symmetries, more terms contribute to the translated and boosted angular-momentum variance (see (A26)).
We now generalize the above analysis to the mixture or impurity system. The notation used hereafter is a straightforward generalization of the notation employed above in the single-species case. Starting from the many-particle wavefunction Ψ(X 1 , Y 1 , X 2 , Y 2 ), we note that each species can be translated and boosted differently, e −i(P X 1 a 1 +P Y 1 b 1 ) e +i(X 1 α 1 +Ŷ 1 β 1 ) × ×e −i(P X 2 a 2 +P Y 2 b 2 ) e +i(X 2 α 2 +Ŷ 2 β 1 ) Ψ(X 1 , Y 1 , X 2 , Y 2 ) ≡ Ψ( a, b; α, β), as is, for instance, the result in the main text for the steered species 1 bosons embedded in un-steered species 2 bosons (see Section 4).
With Ψ( a, b; α, β), it is straightforward to see that the variance of the (mass-weighted) position operator, N M m 1X1 + m 2X2 , total momentum operator,P 1 +P 2 , or any linear combination ofX 1 andX 2 or ofP 1 andP 2 is invariant under translations, boosts, and combinations thereof.
For the variance of the angular-momentum operator of either species 1 or 2, we just use (A26) with the corresponding notations, since 1 . Symmetry can be of assistance as well, yet note that 1 do not vanish at the many-body level (and differences between the mean-field and many-body treatments can thus occur), even for the ground state of a mixture in spherically symmetric traps, due to the coupling introduced by the interspecies interaction λ 12 [121]. All in all, in our case, we have for the angular-momentum variances of the individual species in the driven mixture.
To derive the working expression for the variance of the total angular-momentum operator L Z =L Z 1 +L Z 2 of the mixture, we begin from the transformed quantity e −i(X 2 α 2 +Ŷ 2 β 1 ) e +i(P X 2 a 2 +P Y 2 b 2 ) e −i(X 1 α 1 +Ŷ 1 β 1 ) e +i(P X 1 a 1 +P Y 1 b 1 )L Z × ×e −i(P X 1 a 1 +P Y 1 b 1 ) e +i(X 1 α 1 +Ŷ 1 β 1 ) e −i(P X 2 a 2 +P Y 2 b 2 ) e +i(X 2 α 2 +Ŷ 2 β 1 ) = and its square. The final result is lengthier than (A26), but is otherwise straightforward. For spherically symmetric traps, many terms fall, and the final result can be written as When applied to the driven mixture in the main text, the angular-momentum variance (A30) can be expressed using the variances of time-dependent interspecies-weighted momenta, a 1PY 1 + a 2PY 2 and b 1PX 1 + b 2PX 2 , and position, β 1X1 + β 2X2 and α 1Ŷ1 + α 2Ŷ2 , operators. Alternatively, it can be expressed using the intraspecies momentum and position variances as well as the respective interspecies cross-terms. The contributions of the various intraspecies and interspecies terms to the mixture's angular-momentum variance are different at the many-body and mean-field levels of theory, and are computed and discussed in the main text. In particular, the momentum and position interspecies cross-terms in (A30) vanish at the mean-field treatment of the (driven) mixture.