Monte Carlo calculation of the average neutron depolarization for the NPDGamma experiment

TL;DR

This paper uses Monte Carlo simulations to calculate neutron depolarization effects in the NPDGamma experiment, improving the accuracy of polarization measurements by accounting for spin-flip scattering in liquid hydrogen.

Contribution

It develops methods for tracking neutron spin in MCNPX to accurately model depolarization effects in the experiment.

Findings

Monte Carlo methods quantify neutron depolarization due to spin-flip scattering.

The approach provides a correction factor for neutron polarization measurements.

Results enhance the precision of asymmetry measurements in the NPDGamma experiment.

Abstract

The NPDGamma experiment measures the asymmetry in $\gamma$-ray emission in the capture of polarized neutrons on liquid parahydrogen. The beam polarization is measured using $^3$He spin analysis, but this measurement does not account for the contribution of depolarization from spin-flip scattering primarily due to orthohydrogen in the bulk liquid. This is a systematic effect that dilutes the experimental asymmetry and is modeled using Monte Carlo. Methods for tracking neutron spin in MCNPX were developed in order to calculate the average neutron polarization upon capture for use as a multiplicative correction to the measured beam polarization for the NPDGamma experiment.

Tables

Table 1. Table 1: Molecular hydrogen scattering interactions with spin-flip probabilities as implemented in our MCNPX model using ENDF/B-VII cross sections  [ 15 ] . The uncertainty for H → ortho {}_{\textrm{ortho}}\to H ortho ortho {}_{\textrm{ortho}} is due to the neutron energy spectrum.

Transition	$P_{\text{flip}}$
H${}_{\text{ortho}}{}^{}\to$H${}_{\text{ortho}}$	0.659(1)
H${}_{\text{ortho}}{}^{}\to$H${}_{\text{para}}$	1.0
H${}_{\text{para}}{}^{}\to$H${}_{\text{para}}$	0.0
H${}_{\text{para}}{}^{}\to$H${}_{\text{ortho}}$	1.0

Table 2. Table 2: Isotopes with spin flip probabilities determined using Sears cross sections  [ 26 ] as implemented in the model along with uncertainties due to averaging over the neutron energy spectrum. Zinc is applied over all natural isotopes rather than specific isotopes.

Isotope	$P_{\text{flip}}$
H${}_{\text{bound}}$	0.659(1)
6Li	0.0002(1)
7Li	0.43(5)
14N	0.022(3)
27Al	0.0033(3)
35Cl	0.0326(7)
37Cl	0.00034(4)
55Mn	0.012(5)
63Cu	0.00041(9)
65Cu	0.016(1)
Zn	0.0011(5)

Table 3. Table 3: Polarization correction for the apparatus aluminum at the lower and upper limit orthohydrogen concentrations for each ring and averaged over the detector array.

Ring	${\Delta }_{\text{dep}}^r\left(x_{\text{o,lower}}\right)$	${\Delta }_{\text{dep}}^r\left(x_{\text{o,upper}}\right)$
1	0.9961(3)	0.9799(8)
2	0.9955(3)	0.9722(8)
3	0.9939(5)	0.9581(13)
4	0.9941(6)	0.9556(19)
Total	0.9953(2)	0.9724(5)

Table 4. Table 4: Ring dependent polarization correction for liquid hydrogen at the lower and upper limit orthohydrogen concentrations for each ring and averaged over the detector array.

Ring	${\Delta }_{\text{dep}}^r\left(x_{\text{o,lower}}\right)$	${\Delta }_{\text{dep}}^r\left(x_{\text{o,upper}}\right)$
1	$0.99275\left(11\right)$	$0.94712\left(28\right)$
2	$0.99233\left(9\right)$	$0.94304\left(24\right)$
3	$0.99144\left(10\right)$	$0.93735\left(25\right)$
4	$0.99083\left(13\right)$	$0.93417\left(33\right)$
Total	$0.99187\left(21\right)$	$0.94053\left(56\right)$

Equations

$$\vec{\text{n}}+\text{p}\to \text{d}+\gamma .$$

$$\frac{d\sigma }{d\Omega }\propto \frac{1}{4\pi }(1+A_{\gamma }{\hat{\sigma }}_{\text{n}}\cdot {\hat{k}}_{\gamma }).$$

$$A_{\text{raw}}^i=\frac{\sqrt{{\alpha }_i}-1}{\sqrt{{\alpha }_i}+1},$$

