Self-ordering and cavity cooling using a train of ultrashort pulses

A thin atomic gas in an optical resonator exhibits a phase transition from a homogeneous density to crystalline order when laser illuminated orthogonal to the resonator axis. We study this well-known self-organization phenomenon for a generalized pumping scheme using a femtosecond pulse train with a frequency spectrum spanning a large bandwidth covering many cavity modes. We show that due to simultaneous scattering into adjacent longitudinal cavity modes the induced light forces and the atomic dynamics becomes nearly translation-invariant along the cavity axis. In addition the laser bandwidth introduces a new correlation length scale within which clustering of the atoms is energetically favorable. Numerical simulations allow us to determine the self-consistent ordering threshold power as function of bandwidth and atomic cloud size. We find strong evidence for a change from a second order to a first order self-ordering phase transition with growing laser bandwidth when the size of the atomic cloud gets bigger than the clustering length. An analysis of the cavity output reveals a corresponding transition from a single to a double pulse traveling within the cavity. This doubles the output pulse repetition rate and creates a new time crystal structure in the cavity output. Simulations also show that multi-mode operation significantly improves cavity cooling generating lower kinetic temperatures at a much faster cooling rate.


Introduction
Laser light forces on a dilute atomic gas are enhanced and qualitatively modified within optical resonators as the back-action of the atoms onto the field dynamics cannot be ignored [1]. In particular, even at large detuning, the dipole force does not simply induce a conservative optical potential, but the dynamic nonlinear atom-field coupling opens a wealth of new phenomena [2]. An important consequence is the possibility for opto-mechanical cooling without spontaneous emission via the delayed response of the cavity field to the atomic motion which can be applied to any kind of polarizable particles [3,4,5,6,7,8].
As another surprising phenomenon one finds phase transitions towards regular order and improved collective cooling if the atoms are transversely illuminated by a sufficiently strong pump laser tuned just below the eigenfrequency of a cavity mode [4,9]. This spontaneous ordering process can be traced back to collective (superradiant) coherent Bragg scattering from an emerging periodic atomic density distribution forming selfconsistently at the interference maxima of pump and cavity field [4]. Such spontaneous crystallization can still be observed at zero temperature as a quantum phase transition from a superfluid to a supersolid order [10,11,2,12].
When several laser frequencies and field modes are simultaneously applied, the corresponding atom-field dynamics gets much more complex. Recently ultracold atomic gases in cavities have been used to simulate or emulate a wide class of other exotic solid state phenomena like edge states, topological insulators or quasi-crystal formation [13,14] or general Hamiltonians for quantum enhanced simulated annealing [15,16]. Using atoms with internal spin, cavity mediated interactions have the potential to implement a large class of lattice spin models [17,18,19,20].
In several theoretical proposals it was predicted that cavity cooling is strongly enhanced when the number of modes involved is enlarged [21,22]. Interestingly the use of only two distinct cavity modes and laser frequencies already strongly enlarges the complexity of theoretical modeling. While a single cavity mode creates a virtually infinite range all-to-all interaction throughout the whole cavity volume, by help of a large number of longitudinal modes with corresponding driving frequencies, one should be able to engineer interaction forces of a very general spatial shape [16]. While this sounds extremely complex in practice, one can make use of the fact that cavity modes like comb lines of ultra short pulse trains are equidistant. Hence frequency matching of only two lines will automatically align a huge manifold of frequency components. In general the frequency distribution of the pump light then gives a direct handle on the spatial shape of the interaction forces [23].
Here we generalize the pump light to frequency combs (FCs) with a large bandwidth as generated by a pulse train coming from a mode-locked femtosecond laser. The locking of two comb lines to two different longitudinal modes is then sufficient to guarantee resonant overlap of all modes over virtually the whole pump bandwidth. This way the atoms can simultaneously coherently scatter thousands of different frequencies into the cavity, which can form a newly shaped pulse train between the mirrors. The average optical potential over this pulse train then determines the spatial order of the atoms in a self-consistent way. As in the case of a single-mode the atoms mostly tend to order in such a way that they maximize the total intracavity intensity. This can lead to local periodic order together with long-range clustering of the atoms. Interestingly, the temporal shape of the output pulses and their repetition rate provides for nondestructive observation of the atomic ordering process. As the output pulse train dynamically acquires a different shape and a higher repetition rate than the input pulse, the ordering phenomenon should be closely related to a time crystal formation. This paper is organized as follows: We first introduce the proposed set-up and a semi-classical mathematical model, which includes a discussion about basic physical properties in Sec. 2. Next, in Sec. 3 we identify the lowest energy stationary states of the system and numerically calculate the ordering threshold using a self consistent mean-field approach. The dynamics for a wide range of operating parameters and initial conditions are discussed in Sec. 4. A special emphasis is put on the form of the ordered patterns and their signature in the cavity output light. Finally we discuss the extra benefits for cavity cooling in such a wide bandwidth multi-mode system. We conclude our work with an outlook on future perspectives in Sec. 5.

