First detection of collective oscillations of a stored deuteron beam with an amplitude close to the quantum limit

We investigated coherent betatron oscillations of a deuteron beam in the storage ring COSY, excited by a detuned radio-frequency Wien filter. The beam oscillations were detected by conventional beam position monitors. With the currently available apparatus, we show that oscillation amplitudes down to \SI{1}{\micro \meter} can be detected. The interpretation of the response of the stored beam to the detuned radio-frequency Wien filter is based on simulations of the beam evolution in the lattice of the ring and realistic time-dependent 3D field maps of the Wien filter. Future measurements of the electric dipole moment of protons will, however, require control of the relative position of counter-propagating beams in the sub-picometer range. Since here the stored beam can be considered as a rarefied gas of uncorrelated particles, we moreover demonstrate that the amplitudes of the zero-point (ground state) betatron oscillations of individual particles are only a factor of about 10 larger than the Heisenberg uncertainty limit. As a consequence of this, we conclude that quantum mechanics does not preclude the control of the beam centroids to sub-picometer accuracy. The smallest Lorentz force exerted on a single particle that we have been able to determine is \SI{10}{aN}.

We investigated coherent betatron oscillations of a deuteron beam in the storage ring COSY, excited by a detuned radio-frequency Wien filter. The beam oscillations were detected by conventional beam position monitors. With the currently available apparatus, we show that oscillation amplitudes down to 1 µm can be detected. The interpretation of the response of the stored beam to the detuned radio-frequency Wien filter is based on simulations of the beam evolution in the lattice of the ring and realistic time-dependent 3D field maps of the Wien filter. Future measurements of the electric dipole moment of protons will, however, require control of the relative position of counter-propagating beams in the sub-picometer range. Since here the stored beam can be considered as a rarefied gas of uncorrelated particles, we moreover demonstrate that the amplitudes of the zero-point (ground state) betatron oscillations of individual particles are only a factor of about 10 larger than the Heisenberg uncertainty limit. As a consequence of this, we conclude that quantum mechanics does not preclude the control of the beam centroids to sub-picometer accuracy. The smallest Lorentz force exerted on a single particle that we have been able to determine is 10 aN.