$${\alpha }_i=\frac{N_{\uparrow }^l}{N_{\downarrow }^l}\frac{N_{\downarrow }^m}{N_{\uparrow }^m},$$

$$A_{\gamma }=\frac{A_{\text{raw}}^i-A_{\text{app}}^i}{P_{\text{n}}{\Delta }_{\text{dep}}^i(\lambda ){\Delta }_{\text{sf}}},$$

$$E_{\text{rot}}(J)=CJ(J+1),$$

$$\begin{array}[]{rcl}\frac{d^{2}\sigma^{ss^{\prime}}}{d\Omega^{\prime}dE^{\prime}}&=&\sum_{q}P_{q}\sum_{q^{\prime}}\frac{k^{\prime}}{k}\|\langle q^{\prime}|\sum_{i}e^{\imath\vec{K}\cdot\vec{r_{i}}}U_{i}^{ss^{\prime}}|q\rangle\|^{2}\\ &&\times\delta(\Delta E_{\textrm{n}}+\Delta E_{q}),\\ \end{array}$$

$$U_i^{s{s}^{\prime }}=⟨s^{\prime }|a_{\text{coh},i}-p_i\overrightarrow{S_{\perp i}}\cdot \vec{\sigma }+a_{\text{inc},i}\overrightarrow{I_i}\cdot \sigma |s⟩,$$

$$\begin{array}[]{rcl}U^{++}&=&a_{\textrm{coh},i}-p_{i}S_{\perp zi}+a_{\textrm{inc},i}I_{zi},\\ U^{--}&=&a_{\textrm{coh},i}+p_{i}S_{\perp zi}-a_{\textrm{inc},i}I_{zi},\\ U^{+-}&=&-p_{i}(S_{\perp xi}+\imath S_{\perp yi})+a_{\textrm{inc},i}(I_{\perp xi}+\imath I_{\perp yi}),\\ U^{-+}&=&-p_{i}(S_{\perp xi}-\imath S_{\perp yi})+a_{\textrm{inc},i}(I_{\perp xi}-\imath I_{\perp yi}),\\ \end{array}$$

$$\frac{d{\sigma }^{+-}}{d\Omega }=\frac{d{\sigma }^{-+}}{d\Omega }=\frac{2}{3}{a}_{\text{inc}}^{2}I(I+1),$$

$$\frac{d{\sigma }^{++}}{d\Omega }=\frac{d{\sigma }^{--}}{d\Omega }=\frac{1}{3}{a}_{\text{inc}}^{2}I(I+1),$$

$$P_{\text{flip}}=\frac{2}{3}\frac{{\sigma }_{\text{incoh}}}{{\sigma }_{\text{incoh}}+{\sigma }_{\text{coh}}+{\sigma }_{\text{abs}}},$$

$$E_{\textrm{rec}}=\frac{m_{\textrm{n}}}{m_{\textrm{m}}}\big{(}E_{\textrm{n}}-2\sqrt{E_{\textrm{n}}E_{\textrm{n}}^{\prime}}(\hat{k}\cdot\hat{k}^{\prime})+E_{\textrm{n}}^{\prime}\big{)}=\frac{K^{2}\hbar^{2}}{4m_{\textrm{n}}},$$

$$\Delta E_{\text{n}}=E_{\text{n}}+E_{\text{rot}}(j_{\text{i}})-(E_{\text{n}}^{\prime }+E_{\text{vib}}(n)+E_{\text{rot}}(j_{\text{f}})+E_{\text{rec}}),$$

$$\begin{array}[]{ll}\frac{d^{2}\sigma}{d\Omega dE^{\prime}}|_{o\to o}=&2C(a_{\textrm{coh}}^{2}+\frac{2}{3}a_{\textrm{inc}}^{2})(J_{0}(\alpha)^{2}+2J_{2}(\alpha)^{2}),\\ \frac{d^{2}\sigma}{d\Omega dE^{\prime}}|_{o\to p}=&-6Ca_{\textrm{inc}}^{2}J_{1}(\alpha)^{2},\\ \frac{d^{2}\sigma}{d\Omega dE^{\prime}}|_{p\to p}=&2Ca_{\textrm{coh}}^{2}J_{0}(\alpha)^{2},\\ \frac{d^{2}\sigma}{d\Omega dE^{\prime}}|_{p\to o}=&-\frac{2}{3}Ca_{\textrm{inc}}^{2}J_{1}(\alpha)^{2},\\ \end{array}$$

