Numerical study of the “two components” model and background effects in muonic hydrogen experiments

In this study, we explore the X-ray time spectra of muonic oxygen in various H 2 + O 2 gaseous mixtures by means of numerical simulations based on a sophisticated kineticrate equations model. Our primary focus is on the “two components” hypothesis, which was put forward to account for the experimentally observed two successive exponential behaviors in the radiation time-distribution. To ensure an accurate representation of the physical processes in the experiments, we incorporate precise transition rate coefficients for the elastic and inelastic scattering between the muonic hydrogen atoms and other gas constituents, in our simulation code. Specifically, a crucial aspect in our approach is the utilization of the energy dependence of the muon transfer rate to oxygen, recently extracted from experimental data by the FAMU collaboration. The accurate knowledge of such muonic oxygen X-ray time spectra is vital for cutting-edge experiments like FAMU, where the hyperfine splitting in muonic hydrogen atoms is to be determined with high accuracy, as this time spectra serves as a significant part of the background for the desired signal.


Introduction
The study of light exotic atoms and molecules offers a fresh perspective on fundamental physics.Specifically, muonic hydrogen (pµ) provides valuable insights into the electromagnetic structure of the proton.This is because the muon's significantly larger mass, compared to the electron, results in a greater overlap of the particles' wave functions.This property was leveraged by the research group led by Randolf Pohl at PSI, which conducted measurements revealing a slightly smaller proton charge radius than initially anticipated [1].This discovery led to the so-called 'proton size puzzle', a challenge that has been mostly resolved.
However, our understanding of another crucial proton characteristic, the convolution of its electric charge and magnetic moment distributions known as the Zemach radius, remains limited.Furthermore, substantial disparities exist among various experiments attempting to determine its precise value [2,3], and currently, several ongoing experiments are dedicated to resolving this discrepancy [4,5,6].
All of these experiments rely on precise measurements of subtle effects superimposed on combined radiation resulting from various processes.Thus, a comprehensive understanding of the time-dependent background signal patterns is crucial.This may include the characteristic radiation emitted by higher-Z muonic atoms, which is an intriguing problem in its own right.
Such µZ atoms and charge exchange reactions involving the transfer of muons from pµ atoms to Z-nuclei have been studied by a number of authors [7,8,9] and form the foundation of experiments conducted by the FAMU [4] and CREMA [5] collaborations.
In this study, we investigate the time distribution of X-ray radiation produced by muonic oxygen (µO) effectively immediately after their formation in various gas mixtures of H 2 and O 2 .These investigations are carried out in the context of current or planned pµ experiments, where different configurations of experimental parameters will be explored.Our primary focus is to delve into certain aspects of the 'Two Component' (TC) model proposed by Schneuwly [10].This model posits the existence of a two-component energy distribution of muonic hydrogen immediately following the muon spill -one thermal and the other hot in nature.
The manifestation of the TC model becomes evident in experiments through the structure of the X-ray time spectrum, which can be represented as a weighted sum of two exponential functions.This behavior arises from several key factors.Initially, the 'hot' muonic hydrogen component relaxes rapidly toward lower energy values, thus becoming 'epithermal'.They are added to the epithermal muonic hydrogen pµ produced through another mechanism consisting of transition of pµ to a state with unit total spin F = 1 followed by deexcitation to F = 0 state accompanied with an increase of the kinetic energy [11].During the initial stages of the radiation process, there is a dominant muon transfer to oxygen from the epithermal pµ atoms due to the higher transfer rate λ pO (E) for such intermediate energies.However, over longer time intervals, the epithermal pµ population is depleted, either through the aforementioned transfers or thermalization, and the transfer to oxygen from the thermally distributed muonic hydrogen takes precedence.
In a subsequent study, Werthmuller et al. [7] supported the essence of the TC model by analyzing a series of experiments using gaseous mixtures with varying oxygen concentrations, pressures, and temperatures.They proposed a mathematical representation using a sum of three weighted exponential functions, with the third function serving to correct the behavior at the beginning of the muonic oxygen X-ray radiation.Additionally, through Monte Carlo simulations, they obtained a two-step function for the transfer rate of muons from pµ to oxygen.
However, these findings need revision due to the sensitivity of the observable µO X-ray time spectra to the transfer rate λ pO (E), especially in light of recent results obtained by the FAMU collaboration [12].These results highlight significant restrictions on the behavior of λ pO (E) for µp energies up to 0.1 eV.
To address these developments, we have developed a comprehensive model and conducted numerical simulations to provide a more detailed understanding of the physical processes underpinning the experimental outcomes.Our program code incorporates precise transition rate coefficients for elastic and inelastic scattering among muonic hydrogen atoms, as well as their interactions with other species.Furthermore, we explore various functional dependencies of the transfer rate λ pO (E) while adhering to the constraints outlined in [12].
The paper is organized as follows.In Sec. 2, the kinetic-rate model we use in the computations of the time distribution of muon transfer to oxygen is briefly reviewed.The assumptions and approximations employed in our study, together with the simulation results are presented in Sec. 3. In Sec. 4 the obtained results regarding the behavior of the X-rays emitted in the process of muon transfer are analysed and discussed.Final remarks are given in Sec. 5.

