Monte Carlo simulations of a neutron lifetime experiment with a big gravitational trap

A Monte Carlo model was developed for our current neutron lifetime experiment using storage of ultracold neutrons (UCN) in a big gravitational trap. The model allows us to calculate neutron trajectories in the given geometry, taking into account gravity, and has been used in simulations that reproduce all stages of the experiment. The simulation was proven able to reproduce the time-dependence of the detected neutron count rates. For investigation of systematic effects, the value of the neutron lifetime resulting from the simulation of the experimental procedure was compared to the value entered as an input parameter. In particular, the systematic uncertainty associated with the method of calculating the effective UCN collision frequency in the trap was found not to exceed 0.1 s. In addition, calculations showed that up to 1% uncoated area of the copper trap has a negligible influence on the neutron lifetime result.


Introduction
The neutron lifetime is among the most important fundamental constants for the theory of weak interactions and cosmology. The history of neutron lifetime measurements [1] covers a significant period, starting from first experiments using neutron beams in the 1950ies [2]. Significant progress was achieved using ultracold neutrons (UCNs). One experimental approach consists in storing them in a material trap, with gravity providing vertical confinement. The most precise experiment of this type to date was carried out in 2004 at the UCN facility PF2 at Institut Laue-Langevin (ILL) by a PNPI-JINR-ILL collaboration, which delivered a neutron lifetime value of (878.5 ± 0.7 ± 0.3) s [3,4]. Here we focus on the follow-up project started in 2017 by a PNPI-RAL-ILL collaboration, which features a big gravitational trap at low temperature. The result of first measurements carried out again at PF2 is (881.5 ± 0.7 ± 0.6) s [5,6]. Independent experimental efforts using magnetic trapping have recently been published [7,8]. While the values measured in recent UCN experiments agree reasonably well, they also confirm a discrepancy with values of the neutron lifetime extracted from beam experiments [9][10][11].
Here we review results of simulations performed with a Monte Carlo model of our experiment with a big gravitational trap. The model was designed to reproduce all experimental procedures with the goal to confirm the validity of concepts and to investigate possible systematic errors. Results of complementary investigations that address specific steps of the experiment and related systematic uncertainties are also presented. a e-mail: fomin ak@pnpi.nrcki.ru