$$R=\frac{{\sigma }_{\text{coh}}-\frac{2}{9}{\sigma }_{\text{incoh}}}{{\sigma }_{\text{coh}}+\frac{2}{3}{\sigma }_{\text{incoh}}}=-\mathrm{0.29.}$$

$$P_{\text{ave}}^l=\frac{E_{\uparrow }^l-E_{\downarrow }^l}{E_{\uparrow }^l+E_{\downarrow }^l},$$

$${\Delta }_{\text{dep,ave}}^r=\frac{1}{2}\left[{\Delta }_{\text{dep}}^r(x_{\text{o,lower}})+{\Delta }_{\text{dep}}^r(x_{\text{o,upper}})\right],$$

$${\sigma }_{\text{dep,ave}}^r=\sqrt{\frac{1}{12}}\left[{\Delta }_{\text{dep}}^r(x_{\text{o,lower}})-{\Delta }_{\text{dep}}^r(x_{\text{o,upper}})\right].$$

Full text

1 Introduction

The NPDGamma experiment [1] uses a polarized neutron beam in order to measure the parity violating asymmetry in the angular distribution of the emitted 2.2 MeV $\gamma$-rays from neutron capture on liquid parahydrogen [2]

[TABLE]

The physics asymmetry, $A_{\gamma}$, is defined by

[TABLE]

where $\hat{\sigma}_{\textrm{n}}$ is the neutron spin, and $\hat{k}_{\gamma}$ is the $\gamma$-ray momentum. Neutrons exit the $10\times 12$ cm neutron guide on the Fundamental Neutron Physics Beamline (FNPB) [3] and enter a supermirror polarizer [4], after which is the NPDGamma apparatus [5, 6]. Polarized neutrons then enter the resonant frequency spin rotator (RFSR) [7] and are rotated according to a pattern ($\uparrow\downarrow\downarrow\uparrow\downarrow\uparrow\uparrow\downarrow$) that cancels beam power fluctuations to second order and the signals from opposite spin states are used to isolate the asymmetry signal. The neutrons are incident on the 16-liter liquid hydrogen target [8] that is surrounded by 48 CsI(Tl) detectors arranged in 4 rings of 12 detectors each [9] that detect the $\gamma$-rays from neutron capture. The average neutron beam polarization is measured with an empty liquid hydrogen target vessel, and a 3He spin filter [10] followed by a 3He transmission monitor [11] downstream of the vessel are used to measure the beam polarization. A uniform 9.5 gauss magnetic field aligned within 3 mrad to the vertical axis ($\hat{z}$) transports polarized neutrons to the target vessel [1].

The raw asymmetry for a given spin sequence is determined using the super ratio method for each pair of detectors, $i$, where a pair of opposed detectors, $l=0\ldots 5$ and $m=l+6$, is separated by $180^{\circ}$,

[TABLE]

where $\alpha_{i}$ is the super ratio given by

[TABLE]

and $N_{\uparrow}$ and $N_{\downarrow}$ are the spin up and spin down signals, respectively. The physics asymmetry [2] for a given detector is given by

[TABLE]

where $A_{\textrm{app}}^{i}$ is a term encompassing apparatus asymmetries including noise and beam fluctuations, $P_{\textrm{n}}$ is the measured neutron beam polarization, $\Delta_{\textrm{sf}}$ is the spin flip efficiency, and $\Delta^{i}_{\textrm{dep}}(\lambda)$ is the detector-dependent neutron depolarization correction. The depolarization correction depends on the geometry of the target, detector array, and neutron beam and on the neutron wavelength. Since the detector rings are azimuthally symmetric about the target, the polarization correction can be expressed as a function of detector rings ($\Delta^{r}_{\textrm{dep}}(\lambda)$) rather than a function of detectors or detector pairs.

There are two spin isomers of the hydrogen molecule: parahydrogen $(J=0,2,4\ldots)$ and orthohydrogen $(J=1,3,5\ldots)$. Room temperature hydrogen gas has approximately a 3:1 ortho-to-para ratio, with the equilibrium parahydrogen concentration increasing as the temperature decreases [12]. The energy separation between the $J=0$ ground state of parahydrogen and the $J=1$ lowest energy orthohydrogen state is 14.7 meV and higher order rotational energy levels in the hydrogen molecule are given by

