domain: spi phys

text:
Jouke Hijlkema
email: jouke.hijlkema@onera.fr

Marco De Rosa
email: marco.de.rosa@esa.int

Experimental determination of thermodynamic properties of Hydrazine, MMH and NTO along the saturation line up to the critical point

Thermodynamic properties as density and the speed of sound as a function of pressure and temperature are not very well-known for higher temperatures for typical storable propellants as Hydrazine, MMH or NTO. Most data is based on low temperature experiments that have been extrapolated or calculated from theoretical correlations. To our knowledge no experimental data exists up to the critical temperature except for a single experiment with Hydrazine that was carried out by Lobry de Bruyn [4] in 1896 and suffered from a catastrophic decomposition while trying to reproduce the results.

This paper deals with the design, development and use of an experimental installation capable of measuring the density and the speed of sound of the liquid phase as well as the density of the gaseous phase for these storable propellants from (close to) the triple point up to the critical point. Utmost care has been given to the safety aspects of the remotely operated installation. Low cost solutions have been devised to respect a rather tight budget. The installation is currently in its trial and calibration phase using pure water, this paper contains some results of these tests. The real propellant tests are bound to start in the coming months.

Glossary

moment of inertia Quantity expressing a body's tendency to resist angular acceleration I(T ) 5 Density Mass per unit volume ρ 5 density measurement Pixel correspoding to the position of the blade at temperature T Px(T ) 6 gravitational constant Acceleration due to the atraction of the earth g 5 speed of sound propagation speed of sound in a medium c 1, 8, 9 thermal expansion coefficient Measure of the expansion of a material due to temperature changes α 5 Young modulus Measure of the stiffness of a material E(T ) 5 zero position for the density measurement The pixel corresponding to the position of the unloaded blade

Introduction

The determination of the density and the speed of sound for storable propellants like hydrazine or NTO is a complicated affair. Due to their inherent instability and toxicity, extreme care needs to be taken to avoid catastrophic events and protect personnel. This is the main reason why experimental data of this order is extremely rare.

ESA issued an invitation for tender in 2016 entitled "EXPERIMENTAL DETERMINATION OF STORABLE PROPELLANT PROPERTIES FOR MODELLING PURPOSES".

This article reflects the actual state of the research done at ONERA in this framework. A quick overview of the hardware components, their main functional principles and error estimations is followed by a section dealing with the calibration of the machinery. Some preliminary results from this calibration campaign will be presented as well.

Hardware

The hardware consists primarily of a high pressure vessel with 2 double-glazed sapphire windows. The external walls of the vessel are equipped with 1 kW electrical heaters. 2 isolation valves separate the test vessel from the rest of the hardware which comprises a buffer reservoir, a piston pump feed system, a neutralisation vessel, a vacuum pump, a purge system and a pressurisation system. The whole is operated by an in-house build, Raspberry Pi based, remote control and data acquisition system. The test vessel is supported by an iso-static positioning device that assures that the origin of vessel does not move due to thermal dilatation. This is important for the density measurements (see section 4.2 for details).

Pressure vessel

The pressure vessel is made of a titanium [START_REF] Harold | Theoretical prediction of the thermodynamic behavior of aqueous electrolytes at high pressures and temperatures; II, Debye-Huckel parameters for activity coefficients and relative partial molal properties[END_REF] , cylindrical hollow body with matching titanium top and bottom covers. The outside is covered with 4 125 W lateral heater panels and each cover has a 250 W heater around its edge. Optical access is provided by 2 double surface sapphire windows that allow for a direct view of the internal cavity from 2 opposite sides (See figure 1b for an impression). Because of the nature of the test fluids, it is not possible to plunge a thermocouple directly inside the fluid. Therefore, 2mm diameter holes have been drilled in the vessel wall up to 1 mm from the surface of the internal cavity that allow for the placement of 4 thermocouples as close to the test fluids as possible. The vessel is designed with a burst pressure of 450 bars and has been successfully subjected to a hydrostatic proof pressure test up to 385 bars. To reduce the risks for personnel to a maximum, the test set-up is remote controlled for all but one operation (the filling of the buffer reservoir with the test fluid). This is done through an in-house built remote control and data acquisition system depicted in figure 2. A supervisor (a Python program running on a laptop) is connected to our local network as well as a NAS, the main pilot and the pump pilot.

Remote control and data acquisition system

The last 2 are Raspberry Pi mini-computers situated next to the experimental set-up. Network communication is based on the excellent ZeroMQ library [START_REF] Julien | Characterization and modeling of forged Ti-6Al-4V Titanium alloy with microstructural considerations during quenching process[END_REF] .

The main pilot talks to a DAC3 over a USB link and operates the valves, heater and lights by means of a in-house built controller card (see figure 1e). The pump pilot controls the stepper motor control card that drives the feed system.

A separate (also in-house built) video control system is running on a pc in the remote control room (my desk) at a safe distance from the experiment. This machine controls the live streams from cameras around the set-up as well as the density camera (see section 4.2)

Feed system

The feed system consists of a stepper motor operated piston pump, a 3 way valve and a buffer reservoir. All wetted surfaces are in Teflon to ensure material compatibility. Table top experiments have shown a precision better than 0.01 ml. Figure 1d shows the entirety of the feed system. Before the experiment we fill the buffer reservoir with argon to chase all the air up to V1 (see figure 1a). Then, valves V5, V4, V3, V2 and V1 allow us to create Once the desired pressure drop has been achieved, we close V3 and toggle V4 so that in case of an emergency opening V3 will evacuate the product through the neutralisation vessel. We zero the density measurement system (see section 3), fill the buffer reservoir with the desired quantity of the test fluid (approximately 100 ml), close the lid and leave the test room. This is the only human intervention, the rest of the operations are piloted from a safe distance (200m).

Neutralisation vessel and purge system

Once the experiment is completed the vessel needs to be decontaminated and purged before any human intervention. The entrance of the neutralisation vessel has a diameter of 1 mm to make sure the flow is choked for high pressure conditions to limit the mass flow. The vessel is filled by 2/3 with silicon grit. Before the experiment we top off the grit with a neutralisation fluid (water for NTO or bleach for hydrazine and MMH). With the vessel at a sufficiently high temperature we open the upper isolation valve V3 and route the exit through the neutralisation vessel by using valve V4 (see 1a). The rest heat in the 15 kg heavy test vessel is sufficient to evaporate all the liquid in the system and the high pressure will chase the major part of the product through the neutralisation fluid. The silicon grit breaks the bubbles and extends the transit time and hence the exchange with the neutralisation fluid. The exit of the neutralisation vessel is connected to a closed container that collects the resulting run-off liquid.

Once the system has cooled down sufficiently we operate valves V5, V4, V3, V2 and V1 as well as the feed system pump to rinse all of the interior. The reservoir has an overflow port that connects to the same run-off container. Cameras at both the (transparent) exit of the neutralisation vessel and the reservoir overflow port allow us to verify that both extremities of the system have been rinsed.

Measurement systems

3.1 Data acquisition, pressure and temperature measurements

DAQ