I. INTRODUCTION
The approach to the quantum ground state, the observation of quantum effects in macroscopic systems, and the possibility to detect displacements of macroscopic bodies on the nanometer scale, are the subject of intense theoretical and experimental efforts [1][2][3][4][5]. A notable example is the detection of gravitational waves using an interferometric detector with mirrors in the kilogram range [6]. In all-electric proton storage rings, coherent beam displacements down to the picometer range that are caused by Earth's gravity pull are in principle accessible using the spin rotations of the proton as a detector [7]. We also mention here the ongoing discussions of the possibility to detect gravitational waves via perturbations of the beam orbit in high-energy storage rings [8].
Here we report the first detection of collective oscillations of an intense beam of deuterons in a storage ring with an amplitude close to the quantum limit. The present study is part of an international effort to prepare for the search for the permanent electric dipole moment (EDM) of charged particles. The focus of these studies has been on systematic effects, e.g., imperfection magnetic fields in storage rings [9], and orbit improvements in a machine using beam-based alignment [10], thereby advancing the high-precision frontier in spin dynamics in storage rings. A comprehensive description of this activity and of the proposed stepwise approach leading to a dedicated proton EDM storage ring is presented in Ref. [7].
Experiments searching for electric dipole moments of charged particles using storage rings are at the forefront of the incessant quest to find new physics beyond the Standard Model of particle physics. These investigations bear the potential to shed light on the origin of the anomalously large matter-antimatter asymmetry in the Universe [11], for which the combined predictions of the Standard Models of particle physics and of cosmology fall short of the experimentally observed asymmetry by about seven to eight orders of magnitude [12].
The signal for an EDM is the spin precession in electric fields, where it should be noted that the spins of charged particles can be subjected to large electric fields only in storage rings. The need to eliminate the overwhelmingly stronger spin rotations driven by the magnetic moment in magnetic fields brings to the front an all-electric, so-called frozen spin proton storage ring [7,13]. An important advantage of such a machine is the ability to simultaneously store two counter-propagating proton beams. The concurrent measurement of the EDM-driven spin rotations of the counter-propagating beams would allow to cancel major systematic effects. To this end, to reach an ambitious sensitivity to the proton EDM of d p ≈ 10 −29 e cm, it is imperative to control the relative vertical displacement of the centers of gravity of the two beams to an accuracy of about 5 pm [7]. One may wonder whether such an enormously demanding accuracy is not prohibited by the Heisenberg uncertainty principle. Towards an ultimate precision search for EDMs of charged particles, this particular aspect of the systematics of such measurements had not been investigated so far, and our experiment constitutes the first step in this direction.
Here we report on the measurement of the amplitude of collectively excited vertical oscillations of a deuteron beam orbiting in the magnetic storage ring COSY at a momentum of about 970 MeV/c [14]. The data were taken in 2018 in the course of a dedicated experiment in the framework of systematic beam and spin dynamics studies for the deuteron EDM experiment (so-called precursor experiment), presently carried out by the JEDI collaboration at COSY [15][16][17]. One of the central devices in the precursor experiment is the radio-frequency (RF) Wien filter (WF), shown in Fig. 1, which was designed to provide a cancellation of the electric and magnetic forces acting on the particle. In this operation mode, the Wien filter affects only the particle spins, but does not perturb the beam orbit [18][19][20]. A slightly detuned Wien filter, however, exerts a non-vanishing Lorentz force on the orbiting beam particles. It is shown collective beam oscillation excited by the WF with amplitudes down to 1 µm can be detected with the currently available equipment. Our approach to measuring ultrasmall displacements complements other measurements of ultra-small forces using different techniques [1][2][3][4][5].
In our experiment, the measurement cycles were much shorter than the intrabeam interaction time, and the beam attenuation rate was negligibly weak (see discussion in [21]), therefore we treat the beam as a rarefied gas of uncorrelated particles. Individual particles undergo stable betatron oscillations around the equilibrium orbit in the horizontal and vertical planes, driven by focusing magnetic fields. Apart from their conventional individual betatron motions, all the particles in a bunch participate in one and the same collective and coherent oscillation that is driven by the Wien filter. Therefore, the upper bound of the amplitude of the collective oscillation of the entire beam corresponds to the upper bound of the oscillation amplitude of a single particle. In our approach, access to ultrasmall oscillation amplitudes results from the fact that the measured signal corresponds to a collective response of the electric charge of about N = 10 9 deuterons in the bunch.
The beam tracking simulations were carried out to predict the response of the stored beam to the detuned radiofrequency Wien filter. The simulations use the elements of the ring lattice and realistic time-dependent 3D field maps of the Wien filter. These field maps describe the spatial variation of complex electric and magnetic fields, including the fringe field areas. Furthermore, the tolerances of the elements of the circuit driving the Wien filter are also taken into consideration.
As a reference value for the Heisenberg uncertainty relation, we take an estimate of the amplitude of the single-particle zero-point betatron oscillation amplitude Q. Then, our result for the smallest measured amplitude of the Wien filter-driven single-particle oscillation is only about a factor of ten larger than the quantum limit of Heisenberg's uncertainty relation for vertical single-particle betatron oscillations. The smallest detected oscillation amplitude is by three orders of magnitude smaller than the beam size.
We demonstrate for the first time that the accuracy with which periodic beam oscillation amplitudes can be measured is vastly higher compared to that of static beam displacements generated by steerers using the same BPM. In a broad context, any new precision tool is of interest per se, and the latter point has important implications, for instance, for all-electric frozen-spin EDM storage rings. Here, one aims to control the interfering radial magnetic fields by measuring the vertical spacing of the counter-propagating beams. Our result complements the potential of using beam oscillations to measure this distance, as discussed in Ref. [22].
One must distinguish the RF-driven collective oscillations above the quantum limit Q from the quantum uncertainty of the center of mass of the bunch circulating in a static ring. Specifically, for a rarefied-gas of N uncorrelated particles, the quantum limit of the centroid of the bunch, detected by the BPMs, amounts to Q/ √ N .
The paper is organized as follows. In Sec. II, the measurement principle is introduced, followed by a description of the operation of the radio-frequency Wien filter in Sec. III. The method to determine the beam oscillations is discussed in Sec. IV, and the evolution of the beam to the combined effect of the ring lattice and the Wien filter fields are presented in Sec. V. The time-dependent field maps of the Wien filter are discussed in Sec.V A, and the evaluation of the uncertainties is elaborated in Sec. V B. Experimental results are presented in Sec. VI, followed by conclusion and outlook in Sec. VII.