Model
Let us consider an atomic cloud of N atoms of mass m a trapped along the central axis within a high finesse Fabry-Pérot cavity of length L supporting a wide range of non-degenerate longitudinal modes. The modes have wave-vectors k m = mδk and frequency distance δk = π/L, where the mode indices m range from m low to m high . The atoms are transversally illuminated with a retro-reflected FC creating a standing wave perpendicular to the cavity as depicted in Fig. 1(a). The laser frequency spectrum consists of equidistant narrow lines with spatial frequencies k p,m with a bandwidth ∆k [see Fig. 1(b)]. The line spacing is chosen as δk, such that the round-trip time of a pulse inside the cavity and the pulse repetition period of the ultrashort pulse train match [see Fig. 1(c)]. Pairs of cavity and FC modes can be tuned to be quasi-resonant up to a constant detuning ∆ c = ω p,m − ω m with the cavity frequencies ω m = ck m and the speed of light c. This detuning is adjusted to a value smaller than zero ∆ c < 0. Due to the quasi-resonance, light is scattered from the FC into the cavity via the atoms, where the coupling between pump and cavity mode pairs η is assumed to be constant in m. It depends on the pump strength, the atom-cavity coupling and the atomic detuning [4]. Finally, the so-created cavity field leaks out of the resonator with the rate κ, which we assume to be independent of m for simplicity.
The atoms are confined in a narrow dipole trap along the cavity axis (x axis), where all pump mode functions have a maximum. The theoretical description is thus reduced to a 1D model. Moreover we assume that the atomic cloud is distributed around the center of the cavity at x = 0 and has a size w, which is much smaller than the cavity length L. Then the cavity mode functions can be approximated by alternating cosines and sines, which we write as G c,m (x) = cos(k m x − φ m ) with φ m = 0 for odd m and φ m = π/2 for even m.