Data acquisition is ensured by an in-house built machine since our existing acquisition system is not suited for experiments that take up to 8 hours. We use a Acquisys USB-1608G4 DAQ. These machines have 16, 16 bit, analogue input ports and are capable of a total acquisition frequency of 250 kHz (multiplexed). Data is extracted over a USB link using the Python uldac package [START_REF] Lyakhovitskii | Experimental investigation of acoustic properties of titanium alloys in the temperature range of 20-1000°C[END_REF] and send to the supervisor for storage every 10 seconds. The supervisor sends a copy to the NAS for safe keeping. To reduce noise, temperature and pressure measurements are repeated 1000 times and averaged for each measurement point. This results in a total maximum data acquisition frequency of around 1 Hz.

Thermocouples

Since all known temperature sensors are incompatible with one or more of the storable propellants, the thermocouples (classical type K with a maximum continuous temperature of around 1100 °C) that will monitor the temperature in the test cavity are placed inside the wall of the test cell in holes drilled to within 1 mm of the inner cavity. These thermocouples need to be calibrated in situ in order to take into account the heating of the section of the wiring that touches (or is close to) the hot surface of the test vessel since the correspondence of the measured potential and the actual temperature difference between the thermocouple and the cold junction of type K thermocouples are given for a homogeneous cable temperature. A secondary objective of these calibration test runs is the determination of the time needed for the system to reach thermal equilibrium. When in equilibrium, the temperatures measured inside the wall are considered to be equal to the temperature of the test fluids. To verify this hypothesis a port in the upper cover of the vessel allows the passage of 2 thermocouples. These are plunged in the liquid and the gaseous phase. Once the calibration phase was over this port has been welded shut to prevent micro leakage as much as possible. Section 4.1 deals with these calibration runs.

Pressure transducers

Pressure is measured using GE UNIK 5000 pressure transducers [START_REF] Peytour | Characterization of ceramic/TA6V titanium alloy brazed joints[END_REF] . For the measurements from ambient to critical pressure a 500 bar max range is used. For the measurements from ambient to the triple point pressure we use a 1 bar max range. These sensors have a combined accuracy of ± 0.04% full scale. This gives a precision of 200 mbar for the 500 bar max range and 0.4 mbar for 1 bar max range. The principle of the density measurements is simple. As shown in figure 3, a floater of volume V = d 3 exerts a force F on the tip of a titanium blade.

If we ignore the influence of the pressure on the elasticity of the titanium blade then the deformation δ is a function of the temperature T and the density of the fluid ρ(T ):

δ(T ) = F(T )L 2 (T ) 3E(T )I(T ) (L(T ) + 0.5l(T )) (1)
with E(T ) the Young modulus of TA6V as a function of temperature and I(T ) the moment of inertia. If we define β T = 1 + α(T -T 0 ) with α = 9.10 -6 m/m°C, the thermal expansion coefficient of TA6V, we can express any dimension D(T ) as D(T ) = D 0 β T . Then equation ( 1) becomes:

F(T ) = δ(T )E(T )WH 3 β T 4L 2 (L + 0.5l) = δ(T )E(T )β T Γ (2)
where W and H are the width and height of the blade and Γ = WH 3 4L 2 (L+0.5l) . We know that F(T ) = g(m f -ρ(T )d 3 β 3 T ) with m f the mass of the floater so we can rewrite equation (2) as:

ρ(T ) = m f g -δ(T )E(T )β T Γ gd 3 β 3 T (3)
Classically the temperature dependency of the Young modulus of a metal is considered linear and can be written as This allows to determine the density of the test liquid knowing the temperature T and the deformation δ(T ).

E(T ) = E 0 (1 -a T T m ) = a + bT

THERMODYNAMIC PROPERTIES OF HYDRAZINE, MMH AND NTO

The experimental set-up measures δ(T ) indirectly as the position of a mark on the blade in an image captured with an IDS UI-5480S7 camera equipped with a SIGMA APO Macro 180 mm 1:3.5 D lens. The correspondence between the recorded pixel and the deformation δ(T ) is related to the size of the pixels and the magnifying factor of the optics on one hand and the zero (vertical position of the camera) corresponding to the unloaded blade. This can be expressed as a linear transformation

δ(T ) = c * (Px(T ) -Z) (5)
where Px(T ) is the pixel recorded for a temperature T and Z is the pixel corresponding to the position of the unloaded blade. The coefficient c is related to the pixel size and the optics and can be supposed constant between experiments, Z, however, is sensible to eventual disturbances (movement of the camera or the vessel during preparations, thermal expansion during the test runs is countered by the iso-kinetic support of the test vessel) and needs to be determined before each test. To determine c we fit equation ( 5) to the results of a calibration run with water. Figure 8 shows the results of this operation giving c ≈ 7.366.10 -6 . Now the only thing left is to determine Z before each test. This is easy, we simply record the pixel Px 0 corresponding to position of the blade after pulling a vacuum and before filling the vessel. We know that this corresponds to the deformation due to the weight of the floater (see table 2) and with equation ( 4) and the hypothesis that E(T ) = E = 108 GPa and β T = 1 we see that the corresponding displacement δ 0 is equal to:

δ 0 = 0.042585 1.11932.01 -10 108.10 9 ≈ 0.003523 m (6) then Z = Px 0 - δ 0 c 3.3.

Error estimation

The error estimation for the density measurements can be divided in two.

• The errors due to non-nominal behaviour (drops on the blade, bubbles on the floater, ...)

• measurement errors

The measurement errors are relatively easy to assess. We use the temperature measurements for which we have an error estimation T ≈ 1°(see section 4.1) and images from the camera that introduce an error Px ≈ 6 pixels (from experience). For the other values the error does not come into play since they are constant and included in the calibration so the effective error is 0. Then, if we write the density as a function (see equation 4) of two independent variables ρ(T, Px) :

ρ(T, Px) = C 1 f (T ) + C 2 g(Px)h(T ) with C 1 ≈ 4.250.10 4 C 2 ≈ 1.141.10 -5 f (T ) = 1 (1 + αT ) 3 g(Px) = c(Px + Z) h(T ) = a + bT (1 + αT ) 2
Then we can express the error in ρ(T, Px) as

∆ρ(T, Px) = ∆T δρ δT 2 + ∆Px δρ δPx 2 so with δρ δT = C 1 f (T ) + C 2 g(Px)h (T ) δρ δPx = C 2 g (Px)h(T )
and

f (T ) = - α 3(1 + αT ) 4 g (Px) = c h (T ) = b 1 + αT - α(a + bT ) 2(1 + αT ) 3 we find ∆ρ(T, Px) = ∆T 2 η(T, Px) + ∆Px 2 ζ(T, Px) with η = -C 1 α 3(1 + αT ) 4 + C 2 g(Px) b 1 + αT - α(a + bT ) 2(1 + αT ) 3 2 ζ = C 2 c a + bT (1 + αT ) 2 2
The displacement of the blade is less than 1 mm so g(Px) is inferior to 10 -3. This permits us to fill the following The error due to non nominal behaviour is a lot harder to assess. However, the calibration tests can give us an idea. If we look at figure 8 and we assume that for the third measurement point the error is due to this type of behaviour we see that we would find ρ ≈ 780 kg/m 3 instead. An error of 15 kg/m 3 , 100 times bigger than the typical error in table 3 but still lest than 2%. Once the liquid phase densities are obtained, the gas phase density can be derived as follows. An insert made out of titanium divides the vessel in a lower and upper volume separated by a narrow channel passing in front of the windows. This channel allows for the precise determination of the liquid volume.