[TABLE]

where $C=7.36\textrm{~{}meV}$ is the rotational constant for the hydrogen molecule rigid rotor and $J$ is the angular momentum of the hydrogen molecule. The energy separation between the lowest rotational energy levels in the hydrogen molecule is used as the criterion for discerning scattering event types in MCNPX. The hydrogen target vessel is initially filled with normal hydrogen gas. The hydrogen gas is then liquefied and operated at 15.6 K and the slow natural conversion process is accelerated by circulating the liquid through an ortho-para converter (OPC) [13] containing 150 ml of hydrous iron (III) oxide 30–50 mesh powder[14] as a catalyst. Production data is taken after the hydrogen in the vessel has been condensed and circulated through the OPC over the course of weeks such that the ortho-para ratio has reached a steady state and is assumed to be constant throughout the bulk liquid.

The orthohydrogen scattering cross section for neutrons with energy 1-10 meV (figure 1) is approximately 2 orders of magnitude larger than both the parahydrogen scattering cross section [15] and the hydrogen absorption cross section [16]. The transmission of neutrons through the target as well as the average polarization of neutrons as they propagate through the target is strongly sensitive to the orthohydrogen concentration of the bulk hydrogen. Neutrons can undergo spin flip scattering in the apparatus components and from scattering on orthohydrogen in the bulk liquid hydrogen. Since there is no reliable means of measuring the spin flip scattering rate, the depolarization effect is modeled using Monte Carlo.

There was no in situ method to determine the relative orthohydrogen to parahydrogen concentration for the NPDGamma experiment. The methods to monitor the orthohydrogen concentration available at the time of the construction of the apparatus required feeding electrical signals out of the cryogenic region and hydrogen safety requirements limited the number of wires in the vessel. However, conversion towards a steady state condition with a time constant on the order of days is observed via measurement of the relative neutron transmission and produced a measurement of the parahydrogen scattering cross section [17]. The cross section measurement also yielded a determination of an upper bound (0.15%) on the orthohydrogen concentration in the bulk liquid while the lower limit (0.015%) is given by thermodynamic equilibrium at 15.6 K. The upper and lower bounds on the orthohydrogen concentration are used to as inputs for the calculation of $\Delta^{r}_{\textrm{dep}}$.

The polarization correction, $\Delta^{r}_{\textrm{dep}}$ in equation 5, is calculated with two methods: the first uses a modified version of MCNPX 2.7 [18] and the second uses the Young and Koppel free gas model [19]. These two models will be described next and the two calculation methods are combined to determine the geometry and wavelength dependent polarization for the NPDGamma analysis.

2 Spin-flip scattering formalism

We follow the spin-flip scattering formalism outlined in Moon et al. [20]. The double differential cross section for a scattering process with momentum transfer $\vec{K}=\vec{k}-\vec{k}^{\prime}$ that changes the neutron spin from $s$ to $s^{\prime}$ and the scattering system state from $q$ to $q^{\prime}$ is given by [20]

[TABLE]

where $P_{q}$ is the probability the system is in state $q$, with a sum over the atomic sites $i$, and the $\delta$ function ensures energy conservation for the change in energy of the neutron ($\Delta E_{\textrm{n}}$) and the scattering system ($\Delta E_{q}$). The scattering amplitude is given by [20],

[TABLE]

which has a coherent nuclear scattering term ($a_{\textrm{coh},i}$), a magnetic scattering term ($p_{i}\vec{S_{\perp i}}\cdot\vec{\sigma}$), and an incoherent nuclear scattering term ($a_{\textrm{inc},i}\vec{I_{i}}\cdot{\sigma})$. There are then four combinations of $s$ and $s^{\prime}$ for neutrons with spin in the $\hat{z}$ direction [20],

[TABLE]

where $a_{\textrm{coh}}$ is the bound coherent nuclear scattering length, $p_{i}$ is the magnetic scattering amplitude, $\vec{S}_{\perp}$ is the projection the neutron spin onto the plane perpendicular to the scattering vector, $\vec{I}$ is the nuclear spin, and $a_{\textrm{inc}}$ is the bound incoherent scattering length.

As described by Moon[20], coherent scattering and random isotopic ordering does not contribute to spin-flip scattering. While magnetic scattering contributes to spin-flip scattering in general, it is not relevant for the nonmagnetic hydrogen target in the NPDGamma experiment. Thus, only the incoherent scattering term contributes to spin-flip scattering. The differential scattering cross section per atom becomes simple for the spin-flip and non-spin-flip cases,