II. MEASUREMENT PRINCIPLE
The Cooler Synchrotron (COSY) [14,23] at Forschungszentrum Jülich is a storage ring with a circumference of approximately 184 m. Its principal elements used for the experiments are indicated in Fig. 2. For the investigations presented here, the two key devices are the RF Wien filter, based on a parallel-plates waveguide [18], and a conventional electrostatic beam position monitor (BPM) that is used to monitor the beam oscillations [24]. The Wien filter generates orthogonal and highly-homogeneous electric and magnetic fields. In the present experiment, the Wien filter was operated in the mode with the electric field pointing vertically upward (y-direction), whereas the magnetic field points radially outward (x-direction), and the beam moves in z-direction (see coordinate system in Fig. 2). The effective length of the Wien filter is ℓ = 1.16 m (see Refs. [18,19] for further technical details).
As a spin rotator for the forthcoming deuteron EDM (precursor) experiment [15][16][17], the Wien filter is designed to operate in resonance with the spin precession of the orbiting deuterons [18][19][20]31], and at a vanishing Lorentz force, given by where q denotes the elementary charge and v represents the velocity of the beam particles. Unlike in conventional DC Wien filters, the crossed electric field and magnetic fields ( E and B) of the RF Wien filter are generated simultaneously by exciting the transverse electromagnetic (TEM) mode. The spin resonance tune mapping technique, developed for the Wien filter operation in the deuteron EDM experiment at COSY, is described in [17]. When the electric and magnetic fields in the Wien filter are mismatched, i.e., when the electric and magnetic forces no longer cancel each other, the RF fields excite collective beam oscillations at the frequency at which the Wien filter is operated. In the present experimental set up, a mismatch between electric and magnetic fields provides a vertically mismatched Lorentz force (see coordinate system in Fig. 2).
With a vanishing Lorentz force, the beam performs idle vertical (and horizontal) betatron oscillations where β y (t) is the betatron amplitude function. With the beam revolution period of T = 2π/ω rev , the betatron phase advance ψ y (t) satisfies ψ y (t + T ) − ψ y (t) = ω y T = 2πν y , where ν y is the vertical betatron tune given by ν y = ω y /ω rev . On the other hand, a mismatched Wien filter exerts stroboscopically, i.e., once per turn, a vertical force F y (n) = F y cos(n ω WF T ) on the stored particle, where n is the turn number and ω WF denotes the angular velocity of the RF in the Wien filter (see discussion of Fig. 12 in Sec. VI). The change of the vertical velocity of the stored particle, accumulated during the time interval ∆t = ℓ/v z the particle spends per turn n inside the Wien filter, is given by where γ and m are the Lorentz-factor and the mass of the particle, respectively. The change ∆y of the vertical position y in the Wien filter can be neglected. The coupling of vertical and radial beam oscillations is negligible (see Sec. IV) and it is sufficient here to treat driven oscillations in a one-dimensional approximation. Due to the very strong disparity of synchrotron and fractional WF frequencies, synchro-betatron coupling can be neglected (see discussion in [32]). Furthermore, beam attenuation either by intrabeam scattering or by interaction with residual gas during the data acquisition cycle is very small (see Appendix B), justifying the rarefied gas approximation. According to Eq. (2), the stroboscopic signal of the betatron motion observed at any point in the ring, follows the harmonic law with angular velocity ω y , and we invoke the familiar description of the oscillatory motion in terms of the complex variable z = y − iv y /ω y . With the initial condition z(0) = 0, summing ∆v y (kT ) after n turns, the solution for z(n) behind the Wien filter reads This expression serves as the initial condition for the idle betatron motion during the subsequent (n + 1) turn and so forth. A similar analytic result holds also for generic AC dipole-driven betatron oscillations, discussed in a very different context of machine diagnostics in Ref. [33] (see also references therein). Driven by the mismatched Wien filter, all beam particles participate in one and the same collective and coherent oscillation, and according to Eq. (4), the beam as a whole exhibits oscillations at the Wien filter frequency ω WF . A lock-in amplifier may be used to selectively measure the corresponding Fourier component of the beam oscillation y = ξ y cos(n ω WF T ) from the output of a beam position monitor. Its amplitude is given by where the vertical betatron tune ν y , and the Wien filter tune ν WF , are given by ν y = ω y /ω rev and ν WF = ω WF /ω rev , respectively. When the Wien filter tune is close to the vertical betatron tune, a resonant enhancement of the beam oscillation amplitude ξ y occurs. Equation (5) describes Hooke's law, F y = k H ξ y , and Hooke's constant is given by We invoke an approximate description of the betatron motion by a harmonic oscillator with constant betatron function and evaluate the Heisenberg uncertainty limit Q for the betatron oscillation amplitude ξ y in terms of the zero-point oscillator energy 1 2 ω y , which yields For the present experiment, we obtain With the actual COSY values for the betatron tune ν y and the Lorentz-factor γ of the beam (see Table I, Appendix A), the quantum limit of the vertical betatron oscillations amounts to The interpretation of the measured oscillation amplitudes in terms of the Wien filter parameters requires numerical simulations of the performance of the Wien filter as an element of the storage ring [18]. The details relevant to the present study are described below; the corresponding beam simulations carried out are consistent with the available experimental results on the properties of COSY [21].

III. WIEN FILTER OPERATION
The control of the Lorentz force of the waveguide RF Wien filter is based on the wave-mismatch principle [20]. An impedance mismatch is introduced at the load part of the device to deliberately create reflections that generate a standing wave pattern inside the Wien filter [31]. These standing waves can be represented by the complex-valued field quotient Z q , defined as the ratio of the total electric to the total magnetic field strength, where the superscripts '+' and '−' refer to the forward and backward direction of propagation, Z w is the wave impedance, Z 0 ≈ 377 Ω is the vacuum wave impedance, d = 100 mm is the distance between the electrodes, W = 182 mm is their width [18], and Γ is the reflection coefficient that controls the amplitude and phase of the reflected wave. During the measurements described here, the Wien filter was typically operated at a net input RF power of 600 W. The field quotient Z q is controlled via a specially designed RF circuit [31]. By altering Γ via two variable vacuum capacitors, called C L and C T , a wide range of Z q values can be covered, and the matching point corresponding to the minimum induced vertical beam oscillation amplitude may be determined.  [26], the WASA [27] and the JEPO [28,29] polarimeters, and the Siberian snake [30] are also shown. The coordinate system used is indicated.

