The effect of areal density asymmetries on scattered neutron spectra in ICF implosions

Low mode drive asymmetries are proposed as a major degradation mechanism in inertial confinement fusion (ICF) implosions. When mode 1 drive asymmetries are present, several experimental signatures are measured; these include hotspot bulk flow velocity,^1–4 anisotropy and differences in DT and DD inferred temperatures,⁵ and asymmetric hotspot shape.⁶ Hydrodynamic simulations suggest that the drive asymmetries, which can be used to explain these observations, will also cause corresponding fuel areal density asymmetries that reduce confinement.^1,7–10 Recent work at the National Ignition Facility (NIF) has added to this evidence by exploring the relationship between hotspot velocity and areal density asymmetries measured using activation diagnostics.¹¹ Activation diagnostics, while sensitive to areal density asymmetries, are also affected by hotspot velocity. Additional spectroscopic measurements are required to approximately remove this degeneracy.¹² Scattered neutron images require sophisticated tomographic reconstruction techniques to extract the fuel density distribution.¹³ Scattered neutron spectroscopy is another method commonly used to measure fuel areal densities.^14,15 It can be used at lower areal densities at which activation measurements are difficult, as is the case at OMEGA. While different areal densities can be measured on different detector lines of sight, the effect of areal density asymmetries on the scattered neutron spectrum has not been thoroughly investigated.

In this paper, we explore the effects of fuel areal density asymmetries on the scattered neutron spectrum. The kinematic scattering relationships suggest that properties of the areal density asymmetries should be contained within a single measured spectrum. A physics-based fitting model for the scattered neutron spectrum in the presence of areal density asymmetries is constructed from the nuclear cross sections. This will allow experimental spectroscopic measurements to leverage their high spectral resolution to directly infer more about the areal density distribution than currently possible. By combining multiple lines of sight, the mode, amplitude, and direction of the asymmetries can be measured.

In Sec. II, we will outline the geometry involved in neutron scattering. This will allow the areal density measured in neutron spectroscopy to be precisely defined. A model for a single scatter neutron spectrum from a DT ICF implosion will be described in Sec. III. This will be used to fit results from a neutron transport calculation in Sec. IV. Finally in Sec. V, how to combine results from multiple lines of sight will be discussed. This will include an error analysis as well as a methodology to evaluate the optimal detector arrangement for general projection measurements, such as the hotspot velocity.² 

The link between areal density asymmetries and scattered neutron spectral shape arises due to kinematic energy-angle relationships. These relationships are one-to-one for elastic processes, such as nD and nT scattering. From classical kinematics, in an elastic collision of a neutron with a stationary ion of mass $Amn$⁠, the incoming neutron energy, $E′$⁠, and outcoming neutron energy, E, can be directly related to the scattering cosine, μ_s,

Conversely, inelastic processes, such as n(D,2n)p and n(T,2n)D, produce a distribution of energies for a given scattering angle and vice versa. A feature that is shared amongst all nuclear interactions of interest is azimuthal symmetry. This allows us to define a “scattering cone” at a given scattering angle (or cosine) around which the differential cross section is only energy dependent.

To understand the geometric interaction of scattering cones and areal density asymmetries, we will consider the simplified case for which the birth spectra of primary fusion neutrons are isotropic. This allows separation of the spatial and spectral dimensions. The spatial effects on the neutron spectrum come through the areal density around the scattering cone. In ICF implosions, a central hotspot acts as an extended source of neutrons. Therefore, we must consider the neutron-averaged line-integrated density, or $⟨ρL⟩$⁠, rather than the more common ρR. Working in the scattering cone geometry, we find $⟨ρL⟩s.c.$ along chords at a fixed scattering angle (with respect to our detector). For simplicity, we will consider the detector to lie along the z-axis, i.e., θ = 0. We then consider the areal density seen by a beam of neutrons born at $r→$⁠, with initial direction, $Ω̂$⁠, before they scatter into the detector line of sight,

To neutron-average, we must sum over the whole reacting volume. The $⟨ρL⟩$ and $⟨ρL⟩s.c.$ are then defined as follows:

where Y_n is the fusion neutron yield, ρ is the mass density, and R_n is the neutron production rate. The geometry of these chord integrals is shown in Fig. 1.

For a point source, $⟨ρL⟩$ will coincide with the ρR. However, this is not a good approximation for a typical ICF implosion. With perfect spherical symmetry, $⟨ρL⟩$ is proportional to ρR with a coefficient, which depends on hydrodynamic profiles.^14,16 In asymmetric conditions, $⟨ρL⟩$ is well-defined and measurable, while ρR is not.

As $⟨ρL⟩$ is defined over the surface of a sphere, it is natural to expand it in terms of spherical harmonics, $Ylm$⁠,

Since the scattering cone integrates over the azimuthal angle, $ϕs$⁠, it is simple to show that

where P_l are the Legendre polynomials and μ_s is the scattering cosine [$cos (θs)$]. Thus, a Legendre expansion in $⟨ρL⟩s.c.$ is general as the spherical harmonic expansion of $⟨ρL⟩$ is complete. The areal density Legendre coefficients, $⟨ρL⟩l$⁠, are directly proportional to the m = 0 spherical harmonic coefficients. Modes with $m≠0$ vanish during azimuthal integration so do not affect $⟨ρL⟩s.c.$⁠. Therefore, scattered spectra are blind to $m≠0$ modes. Note that these spherical harmonic expansions are defined in reference to a coordinate system centered on the detector line of sight. The m modes within the experimental geometry, e.g., those in the equatorial plane at NIF, are measurable with an appropriately orientated detector. The coefficients of order l along one line of sight are a linear combination of the coefficients, also of order l, on a different line of sight. Additionally, due to the extended nature of the hotspot neutron source, $⟨ρL⟩$ is likely to be dominated by low modes.¹⁷ 

For mode l = 1 areal density asymmetries, this is particularly simple,

where $Ω̂det$ is the detector line of sight and $⟨ρL⟩→1$ is a vector that points along the axis of the mode 1 and has a magnitude given by the half the peak to trough areal density asymmetry. The addition of an isotropic and vector component resembles the formula for the DT primary centroid shift.¹⁸ Therefore, a similar 4 detector setup^2–4 is required in order to constrain the mode 1 areal density asymmetry using scattered spectra. For every higher mode, additional $2l+1$ detectors are required in order to back out the physical $⟨ρL⟩lm$ coefficients.

The energy spectrum of the scattered neutrons depends on the spectrum of primary fusion neutrons and the areal density distribution. The nuclear interaction differential cross sections provide the link between these spectral and spatial dependencies. Given an birth neutron energy spectrum, $Qb(E′)$⁠, the scattered neutron spectrum is found as follows:

where $m¯$ is the average ion mass and the double differential cross section is the first term inside the integrand in Eq. (9). The dN/dE can be measured via time-of-flight or magnetic recoil spectrometers. The primary neutron spectrum is also measured and this can be used to separately constrain the birth spectrum; this leaves the areal density distribution as the single unknown. These equations will stand as the basis for a model to describe the singly scattered neutron spectrum in asymmetric ICF implosions.

Using the scattering cone integration, one can construct a full singly scattered neutron spectrum model for asymmetric ICF implosions. The single scatter approximation limits the scope of this model to areal densities $≲$200 mg/cm². In this work, we will focus on this areal density regime. The single scatter approximation also implies that the level of attenuation is low since the total cross section is dominated by scattering rather than absorption. Multiple scattering becomes increasingly important at higher areal densities. One can extend the scattering cone method to multiple scattering, but many of the approximations made are invalid. More accurate models for higher areal densities will be the subject of future work. Working with higher energy scattered neutrons can also be used to reduce the impact of multiple scattering.¹⁷ At these higher energies, both the scattering cross section and the number of multiply scattered neutrons are lower. At 200 mg/cm², the ratio of doubly to singly scattered neutrons is just ∼2% within 10–12 MeV compared to ∼10% below the nT edge, and these rise linearly with increased areal density. Hence, the model limit of 200 mg/cm² ensures suitable accuracy across the whole range of neutron energies considered.