Figure 4 shows a schematic overview. Given the liquid density ρ l measured previously (see section 4.2) and the total mass m t of the product, the volume change observed in the capillary section of the cell is a measure of the gas density.

V g = V t -V l ρ g = m t -V l ρ l (T, P) V g

Error estimation

If we state that m t = V l ρ l then we see that the gas density depends on the measurement of the liquid volume, the injected mass as

ρ g (m t , vl) = m t -V l ρ l V g
and we can estimate the measurement error as The speed of sound is measured in the liquid phase by means of an emitter/receiver pair of ultrasonic sensors8 . These sensors are not able to withstand high temperatures and are therefore retractable. We measure the time of flight of the ultrasonic wave through the test vessel and the test liquid (see figure 5). Then we know that

∆ρ g = ∆m t 1 V g 2 +       ∆V m t V 2 g       2 So,
t f light = t ext + t cav = l ext c ext + l cav c liq
where l ext = l 1 + l 2 + l 4 + l 5 and t ext is the flight time of the pulse through l ext . The second term l cav c liq is the flight time through the cavity. If we ignore the influence of the internal pressure on the path lengths (at 345 bars, the maximum deformation of the vessel is 0.082 mm for an initial dimension of the order of 50 mm hence in the order of 0.15%) than the path lengths l ext and l cav are dependent only on the temperature through the thermal expansion coefficient α

l ext (T ) = l ext (T 0 )(1 + α(T -T 0 )) l cav (T ) = l cav (T 0 )(1 + α(T -T 0 ))
The value for α can easily be found in the literature [ [5], [6]]. Now if we do a calibration run with water as the test fluid in which we measure the flight time t H 2 O (T, P) as a function of the temperature and pressure we can write

t ext (T ) = t H 2 O (T, P) - l cav (T ) c H 2 O (T, P)
where c H 2 O (T, P) can be found in the literature [ [3][1]] and all the other terms are known. Once we have t ext (T ) (which is independent of the test fluid and the pressure) we can use it to determine c prod (T, P) for all of our test fluids by measuring the flight time t prod (T, P) c prod (T, P) = l cav (T ) t prod (T, P)t ext (T ) (7)

During the calibration tests we will determine t ext , this is not finished yet.

Error estimation

Equation ( 7) allows us to write:

∆c(t prod ) = ∆t prod C1 (t prod + C2) 2 (8)
With C1 ≈ 0.07 and the values found in table 5 we expect that t prod -C2 ∈ [3.3.10 -5 , 9. If we assume that the error in the measurement of the flight time is of the order of the ultra sound frequency, that is 10 -6 (we wont make an error of more than one cycle) then

∆c(t prod ) ∈ [8.1, 64.3] m/s
this is an error of between 0.96% and 3.1 %.

Gas phase speed-of-sound measurements

Due to the large difference in density between the titanium body of the vessel and the gaseous part of the propellant it is not possible to use the same ultra sound technique as used for the liquid phase (see section 3.5). A completely different approach has been tested to deal with this.

The gas trapped above the liquid phase forms a cavity that we can try to excite acoustically by making the walls vibrate with a well chosen frequency sweep. If we then measure the unsteady component of the gas pressure we should be able to determine the resonance frequency F res and hence the speed of sound inside the gas (We assume the longitudinal dimension L of the cavity known, this is valid since the eventual dilatation effects due to the temperature variation and pressure will be detected during the calibration runs and can be taken into account) C = F res * 2 * L(T, P)

Proof of concept

As a proof of concept the following test setup has been devised (cf. figure 6):

• An aluminium, closed box equipped with a Kistler pressure transducer hooked up to a Kistler charge amplifier (CA).

• A piezo vibrator glued to the side of the box and hooked up to a signal generator.

Both the charge amplifier and the signal generator are connected to and controlled by the main pilot (see section 2.2) through a USB port (this allows us to adjust the range and sensitivity of the system as well as the signal characteristics if need be during operations). The Charge amplifiers output is read by the DAQ and transferred to the main pilot. This set-up does not allow for a precise synchronisation of the recording with the signal generator. This is not a problem if we use a double sweep. The signal generator is programmed to do a double sweep going from 1000 Hz to 5000 Hz and back to 1000 Hz in 20 seconds (10 up, 10 down). The temperature of the air is 25°C so we should measure a speed of sound of 2662 Hz.

Figure 6c show the results of the test. We note the nice symmetries of the resonance peaks and the fact that the recording is not synchronised (the signal is not symmetrical, we miss the tail of the sweep because the recording starts before the sweep. The double sweep avoids the need for a perfect synchronisation of the sub systems).

The analysing software finds the biggest 2 peaks and determines the center of the sweep (the middle between these peaks). This is the gray line in the figure. Then we find the lesser peaks and since we now know the time difference D between the center and the peaks, the start frequency F s and the center frequency F c of the sweep and the half sweep duration S d we can calculate the peak frequencies as:

F1 = F c - D(F c -F s ) S d F2 = F c + D(F c -F s ) S d
We look for peaks with decreasing amplitudes until we find one that is close to our initial guess of the speed of sound (given the graph this does not have to be very precise)

The red lines in the graph represent the values found and we note that the algorithm detects the extremes very nicely. The blue lines in the graph that correspond to the theoretical value of 2662 Hz are invisible since they are nearly fully covered by the red lines (2676.8 Hz). We're off by less than 0.6%! If need be we can adapt the sweep characteristics to be closer to the resonance frequency, this should further increase the acuracy.

Actual experimental set-up

In the actual set-up the walls are much to thick to directly transfer the movement of the vibrator and to induce a P at the surface of the cavity. Since the prime interest of the vessel walls is to withstand the pressure inside the cavity, a very limited number of options are feasible. We've opted for a tunnel drilled into the wall that stops 5 mm from the inner surface. A free floating cylinder acts as a plunger to transmit the movement of the vibrator to the surface of the vessel. (see 6d to get the idea). Sadly this does not work, we tried different vibrators but there is no resonance signal measurable. The most plausible explanation is that the surface of the inner wall that transmit the vibrations is just not big enough to provoke a resonance phenomena. We will therefore refrain from the measurement of the speed of sound in the gas phase. The Kistler pressure transducer was supposed to be mounted on the top port that allowed for the thermocouples used during the calibration runs (see section 3.2.1).

Calibration

Pressure transducers and thermocouples

The pressure sensors are calibrated using a calibration balance present at ONERA. Before every test run the sensor is zeroed against the atmospheric pressure (given to be precise within a range of ± 10 Pa).

During these calibration tests 2 thermocouples are placed inside the cavity (on in the liquid and one in the gaseous phase, see section 3.2.1). These thermocouples permit to quantify the temperature difference between the embedded thermocouples and the actual temperature in the cavity as well as the eventual temperature in-homogeneity between the liquid and the gaseous phase.

The reference temperature is obtained through the pressure measurements. For pure water, injected under vacuum conditions, the relation between pressure and temperature is known perfectly well. A typical test run is given in figure 7a. At regular intervals (150°C,200°C,250°C ...) the system is left to go to equilibrium and a calibration point is obtained. 7b shows the perfectly linear behaviour of the thermocouples and table 6 gives the corresponding coefficients of the best fit for this calibration run. Figures 7c and7d show the results of a verification run using the calibration coefficients found earlier. We see from figure 7c that the measured temperatures join the target temperature of 250 °C around 13:26. Figure 7d show that this corresponds to time when the pressure derivative becomes zero. The same figure show that the derivatives of the temperature measurements inside the cavity join the pressure derivative at 0 around From figure 7c we extract the different temperatures (only the rear thermocouples were present during the test since the front thermocouples needed to be rerouted to make place for a camera). The values are given in table 7 and show a maximum difference of 0.77 °C.

Liquid phase density measurements

Once the temperature calibrations done, the port that was used to place two thermocouples inside the cavity was welded shut. The floater and balance system (see 4.2 for details) were installed and several test runs with pure water have been carried out. The accuracy is really nice for such a rudimentary system. Especially considering that only one point of the curve is fitted.

Conclusion

The work presented in this paper is not finished. During the two years of activities a lot of things went wrong or not as expected and Murphy invited himself to the party a lot more than we had hoped. This is why we are not as far as we had initially estimated and this paper does not contain actual results of the real test fluids as we had hoped to be able to deliver at least partially.

The lessons learned are numerous and we continue to find ways to improve or circumvent the obstacles we constantly find on out way. At the moment I am writing this, the installation is again in hold since we blew a window at a pressure far below the test pressure of 385 bars. We are currently looking for an explanation and are replacing the broken parts.

There is only one more calibration run needed to calibrate the ultra sound measurements before we can start to tackle the real products, beginning with NTO.

Once the real propellant test are completed, we are looking forward to writing the sequel to this adventure. It is more and more clear to us why this has not been done since 1896 but we remain optimistic (most of the time) and persevere.

Figure 2 :

 2 Figure 2: Remote control system

Figure 1 :

 1 Figure 1: Hardware overview

  Figure 3: Principle

  Figure 4: Principle

Figure 5 :

 5 Figure 5: Principle

Figure 6 :

 6 Figure 6: test setup

Figure 7 :

 7 Figure 7: Temperature verification run

Figure 8 :

 8 Figure 8: Calibration density measurements

Table 1

 1 shows the dependency of the Young modulus on temperature as determined by Julien et al.[2] which allows to coherently determine a = 109.954 and b = -0.0698

	E (GPa) T (°C)
	108	28
	30	1031

Table 1 :

 1 Young's modulus as a function of temperature Using the dimensions in table 2 we find that:

	ρ(T ) =	0.042585 -1.11932.01 -10 δ(T )E(T )β T T 9.81 -6 β 3	(4)
		W	0.01 m
		H	0.0002 m
		d	0.01 m
		L	0.0556 m
		l	0.0044 m
		m f 0.004341 kg

Table 2 :

 2 dimensions of the system

table T

 T

		η	ζ	∆ρ (kg/m 3 )
	0	0.016256250 8.5399786e-17 0.1275
	100 0.016168742 7.4364279e-17 0.12715637
	200 0.016081783 6.4161525e-17 0.12681397
	300 0.015995369 5.4775107e-17 0.12647280
	400 0.015909497 4.6188863e-17 0.12613285

Table 3 :

 3 error estimations

  with V g ≈ 2500 mm 3 , ∆V ≈ 5 mm 3 we get the error estimation given in table 4. It should be noted that this does not take into account the error in the liquid density (see section 3.5.2)

	Product	ρ l	Error
	NTO	550 7.2728373 -4
		1500 2.6669666 -4
	MMH	290 1.3793161 -3
		950 4.2107163 -4
	Hyd	230 1.7391350 -3
		1030 3.8837011 -4
	Table 4: gas density error estimation
	3.5 Liquid speed-of-sound measurements
	3.5.1 Principle
			l1 l2	l3	l4	l5
			Emitter	Isolation	Receiver

Table 5 :

 5 3.10 -5 ]. speed of sound expectations

	Product	c min	c max
	NTO	840 1100
	MMH	1200 1700
	Hydrazine 1800 2100