Kinetic-rate model of muonic hydrogen decay
The integral decay of muonic hydrogen arises as an overlay of various decay effects.Despite its complexity, one can still approach its description in a relatively straightforward manner by leveraging the available data for the pertinent averaged scattering rates and decay rates within individual channels.By doing so, as explained bellow, it is necessary to track the kinetic energy distribution of muonic hydrogen which we simplify by discretizating it into n bins.Hence, we model the process using a system of linear differential equations of the kinetic-rate type, characterized by time-independent coefficients.These equations govern the dynamics of the probabilities N (iα) (t) = N (E i , F α , t) to find a muonic hydrogen atom in spin state F α = 0, 1 with kinetic energy E i in the i th bin relative to the laboratory frame.The evolution of the initial probability distribution 2n-dimensional vector N (0) is therefore determined via where the exponentiated 2n × 2n matrix L comprises the rates of the various processes.We can formaly decompose L into two components L = L d + L s .The first component, L d , accounts for the disappearance of µp with negligible kinetic effects related to the changes in the particle's kinetic energy.It is composed of a complex constant λ 0 representing muon decay and nuclear capture rates, as well as rates λ dµ and λ ppµ describing decays due to inelastic collisions in the gas leading to the formation of molecular ions ppµ, muonic deuterium µd, and muon transfer to oxygen λ pO (E) Here, I is the 2n × 2n identity matrix and λ α,β pO (α, β = 0, 1) is an n-dimensional vector with entries λ pO (E i ), i = 1, ...n, and it is assumed that the muon transfer to oxygen λ pO (E) is the same for both possible values of the total spin F = 0, 1.The accurate knowledge of energy dependence of λ pO and its inclusion in (2) are particularly important, as this energy dependence provides an opportunity to employ experiments for precise measurements of hyperfine energy splitting of the ground state of muonic hydrogen.The second component, L s , is related to both elastic and inelastic scattering of the muonic hydrogen on the protons of molecular hydrogen.It can be written in a block-matrix form as where the scattering rates are encoded by the n×n matrices λ α,β , (α, β = 0, 1).All four matrices λ α,β are build by the muonic hydrogen transfer rates describing the transition from a state with spin F = α and kinetic energy E i to a state with F = β and E j .
In the above definitions all the rates but λ 0 are normalized to the liquid hydrogen density (N LHD = 4.25 10 22 atoms/cm 3 ), ϕ is the number density of the atoms of the gaseous mixture in LHD units, and the number concentrations of the x component is denoted by c x .

Simulations of muonic hydrogen experiments
Real-world simulations necessitate the incorporation of specific experimental data, and often prompts some simplifying approximations.In this section we provide further details for the additional assumptions and approximations used in our model, along with the sources of the data for the physical quantities we employ.Furthermore, we present the results of our numerical investigation concerning the "Two Component Model" and compare them with experimental data.
3.1.Approximations and resources 3.1.1.Fundamental constants and parameters In our simulations, we have incorporated the most up-to-date values for fundamental constants and various particle characteristics, including their masses and decay rates, by referencing the latest recommended values from CODATA [13] and the Particle Data Group [14].