Experimental method
The neutron β-decay lifetime is evaluated from the UCN storage lifetime τ st , which can be deduced from the measured number N of neutrons remaining in a trap after different holding times t i . Using two holding times it follows from (see, e.g., Ref. [12]): A rate equation on the other hand relates τ st with the probability rates for neutron β-decay, τ −1 n , and for UCN loss due to other processes, τ −1 loss , respectively: In a high vacuum, losses occurring at wall collisions are dominant. These can be parameterized as where ν (E) is the frequency of wall collisions, and µ (T, E) is the loss probability per neutron and per wall collision, which depends on the wall temperature T and the neutron energy E. Assuming a square barrier wall potential with a real part U 0 and an imaginary part W , and including angular averaging, it can be represented as [13]: (4) Here, η (T ) = W/U 0 is the UCN energy independent loss coefficient, and y = √ E/U 0 = v/v b is the ratio of the actual neutron velocity to the maximum velocity v b at which full reflection from the potential barrier is still possible. Note that in all our simulations is calculated including the angle of impact at each neutron collision with a surface, i.e., applying the formula prior to angular averaging. Capture and inelastic scattering cross sections are both proportional to the neutron wavelength, hence η does not depend on neutron energy. Averaging Eq. (4) over all angles, Eq. (3) can be written as where γ (E) is the calculated effective collision frequency which depends on the UCN energy and on the trap size. The neutron lifetime is deduced from a linear extrapolation of τ −1 st to γ (E) → 0. Experimentally, γ (E) can be varied by using traps with different geometry (referred to as "geometry extrapolation") and/or using a different UCN energy spectrum (this is called "energy extrapolation"). For known trap geometry the value of γ (E) can be calculated analytically. Combining Eqs. (3) and (5) one obtains γ = µν η .
The total collision frequency ν (E) can be calculated from the UCN flux directed to a surface element dS. It is given by ρvdS/4, including the UCN density ρ. Since both the energy and density of the UCNs depend on the height in the trap, Eq. (6) must be integrated over the height and normalized, i.e., H is the maximum height that a neutron with a given energy E can reach. In the gravity field, the UCN density is given by ρ = √ (ε − mgh) /ε, where ε is the UCN energy at the lowest point of the trap, and h is the height above this point. The UCN velocity is, in the same reference system, given by where v 0 is the neutron speed at the trap bottom. Finally, for obtaining the neutron lifetime from the extrapolations mentioned above, γ (E) is to be weighted by the UCN emptying spectra given by the Monte Carlo simulation.
It is a special feature of our setup (a schematic is shown in Fig. 1), that the UCN trap (1) can be cooled down to cryogenic temperatures (∼100 K). It is made of copper in shape of a horizontal half cylinder, with a length of 2 m and a diameter of 1.4 m (hence a trap height of 0.7 m). A copper insert (2) can be rotated into the trap for variation of the UCN collision frequency. It has a diameter of 1.2 m and a length of 1.8 m. The height of the insert is 0.66 m. Two holes at the bottom enable passage of UCNs between the internal volume of the insert and the gap between the insert and the trap. All surfaces are coated as described further below. The trap and the insert can be rotated independently about a common pipe-in-pipe-type axis, using actuators located on opposite sides. It takes 45 s to rotate the trap by an angle of 90 • . A titanium plate (10) serves as an absorber to clean the spectrum from neutrons with energies exceeding the gravitational barrier of the trap. It is rigidly attached to the insert and hence always horizontal during the UCN storage measurements done with the insert fixed either in the upmost or the lowest position. The latter situation is shown in Fig. 1, corresponding to the angle of the absorber/insert set as θ abs = 0 • . The rotation angle of the trap is defined accordingly, i.e., θ = 0 • for its lowest position.
A cylindrical, internal vacuum vessel (3) housing the trap has a diameter of 1.62 m and a length of 2.5 m. The outer vacuum chamber housing the whole setup is not shown in Fig. 1. A vertical neutron guide (4) with a diameter of 0.27 m and a length of 1.19 m connects the bottom part of the vessel with a UCN detector (6). A shutter (5) controls the release of neutrons to the detector. The setup is filled with neutrons from the source through an input guide (diameter 0.14 m, length 2.1 m) (7) via a second shutter (8). A thermal shield (9) is protecting the trap from thermal radiation of the neutron guides.

Monte Carlo model
Our model simulates the whole experimental procedure together with the data processing for the geometry and energy extrapolations. Neutron β-decay is implemented with a preset lifetime value τ n0 . The lifetime value obtained PPNS 2018 from a complete simulation is finally compared with τ n0 . Performing such simulations for various settings of model parameters allows us to identify important effects and draw conclusions about systematic errors. This method was already used for simulation of the earlier experiment including a successful benchmark against data [3,4].
The model simulates the fate of each neutron having entered the neutron guide system from the source. Initial neutron velocities are defined to follow a Maxwellian spectrum, with a maximum neutron energy cutoff set 30 neV above the gravitational holding barrier of the trap. Neutron position and velocity evolve according to the equation of motion of a free projectile in the Earth's gravitational field, with boundaries set by the walls of the guides and the trap. The history ends if one of the following conditions is met: 1) The neutron experiences a β-decay; 2) The neutron returns to the UCN source; 3) The neutron is lost due to exceeding the boundary velocity of the wall material; 4) The neutron is lost in a collision with the wall; 5) The neutron is absorbed by the titanium absorber; 6) The neutron enters the UCN detector.
A typical measurement procedure (with parameters summarized in the Tables 1, 2 and 3 for the different types of simulated measurements discussed here) consists in the following steps. First, the internal vacuum volume (3) is filled for 200-300 s with neutrons from the source via the neutron guide (7) and through the open shutter (8) while the other shutter (5) is closed. During filling, the UCN trap is in its upmost position (θ = 180 • ). At the end of the filling process the trap angle is set to θ = 15 • , thus confining the UCNs. Then, the shutter (8) is closed and the shutter (5) is opened. This step, called "monitoring", is maintained for 300-500 s, for removing those neutrons with energies exceeding the gravitational holding barrier of the trap. Next, the trap is turned into the holding position (θ = 0 • ) in which UCNs are stored for up to 2000 s. The holding period is followed by few (typically 3 to 4) "emptying" positions attained by subsequent turns of the trap, wherein the angle θ is stepwise increased to finally reach 180 • . The counting time of UCNs leaving the trap in each individual tilted position lasts between 200 s and 400 s. The insert is kept at its upmost position (θ abs = 180 • ) all the time when measurements are made "without the insert" and in its lowest position (θ abs = 0 • ) for the measurements "with the insert".
The surfaces of the trap and the insert are coated with a perfluorinated grease (fomblin UT-18) with a boundary velocity of 4.85 m/s and a loss coefficient of 8 · 10 −6 [6]. At the trap walls, the insert, and the vacuum volume, we assume a probability for diffuse neutron reflection of 10% [4]. At the surface of the neutron guides on the other hand, reflections are mainly specular and here we assume a 0.7% probability for diffuse reflection. This value has been evaluated from UCN guide transmission data [14]. In the case of a diffuse reflection, Lambert's law is used; i.e., the probability of reflection at angle α to the normal to the surface is proportional to cos α.