Table 6 :

 6 TC transfer function coefficientsFigure

	TC Id	a	b
	Bottom rear -0.0021 0.9880
	Top rear	0.034	0.9817
	Bottom front 0.0028 1.0087
	Top front	-0.0013 0.9929

Table 7 :

 7 temperature readings at 13:33

	TC id	Measured T (°C)
	Bottom rear	246.40
	Top rear	246.29
	Gas	246.60
	Water	247.03

Table

8

shows the results as well the theoretical density ρ th (P, T ). So if we use T, P, Px and ρ th and equation 6 and fit the data using equation 5 (c = 7.366.10 -6 and Z = -468.14) we find the results in figure

8

.

  THERMODYNAMIC PROPERTIES OF HYDRAZINE, MMH AND NTOT (°C) P (bar) Px (-) ρ th (kg/m 3 )

				318.0	110.5	968	672.0
				294.0	79.4	950	724.2
				269.6	55.4	926	768.1
				245.13	37.2	916	806.0
				196.11 15.126	896	869.2
				147.33	5.346	876	919.5
				13.77	0.972	836	999.2
				Table 8: Density measurement calibration
					Displacement
		0.0036
		0.0034
	δ(m)	0.0032	Theory Measurement
		0.0030
		0.0028
			c=7.36588e-06, Z=468.14
			700	750	800	850	900	950	1000
						ρ(kg/m 3 )

Ta6V because of the material compatibility constraints

http://zeromq.org/

https://www.acquisys.fr/product/usb-1608g/

https://www.acquisys.fr/product/usb-1608g/

https://github.com/mccdaq

https://www.gemeasurement.com/sites/gemc.dev/files/unik_5000_datasheet_english.pdf

https://en.ids-imaging.com/store/ui-5480se.html

DOI: 10.13009/EUCASS2019-297

http://www.epsilon-ndt.com/upload/file/problar-ve-aksesuarlar-.pdf


Jean Brion

E Brion

Floquet thermalisation in a Rydberg-blockaded atomic chain subject to a frequency-modulated drive

Keywords: Rydberg atoms, Rydberg blockade, Floquet thermalisation in closed quantum systems

may come

Floquet thermalisation in a Rydberg-blockaded atomic chain subject to a frequency-modulated drive

Introduction

For the last two decades, Rydberg atoms [START_REF] Gallagher | Rydberg atoms[END_REF] have been the subject of an intense research activity, mainly due to their potential use in quantum technologies [START_REF] Saman | [END_REF]3]. This interest is motivated by the so-called Rydberg blockade which results from the strong dipole-dipole interactions that two atoms exhibit when prepared in a highly excited level [4,5]. Since its discovery, Rydberg blockade was put forward as a powerful way to create entanglement among information-carrying atoms [6,7] and therefore implement quantum computation [5,8,9], communication [1012] and simulation protocols [13,14]. The highly non-linear optical response of Rydberg-excited atomic media and its applications were also addressed experimentally and theoretically, either in free-space or in cavity qed setups [1526].