The set of reactions we will consider are D(T,n)α, D(D,n)³He, T(T,2n)α, nT elastic, nD elastic, n(D,2n)p, and n(T,2n)D. Cross-sectional data for these reactions were taken from the ENDF¹⁹ and CENDL²⁰ nuclear data libraries. The Bosch–Hale²¹ DT and DD fusion reactivities and the Appelbe²² calculation of the TT temperature dependent spectral shape were used. The DT and DD primary moments were calculated using relativistically correct expressions derived by Ballabio.²³ For stationary target ions, the differential cross sections can simply be taken from nuclear data libraries, as will be done in this model. The effects of scattering ion velocities can be included for elastic processes using the method outlined in Crilly et al.²⁴ These ion velocity effects are particularly important for the spectral shape of the backscatter edges. It is worth noting that most of the scattering interactions have been well measured,^25,26 but this is not the case for n(T,2n)D. Models disagree on the cross section at 14 MeV (Refs. 19 and 27) and experimental data are limited²⁸—this prompts further theoretical and experimental effort to reduce these uncertainties.

In current indirect drive implosions, remaining ablator areal densities²⁹ can be several hundreds of mg/cm². The higher mass nuclei mean that these ablator areal densities are generally less effective scatterers than equivalent fuel areal densities. However, scattering from the remaining ablator can prove a non-negligible contribution to the scattered neutron spectrum.¹⁷ A minimal extension to the spectral model, described above, could include the cross sections for ablator scattering interactions in a similar manner to those of the fuel. Carbon γ-ray measurements^29–31 could be used to infer the carbon areal density as an input or constraint to this extended model. In this work, we will focus on the dominant fuel contribution to the scattered neutron spectrum.

First, we will look at the form of the scattering cone neutron spectrum as given in Eq. (9). The double differential cross section term can be calculated independently of the areal density distribution. The birth neutron spectrum of the primary fusion reactions can be calculated using the models listed above. For a stationary uniform plasma, the spectral shapes and yields of the primaries are only dependent on a single parameter, the ion temperature. More complex source conditions will alter the birth spectrum, but the model will aim to describe the birth spectra using only a single temperature, taken as the burn averaged DT ion temperature, $⟨Ti⟩DT$⁠. Figure 2 shows the differential cross section terms for a 4 keV ion temperature case.

While limiting the energy range of a spectral measurement might be dominated by a particular angular range, there are still contributions from other scattering angles. This becomes increasingly important at lower energies where many scattering sources have similar amplitude. Areal density measurements found using a down-scatter ratio (DSR) technique will, therefore, find the following average measurement of the areal density distribution:

where the notation on the double differential cross sections indicates that they have been integrated over the DT peak only. Note we convert the DSR to a $⟨ρL⟩$⁠, not a ρR. The conversion to ρR is sensitive to hydrodynamic profiles and, therefore, varies between implosions. Estimates for this geometry effect exist in the literature.^14,16Figure 3 shows the angular weighting function for various DSR ranges. Inferring areal density asymmetries using DSR techniques from multiple lines of sight will, therefore, have to account for this averaging effect as well as the line of sight projection effect of the $⟨ρL⟩$ distribution. Therefore, it may be advantageous to perform a fit to scattered spectrum over a large energy range to infer the $⟨ρL⟩l$ coefficients directly. Fitting the spectrum also has the additional advantage that a single line of sight measurement can be used to quantify areal density asymmetries, whereas DSRs require multiple lines of sight to do this. This will be investigated further in the following sections of this work. Finally, we will consider the effect of cross section uncertainty on areal density inference. This can be estimated by calculating the fractional uncertainty in the average scattering cross section, $⟨σ⟩$⁠. The uncertainty on the inelastic processes is likely the dominant source of error, leading to

where Δ indicates a fractional error and subscript i indicates the inelastic component of the cross section. The ratio $⟨σi⟩/⟨σ⟩$ quantifies the fraction of scattered neutrons within the DSR range, which have undergone inelastic processes. This fraction is 2%, 11%, and 38% for 10–12, 6–12, and 3.5–4 MeV, respectively. We can, therefore, see that 3.5–4 MeV is the most sensitive to inelastic processes. For this range, fractional errors in inelastic cross sections of ∼10% (Ref. 26) produce ∼4% error in areal density inference.