Number density of atoms ϕ
To obtain the number density of atoms ϕ, normalized to LHD, we adopt the ideal gas approximation.This simplifies the calculation and reduces the number of thermodynamic parameters influencing the model to only pressure and temperature.Strictly speaking, we focus on describing the non-stationary processes that lead to X-ray radiation from the µO atoms.However, it is acknowledged that the percentage of gas constituents involved in these relaxation processes is small.
3.1.3.Muonic hydrogen scattering rates λ α,β We account for the muonic hydrogen scattering rates λ α,β , (α, β = 0, 1) on hydrogen molecules by making use of the corresponding values for T = 300K provided in Adamczak et al. [15].These values are obtained through intricate calculations involving numerous approximations, which we deem appropriate for the description of considered experiments.
3.1.4.Muons to oxygen transfer rate λ pO (E) The energy-dependent transfer rate of muons to oxygen λ pO (E) plays a crucial role in shaping the observed X-ray radiation patterns.Here, we make the assumption that the muon transfer rate to oxygen is identical for both possible values of the total spin λ 0 pO (E) = λ 1 pO (E).We generate and evaluate the performance of various rates that comply with the constraints derived in [12], contrasting them with experimental data.
3.1.5.X-ray radiative decay of muonic oxygen dN O /dt The measured quantity in the experiment is the X-ray radiative decay of muonic oxygen atom, denoted as dN O /dt.In our calculations, we make the assumption that the number of X-ray events per second is equal to the number of muon transfers to oxygen.This assumption is grounded in the rapid radiative transition of the initially formed µO atom from its excited state to its ground state.
3.1.6.Initial distribution of muonic hydrogen Inline with the TC-model the probability distribution of muonic hydrogen pµ is selected to account for a portion of the population following Maxwell-Boltzmann statistics, while the rest is concentrated in the high energy bins.Various choices for the description of the components are considered to both validate the model and describe the experimental scenarios.However, we have inferred that the dependence on slight variations in factors such as relative weight of the components, number of bins that spans hot component and its shape, is relatively weak.In our simulations we choose to assume an equally shared population between both components.

Numerical results
In this work we study the time-distribution of X-rays, emitted due to the transfer of muons to oxygen, by numerical simulations with the code described above.To assess the goodness of our results we need a reliable benchmark.As such, we use the three weighted exponential function given by Werthmuller et al. [7]: For the particular combination of initial conditions: oxygen concentration c O = 0.396 %, gas pressure P = 15.06 bar, and temperature T = 300 K, the parameters λ 1 = 0.006329, λ 2 = 0.02941, λ 3 = 0.08333 are given in [7].The coefficients a = 1212 and b = 1868 are found by fitting the experimental data, digitized from the graphical results presented in the same work, by Eq. ( 4).The obtained curve is shown on Fig. 1 with solid black line and the blue crosses depict the digitized experimental data points.We will note that the first peak in the experimental data corresponds to the initial capture of muons directly by the oxygen atoms.As we focus on muon transition from pµ to oxygen, the former process will not be investigated.With light-green is marked the area of λ pO allowed by the latest FAMU results [12].For E > 0.1 eV there are no restrictions imposed by FAMU.(Right) Time distribution of the characteristic oxygen X-rays following the transition of muons from hydrogen to oxygen.Three distinct simulation curves are depicted, each corresponding to λ pO (E) given on the left picture with the same color.The solid black curve is the fit function from Fig. 1.In both figures, the dashed green, dotted blue, and dot-dashed red lines correspond to λ S pO (E) (taken from [12]) and the transfer rates given by and Eqs. ( 5) and ( 6).
Our investigations show that the quantity which has the greatest impact on the X-rays time distribution behavior, is the transfer rate of muons to oxygen λ pO .To demonstrate it, we have performed simulations of the characteristic oxygen X-ray events with a few different λ pO (E) and the obtained results have been compared with our benchmark curve shown on The behavior of λ A pO (E) and λ B pO (E) is subject to stringent constraints defined by two key factors.Firstly, the latest FAMU results [12] determine the energy dependence of the muon transfer rate with a relatively small error in the range 0 ≤ E ≤ 0.1 eV.Thus we have to use functions that comply with this result, i.e. they have to be within the light-green area shown on Fig. 2(Left) that represents the uncertainty band of λ S pO (E).Secondly, the timing of Xray events in relation to the moment of injecting muons in the target imposes an additional limitation.On Figure 1, muons are emitted at t = 0 ns while the first distinctive µO X-rays are only detected approximately 22 ns later.This temporal relationship serves as a defining parameter that exclusively establishes the upper energy threshold for the proposed λ X pO (E) functions.
On Fig. 2(Right) are given the time distributions of the characteristic oxygen X-rays obtained with the three distinct λ pO given on the left picture with the same colors.The solid black line is our benchmark shown on Fig. 1.At considerably large time all distribution curves become close together with similar slopes that correspond to the transfer of thermalized muonic hydrogen atoms.
We have performed numerical simulations of the muon transfer from hydrogen to oxygen under five more combinations of initial conditions.Our code have been run with different values for the oxygen concentration and gas pressure and fixed temperature T = 300 K.The initial parameters and the corresponding results are summarized in Table 1.In columns 4 to 13 are given the coefficients (λ 1 and λ 2 ) corresponding to the slopes in the characteristic X-ray time distribution.The transition rates' coefficients without a superscript (λ 1 and λ 2 ) are derived by Werthmuller et al. [7] for concentration, pressure, and temperature given in the first three columns of the table.While λ 1 and λ 2 with superscripts are extracted from our simulation results for different λ pO : those indexed by "S", "A", and "B" are obtained with the muon transfer rates λ S pO (E) and the ones given by Eqs. ( 5) and ( 6) respectively.
Table 1.Numerical values of λ 1 and λ 2 , representing the two characteristic slopes of the oxygen X-ray time spectra curve.Each row corresponds to a certain combination of oxygen concentration (as a percentage of the total particles) and pressure of the gas mixture.The transition rates's coefficients without a superscript are taken from [7].The others are obtained with our simulation code for different λ pO : those with indices "S", "A" and "B" correspond to λ S pO (E) and the ones given by Eqs. ( 5) and ( 6) respectively.