Simulations of the experiment
We performed three different main simulations to derive neutron lifetime values: one with the monitoring at trap orientation θ = 15 • but without the absorber, one without monitoring but with the absorber, and a last one with both, monitoring and the absorber. Parameters and results are quoted in three corresponding subsections. The simulations in Sect. 4.1 were carried out at the design stage of the setup. Therefore, the loss coefficient and the heights of the trap and insert are not the same as in Sects. 4.2 and 4.3. For each case we report differences of simulated neutron lifetime values from the model preset value τ n0 , for a measurement without the insert, with the insert, and for the joint extrapolation of these two measurements. Each simulation consisted of 4 runs with the same number of starting neutrons: for short and long holding times, with and without the insert. The number of neutrons was chosen to obtain a final accuracy of 0.1 s for each joined extrapolation. The total number (i.e., 4 times the number per run) is reported together with each group of results. Note that it includes the neutrons with energy up to 30 neV above the gravitational holding barrier of the trap.

Simulation with monitoring but without an absorber
This simulation is performed for measurements using the two holding times t 1 = 300 s and t 2 = 2000 s and four emptying steps with periods and tilt angles θ as listed in Table 1. The value of loss coefficient for fomblin coating is taken as 2 · 10 −6 (i.e., the value attained in our earlier experiment [4], using a low-temperature fomblin oil). In    Figure 2 shows the simulated time diagram of counted neutrons. From this, we extrapolated the storage lifetimes to γ (E) = 0. The resulting value for the simulated neutron lifetime deviates from the model preset value τ n0 by (0.05 ± 0.54) s for the measurements without the insert, by (-0.16 ± 0.32) s for the measurements with the insert, and by (-0.05 ± 0.10) s for the joint extrapolation [15]. The total number of neutrons injected into the 4 simulation runs was 6.5 · 10 11 .