IV. BEAM OSCILLATIONS
In this experiment, the electric field of the Wien filter is oriented vertically and the magnetic field horizontally. This implies that the oscillations mainly take place along the y-axis [see Eq. (1)]. For the detection of the vertical beam oscillations, a conventional beam position monitor has been employed. In order to be most sensitive, BPM 17 located in the straight section opposite to the Wien filter (see Fig. 2) with a large vertical β function was used, β BPM y ≈ 15.3049 m, while at the Wien filter location, β WF y ≈ 2.6784 m, as shown in Fig. 3. The arguments to pick BPM 17 are further discussed below in Sec. V.
In order to measure small beam oscillations, a technique based on lock-in amplifiers 1 was developed [34]. These devices operate in the frequency domain and lock onto a signal whose frequency is set as a reference, which  (2) and as far as the Fourier spectrum of the beam oscillations is concerned, without loss of generality, the BPM can be considered to be located right behind the Wien filter, and the induced voltages can be represented by where the index 't' refers to the + sign and the index 'b' to the − sign, respectively. The harmonic factor cos(ω rev t) emphasizes the pulse repetition frequency, although cos(ω rev t) = 1 for t = nT . The voltage U t, b is non-zero only at the time the beam passes through the BPM. Here, U 0 denotes the voltage proportional to the beam current that is induced when the beam passes exactly through the center of the BPM, and ∆U (∆y) represents the voltage variation induced by a beam that is vertically displaced by ∆y.
For small beam displacements, the beam position monitor operates in its linear regime, which implies that the induced voltages take the form where κ is a calibration factor that needs to be determined. At a momentum of 970 MeV/c, the revolution frequency of deuterons orbiting in COSY is f rev ≈ 750 kHz. The Wien filter is operated at the k th sideband of the spin-precession frequency f s , given by Here ν s = Gγ denotes the spin tune, i.e., the number of spin precessions per revolution, G ≈ −0.1430 is the magnetic anomaly of the deuteron, and the spin precession frequency f s = ν s f rev [18]. It should be noted that in view of ω rev T = 2π, trains of beam oscillation pulses do not depend on the actual choice of the sideband. In the present experiment the Wien filter was operated at k = −1, which corresponds to f WF ≈ 871 kHz 2 . The induced oscillations of amplitude ξ y contribute to Eq. (11) the harmonic voltage variation ∆U (y(t)). The BPM in conjunction with the lock-in amplifiers is used to measure at times t = nT the beam positions at the reference frequencies, i.e., at f WF and at f rev . Given that y(t) can be evaluated at the spin precession frequency, the BPM signals of the upper and lower electrodes can be written as follows 2 For the considerations presented in this paper, negative and positive frequencies are considered equivalent.
Here ω ∆ and ω Σ represent sidebands of the Wien filter frequency at (a) Magnitude of the field quotient |Zq|, evaluated integrally, where |Zq| int = |Zq|dℓ. Ideally, with |Zq| close to 176 Ω, the electric and magnetic forces are equal.
(b) Phase of the field quotient ∠Zq evaluated integrally, where ∠Z int q = ∠Zqdℓ. A non-vanishing ∠Zq implies a phase shift between the electric and magnetic fields. In order to measure the beam oscillations, four lockin amplifiers [34] were used, two for the horizontal and two for the vertical direction. For each direction, one lock-in amplifier detects the Fourier amplitudes at f rev ≈ 750 kHz and a second one at f Σ = f rev + f s ≈ 871 kHz. The lock-in amplifiers receive reference frequencies from the signal generator of the Wien filter and from the master oscillator of COSY. The four Fourier amplitudes of the top and bottom electrodes are determined practically in real-time, yielding The amplitude of the vertical oscillation ξ y can then be determined from The uncalibrated raw asymmetry of the four Fourier amplitudes is denoted byξ y . The readout scheme, shown in Fig. 4, was used to concurrently record radial beam oscillations, and the above analysis has been repeated for the correspondingξ x . The main result is that the coupling of vertical and radial betatron oscillations is negligibly weak, |ξ x /ξ y | < 2 × 10 −2 , which justifies treating the RF-driven beam oscillations as one-dimensional.
The determination of the calibration constant κ, required to calibrate the vertical oscillation amplitude, is described in detail in Appendix B. It amounts to κ = (5.82 ± 0.43) · 10 −6 µm −1 .
During the experiments, the vertical betatron tune of the machine amounted to about ν y ≈ 3.6040. 3 The frequency f Σ ≈ 871 kHz, at which the Wien filter is operated, is well separated from the lowest intrinsic spin resonances 4 at 297 kHz, 453 kHz, 1048 kHz, and 1204 kHz.
The two variable and highly accurate capacitors, C L and C T , are driven by stepper motors. They constitute the main dynamical elements of the driving circuit. Each pair of capacitor values yields a well-defined field quotient |Z q |, as shown in Fig. 5 (a). Away from the matching point, a phase shift ∠Z q occurs between electric and magnetic fields, as shown in Fig. 5 (b). The corresponding Lorentz force leads to the measured beam oscillations, i.e., the function ξ y = f (C L , C T ), which can be visualized in the form of a 2D map, as shown in Fig. 6. The experimental data were taken on a grid of (7 × 6) points of C L and C T , with corresponding grid spacings of (94.5 ± 1.0) pF for C L and (95.8 ± 1.0) pF for C T . Each grid spacing corresponds to 1000 steps of the corresponding stepper motors. The calibration of the capacitances C L and C T as a function of step number is discussed in detail in [31]. The grid spans over C L ∈ [318.88, 885.58] pF and C T ∈ [428.99, 907.79] pF. The uncertainties of the grid spacings are systematic ones. 5 The map of the measured and calibrated vertical beam 3 The numerical values used for the simulation calculations are listed in Table I of Appendix A. 4 An intrinsic depolarizing resonance is encountered, when the betatron motion of the particles is in sync with the spin motion, hence, when the condition fs = νsfrev = fy = (nP ± ν ′ y )frev is fulfilled [35], where n ∈ N, P denotes the superperiodicity of the lattice, and ν ′ y the fractional tune. During the experiments described here, P = 1 (see also Fig. 3). 5 The individually measured uncertainties of the capacitors are ac-oscillations ξ y is shown in Fig. 6. The parameters of the matching point are given by C L = (697.1 ± 1.0) pF , and C T = (503.0 ± 1.0) pF , (19) and the corresponding minimal detected beam oscillation amplitude at the location of BPM 17 amounts to The above accuracy of δξ min y BPM = 0.52 µm can be compared to the accuracy of measurements of static distortions of the beam orbit, which is about 20 µm (see Table II). As expected, the accuracy of the beam oscillation amplitudes is better by a factor of about 40 compared to the static amplitudes measured using the same BPM.
Upon rescaling the oscillation amplitudes using Eq. (2) with β functions listed in Table I to the WF location, we obtain which should be compared to the value of Q ≈ 41 nm, given in Eq. (9). The largest measured amplitude of driven beam oscillations at a strongly mismatched point with 600 W of input RF power amounts to Here one must bear in mind that the sensitivity to a periodic signal scales inversely with the square root of the observation time. 6 The frozen-spin proton EDM experiment aims at the accumulation of the EDM signal for a duration of about 10 7 s [7]. The accuracy δξ min In Fig. 7 (a), the data measured at the matching point [Eq. (19)] are shown. Each sample was recorded by the lock-in amplifiers with an integration time set to 0.5 s, corresponding to an average of 5000 measurements. A Monte Carlo error propagation model was applied to treat the uncertainties of the still uncalibrated raw position asymmetriesξ y and the calibration coefficient κ [37]. The results are fitted with a normal distribution, as shown in Fig. 7 (b), from which the mean value µ ξy and tually much smaller than the stated uncertainty of 1.0 pF. However, other factors, such as the capacitances and inductances of the connectors and cables and their power dependencies, also contribute to the aforementioned uncertainties. 6 The relevant discussion is found in Ref. [36]. See also the observation of white noise suppression by two orders in magnitude when using 5 h signal averaging in a test bench experiment with SQUID BPMs [22]. the error of the measured beam oscillations σ ξy are estimated. The latter represents the systematic error of the measurement. It should be noted that the map shown in Fig. 6 is actually a function of all the circuit elements. The uncertainties of ξ y are influenced by the uncertainties of all circuit elements and also by the ones of the BPM itself, which include its readout electronics, i.e., the lock-in amplifiers.
To appreciate the result given in Eq. (20), one can compare the oscillation amplitude to the 1σ vertical beam size. The latter has been deduced from the 1σ beam emittance ǫ y and the amplitude of the β function at the position of the BPM, yielding In the present experiment, the beam emittance was not monitored. The above numerical estimate of σ BPM y is based on rescaling the experimental result for the 2σ beam emittance of 49.3 MeV protons in COSY of ǫ y = (0.92 ± 0.15) µm [21] to the conditions of the present experiment. It is noteworthy that with the present equipment it is possible to access coherent beam oscillations with amplitudes that are more than three orders of magnitude smaller than the beam size.