Mathematical model and its properties
When the temperature of the atomic cloud is larger than the recoil energy, i.e.
for all m, the coupled dynamics of the dilute gas and the In order to match to the cavity modes we choose δk ≡ π/L. (c) This line spectrum corresponds to a train of phase-stabilized ultrashort pulses in the time domain, which translates to a regular train of pulses with width 2π/∆k and distance 2π/δk ≡ 2L in space, along the FC propagation axis.
field is well described by the semi-classical stochastic differential equations [4] x j = p j m a (1a) coupling atomic positions x j , momenta p j and cavity fields α m . The diffusion due to photon loss is taken into account by the Wiener process ξ m , which obeys the noise correlations ξ m (t)ξ * n (t ) = δ(t − t )δ mn . The crucial quantities describing the amount of scattering into the m-th cavity mode are the atomic order parameters We neglect here the dynamical Stark shift of the cavity field and assume that the mode spacing is large compared to the coupling strength and cavity detuning such that there is no cross-scattering between different cavity modes [24].
In the single-mode case, it is well known that the described set-up and Eqs. (1) lead to a stationary state for ∆ c < 0 [25]. Below a threshold of the laser power, the stationary spatial distribution of the atoms is homogeneous, while above this threshold it exhibits a density modulation with the period of one wavelength [4]. With the rescaled pump intensity [25,26] the threshold is at ζ = 1 for an atomic cloud at the self-consistent stationary temperature The stationary momentum distribution is well approximated by a thermal state with the self-consistent stationary temperature T st as long as |∆ c |/ω R 1 and |∆ c |/ω R 4ζ [25], i.e. when the pump intensity is not too high. This sometimes allows for employing equilibrium thermodynamics to describe the stationary properties of the system, even though the system is open [26]. When the cavity parameters are tuned such that the stationary temperature is smaller than the initial temperature, the dynamics cools the atomic cloud -a phenomenon called cavity cooling, which allows for sub-Doppler cooling without relying on the internal structure of the atoms [3].
Within the thermal limit, the stationary momentum distribution is not altered when pumping two modes as long as κ and ∆ c are uniform [27]. In such a case the definition of Eq. (4) provides a unique stationary temperature. The thresholds for self-ordering and the spatial form of the stationary states, however, strongly depend on the mode structure [28]. Moreover, the type of the phase transition can change from second to first order, resulting in quasi-stationary states depending on the initial spatial condition and temperature [24].
In this work we present a further generalization to the two-mode case by pumping many adjacent modes, which have a similar frequency. Since there is no cross-scattering between the modes, the extension from the two-mode case is theoretically straightforward, the sums simply run over larger range. Thus the results about the stationary temperature can be generalized to many modes for uniform κ and ∆ c . The mode structure, however, is obviously changed leading to an alteration of the thresholds and the spatial stationary distribution.