Simulation without monitoring but with the absorber
An issue of the monitoring step is that it reduces the number of detected neutrons. Indeed, while at θ = 15 • the spectrum is cleaned from neutrons with energy exceeding the gravitational holding barrier (E trap = 71.7 neV for the trap used in the experiment [5,6], see Table 2), neutrons with smaller energies in the range E trap = (53.2−71.7) neV also leak out from the trap. A further significant reduction of trapped UCNs during a monitoring period as long as the 500 s in the previous set of simulations (Sect 4.1) is due to β-decay and losses at wall collisions.
It is therefore useful to study a measurement procedure without monitoring but using the fixed horizontal absorber described above. In experiments using UCNs, a wide range of materials can be used as absorbers [16]. Titanium was chosen from a constructive point of view. Figure 3 shows the energy of those UCNs being absorbed in the course of time (without taking quantum reflections into account). One can see that neutrons with energies slightly exceeding the threshold required to reach the absorber surface can remain in the trap for a rather long time. This can cause a systematic error if they stay exposed to the absorber during the holding period.
A solution of this problem would be to lift the absorber for the holding and emptying period. However, this is technically not possible in the current setup. On the other hand, it is possible to rotate the absorber (together with the insert) by some angle, thus exhibiting some absorbing surface at a somewhat lower height in the trap. After cleaning the spectrum that way, before the holding period, it has to be rotated back to the horizontal position to bring the whole absorbing surface to maximum height. Figure 4 shows the results of a simulation of neutron absorption events at the absorber permanently tilted by θ abs = 4.78 • such that its lower edge was 5 cm below the trap axis (defined as z = 0). One can see that after cleaning no more absorptions occur at z > 0, so that no more neutrons will hit the absorber rotated back to the horizontal position (z = 0).
Simulations of a complete set of neutron lifetime measurements were performed with the set of parameters shown in Table 2. The following discrepancies of the deduced neutron lifetime from the preset value τ n0 were obtained: (-0.17 ± 0.53) s for the measurements without the insert, (-0.03 ± 0.33) s for the measurements with the insert, and (-0.07 ± 0.10) s for the joint extrapolation [17]. The total number of neutrons injected into the 4 simulation runs was 3 · 10 11 .
We also investigated the influence of uncoated areas of the trap, which in the real experiment do not exceed 0.1% of the total surface. Since the substrate material of the trap and the insert is copper with a boundary velocity of 5.6 m/s and a loss coefficient of 3 · 10 −4 [13], one may expect the impact of such areas to be negligibly small. The simulations shown in Fig. 5 confirm this: for all the simulated variants (0%, 0.1% and 1% of uncoated  area), the model neutron lifetime τ n0 is reproduced by the extrapolations of τ −1 st with an accuracy of 0.1 s: deviations were found to be (0.08 ± 0.10) s for 0.1% uncoated area, and (0.05 ± 0.10) s for 1% uncoated area.

Simulation with monitoring and with the absorber
Despite the good agreement of neutron lifetime values obtained in the previous case, it is nonetheless reasonable to do some monitoring in the real experiment. Online detection of those neutrons with energies exceeding the gravitational holding barrier and leaving the trap during the monitoring provides a valuable experimental control of the cleaning process of the trapped UCN sample. For this, the monitoring angle can be chosen smaller than 15 • . According simulations were performed for measurements with the parameters listed in Table 3. The calculated   data have been compared with the experimental ones (Fig. 6). Best agreement in the region of the monitoring was obtained for the case when 50% of the surface of the titanium absorber was oxidized and assumed not to absorb neutrons. The agreement with the experiment in the region of the emptyings does not depend on this parameter.
The discrepancies of the deduced neutron lifetime from the preset value τ n0 of the model are (0.73 ± 0.50) s for the measurements without the insert, (-0.07 ± 0.31) s for the measurements with the insert, and (-0.06 ± 0.10) s for the joint extrapolation [18]. The total number of neutrons injected into the 4 simulation runs was 4.8 · 10 11 .
As additional information, Fig. 7 shows the influence of a variation of the tilt angle of the trap for the first UCN extraction. It shows the sensitivity to a correct setting of this angle. Corresponding UCN spectra are shown in Fig. 8.

Conclusion
The developed model proved to be an effective tool for analysing neutron lifetime experiments with gravitational UCN traps. The simulations presented here have shown that even an unrealistically high fraction of 1% of area of the copper trap uncoated by fomblin has a negligible influence on the neutron lifetime result. The calculated and experimental time dependence of the detector counts coincided within the error bars. The difference between the preset and the extrapolated value is considered to be a systematic error associated with the method for calculating the effective collision frequency of UCNs in the trap. The error obtained by this procedure does not exceed 0.1 s [5,6].
Calculations were performed using the computer cluster of the Information Technologies and Automation Division as well as the Data Center of reactor PIK in NRC "Kurchatov Institute" -PNPI.