Rydberg atoms furthermore appeared as versatile systems to investigate fundamental paradigmatic problems, such as quantum dynamics after a quench. Atomic ensembles suddenly submitted to a resonant laser excitation towards a Rydberg level were indeed numerically shown to thermalise [27,28]: the expectation values of certain observables specically single-particle ones such as the populations exhibit rapidly damped oscillations before entering a thermal regime, in which they slightly uctuate around an average thermal-like value [2936]. This eect appears in simulations though the model does not include any coupling to the environment: Rydberg levels being long-lived [START_REF] Gallagher | Rydberg atoms[END_REF], their spontaneous emission may indeed be safely neglected on the typical timescales of experiments and the atomic system may therefore satisfactorily be considered closed. Here, equilibration is fundamentally due to the dephasing process resulting from the combined action of single-particle laser excitation and long-range two-particle dipole-dipole interactions. To be more explicit, the (typically factorized) state of the ensemble just before the quench decomposes on a large number of eigenstates of the quenched Hamiltonian: if the associated eigenenergies are distributed in a suciently wide range, the dierent eigencomponents in the system state vector dephase, hence the observed eective thermalisation. According to the Eigenstate thermalisation Hypothesis -ETH [35], the reached thermalised state can be described by a canonical thermal ensemble.

Periodically driven (so-called Floquet) systems, may also exhibit stroboscopic thermalisation. In other words, the discrete-time expectation values of observables measured every period may converge towards a thermalised value, thus implying a synchronization of the system with the drive [37]. Floquet thermalisation was investigated in various systems and models, including interacting fermions [38], driven Hubbard model [39], and the O (N ) model [40]. Other behaviours may arise in Floquet systems, including the emergence of Floquet time crystals [41], which exhibit macroscopic oscillations at integer multiples of the driving period. Floquet stroboscopic evolution can also be specically engineered in order to emulate exotic phases of matter [42,43] or light [44]. Recently, a protocol for a Rydberg-based Floquet time crystal has been theoretically suggested [45] and analyzed and the implementation of interacting Floquet Rydberg cavity polaritons has been demonstrated [46].

The present work focuses on Floquet thermalisation in a chain of atoms excited towards a Rydberg level by a frequency-modulated laser beam. Here, we simulate the stroboscopic dynamics of the population in the Rydberg state, analyze the conditions for its equilibration and numerically investigate the inuence of the dierent parameters of the system. We also present a potential implementation of our scheme in a system of atoms trapped in optical tweezers.

The article is structured as follows. In Sec. II, we briey recall basic notions of Floquet systems and thermalisation, and introduce useful notations. In Sec. III, we present our model, provide numerical simulations and interpret our results. In Sec. IV we conclude and give perspectives of our work.

2. Stroboscopic thermalisation of a closed system submitted to a periodic drive In this section, we introduce the general notions and notations which will be used in the remainder of the article. To this end, we consider a generic many-body closed system of natural Hamiltonian H 0 , initially prepared in the state |Ψ 0 . At time t = 0, a T -periodic time-dependent drive V (t) is switched on: the full Hamiltonian of the system H = H 0 + V is therefore itself T -periodic, i.e.

H (t + T ) = H (t)

After one period, the system is in the state

|Ψ 1 ≡ U F |Ψ 0 where U F ≡ T exp -i ¡ T 0 dt H (t)
is the evolution operator over one period and T is the chronological product. As a unitary operator, U F can be put under the diagonal form U F = α exp (-iω α T ) |α α| where exp (-iω α T ) and |α denote the eigenvalues and eigenstates of U F , respectively, with ω α T uniquely determined in the interval [-π, π]. The Floquet Hamiltonian H F can then be uniquely and unambiguously dened by H F ≡ α ω α |α α|. The evolutions induced by H and H F at an arbitrary time t are a priori dierent but they do coincide at times equal to multiples of the period T , i.e. at t n = nT where n ∈ N.

For an arbitrary observable O of the system, expectation values at times nT (for n ∈ N), denoted by O n ≡ O (nT ), therefore simply write

O n = Ψ 0 U † F n O (U F ) n Ψ 0 = αβ e inT (ωα-ω β ) c * α c β O αβ (1)
where O αβ ≡ α |O| β and c α ≡ α| Ψ 0 is the projection onto the eigenvector |α of the initial state of the system just before switching the drive on. [36].

After [36], we introduce the frequency gaps

G αβ ≡ (ω α -ω β ) (2)
their respective amplitudes

v O,αβ ≡ c * α c β O αβ ∆ O (3)
where ∆ O is the dierence between the highest and lowest eigenvalues of O, and their respective relevances