Frequency comb and mode regimes
We now aim to identify the crucial properties and parameters of the mode structure of a frequency comb pumped cavity. The spectrum of the frequency comb can be formally written as a convolution of an envelope function η env (k) and a Dirac comb. In units of the pump-cavity coupling it takes the general form For our choice of uniform coupling η, the envelope function is a rectangle with height η and width ∆k centered around k c . The spatial shape of the laser pulses is related to the envelope via a Fourier transform and thus has a width of about 2π/∆k, a pulse-to-pulse distance 2π/δk and a characteristic sinc form [see Fig. 1 In this work the spectral lines of the FC match all longitudinal cavity modes within a certain range, leading to the mode spacing δk = π/L and a pulse distance of one cavity round trip length 2L. Put differentily, this means that the (spatial) shapes of adjacent modes dephase only on a length scale of the cavity length L, and since w L they do not differ within the atomic cloud. Thus adjacent sine-cosine mode pairs can be assumed to effectively have the same frequency. Within this "dense" mode parameter regime wδk/π 1, the behavior is left to depend on the pulse width in relation to the atomic cloud size, i.e. on the parameter For χ 1, the spatial pulse width is large compared to the cloud. In this regime all modes have the same shape within the extent of the atomic cloud, and thus can be effectively assumed to have the same frequency [see Fig. 2(a)]. In the center of the cavity there are thus two groups of spatially coinciding modes, M/2 sine modes and M/2 cosine modes. Since the physics can be effectively described by two modes, we refer to this regime as the effective two-mode regime.
For χ 1 the spatial pulse width becomes smaller than the atomic cloud size and thus introduces a new length scale to the system. Modes which are far apart in the spectrum spatially dephase within the atomic cloud [see Fig. 2(b)], and thus we refer to this regime as the dephasing regime.
In the following we will investigate these two regimes and the transition from one to the other by tuning χ. Practically, this is done by increasing the bandwidth while keeping the mode density and the atomic trap size constant, or by increasing the size of the atomic cloud keeping the bandwidth constant.

Stationary states
In this section we consider the form of the stationary states and aim to find the threshold for different χ. For this we first introduce the stationary potentials. The stationary The smallest length scale is given by the comb line wavelengths as e.g. the central comb wavelength λ c = 2π/k c . The beat length between neighboring comb lines is π/δk ≡ L, which typically is much larger than the atomic cloud of size w. In the effective two-mode regime (a) 2π/∆k is larger than the atomic cloud w and thus none of the involved modes differs significantly within the atomic cloud. The blue curve shows the fields of 8 adjacent cosine modes over the interval w = 20λ c with 2π/∆k ≈ 45.9λ c . (b) Dephasing regime: Distant modes run out of phase within the atomic cloud size and the spatial laser pulse width is smaller than the atomic cloud. The green curves show 50 cosine modes with 2π/∆k ≈ 9.2λ c , which clearly dephase within interval w = 20λ c . The mode spacing δk is the same in both cases. For visual clarity we depict only the cosine modes, omitting sine modes.
are attained after a time κ −1 and can be calculated from Eq. (1c) by settingα m = 0 and omitting the noise. Replugging this expression into Eq. (1b) yields an adiabatic force, which can be integrated to the single-particle potential This is a mean-field potential for a single particle at position x created by all particles via Θ m (x 1 , ..., x N ). Moreover, summing over all particles reveals the total potential energy [24] U where the square of the total order parameter is the crucial quantity describing ordering to all modes. A minimization of the potential energy thus corresponds to a maximization of the square of the total order parameter.

Low energy states far above threshold
For a high pump power ζ far above the self-ordering threshold, the potential depth is large compared to the stationary temperature T st . The stationary spatial distribution then minimizes the potential energy U and maximizes the square of the total order parameterΘ 2 . It is thus worthwhile taking a closer look at this quantity.
As discussed in the previous section, adjacent cosine-sine pairs have the same effective frequency within the atomic cloud. In the expression in Eq. (10) we can thus summarize these adjacent pairs using the identity cos(α) cos(β) + sin(α) sin(β) = cos(α − β) for each pair, and obtainΘ where the factor 1/2 emerges since we still sum over all M modes. Rewriting this expression as a convolution of a rectangular envelope and a Dirac comb in (spatial) frequency space as in Eq. (5), further yields the expression Here we made an approximation by considering only a single sinc, which is very accurate since the next sinc is shifted by 2π/δk = 2L and is thus far outside of the atomic cloud. Note that in the dense mode regime the total order parameter is bound byΘ ≤ 1/ √ 2. This is a consequence of the fact that the atoms cannot simultaneously order to a sine and a cosine mode with the same wavelength.
The multi-mode behaviour enters via the sinc depending on the parameter χ. In the effective two-mode regime χ 1, the sinc is always approximately 1. A low energy distribution is thus left to maximize the cosine, which favors ordering of the atoms in a λ c -periodic grating, as in the single-mode case. The absolute position in space, however, is not fixed since only the distances of particles enter. This translational invariance comes from the coupling to sine and cosine mode pairs with the same (approximate) wavelength. Considering one mode pair only, the translation over one wavelength of an atomic distribution would then shuffle photons from the sine mode to the cosine mode and back, while leaving the potential energy unchanged. This resembles the U(1) symmetry described for crossed cavities [29] and atoms in ring cavities [30]. Let us emphasize that in the present case this only valid close to the center of the cavity and when the dynamical Stark shift of the individual modes can be neglected.
In addition to the λ c -periodic grating, a localization of the atomic cloud to some width smaller than the pulse width 2π/∆k becomes energetically favorable in the dephasing regime (χ > 1). In such a way the atoms stay in the in-phase region of the modes [see Fig. 2(b)] and the sinc takes a value close to one. A wider distribution could scatter only a narrower bandwidth and thus less power.

Self-ordering threshold from a mean-field model
We now aim to find the threshold of the total laser pump power ζ tot = M ζ for selfordering at the self-consistent temperature T st as a function of χ. Remembering that the stationary state can be approximated by a thermal state in typical parameter regimes [26], we assume that the stationary state is a Boltzmann distribution with the temperature T st and the mean-field potential given in Eq. (8). To obtain the self-consistent cavity fields creating this potential, we employ the method introduced in Ref. [31] for the single-mode case. After choosing an initial spatial atomic distribution P 0 (x), one calculates its order parameters Θ m [P 0 ], which can be written as analogously to Eq. (2). These specify a mean-field potential U MF (x; Θ 1 [P 0 ], ..., Θ M [P 0 ]) created by the distribution P 0 . The updated atomic configuration P 1 in the next iteration is the normalized Boltzmann distribution in this potential with the temperature T st [Eq.
(4)]. In this vein the n-th step of the algorithm can be formulated as The iteration is then repeated for N iter times, until the distribution does not change anymore. Slow ramping of the laser power over the threshold can be emulated by the following procedure. For a fixed χ we start with a low pump strength below the threshold and a homogeneous distribution, and calculate the stationary state using the above described algorithm [Eq. (14)]. In the next step, we use this stationary state as initial state for the iteration with a slightly higher pump strength. Continuing this way we simulate an up-sweep of the laser power. Analogously, for the down-sweep we start with the stationary state obtained from the up-sweep far above the threshold and then reduce the pump power in small steps.
The atomic distribution is confined to a finite interval with the width w = 300λ c centered around the cavity center. The parameter χ is varied by changing the bandwidth ∆k while keeping the mode spacing (or mode density) δk constant. The maximal bandwidth we consider is ∆k max /k c ≈ 0.0213 = 2.13% with M = 998. Correspondingly,  In the two-mode case (χ ≈ 0), these curves coincide and we get the a typical second order self-organization phase transition. In the multi-mode regime for larger χ, jumps in the order parameter and a hysteresis as signs of a first order phase transition appear. (b) The corresponding phase diagram as a function of χ and the laser power ζ tot . The red curves depict the thresholds for an upward-sweep (solid line) and a downward-sweep (dotted line) of the laser power. The black curve is a theoretical prediction of the phase boundary depending on the mode structure given in Eq. (17). In the large χ-regime there is a metastable region and hysteresis. The insets depict atomic distributions at different positions of the phase diagram above the threshold. Apart from the wavelength grating, the distribution is flat for small χ (blue) and is localized for large χ (green). the mode spacing is given by δk/k c = ∆k max /k c /998 = 2.13 · 10 −5 , which is deep in the dense mode regime (wδk/π ≈ 0.013 1). Note that in theory the maximum number of modes (here 998) can be chosen quite arbitrarily as long as the dense mode condition wδk/π 1 is fulfilled, while in the experiment it is fixed by the comb-cavity set-up. The cavity parameters are κ = 400ω kc R and ∆ c = −κ and the stationary temperature hence has the simple form k B T st = κ/2.
As shown in the insets of Fig. 3(b), above the threshold the atomic distributions have a form similar to the low energy states maximizing Eq. (10) discussed in the previous section 3.1. In the effective two-mode regime (χ 1), the λ c -periodic density modulation extents over the whole interval and the peaks get narrower for larger pump intensities (blue). In the dephasing regime (χ > 1), the atoms additionally localize within an interval smaller than the pulse width, which is 46.9λ c for χ ≈ 6.4 (green).
Interestingly, for small χ the order parameter continuously grows at the transition point, and the results from up-and down-sweep coincide [see Fig. 3(a)]. For larger χ, in contrast, the order parameter exhibits a jump at the transition point and the threshold depends on the sweep direction. This hysteresis hints that the type of the phase transition changes from second order to first order upon increasing χ.
The threshold ζ c tot increases as a function of χ (see Fig. 3), signifying that more total input power is needed to create a sufficiently deep potential to stabilize an ordered pattern for larger χ. The reason for this is that upon increasing the bandwidth, less modes are being scattered into by an atomic distribution with a fixed width and less power is transmitted to the cavity. In other words, at some χ not the whole spectrum can be transmitted anymore.
In this vein, in order to find an expression for the multi-mode threshold, we rescale the single-mode threshold ζ c = 1 by a factor representing the efficiency of power transmission into the cavity compared to the single-mode case. This efficiency factor is given by the square of the total order parameterΘ 2 [P c ], where P c is an atomic distribution which optimally couples to a single-mode, i.e. maximizes the cosine in Eq. (12). That is, P c is a λ c -periodic grating with some arbitrary width. In this quite general sense the multi-mode threshold is modeled by A rather simple expression can be derived for the up-sweep. Coming from below the threshold, the atoms are distributed over the whole interval, and thus P c has the width w. Explicitly, we write with N λ = w/λ c . With this distribution, the estimated multi-mode threshold becomes This quantity is depicted as a black line in Fig. 3(b) and very well resembles the results from the self-consistent iteration. The situation for the down-sweep is different, since by coming from above the threshold the width of P c depends on χ. Figure 3 shows the numerical results for the down-sweep boundary. Again, for larger χ, the up-and down-sweep phase boundaries depart from each other, spanning a region where several metastable states exist.

Dynamics
In addition to the analytical and stationary considerations above, we now simulate the dynamics of the frequency comb set-up by numerical integration of the stochastic differential equations of motion given in Eqs. (1).

Self-ordering
First we compare the dynamics of N = 100 atoms for different pump strengths for the initial temperature k B T = κ/2. The initial spatial distribution is uniform on an interval of the width w ini , which we choose as either 20λ c or 300λ c . Note that the situation in the dynamical case here is different to the self-consistent calculations above since the particles can move out of the interval and due to the relatively small particle number.
Again we choose the cavity parameters as κ = 400ω kc R and ∆ c = −κ. The bandwidth ∆k max /k c ≈ 2.13% is spanned by 50 modes now, leading to the mode spacing δk/k c = ∆k max /k c /50 = 4.26 · 10 −4 . That means that we are still in the high mode density regime where adjacent modes do not dephase significantly. Figure 4 depicts the trajectory averages of the total order parameter Θ , the total intra-cavity photon number n = m |α m | 2 and the kinetic energy E kin = p 2 /(2m a ). For relatively free particles it is related to the temperature by E kin = k B T /2 due to the equipartition theorem. In Fig. 5 the position of the particles is plotted for an example trajectory. Moreover, we depict the averaged field distribution inside the cavity with the normalized spatial form of the m-th cavity mode field Θ m cos(k m x − φ m ). The envelope of this quantity gives insight about where in space the particles are ordered and can be extracted as the smoothed field distribution with the rectangle function rect(x) = 1 for −1/2 < x < 1/2 and 0 else. From our simulations we obtain the time evolution of the trajectory averaged order parameters Θ m (t), which enables us to consider the time evolution of the function F(x). Note that also experimentally the order parameters can be indirectly measured via a measurement of the cavity output fields. When the dynamics of the cavity fields α m is much faster than the atomic motion, we can use the stationary approximation of the fields from Eq. (7) to obtain an approximate expression for the order parameters 20) in terms of the measurable field quadratures Finally, Fig. 6 depicts the laser output intensity for the trajectories in Fig. 5 at the final time κt = 10 6 . We start with the case of two modes for comparison, labeled (A) in Figs. 4, 5 and 6. The total and single-mode pump strength are ζ tot = 15 and ζ = 7.5, respectively, where the latter is much larger than one and thus far above threshold. The dynamics can be split into two time regions: The fast initial increase of the order parameter Θ and the photon number n is accompanied with an increase of the kinetic energy. Here the atoms fall into their self-created potential. On a much longer time scale, the atomic cloud is cooled down to the approximate stationary temperature (dashed line), leading to further localization in the potential and thus a further increase of the order parameter. This characteristic dynamics resembles the single mode case and has been discussed in detail e.g. in Ref. [32]. In contrast to the single-mode case, the potential is nearly translationally invariant, which allows for a collective translational motion of the ordered atoms [see Fig. 5(A)]. The field distribution F(x) is spatially uniform for two modes (cosine-sine pair), since the light fields do not contain any information about the absolute position. Just as the cavity input light, the output light is continuous for two modes and two comb lines (see Fig. 6).
Case (B) shows the multi-mode scenario (M = 50) for the same total pump intensity ζ tot = 15 (ζ = 0.3) and a small initial width w ini = 20λ c . Since χ ini = w ini ∆k/(2π) ≈ 0.426 is well below one, the initial distribution is narrower than the pulse width and the atoms are in the in-phase region. Due to this and since the pump strength is above the multi-mode threshold, the atoms immediately order into a (approximately) λ c -periodic grating, which hinders the atomic cloud from expanding. The atoms thus stay in the in-phase region and the curves of the collective quantities in Fig. 4 closely resemble the two-mode case. When an atom leaves the in-phase region, it is not trapped anymore, in contrast to the two-mode case. Also the field distribution F(x) is not uniform anymore, but is localized on the dephasing length scale. The cavity output field phases thus include information about the absolute position of the ordered atoms. As depicted in Fig. 6, the cavity output is pulsed and the shape of the input pulse is not considerably changed by the Bragg reflection on the atoms. The situation drastically changes for an initial distribution of w ini = 300λ c , which is much wider than the pulse width (χ ini = w ini ∆k/(2π) = 6.39 > 1). Case (C) depicts this situation with the same parameters as in (B). Initially, there is only a small increase in the order parameter since the atomic cloud is too wide to order to all modes. This comes hand in hand with a smaller increase of the kinetic energy, and faster cooling, which will be examined in the next section. The potential depth is too small to permanently trap the atoms, and does not hinder the expansion of the atomic cloud (see Fig. 5). With increasing (dynamical) width of the cloud, less and less modes can couple to the cloud and the total order parameter starts to decrease again. At intermediate times the field distribution shows several regions of ordering, which later disappear.
Finally in (D) we consider a very large pump intensity ζ tot = 60 above the singlemode threshold (ζ = 1.2). In spite of the wide initial distribution w ini = 300λ c , the cloud organizes into some stable pattern, since even a single-mode can sustain ordering here. This initial freezing does not only impede the expansion, but also further localization of the atoms. As depicted in Fig. 5 and the inset in Fig. 6, the atomic cloud orders to different modes in different points in space. This leads to modified output pulse shapes, as depicted in Fig. 6. Summarizing, we find that the initial width of the distribution is the main factor deciding whether strong ordering to many modes occurs. Initial widths smaller than the dephasing length often result in ordering which can be stable on large time scales even below the multi-mode threshold obtained using the self-consistent iteration. For initial widths larger than the dephasing length instead, the atoms do not find the multi-mode stationary state in the cases we considered and order only to few modes simultaneously. Even for pump strengths slightly above the multi-mode threshold obtained before, the trapping is not strong enough to stabilize the atoms. Finally, above the single-mode In the case of an atomic distribution localized in the cavity center, these pulses are transmitted to the cavity without changing their shape significantly. The black curve shows the input pulse shape for reference. (C) The weakly ordered case leads to a weak cavity output with irregular aperiodic structure. (D) A self ordered atomic distribution results in periodic pulses of different shape representing the intracavity order. Since ordering happens at several locations the output pulses have multiple peaks breaking the pump pulse symmetry. Insets show the order parameters for a range of cosine (magenta) and sine (cyan) modes. threshold, where one mode is sufficient to sustain ordering, the atoms freeze close to their current position. A better ordering also in the case of wide initial distributions might be obtained by slowly ramping up the pump strength to its final value or by using higher initial temperatures [24].