V. BEAM DYNAMICS SIMULATIONS
To improve our understanding of the measured results, a computer code was developed to model the beam dynamics in the COSY storage ring. The modeled storage ring consists of a sequence of drift regions, quadrupole and dipole magnets, the Wien filter, and beam position monitors. These elements are represented by transfer  matrices, which are well understood and documented in the literature [38]. In the model of the ring, the actual settings of the beam optics elements of COSY were those used at the time when the experiment took place. Simulations are based on the Hamiltonian formulation as presented in Ref. [38]. The Wien filter is modeled by a timedependent matrix that also takes into account the arrival time of the particles.

A. Time-dependent Wien filter field maps
In order to be able to perform reliable beam simulations, we have placed great emphasis on good spatial resolution and the accuracy of the 3D field maps inside the Wien filter 7 , computed using a 3D electromagnetic simulation tool 8 . The fringe fields of the Wien filter are included, because they are of particular importance for the beam oscillations, as will be discussed later. An example of the computed 3D fields of the Wien filter at the experimentally determined matching point is shown in Fig. 8. The beam-tracking simulations use the three vector components of the electric and magnetic fields. The Wien filter is implemented as an RF kicker, as described by Eq. (3). Inside COSY there are 32 BPMs available to control the horizontal and vertical beam position during operation. In order to select one of them with a good sensitivity to determine the beam oscillations induced by the Wien filter, a number of particles were tracked, as described above, and the orbit response induced by a field change at the location of the Wien filter was calculated at each BPM location 9 . As a result, BPM 17, located about 70 m downstream of the Wien filter (see Fig. 2), was chosen because it offered good sensitivity to both radial and vertical beam oscillations. Figure 8 shows the main field components E y and B x . When mismatched, the Wien filter generates periodic transverse perturbations of the trajectory. Switching off the Wien filter eliminates such oscillations. The maximum amplitude of the observed oscillations in the simulation is then taken for ξ y . The two simulated vertical beam oscillation amplitudes of BPM 17 and Wien filter read which agree well with the experimentally measured results, given in Eqs. (20) and (21). A detailed description of the determination of the uncertainties of the beam simulations is discussed in the next section. For each and every measured point on the C L versus C T grid, a beam dynamics simulation was carried out. For each of these points, a 3D field map of the Wien filter was generated and then used for the beam tracking simulations. The results of these simulations are shown in Fig. 9, and are later compared with the results of the measurements.