Discussion
The monitored quantities in the experiments under consideration are the time spectra of X-ray radiation emitted by muonic oxygen atoms immediately after their formation.As previously discussed in Sec.1, the observed behavior of these spectra characterized by the successive exponential descents stems from some peculiarities related to the energy dependence of the muon transfer rate to oxygen.Namely, the initial steeper descent of the function primarily arises from the predominant transfer from epithermal muonic hydrogen, and the subsequent shallower descent is due to the prevalent transfer from the thermally distributed pµ.We have established that the initial steep descent is indeed modified due to the presence of the hot pµ component at the beginning and strongly dependent on the specific choice of the transition rate λ pO (E).Therefore, we have considered various functional dependencies for λ pO (E) that align with the latest developments on the subject [12].Three of them -the one λ S pO (E) proposed in [12], and two of its modifications λ A pO (E) and λ B pO (E), are depicted in Fig. 2(Left).We demonstrate on a specific example, see Fig. 2(right) and the accompanying explanations, that it is possible to reproduce only the slope of the slower descent described by λ 1 when using the rate λ S pO (E) from [12].However, the prolonged in energy modifications made to λ S pO (E) as seen in λ A pO (E) and λ B pO (E) show some improvement in matching the initial stage in the benchmark curve.This is obviously due to the fact that more muons are transferred to oxygen at these intermediate energies.The improved performance of these modified transition rates at the initial stage of the process is further validated through a comparison of pairs of slopes obtained with their use and pairs of slopes extracted from experimental data [7], see Table 1.

Conclusions
We have developed a physical model and conducted numerical simulations to describe the X-ray radiation produced by µO atoms in various gaseous mixtures H 2 + O 2 .Given the sensitivity of the initial descent to the muon transfer rate to oxygen Λ pO (E), we tested several modifications that adhere to the stringent restrictions imposed by recent analysis of the experimental data [12].Our computations consistently recover the slope of the subsequent slower descent, but we found it challenging to match the slope of the initial steeper descent with the experimental data.This suggests that our model, in combination with the additional approximations, accurately describes the µO X-ray time spectra only after time interval of the order of a few tens of nanoseconds.
Enhancing the considered model and obtaining more precise knowledge of the physical quantities is essential to comprehensively describe the experiments.Furthermore, additional experimental data is needed to assess the reliability of the presented approach and to describe better the behavior of the muon transfer rate function λ pO (E) for muonic hydrogen energies greater than 0.1 eV.The results presented here are of potential use in a number of planned and ongoing muonic experiments.