Cavity cooling
A central application of the coupled atom field dynamics is motional cooling without spontaneous emission [1], which is known to have some important restrictions for larger ensembles [7]. Here we finally focus on the cooling capabilities in the extreme broadband case for different bandwidths and pump strengths. For this we simulate the dynamics of N = 1000 atoms, starting with an initial temperature k B T ini = 10 κ/2 which is ten times higher than the expected stationary temperature. Moreover we add a harmonic trap with the frequency ω T = 0.046ω kc R preventing the atomic cloud from expanding when they are not trapped by the self-consistent cavity potential. This is done by adding a restoring force −m a ω 2 T x in Eq. (1b). The resulting dynamics is depicted in Fig. 7 for different pump strengths and  Figure 7. Dynamics of (a) the total order parameter Θ and (b) the kinetic energy E kin for N = 1000 with an initial temperature k B T ini = 10 κ/2 for different mode numbers (or bandwidths) M = 1 (solid), M = 2 (dashed), M = 10 (dotted) and M = 50 (dashdot). The depicted quantities are averages over 25 trajectories. While for ζ tot = 0 (blue) there is obviously no cooling, a non-zero pump strength leads to cooling. When pumping below the multi-and single-mode threshold ζ tot = 0.5 (orange) there is no ordering and no difference in the kinetic energy curves for different mode numbers (or bandwidths). For ζ tot = 1.5 (green), below the multi-mode threshold, but above the single-mode threshold we observe slower cooling for the M = 1 and concurrent ordering. This difference becomes large for ζ tot = 3 (red).
bandwidths. We ask how the cooling speed differs depending on the bandwidth for a certain ζ tot . When pumping below the multi-and single-mode threshold ζ tot = 0.5 there is no ordering and the kinetic energy curves for different bandwidths (or mode numbers) coincide (see the orange curve). Below the multi-mode threshold, but above the single-mode threshold (for ζ tot = 1.5 < 2) there is ordering for M = 1 only. This leads to slower cooling for M = 1 [25] compared to M = 2 or M = 50, as can be seen from the zoomed inset. For ζ tot = 3 (red curve) this anti-correlation between ordering and fast cooling becomes even more evident. It is made explicit in Fig. 8, where we plot the final values as a function of the bandwidth. The kinetic energy at κt = 10 6 is decreased by a factor of about 0.49 for the bandwidth ∆k max /k c ≈ 0.0213 = 2.13% at M = 50 compared to the single-mode case. Note that the kinetic energy follows a damped oscillation. When the cloud does not organize fast enough, it expands freely and bounces back from the fixed harmonic trap, leading to the oscillations in the beginning.
As already known from previous work [25] the cooling speed increases with larger pump strength until the self-organization threshold is traversed. Then the ordering of the atoms slows down cooling and, for high pump strengths, also increases the stationary temperature. Using a train of short pulses with a large bandwidth shifts up the (multimode) self-organization threshold, avoiding self-ordering up to large pump powers. This allows for implementing faster cavity cooling.