[TABLE]

[TABLE]

and which should be compared to equation 15 below. The spin flip terms for a neutron with spin along the $\hat{z}$ axis, the vertical polarization axis in the NPDGamma experiment, contain only $S_{\perp xi}$ and $S_{\perp yi}$ terms, such that only nuclear spins perpendicular to the neutron spin can cause a spin flip and this gives rise to the $\frac{2}{3}$ factor in equation 10. The probability that any given particle interaction leads to a spin flip is given by

[TABLE]

where $\sigma_{\textrm{incoh}}$ is the incoherent scattering cross section, $\sigma_{\textrm{coh}}$ is the coherent scattering cross section, and $\sigma_{\textrm{abs}}$ is the absorption cross section.

2.1 Neutron-hydrogen spin-flip scattering

There are four possibilities of initial and final states for neutron scattering from a hydrogen molecule (see table 1). H${}_{\textrm{para}}\to$H${}_{\textrm{para}}$ scattering events do not change the spin of the neutron for the neutron energies relevant to the NPDGamma experiment. All H${}_{\textrm{para}}\to$H${}_{\textrm{ortho}}$ and H${}_{\textrm{ortho}}\to$H${}_{\textrm{para}}$ scattering events are inelastic in nature and necessarily cause a spin flip. A neutron with energy well below 14.7 meV does not have sufficient energy to undergo an inelastic spin-flip scatter from a parahydrogen molecule such that all interactions with parahydrogen do not change the neutron spin, unless the neutron has been upscattered above the 14.7 meV threshold. An H${}_{\textrm{ortho}}\to$H${}_{\textrm{ortho}}$ scattering event may cause a spin-flip if the molecular spin is orthogonal to the neutron spin, which follows the spin flip probability in equation 12. The Young-Koppel gas model as well as the Keinert-Sax liquid model were used separately in order to estimate the contributions of these scattering processes to the depolarization of neutrons for the NPDGamma experiment.

2.1.1 Young-Koppel free gas model treatment

Young and Koppel [19] give analytical expressions for the differential scattering cross sections for orthohydrogen and parahydrogen gases with the inclusion of spin correlations, vibrations, and rotations. The momentum transfer is given by $\vec{K}=\vec{k}-\vec{k}^{\prime}$, and the recoil energy of the molecule is given by,

[TABLE]

where $E_{\textrm{n}}$ and $E_{\textrm{n}}^{\prime}$ are the initial and final neutron energies, $\hat{k}$ and $\hat{k}^{\prime}$ are the initial and final neutron velocity vectors, $m_{\textrm{n}}$ is the neutron mass, and $m_{\textrm{m}}$ is the mass of the molecule. The energy balance for the scattering event is given by,

[TABLE]

where $E_{\textrm{vib}}(n)$ is the vibrational energy of the molecule in the $n$th mode. At low temperatures, the molecular vibrational states can be ignored because molecules are in the vibrational ground state and cold neutrons do not have sufficient energy to excite higher order states. The Y-K model ignores intermolecular interactions and the cross section is given by the molecular “self”-term. The simplified cross sections (for $J=0,1$ )for the spin isomer transitions in hydrogen can be expressed analytically as follows from Zoppi [21][22],

[TABLE]

where C=\sqrt{\frac{E_{\textrm{n}}^{\prime}}{\pi E_{\textrm{n}}E_{\textrm{rec}}k_{B}T}}\exp{\big{[}-\frac{(\Delta E_{\textrm{n}})^{2}}{4E_{\textrm{rec}}k_{B}T}\big{]}} is a factor common to each, $a_{\textrm{coh}}$ and $a_{\textrm{inc}}$ are the bound coherent and incoherent nuclear scattering lengths, $J_{l}(\alpha)$ is the spherical Bessel function of order $l$, $d$ is the equilibrium separation between the atoms in the molecule, and $\alpha=\frac{1}{2}Kd$. The $3J$ symbol and integral terms have been evaluated in order to supply concise expressions for the $J=0,1$ transitions, which are the only terms that contribute for cold neutrons and liquid hydrogen in the vicinity of 15 K.