Given the scattering cone neutron spectrum, one can compute the resultant total neutron spectrum by integration over scattering cosine [as shown in Eq. (10)]. Figure 4 shows the spectrum for various l = 1 and l = 2 areal density asymmetry amplitudes. Deviations from the symmetric case are notable and different behavior is observed for l = 1 and l = 2. A simple interpretation based on the elastic energy-angle relationship in Eq. (1) fails to predict some of the behavior. For example, for l = 1, one might expect a positive mode 1 to be lower amplitude compared to the symmetric case in the backscatter regions. This is not the case due to the (n,2n) reactions, which have a very forwarded peaked differential cross section (c.f., Fig. 2). Their contribution to the elastic backscatter regions samples areal density from the forward direction—opposite to that sampled by the elastic interactions (c.f., the 3.5–4 MeV result in Fig. 3). A diagrammatic representation of this effect is given in Fig. 5. This demonstrates the need for the model outlined in this section in order to understand the effects of areal density asymmetries on neutron spectra.

The model laid out in Sec. III will be compared to results from a more complete neutron transport description of an ICF implosion. The inverse ray trace method described in Crilly et al.¹⁷ allows the calculation of neutron spectra given 3D grids of hydrodynamic conditions. As a suitable test problem, we consider an ice block model for a mode 1 asymmetry where a stationary spherical hotspot is displaced within a dense cold fuel layer, as in Fig. 5. The hotspot was chosen to have a parabolic temperature profile to test the spectral model's reliance on a single average temperature. The ice block configuration was then chosen such that $⟨Ti⟩DT$ = 3.9 keV, $⟨ρL⟩0$ = 150 mg/cm², and $⟨ρL⟩1/⟨ρL⟩0$ = 30%. The detector line of sight was then rotated with respect to the mode 1 axis to explore the relationship derived in Eq. (8). Neutron spectra were calculated for these different detectors using the inverse ray trace method. These spectra were fitted using the model outlined in Sec. III to obtain $⟨ρL⟩0$ and $⟨ρL⟩1$ for each detector.

The results of this analysis are shown in Fig. 6. Every fit retrieved a value of $⟨ρL⟩0$ within 4 mg/cm² of the true value of 150 mg/cm². This error originates in the approximations of the model, namely, single temperature primary spectra and no attenuation. For $⟨ρL⟩1$⁠, it was found that a fit to the spectrum between 3 and 6 MeV produced the best results. This range of energies is sensitive to all scattering angles (c.f., Fig. 2) but also narrow enough that the effect of differential attenuation of the scattered neutrons is small. Narrow high energy ranges are more robust to attenuation effects but have lower angular sensitivity.

It is common for the primary spectra to exhibit anisotropy due to bulk fluid flow of the fusing plasma. For a mode 1 drive asymmetry, the hotspot velocity is aligned with the areal density asymmetry.^1,7,8,11 According to the experimental data analysis by Rinderknecht et al.,¹¹ the following linear relationship between mode 1 areal density asymmetry and hotspot velocity matches experimental data:

Thus, the 30% asymmetry presented here corresponds to a hotspot velocity ∼77 km/s. To see the effect of a spectral anisotropy of this magnitude, a uniform hotspot velocity (of 77 km/s) is introduced to the ice block model. For the fitting method, the spectral anisotropy was included through the birth spectrum term, $Qb(E′,μs,vHS)$⁠, in Eq. (9) with the assumption that the mode 1 areal density asymmetry and hotspot flow velocity were aligned. It was found that the same results as shown for the stationary hotspot case were obtained when the hotspot flow velocity was properly accounted for. If the anisotropy of the birth spectrum was not included, then an additional absolute $∼1%$ error was introduced to the inferred $⟨ρL⟩1/⟨ρL⟩0$⁠.

Finally, an areal density inference based on DSRs was performed. For a pure mode 1, the $⟨ρL⟩DSR$ as defined in Eq. (12) can be simplified as