q O,αβ ≡ |v O,αβ | 2 ab |v O,ab | 2 (4)
When the gaps G αβ which contribute to the sum in equation ( 1

σ2 O ≡ σ 2 O (n → +∞) = α =β (c * α c β O αβ ) 2 (6)
In the next section, we apply the formalism introduced above to a linear chain of atoms submitted to a laser eld whose frequency is periodically modulated close to resonance with a transition towards a Rydberg state.

3. Floquet thermalisation in a Rydberg-blockaded atomic chain

Model and assumptions

In this section, we consider a system comprising N atoms with a two-level structure involving a ground and a Rydberg states, respectively denoted by |g and |r (see gure 1 a). The atoms are regularly arranged along a line with the step a. They are submitted to a laser beam which near-resonantly drives the transition |g ↔ |r with the Rabi frequency Ω and the time-dependent T -periodic detuning ∆ (t). Atoms excited in the Rydberg state moreover exhibit strong dipoledipole interactions leading to the so-called Rydberg blockade [25]. In the usual rotating wave approximation, the full time-dependent Hamiltonian takes the form H (t) = H 0 + V (t) with

H 0 ≡ C 6 a 6 m =n σ (m) rr σ (n) rr |m -n| 6 (7)
V (t) ≡ n ∆ (t) σ (n) rr + Ω 2 σ (n) rg + σ (n) gr (8)
where σ

(m) kl ≡ I ⊗ • • • ⊗ |k m l| ⊗ • • • ⊗ I is the operator σ kl ≡ |k
l| acting on the m th atom in the chain, and C 6 denotes the coecient of the van der Waals interaction between two atoms in the state |r .

Emission from the Rydberg state can usually be neglected provided experimental runs are short enough. For instance, for n > 50, the s, p, d, and f states of rubidium have typical lifetimes (including radiative and blackbody eects) greater than about 50µs at room temperature while experiment sequences can be achieved on the MHz-scale [START_REF] Saman | [END_REF]. We therefore shall not consider spontaneous emission in our model.

Numerical results

In this section we present the numerical results we obtained for a system of N = 10 atoms submitted to a frequency-modulated laser beam. We set ΩT = 2π and consider the following two types of modulated detunings (see gure 1 b)

∆ (t) = ∆ 0 × 1 -4 t T -n , for n ≤ t T ≤ n + 1 2 4 t T -n -3, for n + 1 2 ≤ t T ≤ n + 1 (9)
and

∆ (t) = ∆ 0 cos 2πt T (10)
For a given set of parameters {∆ 0 T, C6T / a 6 } we numerically compute the Floquet evolution U F by splitting the period [0, T ] into K short time steps τ ≡ T K and resorting to the approximate discretized form

U F ≈ e -iτ H(Kτ ) × • • • × e -iτ H(2τ ) × e -iτ H(τ ) (11)
For the parameters we use in our simulations, we numerically check that K = 100 steps are enough to obtain a converged value of U F ; in other words, increasing K beyond 100 does not signicantly change U F nor the results obtained, e.g., for the dynamics of populations. The eigenenergies ω α , eigenstates |α and the Floquet Hamiltonian H F are derived as described in section 2. Assuming the system is initially in the factorized ground state |Ψ 0 ≡ |g • • • g , we can straightforwardly compute the stroboscopic evolution of the expectation value of the observable of interest, i.e. the population of the Rydberg level in the ensemble We can also derive the respective relevances q αβ of the gaps G αβ for the observable of interest, as well as the thermalised value Pr and associated standard deviation σr . Before presenting and analysing our results, let us note that our choice for |Ψ 0 ≡ |g • • • g corresponds to the typical initial state in Rydberg experiments. This conguration was considered in the rst articles on thermalisation in Rydberg atomic systems (see, e.g., [27,28]). Other choices for |Ψ 0 = α c α |α should a priori lead to dierent Floquet thermalised states ρ = α |c α | 2 |α α|. It was, however, shown that a nonintegrable periodically driven model generically leads to a Floquet thermalised state independent of the initial condition [47]. By contrast to this generic scenario, recent works pointed out the existence of quantum scars in certain Rydberg atomic systems, i.e. eigenstates leading to persistent coherent oscillatory dynamics [48]. According to the preliminary numerical checks we performed on systems with moderate numbers of atoms, starting from randomly picked initial states |Ψ 0 (with one or several Rydberg excitations), our system seems to follow the generic case. But, in the absence of a systematic joint numerical and theoretical study, we cannot claim quantum scars do not exist in our system. This will be the subject of future work.

P r = N n=1 σ (n) rr . We set P r,n ≡ Ψ 0 U † F n P r (U F ) n Ψ 0 .
We start with the trivial case C6T a 6 = ∆ 0 T = 0 which corresponds to non-interacting atoms undergoing independent and identical Rabi oscillations at the frequency 2π /Ω = T induced by the resonant laser beam. Figure 2 displays the continuous-time evolution of the population P r (t) over the rst period of the time-dependent detuning (though in this case, the amplitude of modulation is zero): one full oscillation can be clearly identied. Measuring the population P r (t) at the beginning of each period therefore yields zero and accordingly P r,n remains constantly zero.

By contrast, gure 3 shows P r,n over 800 periods for ∆ 0 T = 5 and C6T / a 6 = 8 in the case of a triangular T -periodic time-dependent detuning ∆ (t) dened in equation ( 9). Now, one observes a quick equilibration of P r,n around the average value Pr ≈ 4.62 with reduced long-term uctuations σr ≈ 0.10. With this choice of parameters, the simultaneous excitation of two neighbouring atoms is highly prevented, i.e. one imposes nearest-neighbour Rydberg blockade, as conrmed by the average number of excited along the chain, Pr ≈ N /2. The dynamical behaviour of P r,n may be qualitatively accounted for through analyzing the distributions of energy gaps and associated relevances of the system, represented in gures 4, 5 and 6. We rst make two remarks: A. For sake of readibility of gures 4, 5 and 6, we chose not to display : i) degenerate state pairs (|α , |α ) of zero energy gaps which play no role in the time-dependent part of equation (1) and ii) the state pairs (|α , |β ) of negligible relevances q r,αβ , i.e. for which q r,αβ ≤ 10 -6 . This is why, e.g., only ≈ 3500 state pairs appear in gure 4, while the actual number of state pairs in the system is 2 20 ≈ 10 6 . B. The relevance distribution is symmetric with respect to the zero energy gap. This results from the reality of the expectation value in equation [START_REF] Gallagher | Rydberg atoms[END_REF].

We now turn to the interpretation of gures 4, 5 and 6. The relevance distribution on gure 4 can be roughly split into a (dominating) low-frequency and a high-frequency parts, respectively displayed on gures 5 (a) and 6 (a). The boundary between low-and high-frequency contributions was xed at |G αβ | ≈ 0.1 × Ω which, in order of magnitude, corresponds to 0.1 × T -1 for the parameters used in our simulations. Following the (oversimplistic but intuitive) picture of populations performing slowly damped Rabi oscillations, the separation between low and high frequencies thus dened therefore allows us to isolate the slow damping from rapid Rabi oscillations. The associated reduced low-frequency and high-frequency dynamics are represented on gures 5 (b) and 6 (b), respectively, and compared to the exact full dynamics on gures 5 (c) and 6 (c), respectively. Note that the zero-frequency component (associated to degenerate state pairs (|α , |α )) was implicitly included in the calculation of both low-frequency and high-frequency reduced dynamics to make the comparison with the full dynamics easier.

The low-frequency part is dominated by two symmetric peaks denoted by A on gure 5 (a), of The high-frequency part is dominated by a pair of symmetric peaks, denoted by A on gure 6 (a), of gaps ±G A ≈ ±0.20 × T -1 associated to relatively slow oscillations of period 2π /G A ≈ 31 × T which can be qualitatively identied on gures 6 (b) and (c). Because of its complexity, the rest of the distribution does not lend itself to direct interpretation : it comprises many high-frequency peaks which induce fast periodic dynamics of respective periods ≈ 1, 2, 6, 10, 20, • • • × T (among others) and result in the noisy uctuations observed on gures 6 (b) and (c).

Floquet thermalisation is also accompanied by a synchronization phenomenon: in other words, the continuous-time function P r (t) tends to synchronize with the drive and becomes itself Tperiodic. This can be seen on gure 7 which shows superimposed plots of the function P r (t) on the rst (top) and last (bottom) ten periods of the range considered in gure 3: while the behaviour of P r (t) is quite dierent on the rst ten periods, the traces of P r (t) on the last ten periods almost coincide.

Figures 8 and9 analyze the dependence on the parameters ∆ 0 T and C6T / a 6 , respectively, of the quantities Pr , σr and σr / Pr. The latter is particularly signicant since it quanties how eective Floquet thermalisation is. As shown on gures 8 and 9, thermalisation is robust with , and associated relevances q r,αβ are by a full (black) and a dashed (red) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. (b) The reduced stroboscopic dynamics of the Rydberg population P red r,n induced by the low-frequency gaps only including zero-frequency components Gαα not displayed on gure (a), is represented by a full (blue) line as a function of the number of periods n, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dasheddotted (black) lines. (c) The exact and reduced stroboscopic dynamics of the Rydberg population, Pr,n and P red r,n , are represented by a full (blue) and a dashed (red) lines, respectively, as functions of the number of periods n. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5. respect to variations in the system's parameters: in the range of ∆ 0 T values investigated around ∆ 0 T = 5, σr / Pr ≤ 0.022, while, in the range of C6T / a 6 values investigated around C6T / a 6 = 8, σr / Pr ≤ 0.022.

In the same spirit, one can investigate the robustness of the phenomenon against a slight distortion of the modulation scheme. Our ambition here is not to provide a quantitative and systematic analysis of the inuence of the frequency modulation shape on, e.g., the spectral properties of H F and/or the population dynamics, this is much beyond the scope of the present work. More modestly, on the specic example presented below, we merely want to point out that the observed thermalisation phenomenon seems to remain qualitatively unchanged when slightly modifying the modulation scheme.