For the Y-K gas model, the probability of a spin flip interaction for H${}_{\textrm{para}}\to$H${}_{\textrm{ortho}}$ and H${}_{\textrm{ortho}}\to$H${}_{\textrm{para}}$ transitions is unity. Cold neutrons do not have sufficient energy for spin flip interactions from H${}_{\textrm{para}}\to$H${}_{\textrm{para}}$ scattering events. The depolarization factor for H${}_{\textrm{ortho}}\to$H${}_{\textrm{ortho}}$ transitions and 5 meV neutrons is determined from the cross section [21] in equation 15 and is given by

[TABLE]

The Y-K model only includes hydrogen transitions and does not incorporate spin-flip scattering from apparatus components.

2.1.2 Keinert-Sax kernel model treatment

The basis for hydrogen scattering kernels is the Young and Koppel model, which accurately describes rotational, vibrational, spin correlations, and free translations for gaseous hydrogen. However, translational modes in liquid hydrogen are not free and the Young and Koppel model is not sufficient to describe the scattering of cold neutrons from liquid hydrogen. The Keinert and Sax model [23] improves the description of the translational modes compared to the Young and Koppel model in order to more accurately describe scattering from liquid hydrogen. The scattering kernels compile for MCNPX from MacFarlane [24] are based on the hindered translation model from Keinert and Sax and also incorporate interference between scattered waves from different molecules. The hydrogen molecules are represented as if there are solid-like clusters of approximately 20 molecules that diffuse through the liquid. The MacFarlane model better describes the steep drop in the parahydrogen cross section above 20 meV due to intermolecular spin correlation and the drop around 3-4 meV due to intramolecular interference [24] when compared to the Young and Koppel gas model.

Despite the fact that MCNPX does not track spin or determine when a spin exchange scattering event takes place, the kinematics of a neutron interacting with orthohydrogen or parahydrogen carries the signature of spin exchange and this is used to determine whether a spin-flip interaction has taken place. The first excited state of the hydrogen molecule with $J=1$ is 14.7 meV above the ground state. Therefore, the signature of the H${}_{\textrm{ortho}}\to$H${}_{\textrm{para}}$ transition is a substantial increase in the energy of the neutron with an expected value near 14.7 meV, while the H${}_{\textrm{para}}\to$H${}_{\textrm{ortho}}$ transition signature is an energy decrease of the same magnitude.

A broad and flat input neutron source that does not match the FNPB energy spectrum was used to generate a map of the kinematics in $\Delta E_{\textrm{n}}=E_{\textrm{n}}^{\prime}-E_{\textrm{n}}$ and $E_{\textrm{rec}}=\frac{K^{2}\hbar^{2}}{4m_{\textrm{n}}}$ space for neutron scattering from parahydrogen and orthohydrogen. This was done for illustration and to determine whether the kinematics could be used as a marker for spin-flip scattering interactions by isolating a region around $\Delta E_{\textrm{n}}=\pm 15$ meV. Figures 2 and 3, parahydrogen and orthohydrogen respectively, illustrate two dimensional histograms of the kinematics as calculated in MCNPX 2.7 using ENDF-VII [15] cross sections.

The regions near $\pm$14.7 meV for each were fit to a 2 dimensional Gaussian in $\Delta E_{\textrm{n}}$ and $E_{\textrm{rec}}$. A 3$\sigma$ cutoff from these fits was used in the depolarization calculation such that any scattering event with $\Delta E_{\textrm{n}}$ and $E_{\textrm{rec}}$ inside 3$\sigma$ is labeled as a H${}_{\textrm{ortho}}\to$H${}_{\textrm{para}}$ or H${}_{\textrm{para}}\to$H${}_{\textrm{ortho}}$ and the neutron spin is flipped with unit probability. If the event is outside of this ellipse, the interacting species and table 1 are used to determine whether the interaction is H${}_{\textrm{para}}\to$H${}_{\textrm{para}}$ or H${}_{\textrm{ortho}}\to$H${}_{\textrm{ortho}}$.

3 Monte Carlo calculations of spin-flip scattering

Two separate Monte Carlo models have been developed in order to estimate the depolarization correction for the NPDGamma experiment. The first model uses the Y-K gas cross sections in order to estimate the contribution from orthohydrogen and parahydrogen. For the second, MCNPX 2.7 was modified in order to provide neutron spin state tracking throughout particle histories such that the neutron spin could be flipped probabilistically and the signal contributions from each spin state recorded for each detector.

3.1 Gaseous cross-section based simulation