Conclusions
In this work we theoretically studied the dynamics of a cold cloud of atoms in a multi-mode cavity which is illuminated by a train of phase-stabilized ultrashort pulses with a wide frequency spectrum. For red cavity and atomic detuning (high field seeking particles), the lowest energy state of the system is connected to a spatial atomic distribution, which optimally scatters the majority of the pulse frequency components into the cavity. In order to fulfill this requirement, in addition to a periodic Bragg grating with the period of the wavelength known from single-mode self-organization, the atoms need to be localized on a length scale smaller than the pulse length to adapt optimally to the mode structure. We have shown that the threshold power required for self-ordering increases when the spatial pulse width is smaller than the atomic distribution as it is not possible for all the particles to commonly scatter in phase and ordering to all modes becomes increasingly difficult. Thus even at optimal single-mode order, less power is transmitted to the cavity and a fracturing of the atomic density into several lobes appears. Moreover, in the corresponding simulations we identify a hysteresis hinting that the phase transition switches from second to first order in the extreme multi-mode regime.
Dynamically, the atoms do not attain this localization themselves. For high pump powers, they freeze in the wavelength Bragg grating, while for weak pump powers the atomic cloud expands out of the in-phase region of the modes. On the other side such a spatial dephasing of the modes generally has a positive effect on the cooling time and temperature. Using a large bandwidth avoids self-ordering and the concurrent increase of kinetic energy at too low pump power, while cavity cooling itself is hardly affected by the bandwidth.
Since the transmission of the pulse into the cavity depends on the atomic pattern, the cavity output intensity as a function of time can be used to gain information about the longitudinal form and location of the atomic cloud. This is similar to an atomic kaleidoscope [33] but operated in the longitudinal direction in the time domain. The formation of such an output pulse train with different periodicity and shape than the input pulse can be viewed as a spontaneous formation of a time crystal. While many of the discussed phenomena might survive in the zero temperature regime of quantum gas cavity QED, new phenomena could be expected from direct inter-particle interactions.