B. Uncertainty evaluation
The accuracy with which the Lorentz force and the resulting amplitudes of the beam oscillations can be tuned depends on the accuracy with which the field quotient Z q can be integrally set to the desired value. Z q depends on the hardware elements in the driving circuit. In order to evaluate the effects of uncertainties of these elements, extensive coupled circuit electromagnetic simulations have been conducted, as discussed in Ref. [31]. The uncertainties involved are listed in Table 6 and shown in Fig. 16 of Ref. [31]. As far as the Lorentz force is concerned, most important are the uncertainties of the fixed inductance L f and the fixed resistance R f . Once these uncertainties are known, one can compute the electric and magnetic fields and the corresponding Lorentz force, including their corresponding errors. Figure 10 shows a few examples of the main components of the electric and magnetic fields, computed with the above mentioned circuit uncertainties. As will be explained below, these 3D fields, together with their uncertainties, are subsequently used as input to the beam simulations.
The algorithm used to compute the uncertainties of the beam simulations is the polynomial chaos expansion (PCE), as explained in Refs. [19,31] and in Appendix C. The PCE has been proven in many applications in science and engineering to be just as accurate as the computationally much more expensive Monte-Carlo counterpart [19,[40][41][42].
To compute the uncertainties σ ξy , the PCE algorithm requires a random set of the simulated ξ y , alongside a set of randomized input parameters according to their uncertainties to generate the output. The set of ξ y is produced using a number of beam-tracking simulations, where for each instance, a 3D field map of the Wien filter is generated, according to the randomized input parameters. An example of the electric and magnetic fields evaluated at the center of the Wien filter for the matching case [see Eq. (19)] is shown in Fig. 10. The magnitudes of the fields vary as a function of the uncertainties of the driving cir-  Table 6 of [31].
cuit [31]. The numerical tracking of the particles through these fields generates a collection of different ξ y values that the PCE algorithm can use to project the output onto orthogonal polynomial functions. These functions serve as basis functions, from which the expansion coefficients are determined that are used to generate a large sample of outputs to compute the uncertainties of the beam simulations. In Fig. 11 (a), the simulated values of ξ y are shown for the matching case. The detailed steps to achieve this result are discussed in Appendix C. As shown in Fig. 11(b), fitting these data to a Gaussian yields a standard deviation of σ ξy = 0.082 µm. This number is of considerable importance, because, given the uncertainties of the driving circuit, it sets the lower limit that can be achieved by minimizing the amplitude of the vertical beam oscillations when more sharply tuning the driving circuit of the Wien filter. The same procedure is performed on each point of the map shown in Fig. 9.

VI. COMPARISON OF SIMULATION AND EXPERIMENTAL RESULTS
The simulations yield the net Lorentz force exerted by the Wien filter on beam particles and the corresponding oscillation amplitudes for each measured point of the C L versus C T map, shown in Fig. 6. The only variables in this case are the field maps of the Wien filter itself. After 1000 turns, the beam position is computed at the same location in the ring, where the measurement using BPM 17 took place (see Fig. 2). The net Lorentz force is a result of local cancellations between the electric and magnetic field components, as illustrated in Fig. 12 for the matching point given in Eq. (19) with the minimal measured oscillation amplitude.
In Fig. 12(a), the local Lorentz force is shown along the trajectory for 5 randomly chosen passes though the Wien filter. The trajectory of the same particle changes from pass to pass, thereby different Wien filter fields and consequently different values of the Lorentz force F y will be picked up. As shown in Fig. 12, even at the matching point, the matching is still imperfect, and the largest local F y contributions are caused by the fringe fields at the entrance and exit of the Wien filter. Despite the different location of the particle in the vertical and horizontal phase space at the entrance of the Wien filter upon subsequent passes, the integration of these local forces along the particle trajectories exhibits nevertheless a perfectly harmonic time dependence with the frequency f s , as shown in Fig. 12(b). The points encircled in blue correspond to the randomly selected passes through the Wien filter, shown in Fig. 12(a).
In the left panel of Fig. 13, the amplitude of the simulated Lorentz force F sim y is plotted versus the simulated oscillation amplitude ξ sim y at the Wien filter position. As expected, it exactly follows Hooke's law with a spring constant of k H = (151.2±0.2) MeV/m 2 . In Fig. 13(b), the measured amplitudes are compared with the ones simulated for the location of BPM 17. The two sets ξ exp y and ξ sim y are in very good agreement with each other. The horizontal and vertical error bars are derived from the uncertainties of the measurements and simulations, represented by the width of the distributions, as shown in Figs. 7(b) and 11(b). It is important to note that the error bars refer to systematic uncertainties and should not be confused with statistical ones. This implies that repetitions of either the measurements or the simulations will neither reduce the systematic error of the readout electronics of BPM 17, nor will it affect the uncertainties of the elements of the driving circuit.  . Depending on the initial coordinates in the vertical and horizontal phase space, the particle travels along different trajectories, and therefore picks up different field components Fy.
The fit shown in Fig. 13 yields χ 2 /ndf = 45.5/41, very close to unity [43]. The linear fit yields a slope of 0.999 ± 0.018, which is perfectly consistent with unity. The intercept parameter of the fit yields (−0.93± 0.31) µm, and within three standard deviations, it agrees with zero.
The very good agreement between measurements and simulations reflects our good understanding of both the electromagnetic fields generated in the Wien filter and of the underlying beam dynamics in the machine. This point is further substantiated by comparing the simulated amplitudes at the Wien filter and at the positions of the BPMs with the estimated amplitudes expected from rescaling based on the β functions 10 , taking into account the numerical values, listed in Table I   The good agreement between these two numerical values in Eq. (25) indicates that the observation of the oscillation amplitude at one location in the ring can be reliably transferred to some other place in the ring by use of Eq. (2). The above quoted value of ξ WF y sim = 0.435 µm is about a factor of 10 larger than the quantum limit of the vertical oscillation amplitude Q, given in Eq. (9).
In searches for EDMs in dedicated all-electric storage rings, a continuous monitoring of the orbits of the two counter-rotating beams is mandatory during data acquisition within the horizontal spin-coherence time [7]. When intrabeam scattering can be neglected [21], which is arguably justified within the horizontal spin-coherence time, the beam can be described as a rarefied gas of particles, i.e., the zero-point oscillations of individual particles are uncorrelated. We repeat the point from the introduction that in a static regime, the quantum limit of the center of mass of a bunch with N particles can be estimated via Q bunch = Q/ √ N . For a bunch of N = 10 10 stored particles, one obtains Q bunch ≃ 0.4 pm. It follows that Heisenberg's uncertainty relation does not present an obstacle to achieving a sensitivity of 5 pm for the vertical separation of clockwise and counter-clockwise beams -the real challenge is to develop compact BPMs with a sensitivity improved by a factor of about 300 compared to those used here [25].
The theoretical estimate of k th H , calculated using the numerical values listed in Table I

VII. CONCLUSION AND OUTLOOK
As part of several studies to investigate the performance of the waveguide RF Wien filter, exploratory data were taken to provide a benchmark on the sensitivity to very weak collective vertical beam oscillations of deuterons stored in the COSY ring. To a good approximation, the beam can be viewed as a rarefied gas of uncorrelated particles, and the sensitivity limit is applicable to the classical motion of individual particles, propagating along the ring circumference in the confining oscillatory potential. Simulations of the beam dynamics in the COSY ring equipped with an RF Wien filter suggest that with the present apparatus, the sensitivity to collective beam oscillations on the sub-micron level is only a factor of about 10 larger than the amplitude of single-particle zero-point quantum oscillations of the stored deuterons. From the perspective of future EDM experiments, our finding confirms that, as far as the Heisenberg uncertainty relation is concerned, a separation of the centroids of two counter-propagating beams may be determined to sub-picometer accuracy.
The reported excellent agreement between simulated and experimentally observed vertical beam oscillations at COSY suggests that a further increase in sensitivity to collective beam oscillations is possible. Specifically, the simulation on finer capacitor grids indicates that by further optimization of the Wien filter settings to C L = (692.76 ± 1.00) pF and C T = (495.77 ± 1.00) pF, an oscillation amplitude at the Wien filter location of ξ y = (0.077 ± 0.032) µm may be achieved. Thus in that case, the vertical oscillation amplitude would only be about a factor of 2 away from the quantum limit, with a corresponding Lorentz force of F y ∼ 3 aN.

ACKNOWLEDGMENTS
We would like to thank I. Bekman, B. Breitkreutz, J. Hetzel and R. Stassen for their support in setting up the COSY accelerator for the experiment. The work presented here has been performed in the framework of the JEDI collaboration and is supported by an ERC Advanced Grant of the European Union (proposal number 694340). In addition, it was supported by the Russian Fund for Basic Research (Grant No. 18-02-40092 MEGA) and by the Shota Rustaveli National Science Foundation of the Republic of Georgia (SRNSFG Grant No. DI-18-298: High precision polarimetry for charged-particle EDM searches in storage rings).

Appendix A: Quantities used in beam simulations
In order to provide a consistent calculation of all effects in the storage ring, the beam simulations were carried out using the set of quantities given in Table I as an input. The vertical machine tune ν y is a result of simulations with the known COSY lattice, reflecting the actual currents of the magnetic elements in the machine at the time when the experiment was conducted. The simulations provide the uncalibrated parameters of the vertical beam oscillations to about per mill accuracy, and giving the kinematic, ring, and Wien filter parameters to four digits appears therefore sufficient. It should be noted that within the simulation calculations carried out in the context of the present work, all quantities have been computed to double precision (machine epsilon of 1.11 × 10 −16 ). Of the physical quantities, the highest sensitivity to the vertical betatron tune is exhibited by the theoretical estimate for Hooke's constant, dk th H /dν y ≈ 2 × 10 3 MeV/m 2 . The largest uncertainty contributing to the error of the detected oscillation amplitudes arises from the calibration factor κ of the beam position monitor, given in Eq. (B2). It amounts to about 7.3% and is considered a systematic scale-factor uncertainty (see Appendix B). and A rev b , given in Eqs. (16), are related to a vertical beam displacement ∆y in the following way, The calibration constant κ is experimentally determined by introducing local vertical beam bumps in the ring at the location of BPM 17. The orbit positions y and the orbit displacements ∆y, listed in Table II, differ by the position of the unperturbed orbit, and are generated by altering the current of a set of vertical steerers. The steerer magnets have well-known conversion factors from current to magnetic field. The calibration factor κ is obtained by fitting the ratio R from Eq. (B1) as a function of the vertical orbit variation ∆y, exhibiting the nearly linear relationship shown in Fig. 15. The slope corresponds to κ = (5.82 ± 0.43) · 10 −6 µm −1 . (B2) can be used to generate an arbitrarily large sample of observables, from which the uncertainties can be computed to the desired statistical accuracy. The PCE algorithm has been compared with the MC method in many applications and has been shown to provide very reliable results [40]. The PCE requires much fewer simulations to converge compared to the MC method. For instance, for the present case, 200 beam tracking simulations per point in the 2D the map of beam oscillations, shown in Fig. 9 were sufficient to reach convergence. In cases where the number of random input variables m is larger than 10, the PCE method offers clear advantages over the MC method. The reason is that the number of basis functions in the PCE method increases enormously as a consequence of the tensor product of the involved polynomials. Therefore, the algorithm has been improved further to allow for a reduction of the number of simulations required. Such an approach is also adopted here, as described in Algorithm 1. The hyperbolic truncation scheme together with the Least-Angle Regression (LAR) method form a sparse version of the original algorithm.
An m-dimensional set is first created, representing N combinations of simultaneous random variables. Many methods can be used to generate such sets, and here the Latin-hypercube sample scheme is adopted [31]. Subsequently, the set is standardized for convergence reasons.
Depending on the distribution of the data, the basis functions, here Hermite polynomials, are determined. The number of basis functions restricts the lower limit of the number of simulations (full-wave and tracking) which are usually computationally expensive. As a rule of thumb, with N basis functions, the PCE algorithm requires at least 1.5 × N (in this case, full-wave) simulations to converge. The number of basis functions itself can, however, be reduced by the hyperbolic truncation scheme that eliminates higher-order terms that do not have a significant impact on the observation objects [44,45]. Furthermore, by applying the LAR algorithm, the number of remaining basis functions can be further reduced substantially, whereby the problem becomes computa-tionally solvable in a very efficient fashion. The matching point, specified in Eq. (19), yields the minimum measurable beam oscillations, as given by Eq. (20). This experimental result can be estimated using the beam-tracking calculations. Subsequently, the concrete steps of the application of the PCE algorithm are discussed.
All the reasonable sources of uncertainties of the circuit are represented by 15 random parameters that are allowed to vary simultaneously. At first, a sample of (200 × 15) entries is generated using the Latin-hypercube sampling scheme. As an example of this sample, the variation of the three circuit elements C L , C T , and the load resistor R f is shown in Fig. 16(a).
All the uncertain parameters in the electromagnetic circuit simulations are used to generate the electric and magnetic fields shown in Fig. 10. These are subsequently used in the beam-tracking calculations. For the matching point of the map [Eq. (19)], N = 200 full-wave simulations were conducted. The import of these field maps into the beam-tracking calculations resulted in a set of N = 200 values of ξ y . This set is not directly used to conduct the statistical analysis. Instead, in conjunction with the input samples, these data are used as input to the sparse PCE algorithm.
The optimum set of basis function is determined using the LAR algorithm, as shown in Fig. 16(b). With an expansion order of p = 6 and a truncation norm q = 0.35, executing the PCE algorithm required 91 basis functions to converge, reflected by the low value of the leave-one-out error LOO err = 1.7 × 10 −4 . Subsequently, the expansion coefficients are computed, qualitatively depicted in Fig. 16(c). It is shown in Fig. 16(d) that the PCE algorithm perfectly reproduces the tracking results using these expansion coefficients. Finally, these coefficients are used to reconstruct a larger sample of ξ y to (c) Expansion coefficients on a semi-log scale. 91 coefficients have been selected after applying the LAR algorithm to the matching point.
(d) Comparison between the tracking results and the PCE with respect to the oscillation amplitude, determined using the expansion coefficients of (c).  Table III. estimate the error σ ξy . Figure 11 shows 10 4 of the 10 6 reconstructed samples. The PCE parameters used are summarized in Table III.
The fitting of these results with a Gaussian, as depicted in panel (b) of Fig. 11, yields a standard deviation of σ ξy = 0.082 µm. The same technique is repeated for each point in the map.