Figures 10111213141516address the case of a T -periodic cosine time-dependent detuning, dened in equation (10) with the same parameters ΩT = 2π, C6T / a 6 = 8 and ∆ 0 T = 5 as above. The results obtained are very much alike previous ones. Equilibration of P r,n is observed around its average Pr ≈ 4.74 with reduced long-time uctuations σr ≈ 0.09 (gure 10), as well as synchronization of P r (t) (gure 14).

As in the triangle-frequency-modulation case, the dynamical behaviour of P r,n is qualitatively accounted for through analyzing the distributions of energy gaps and associated relevances of the system, represented in gure 11. Again, the relevance distribution on gure 11 can be roughly split into a (dominating) low-frequency and a high-frequency parts, respectively displayed on gures 12 (a) and 13 (a). The associated reduced low-frequency and high-frequency dynamics are represented on gures 12 (b) and 13 (b), respectively, and compared to the exact full dynamics on gures 12 (c) and 13 (c), respectively.

The low-frequency part is again dominated by two symmetric peaks denoted by A on gure The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

2 000 4 000 1 000

3 000 500 1 500 2 500 3 500 Finally gures 15 and 16 show the strong robustness of Floquet thermalisation with respect to changes in the system's parameters ∆ 0 T and C6T / a 6 : in the range of ∆ 0 T investigated around ∆ 0 T = 5, σr / Pr ≤ 0.020, while, in the range of C6T / a 6 values investigated around C6T / a 6 = 8, σr / Pr ≤ 0.021.

Experimental implementation

To close this section we suggest a possible experimental implementation to investigate the Floquet thermalisation identied above. We consider the same conguration as in [45], i.e. rubidium atoms 87 Rb with a ground and Rydberg states, respectively |g = 5s1 /2 ; F = 1; m F = 1 and |r = |60s; j = 1 /2; m j = 1 /2 , coupled through a two-photon process via the non-resonant intermediate state |e = 5p1 /2 ; F = 1; m F = 0 . Adiabatic elimination of |e [49] yields an eective model of interacting two-level atoms as considered above. The van der Waals coecient of the Rydberg state is C6 = 2π × 1.4 × 10 11 Hz • µm 6 and its spontaneous decay rate is γ r = 2π × 2.0 kHz.

Setting T ≈ 80ns, ∆ 0 ≈ 2π × 10MHz, Ω = 2π × 12.5MHz and a ≈ 4.5µm which are feasible parameters in current experiments with atoms trapped in optical tweezers [13] we get C6T a 6 ≈ 8 and ∆ 0 T = 5 which allow to observe the thermalisation identied above.

Note that with our choice of parameters, γ r T ≈ 10 -3 . For approximately the rst hundred periods, spontaneous emission will therefore have no noticeable inuence on the dynamics of the system. This is sucient for the system to reach the thermal regime and therefore exhibit Floquet thermalisation and synchronization in both triangle and cosine frequency modulation schemes studied here. The analysis of the impact of spontaneous emission is, however, important since it may compete with the phenomena investigated here and lead to qualitatively dierent phenomena in the long-time limit. It will be addressed in a future work.

Conclusion

In this article we investigated the emergence of Floquet thermalisation and synchronization in a chain of atoms excited by a frequency-modulated laser beam on a transition towards a Rydberg state. We related this phenomenon to properties of the Floquet Hamiltonian spectrum, i.e. the distribution of energy gaps contributing to the Rydberg population dynamics. We found this process is robust against small variations of the parameters, such as the detuning modulation amplitude or the distance between atoms. Finally we suggested a possible experimental situation in which Rydberg-blockade induced Floquet thermalisation could be investigated. , and associated relevances q r,αβ are represented by a full (black) and a dashed (red) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. (b) The reduced stroboscopic dynamics of the Rydberg population P red r,n induced by the high-frequency gaps only also including zero-frequency components Gαα not displayed on gure (a), is represented by a full (blue) line as a function of the number of periods n, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dashed-dotted (black) lines. (c) The exact and reduced stroboscopic dynamics of the Rydberg population, Pr,n and P red r,n , are represented by a full (blue) and a dashed (red) lines, respectively, as functions of the number of periods n. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5. The impact of unavoidable spontaneous emission from the Rydberg level remains to be systematically investigated. The inuence of the spatial arrangement of the atoms on the phenomenon is also worth studying, in particular the choice of other geometries than the line or the introduction of disorder in the atomic positions. Richer level congurations involving several Rydberg states could also be considered, such as, e.g., the excitation of a Förster resonance. These aspects shall be considered in future works.

Figure 1 .

 1 Figure 1. (a) The model system : N atoms with two levels {|g , |r } are regularly arranged on a line with the step a, and submitted to a near-resonant laser eld of Rabi frequency Ω and time-dependent detuning ∆ (t). (b) Typical shapes used for ∆ (t) in simulations (see equations 9,10) in the text).

Figure 2 .

 2 Figure 2. Trivial non-interacting and non-modulated case The Rydberg population Pr is represented as a function of time on the rst period of the drive. The parameters used in this simulation are ΩT = 2π and C 6 T / a 6 = ∆ 0 T = 0.

Figure 3 .

 3 Figure 3. Triangular frequency modulation : Floquet thermalisation of the Rydberg state population The Rydberg population Pr,n is represented by a full (blue) line, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dashed-dotted (black) lines. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

gaps

  ±G A ≈ ±0.02 × T -1 , surrounded by a band of gaps of width ∆G ≈ 0.1 × T -1 . The associated dynamics consists of very slow damped oscillations of period 2π /G A ≈ 292 × T and typical damping timescale τ d = 2π /∆G ≈ 65 × T . The two pairs of secondary peaks outside the main band, denoted by B and C on 5 (a), induce surviving oscillations of typical period is given by 4π /(G B -G C ) ≈ 75 × T . This simple analysis is in good qualitative agreement with the dynamics observed on gures 5 (b) and (c).

Figure 5 .

 5 Figure 5. Triangular frequency modulation : low-frequency dynamics (a) Lowfrequency energy gaps G αβ , i.e. for which G αβ ≤ 0.1 × T -1, and associated relevances q r,αβ are by a full (black) and a dashed (red) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. (b) The reduced stroboscopic dynamics of the Rydberg population P red r,n induced by the low-frequency gaps only including zero-frequency components Gαα not displayed on gure (a), is represented by a full (blue) line as a function of the number of periods n, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dasheddotted (black) lines. (c) The exact and reduced stroboscopic dynamics of the Rydberg population, Pr,n and P red r,n , are represented by a full (blue) and a dashed (red) lines, respectively, as functions of the number of periods n. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

Figure 6 .

 6 Figure 6. Triangular frequency modulation : high-frequency dynamics (a) Highfrequency energy gaps G αβ , i.e. for which G αβ ≥ 0.1 × T -1, and associated relevances q r,αβ are represented by a full (black) and a dashed (red) lines, respectively. For sake of clarity, the x-axis does not show state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. (b) The reduced stroboscopic dynamics of the Rydberg population P red r,n induced by the high-frequency gaps only also including zero-frequency components Gαα not displayed on gure (a), is represented by a full (blue) line as a function of the number of periods n, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dashed-dotted (black) lines. (c) The exact and reduced stroboscopic dynamics of the Rydberg population, Pr,n and P red r,n , are represented by a full (blue) and a dashed (red) lines, respectively, as functions of the number of periods n. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