where $⟨μs⟩DSR$ is the average scattering cosine being sampled in the DSR range. The 10–12 MeV (Ref. 14) and 3.5–4 MeV (Ref. 32) down scatter ranges have $⟨μs⟩DSR$ values of 0.70 and −0.22, respectively (based on the results of Fig. 3). Between 3.5 and 4 MeV, if the nD and (n,2n) backgrounds under the nT backscatter edge can be successfully subtracted, this range becomes much more backward focused with a $⟨μs⟩DSR$ = −0.89. However, the subtraction method must be suitably general as the scattering backgrounds are themselves sensitive to areal density asymmetries.

Unlike the fitting method, values of $⟨ρL⟩0$ and $⟨ρL⟩1$ cannot be inferred separately on a single line of sight when using DSRs. Figure 7 shows the inferred $⟨ρL⟩DSR$ for the ice block spectra presented in Fig. 6. Without taking into account the average scattering cosine, $⟨μs⟩DSR$⁠, the DSRs as a function of viewing angle are a poor representation of the underlying areal density asymmetry. This is due to the fact that the areal density probed by the DSR measurement is not along the detector line of sight. Because of this, interpreting the $⟨ρL⟩$ through measured DSRs requires additional care. If the underlying areal density distribution is assumed to be mode 1, a linear fit can be performed once the measurements have been projected using the appropriate angle, $⟨μs⟩DSR$⁠. For higher modes, Eq. (15) must be expanded to higher order introducing higher moments of μ_s, complicating further the relationship between DSR and underlying physical areal density asymmetries. The results of a mode 1 fit to the data in Fig. 7 are values of $⟨ρL⟩0$ = 148 mg/cm² and $⟨ρL⟩1/⟨ρL⟩0$ = 30%, within a few percent of the true values. Here, we have used the fact we know the direction of the mode 1 a priori to analyze the data, and in the next section, we will consider the general case of a completely unknown $⟨ρL⟩→1$⁠.

As shown in Sec. II, a single measurement of the scattered neutron spectrum contains incomplete information about the physical $⟨ρL⟩(θ,ϕ)$ distribution. Multiple measurements are required in order to infer the underlying spherical harmonic moments. The inference of hotspot flows^2,4 and apparent ion temperature anisotropy³³ from primary neutron spectra requires similar analysis.

In this section, we will consider a more general approach to inferring the mode 1 areal density asymmetry vector, $⟨ρL⟩→1$⁠, and the error analysis. We will approach this as a general projection measurement first. We will consider the case where we have n_d detectors and $nd≥np$⁠, where n_p is the number of physical parameters we are inferring. The measurements can then be related to the underlying physical parameters with the following system of linear equations:

where $m→$ is the length n_d vector of measurements, A is the projection $nd×np$ matrix, and $p→$ is the length n_p vector of physical parameters. The normal equation can be used to find the least squares distance between measurement and inferred physical parameters. This leads to the definition of the pseudo-inverse, $A+$⁠,

We are also interested in how error propagates from measurement to inferred values. By considering the covariance matrix of $p→$⁠, which is denoted Σ_p, one can find a detector arrangement that minimizes this error,

In order to proceed with the analysis, we will assume that the detectors are independent and have the same error, σ_m. Then, the covariance matrix of the measurements $Σm=σm2I$⁠. Using the properties of symmetric matrices for $ATA$⁠, this simplifies Eq. (21),

The determinant of Σ_p is a measure of the total error in $p→$⁠. Therefore, the determinant of $ATA$ must be maximized to minimize this error. This defines a method to evaluate error propagation and optimal detector setup for measurements which can be written as linear projection problems.^2,33

We will investigate the optimal detector arrangement for multiple DSR measurements of a mode 1 areal density asymmetry as a worked example of this methodology. The projection matrix, A, for the DSR method can be found by combining Eqs. (8) and (15),

where the bracketed superscripts denote the detector index. If one sets $⟨μs⟩DSR(i)=1$⁠, then this projection matrix also applies for the measurement of vector hotspot velocity and isotropic Gamow shifts for primary DT and DD neutron spectra.^2,4 Evaluating $ATA$ gives the following 4 × 4 matrix:

In order for the inferred quantities to have no covariances, then this matrix must be diagonal. Setting the off diagonal components to zero will set requirements on the detector lines of sight. As discussed earlier, the minimal error detector arrangement is the one that maximizes the determinant,

where Σ_u is the 3 × 3 covariance matrix of detector projection vectors, $u→$⁠. We will investigate the case where every detector uses the same DSR range such that $⟨μs⟩DSR(i)=⟨μs⟩DSR$⁠. The angular factors can then be factored out,

leaving only dependence on the (unit length) line of sight vectors, $Ω̂(i)$⁠. Therefore, we can now use statistical results for the analysis of spherical data. The normalized orientation matrix, $T=⟨Ω̂Ω̂T⟩$⁠, has positive eigenvalues, τ_i, that sum to unity.³⁴ Its determinant is then maximized when all these eigenvalues are equal. This describes a uniform isotropic distribution for which

leading to

These requirements also lead to a diagonal $ATA$ with elements (n_d,n_d/3,n_d/3,n_d/3). For n_d = 4, the optimum is found for a tetrahedral arrangement. The solution to this problem is analogous to arranging masses on the unit sphere such that the center of gravity is at the origin and all principal moments of inertia are equal. This second constraint is satisfied if there are two axes with n-fold (n $≥3$⁠) rotational symmetry.³⁵ 

Using the analysis outlined in this section, we can calculate new optimal lines of sight for mode 1 projection problems given current experimental detector arrangements. This will apply for hotspot velocity inference and mode 1 areal density measurements, which use the same DSR range on all detectors. Given the current (n_d = 4) NIF nToF arrangement,² the next optimal line of sight would be along a (θ–$ϕ$⁠) of (77–68). At OMEGA, for the current nToF arrangement⁴ (n_d = 5 unique lines of sight), the next optimal line of sight is along (135–301). The analysis is more involved if multiple different DSR ranges are used. This is because the scattering angle terms cannot be factored out [c.f., Eq. (29)]. The optimization problem must then be solved separately for any given combination of DSR ranges. Previous work in optimizing new lines of sight has been carried out for scattered neutron imaging using the NIF neutron imaging system (NIS).³⁶ The optimal new detector location was taken as the one which would maximize coverage of the sky. Using the analysis in this section produces a different optimal line of sight to that which maximizes coverage. This highlights the different considerations required for different diagnostic techniques.

Scattered neutron spectroscopy contains valuable information on areal density asymmetries in ICF implosions. Areal density asymmetries lead to a reduction in confinement, and so accurate measurement of such asymmetries is key in understanding current perturbation sources and improving experimental performance.

In this work, we developed a model for the singly scattered neutron spectrum, which can fully account for the effects of an asymmetric areal density distribution. This required a more general description of the scattering geometry, which defined the measurable neutron averaged areal density, $⟨ρL⟩$⁠. The model outlined in this work will allow fits to spectroscopic data to infer the amplitude and mode number of the areal density asymmetries. The model was tested on neutron transport results for a simple mode 1 test case. Fitting the spectra from multiple lines of sight reproduced the theoretical results, even in the presence of non-uniform temperature and anisotropic birth spectra. This analysis provides the most direct measurement of areal density with little sensitivity to other factors. However, going to higher areal densities presents challenges due to the increasing level of multiple scattering.

Due to the azimuthal symmetry of scattering, a single spectroscopic line of sight is insensitive to variation in $⟨ρL⟩$ along the azimuthal direction. Therefore, combining measurements from a sufficient number of detectors is required to obtain a complete picture of the $⟨ρL⟩$ distribution across the whole sphere. The error propagation for combining the measurements from multiple lines of sight was explored. This led to a methodology for calculating the optimal detector arrangement, which minimized the error in the inferred physical quantities. This error analysis is also relevant to the measurement of primary DT and DD spectra, which involve similar line of sight dependencies.^2,4

This work was supported by the Lawrence Livermore National Laboratory through the Academic Partnership Program.

The data that support the findings of this study are available from the corresponding author upon reasonable request.