While the Y-K model does not adequately describe the translational modes in the liquid phase [25], the spin-correlation assumptions in the Y-K model are valid. The Y-K cross sections in equation 15 are used to directly model the average polarization of a neutron beam as it passes through a parahydrogen/orthohydrogen gaseous mixture.

The source neutron beam in this model is initialized with 100% polarization and is monochromatic with an energy of 5 meV. The mean-free-path of the neutron is calculated and the neutron is propagated to a new location inside the hydrogen volume and the event type is chosen probabilistically at that location according to the value of $R$ in equation 16. If the event is a capture event, the neutron history is terminated, the average polarization is recorded, and a new neutron track is generated at the source. If the event is a scattering event, then the neutron spin is flipped as appropriate, the kinematics of the collision are determined, and the calculation continues with a new determination of the mean-free-path and calculation continues until the neutron track is terminated when it has either escaped the liquid hydrogen volume or has captured on hydrogen.

3.2 MCNPX simulation

MCNPX does not have any internal way of defining and tracking neutron spin. This capability was implemented by adding a spin state tag for neutrons. As neutrons propagate in MCNPX, all information regarding the scattering material, interaction type, outgoing momentum, and possible newly created particles is determined internally. An additional step was incorporated to probabilistically change the neutron spin tag from spin up to spin down at each scattering event. Hydrogen in plastics and rubber components in other sections of the apparatus is denoted as $H_{\textrm{bound}}$.

The MCNPX calculation incorporates implicit capture, such that capture events produce secondary particles and reduce the weight of that neutron and continue propagating it. The materials in the MCNPX model of the NPDGamma experiment contain complicated alloys rather than mono-isotopic materials, and MCNPX handles the bookkeeping for the selection of the isotope involved in a scattering event. The isotopes relevant to the NPDGamma experiment with non-zero nuclear spin and an incoherent cross section from Sears [26] are shown in table 2. In the case of H${}_{\textrm{bound}}$, the spin flip probability corresponds to the probability of a non-absorption event multiplied by a factor of $\frac{2}{3}$ due to the lack of separate coherent and incoherent cross sections in MCNPX. The material definitions in the simulation used ENDF/B-VII cross sections [15] for the MCNPX calculations.

The source neutron beam begins with 100% polarization and has the energy spectrum shown in figure 4. This spectrum is representative of the spectrum determined by the beamline chopper settings [3] without frame overlap. As neutrons propagate in the model and create $\gamma$-rays, the $\gamma$-rays inherit the tag indicating the spin of the neutron as the $\gamma$-ray was created. Additionally, Compton electrons inherit this same tag as they are created through $\gamma$-ray scattering processes.

The average polarization is determined by binning the energy deposition in each detector by the spin state tag, and a user supplied subroutine performs these tallies. The geometrically averaged polarization as measured for detector, $l=0...47$, is given by,

[TABLE]

where $E^{l}_{\uparrow}$ and $E^{l}_{\downarrow}$ are the spin up and spin down energy deposition tallies respectively. The detector rings, $r$, are numbered from 1 to 4 with ring 1 being the most upstream ring. Each detector ring should observe approximately the same average polarization and downstream rings should observe a slightly lower average polarization because further downstream detector rings measure signals from neutrons that have penetrated deeper into the liquid and are therefore more likely to have undergone a spin-flip scattering event. Additionally, the penetration depth for neutrons is dependent on the neutron energy, and thus this method of determining the polarization correction captures both the geometric ring-dependence and energy dependence of $\Delta^{r}_{\textrm{dep}}(\lambda)$.

4 Calculated polarization correction

The MCNPX model provides the versatility to calculate the $\Delta^{r}_{\textrm{dep}}$ correction for multiple target materials in addition to hydrogen. In particular, this correction factor was calculated for auxiliary NPDGamma targets such as the chlorine-35 target and the aluminum 6061 target.

4.1 Chlorine-35

The chlorine-35 parity-violating asymmetry is large and well-known, with a world average of $(2.39\pm 0.136)\times 10^{-5}$ [27][28][29]. A Teflon target filled with carbon tetrachloride was used as a systematic check on the apparatus components. Natural chlorine is 75.8% chlorine-35 and 24.2% chlorine-37[30], and chlorine-35 has a significant spin flip probability (see table 2). The polarization correction for the chlorine target is $0.9840\pm 0.0003$. Since the chlorine target is thin rather than comparable to the length of the detector arrays, there is no ring dependence to the polarization correction for chlorine.