Figure 7 .

 7 Figure 7. Triangular frequency modulation : Synchronization with the drive The plots of the continuous-time function Pr over the rst (top) and last (bottom) ten periods are superimposed. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

Figure 8 .

 8 Figure 8. Triangular frequency modulation : inuence of parameter ∆ 0 T on the thermalisation phenomenon The thermalised Rydberg population Pr, standard deviation σr and ratio σr/ Pr are represented by (blue) stars, (red) triangles and (black) circles, respectively, as functions of ∆ 0 T . The parameters used in this simulation are ΩT = 2π and C 6 T / a 6 = 8.

8

 8

Figure 9 .Figure 10 .

 910 Figure 9. Triangular frequency modulation : inuence of the parameter C 6 T / a 6 on the thermalisation phenomenon The thermalised Rydberg population Pr, standard deviation σr and ratio σr/ Pr are represented by (blue) stars, (red) triangles and (black) circles, respectively, as functions of C 6 T / a 6 . The parameters used in this simulation are ΩT = 2π and ∆ 0 T = 5.

  12 (a), of gaps ±G A ≈ ±0.039 × T -1 , surrounded by a band of gaps of width ∆G ≈ 0.094 × T -1 . The associated dynamics consists of very slow damped oscillations of period 2π /G A ≈ 160 × T and typical damping timescale τ d = 2π /∆G ≈ 70 × T . The pair of secondary peaks outside the main band, denoted by B on 12 (a) (and peaks close to B), of gaps ±G B ≈ ±0.06 × T -1 induce surviving oscillations of typical period by 2π /G B ≈ 105 × T . Here again our simple analysis reproduces the essential features observed on gures 12 (b) and (c).The high-frequency part is dominated by a pair of symmetric peaks, denoted by A on gure 13 (a), of gaps ±G A ≈ ±0.135 × T -1 associated to relatively slow oscillations of period 2π /G A ≈ 46 × T which can be qualitatively identied on gures 13 (b) and (c). The presence of many secondary peaks makes, however, further analysis illusive: the noisy uctuations observed on gures 13 (b) and (c) results from the superposition of oscillations with periods ≈ 1, 2, 6, 10, 13, 17, 23, • • • × T (among others).

Figure 13 .

 13 Figure 13. Cosine frequency modulation : high-frequency dynamics (a) Highfrequency energy gaps G αβ , i.e. for which G αβ ≥ 0.1 × T -1 , and associated relevances q r,αβ are represented by a full (black) and a dashed (red) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. (b) The reduced stroboscopic dynamics of the Rydberg population P red r,n induced by the high-frequency gaps only also including zero-frequency components Gαα not displayed on gure (a), is represented by a full (blue) line as a function of the number of periods n, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dashed-dotted (black) lines. (c) The exact and reduced stroboscopic dynamics of the Rydberg population, Pr,n and P red r,n , are represented by a full (blue) and a dashed (red) lines, respectively, as functions of the number of periods n. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

Figure 14 .

 14 Figure 14. Cosine frequency modulation : Synchronization with the drive The plots of the continuous-time function Pr over the rst (top) and last (bottom) ten periods are superimposed. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

Figure

  Figure Cosine frequency modulation : inuence of the parameter ∆ 0 T on the thermalisation phenomenon The thermalised Rydberg population Pr, standard deviation σr and ratio σr/ Pr are represented by (blue) stars, (red) triangles and (black) circles, respectively, as functions of ∆ 0 T . The parameters used in this simulation are ΩT = 2π and C 6 T / a 6 = 8.

8

 8

Figure 16 .

 16 Figure 16. Cosine frequency modulation : inuence of the parameter C 6 T / a 6 on the thermalisation phenomenon The thermalised Rydberg population Pr, standard deviation σr and ratio σr/ Pr are represented by (blue) stars, (red) triangles and (black) circles, respectively, as functions of C 6 T / a 6 . The parameters used in this simulation are ΩT = 2π and ∆ 0 T = 5.

  ), i.e. whose relevances q O,αβ do not vanish, distribute on a suciently large range ∆G and in a suciently homogeneous manner without

	degeneracy, the dierent spectral components in equation (1) dephase and interfere destructively.
	For t ≥ t ≈ 2π /∆G, i.e. for a number of periods n ≥ n = [ 2π /∆G×T ] where [•] denotes the integer
	part the expectation value O n equilibrates around the value
	Ō ≡	|c α | 2 O αα	(5)
	α
	which can be seen as the mean value Tr [ρO] calculated in the mixed state ρ ≡ α |c α | called the thermalised state. Fluctuations of O n around Ō can be characterized by their variance 2 |α α| also σ 2 O,n ≡ 1 n n j=0 O j -Ō 2 and its long-time limit

  Figure 4. Triangular frequency modulation Energy gaps G αβ and associated relevances

	6					0.0018
	5					0.0016
	4
						0.0014
	3
	2					0.0012
	1
						0.001
	0
	-1					0.0008
	-2					0.0006
	-3
						0.0004
	-4
	-5					0.0002
	-6
	500 500	1 000 1 000	1 500 1 500	2 000 2 000	2 500 2 500	3 000 3 000

q r,αβ are represented by a dashed (black) and a full (blue) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

  Figure 11. Cosine frequency modulation Energy gaps G αβ and associated relevances q r,αβ are represented by a dashed (black) and a full (blue) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and ∆ 0 T = 5.

								0.0014
	5
								0.0012
								0.001
								0.0008
								0.0006
								0.0004
								0.0002
	500	1 000	1 500	2 000	2 500	3 000	3 500	4 000

  Figure 12. Cosine frequency modulation : low-frequency dynamics (a) Low-frequency energy gaps G αβ , i.e. for which G αβ ≤ 0.1×T -1 , and associated relevances q r,αβ are represented by a full (black) and a dashed (red) lines, respectively. For sake of clarity, the x-axis does not show degenerate state pairs (|α , |α ) of zero energy gaps, nor state pairs (|α , |β ) for which q r,αβ≤ 10 -6 . The remaining pairs are arbitrarily numbered in such a way that G αβ 's appear in decreasing order. (b) The reduced stroboscopic dynamics of the Rydberg population P red The parameters used in this simulation are ΩT = 2π, C 6 T / a 6 = 8 and

	induced by the low-frequency gaps only including zero-frequency components	r,n

Gαα not displayed on gure (a), is represented by a full (blue) line as a function of the number of periods n, the thermalised value Pr by a dashed (red) line, the uctuation area is delimited by dashed-dotted (black) lines. (c) The exact and reduced stroboscopic dynamics of the Rydberg population, Pr,n and P red r,n , are represented by a full (blue) and a dashed (red) lines, respectively, as functions of the number of periods n. ∆ 0 T = 5.

Acknowledgments E. B. thanks F. Alet for helpful discussions.

* retired, previously at the Groupe de Spectrométrie Moléculaire et Atmosphérique, UMR CNRS 7331, University of Reims Champagne-Ardenne, F-51687, Reims Cedex 2, France.