4.2 Aluminum-27

The apparatus components, including the hydrogen target vessel and vacuum windows, are constructed of aluminum 6061 and neutron capture on these components can carry a parity violating signal that is a dilution factor in the measurement of the hydrogen asymmetry. The polarization correction for the apparatus aluminum is calculated using the MCNPX model and is dependent on the orthohydrogen concentration in the hydrogen vessel. The results are shown in figure 5 and tabulated for each ring in table 3.

There is a detector ring dependence seen in 5 that is the result of the geometry of the target vessel and the fact that it is distributed over $\approx$ 30 cm. Downstream detectors are more likely than upstream detectors to detect $\gamma$-rays from neutrons that travel further through the liquid hydrogen volume before being captured. Neutrons that travel further through the liquid hydrogen are more likely to have scattered before being captured and therefore are more likely to have experienced a spin-flip scattering event than neutrons that captured a short distance into the target.

4.3 Hydrogen

Similar to the aluminum calculation, the $\Delta^{r}_{\textrm{dep}}$ correction for hydrogen has a detector ring dependence due to being a distributed volume. $\Delta^{r}_{\textrm{dep}}$ is sensitive to the orthohydrogen concentration ($x_{\textrm{o}}$) and calculations were performed at the upper and lower bound orthohydrogen concentrations that were determined along with the measurement of the parahydrogen scattering cross-section [17], which are 0.15% and 0.015%, respectively. Since there was no in situ method of determining the $x_{\textrm{o}}$ for the NPDGamma experiment, the $x_{\textrm{o}}$ distribution is assumed to be flat between the upper and lower bounds, with an uncertainty given by $\frac{1}{\sqrt{12}}(x_{\textrm{o,lower}}-x_{\textrm{o,upper}})$. Similarly, the polarization correction is given by

[TABLE]

with uncertainty equal to

[TABLE]

The neutron energy dependence of $\Delta^{r}_{\textrm{dep}}(\lambda)$ was calculated using a series of calculations with monochromatic neutron source terms and varying the energy. The results of this calculation are shown in figure 6, where the results from the four rings have been averaged together for each neutron energy. The small deviation from 100% polarization at zero orthohydrogen concentration is due to the upscattering of neutrons above the spin-flip threshold as well as the combined kinetic energy of the molecules and neutron above the spin-flip threshold.

Simulations with the MCNPX model were performed at multiple orthohydrogen concentrations within this range and the results are shown in figure 7. The polarization upon capture is a linear function of orthohydrogen concentration within this range. Calculations with monochromatic neutron beams were also performed and show that the polarization correction decreases with increasing energy within the typical NPDGamma neutron spectrum range, which is below the 14.7 meV rotational transition energy. The average contribution to $\Delta^{r}_{\textrm{dep}}$ for hydrogen due to materials other than orthohydrogen and parahydrogen in table 2 is $2\times 10^{-4}$ and is independent of orthohydrogen concentration, since the majority of this contribution is due to spin-flip scattering that takes place upstream of the target vessel in vacuum windows and other apparatus components.

The uncertainty due to the unknown orthohydrogen concentration is much larger than both the uncertainty due to the energy dependence and the ring dependence, thus the polarization correction is averaged over the energy spectrum and over all detector rings to produce a single value that applies to all detectors. The resulting polarization correction from the MCNPX calculation is between 0.99187 and 0.94053 with an average of $0.966\pm 0.015$. By comparison, the Y-K model uses a monochromatic neutron beam at 5 meV and yields an average polarization between $0.97$ and $0.9959$ with an average of $0.983\pm 0.007$. In order to be conservative, we take the maximum and minimum of the results from the Y-K and MCNPX results and assign the bounds on $\Delta_{\textrm{dep}}$ to be between $0.94053$ and $0.9959$ and take the probability distribution between these to be flat. Therefore, the multiplicative polarization correction to the NPDGamma experiment, $\Delta_{\textrm{dep}}$, is determined to be $0.968\pm 0.016$. This uncertainty is comparable to that of other multiplicative corrections for the NPDGamma experiment, namely the geometrical factors calculation [31] and the measured beam polarization [10] and significantly smaller than the statistical uncertainty in the measurement.

5 Acknowledgments

This material is based upon work supported by the U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences under contract number DE-AC05-00OR22725. This work was also supported by DE-FG02-03ER41258.