734
Views
3
CrossRef citations to date
0
Altmetric
Technical Papers

Partial Cross-Section Calculations for PGNAA Based on a Deterministic Neutron Transport Solver

ORCID Icon, &
Pages 1114-1123 | Received 23 Apr 2021, Accepted 06 Dec 2021, Published online: 10 Feb 2022

Abstract

The measurement facility QUANTOM is being developed for the material analysis of radioactive waste packed up in 200-L drums. QUANTOM enables a spatially resolved elemental analysis based on prompt gamma neutron activation analysis. The evaluation of the spatially resolved gamma spectra relies on the calculation of partial (n,γ) cross sections. Hereby, the neutron flux spectrum enters as a parameter, which needs to be simulated in the full three-dimensional geometry of the measurement facility. To ensure that the simulations can be carried out within an acceptable time frame, we use a deterministic neutron transport code specially developed for this purpose based on the SPN approximation of the linear Boltzmann equation. The following question arises: Does the approximation in the neutron transport model still allow a calculation of the partial cross sections at a sufficient level of accuracy. Therefore, in this paper, we study the calculation of partial cross sections in light of the approximation in the neutron transport model in the geometrical setting of the measurement facility. In a simulation study we consider four typical matrix materials and compare cross sections for all elements of the periodic table to reference results obtained by Monte Carlo simulations.

I. INTRODUCTION

For the storage of low-level and medium-level radioactive waste in the German final repository Konrad, a comprehensive documentation of the material composition of the waste is required.Citation1 Strict waste acceptance requirements were defined based on the results of a site-specific safety assessment, which require the quantification of various hazardous substances (e.g., mercury, cadmium, copper, or arsenic).Citation2 The existing documentation is usually not sufficient, especially for so-called historic waste.Citation3

With the measuring system QUANTOM (CitationRef. 4), for the first time, a technology is being developed that is able to determine in a nondestructive manner the material composition of 200-L drums filled with radioactive waste. It utilizes prompt gamma neutron activation analysis (PGNAA) to identify all chemical elements. Neutrons are generated by a deuterium-deuterium (D-D) neutron generator and moderated by surrounding graphite such that a (mostly) thermal neutron flux field is formed inside the waste drum. Atoms inside the measurement facility interact with free neutrons, resulting in the emission of element-specific (even isotope-specific) gamma radiation.Citation5 Spatially resolved gamma spectra are acquired during the measurement process by two symmetrically mounted high-purity germanium detectors. By rotation and vertical translation of the waste drum inside the measurement chamber, the surface is completely and disjunctively scanned by the detectors at discrete positions. A joint evaluation of the spectra is used to reconstruct spatial element distributions within the waste drums. Therefore, the drum content is discretized into smaller subvolumes of the drums, called “sectors” in this paper. The gamma emission at a given characteristic energy emitted from each sector is given by the element mass in this volume, the partial (n,γ) cross section, and the incident neutron flux. By also including photon transport and summing over all disjunct sectors, one obtains a model for the measured signal. The reconstruction of spatially resolved element masses can then be interpreted as an inversion of this model.

The evaluation of the spatially resolved gamma spectra hence relies on the calculation of spatially resolved partial (n,γ) cross sections in all sectors of the waste drum. Hereby, the neutron flux spectrum enters as a parameter, which needs to be simulated in the full three-dimensional (3-D) geometry of the measurement facility.

In contrast to classical instrumental neutron activation analysis at research reactors,Citation5 where the sample mass is small and the neutron flux is approximately homogeneous in the sample volume, the interactions of the sample material with neutrons cannot be neglected or corrected by a simple factor. The complex geometry of the measurement facility suggests the use of Monte Carlo methods, but for the practical application, their performance is insufficient. The inversion method iterates over the material composition, and a full neutron transport simulation is necessary for each measurement position in every iteration, accounting for hours of simulation time per iteration. To ensure that the simulations can be carried out within an acceptable time frame for a large-scale industrial application, we use a deterministic neutron transport code specially developed for this purpose based on the SPN approximation of the linear Boltzmann equation. This code is called SPN Approximation for Remnant Characterization (SPARC). It solves the multigroup SP1 or SP3 approximation to the linear Boltzmann equation using a finite element numerical scheme. For the implementation of the corresponding solver, we rely on the software package FEniCS (CitationRef. 6).

Since the neutron spectra in the measurement facility are not strictly thermal, the diffusive SPN approximations that we used are not perfectly suited to determine neutron transport at the highest order of accuracy. On the other hand, computational efficiency, user-friendly implementation, and the possibility to handle complex 3-D geometries are important factors when choosing the deterministic transport model. This can be achieved using SPARC for the simulation of neutron transport. The following question arises: Does the approximation in the neutron transport model still allow a calculation of the partial cross sections at a sufficient level of accuracy? Therefore, in this paper, we study the applicability of diffusive SPN approximations in partial cross sections for PGNAA in the geometrical setting of the facility. In a simulation study we consider four typical matrix materials and compare cross sections for all elements of the periodic table to reference results obtained by Monte Carlo simulations.

To the authors’ knowledge, the only other approach to use a deterministic multigroup code for PGNAA is described in CitationRef. 7, where the neutron transport simulations are carried out with Attila.Citation8 However, the main focus there lies on a new method for generating multigroup neutron-photon cross sections that is verified for a series of simple benchmark problems. The applicability of an approximate neutron transport model in the calculation of partial (n,γ) cross sections was not a subject of the study.

The structure of this paper is as follows. In Sec. II we define the partial cross sections as the quantity of interest and recall the underlying neutron transport theory. Section III is dedicated to the numerical method, geometrical setup, and implementation details. Computational results regarding the applicability of the approximate neutron transport model in computing the quantity of interest are presented in Sec. IV, for which simulation studies with four typical matrix materials were carried out. This is followed by a conclusion in Sec. V.

II. NEUTRON TRANSPORT AND CROSS SECTIONS

In PGNAA, gamma lines with given energies are detected in order to identify elements. An important quantity for the evaluation of gamma spectra to quantify element masses is partial (n,γ) cross sections. In order to compute these cross sections, it is necessary to determine the neutron flux first. The multigroup form of the neutron transport equation for fixed-source problems is given by

(1) (Ω+Σtg(x,Ω))ϕg(x,Ω)=g=1G4πsgg(x,ΩΩ)ϕg(x,Ω)dΩ+Qg(x,Ω) ,(1)

whereϕg(x,Ω)=groupflux

Σtg(x,Ω)=totalcrosssectionofgroupg

Σsgg (x,ΩΩ )= group-to-group scattering cross section.

The angular discretization is carried out by applying the SPN approximation, originally introduced by GelbardCitation9 in an ad hoc way. A strict derivation of the SP3 equations based on variational analysis can be found in CitationRef. 10. This derivation also yields the appropriate boundary and material-interface conditions. A detailed review of the theory underpinning the SPN equations can be found in CitationRef. 11. The results shown in Sec. IV are obtained using either the SP1 or the SP3 equations. The SP3 equations in the case of an isotropic source and isotropic scattering read as

Δ13Σtgϕ0g(x)+2ϕ2g(x)+(ΣtgΣsgg)ϕ0g(x)=Qg(x)

and

(2) Δ935Σtgϕ2g+Σtgϕ2g=25(ΣtgΣsgg)ϕ0gQg(x) ,(2)

where ϕ0g(x) and ϕ2g(x) are the zeroth and the second Legendre moment of the neutron flux, respectively. The source term Qg(x) is given by

(3) Qg(x)=gg Σsgg ϕ0g(x)+Qextg(x) .(3)

The SP1 equation, also widely known as the neutron diffusion equation, is obtained from EquationEq. (2) by setting the higher-order moments to zero. In order to define the partial (n,γ) cross sections, it is necessary to first introduce two auxiliary quantities: the neutron energy distribution and the isotopic capture cross section. The neutron energy distribution can be defined as

(4) υ(E)=ϕ0(E)0ϕ0(E)dE .(4)

For an incident neutron flux that is not monoenergetic, the isotopic capture cross section is defined as

(5) σi=0σi(E)υ(E)dE,(5)

where σi(E) is the energy-dependent microscopic capture cross section for the given isotope. The partial (n,γ) cross sections are then defined as

(6) σEγ=IEγθσi,(6)

where

σi ==

isotopic capture cross section

θ ==

natural abundance of the given isotope in the element

IEγ ==

so-called intensity of the gamma line.

This can be interpreted as the number of photons with energy Eγ emitted per capture. According to CitationRef. 12, the energy dependence of the intensity in the thermal energy range is negligible. The argumentation found there utilizes previous work on neutron capture models, based on either statistical theoryCitation13 or direct capture.Citation14–16 The partial (n,γ) cross sections can be characterized as the probability of producing a gamma ray with energy Eγ per atom of the given element per neutron capture. Since the gamma lines are characteristic for each isotope, PGNAA can in principle identify single isotopes. In most applications including ours, this is not needed; a reconstruction of the element masses is sufficient. The natural isotope ratio also applies in our case, so the natural abundance θ is incorporated into EquationEq. (6) to make sure that the partial (n,γ) cross sections directly relate the measured values and the masses of the elements. A more detailed introduction of this topic can be found in CitationRef. 17.

III. Methodology and Model Description

In order to solve the SP1 and SP3 equations for this specific application, a special neutron transport code was developed, as mentioned in Sec. I. It is completely based on open source libraries. For the spatial discretization of EquationEq. (2), the finite element method was used since it features great geometrical flexibility. An additional benefit is that the finite element method works very well with Poisson-like equations, in this case the odd-order SPN equations. The FEniCS libraryCitation6,Citation18 with PETSc (CitationRef. 19) as the linear algebra backend is used for the implementation. It offers a large amount of linear solvers and is complemented by hypre,Citation20 a library of high-performance preconditioners. For this work, a gmres solverCitation21 is used in conjunction with an algebraic multigrid preconditioner.Citation22 As a mesh generator, Gmsh (CitationRef. 23) is used. The multigroup cross sections are computed using the OpenMC Monte Carlo code.Citation24 It also yields the neutron flux that is used for benchmarking purposes.

The geometric model used for both SPARC and OpenMC is a simplification of the computer-aided-design representation of the measurement facility. depicts horizontal and vertical cross sections through the finite element mesh used by SPARC. It shows the most important components of the facility, namely, the neutron generator (red), the graphite moderator (black), the drum (yellow), and the outer radiation shielding (gray), as well as some air gaps within the facility. also shows a drum sector (green). As a future improvement it is planned to divide the whole drum into 48 sectors with possibly different filling materials in order to enable inhomogeneous drum fillings. The model has a length of 219 cm, a width of 232 cm, and a height of 273 cm. In reality, the facility is larger; however, the parts outside the radiation shielding have a negligible influence on the neutron flux in the drum. Monte Carlo simulations show that only 0.55% to 0.56% of the neutrons leave the modeled parts of the facility. The shielding consists of highly absorbing borated polyethylene, so possible room shine effects with neutrons that leave the facility and then return are insignificant. The outer parts of the facility can therefore be neglected for this study. The characteristic length of each finite element is not larger than 2.5 cm, resulting in a mesh of the measurement facility that consists of 3.88 mio. finite elements. The exact size of the finite elements varies to adapt better to the geometry, resulting in an unstructured mesh. This ensures that the discontinuous transitions in the cross sections on material interfaces are represented correctly. depicts the picture section in the blue frame in in greater detail, showing how the finite elements vary in size and form to adapt to the material boundaries.

Fig. 1. Slice through the finite element mesh used by SPARC. Depicted are the drum content (yellow), a sector of the drum (green), the neutron generator (red), the moderating graphite (black), the outer shielding (gray) and air gaps within the facility (white). shows the picture section within the blue frame in greater detail.

Fig. 1. Slice through the finite element mesh used by SPARC. Depicted are the drum content (yellow), a sector of the drum (green), the neutron generator (red), the moderating graphite (black), the outer shielding (gray) and air gaps within the facility (white). Figure 2 shows the picture section within the blue frame in greater detail.

Fig. 2. Picture section of in greater detail. The variation of the finite elements in size and form ensures that the mesh preserves the geometric properties of the model. The slice does not divide single finite elements, and since the mesh is unstructured, not all finite elements lie on the same height. The varying brightness of the finite element facets indicates the different angles between the image plane and the facets.

Fig. 2. Picture section of Fig. 1a in greater detail. The variation of the finite elements in size and form ensures that the mesh preserves the geometric properties of the model. The slice does not divide single finite elements, and since the mesh is unstructured, not all finite elements lie on the same height. The varying brightness of the finite element facets indicates the different angles between the image plane and the facets.

Since there is no significant amount of transuranic elements present in the measurement facility, we can assume that there is no fission but that there is a monoenergetic point source with an energy of 2.45 MeV, the energy of the D-D fusion in the neutron generator. The neutrons are emitted isotropically. It is also remarkable that the huge amount of graphite used as a moderator causes strong upscattering even at room temperature; the neutron flux with and without upscattering differs by an order of magnitude for thermal neutron energies.

In order to represent a wide variety of possible drum contents, the simulations were carried out with four different filling materials. The chosen materials were concrete, polyethylene, cast iron, and lead. The exact material composition can be found in CitationRef. 25 under the material numbers 97, 248, 162, and 271. The first two materials have a relatively low density. They were chosen since concrete is a filling that is used very often, while polyethylene represents light mixed wastes quite well. Cast iron and lead represent regular metallic and heavy metallic materials, respectively. In the following discussion, the drum consists mainly of one of these materials. They serve as matrix materials, and small amounts of other materials are embedded in this matrix. It is assumed that the small amounts of other materials have a negligible influence on the neutron flux, so the transport calculations are carried out taking only the matrix material into account.

After the neutron transport simulations are carried out as described in Sec. II, the isotopic and the partial (n,γ) cross sections are calculated for each element using the computed neutron energy distribution densities. In the following discussion, the results for the gamma lines with the highest σEγ for the respective element are shown. It is reasonable to concentrate on these since these are typically the best detectable and therefore the most relevant ones for the practical application.

In this analysis, the isotopic and therefore partial (n,γ) cross sections are not computed for the full energy range but are computed for three distinct energy ranges. The first energy range includes neutrons with energies smaller than 0.25 eV and is named the thermal region in the following discussion. It is the most important energy range since the high microscopic cross sections in this region in combination with a neutron flux that is nearly completely thermalized ensure that these energies have a dominant influence on the final result. The second energy range is named the resonant energy range and includes neutrons with energies between 0.25 eV and 1 MeV. The inclusion of the resonant energy range is necessary in this case since the neutron flux is not fully thermalized. For experiments at research reactors, where the neutron flux is essentially a thermal one,Citation26 resonant neutrons can usually be neglected. Neutrons with energies greater than 1 MeV, here named fast neutrons, have an insignificant influence on the final result. The absorption cross sections in this region are by orders of magnitude smaller than the cross sections in the thermal energy range while the fast neutron flux in the facility is also significantly smaller than the thermal neutron flux. In the following discussion, the integration boundaries in EquationEqs. (4) and (Equation5) are adjusted in order to compute the cross sections for the thermal and resonant energy ranges while the fast energy range is neglected for the reasons mentioned above.

IV. RESULTS AND DISCUSSION

In this section, the results for the partial (n,γ) cross sections will be discussed for the application cases introduced in Sec. III. The necessary neutron energy distributions were computed with the SPN neutron transport code SPARC and with the OpenMC Monte Carlo code. The purpose of this analysis is to show that the SPN approximation, although having a higher model error than Monte Carlo, is nevertheless appropriate for computing σEγ for the evaluation of PGNAA spectra. In the following discussion, the OpenMC solutions serve as reference solutions for SPARC. All so-called errors are to be understood as deviations of SPARC from OpenMC. The results have been divided into two subgroups: materials with lower density (Sec. IV.A) and dense metallic materials (Sec. IV.B). The performance of SPARC and OpenMC is compared in Sec. IV.C.

IV.A. Light Materials

shows the neutron energy density for a drum filled with concrete, computed with both SPARC and OpenMC. In the relative deviation of the SPARC results in comparison to OpenMC is depicted. For thermal energies, the results for both codes are almost identical. The deviation gets larger for higher energies, and it reaches its maximum for fast neutrons where SPARC underestimates the OpenMC flux by over 60%. For the application at hand, the quantity of interest is not the neutron flux but rather the partial (n,γ) cross sections. However, these results show that the SPN approximation is insufficient for a full neutron transport simulation in this application case.

Fig. 3. Neutron energy density υ(E) within a drum filled with concrete computed with OpenMC and SPARC.

Fig. 3. Neutron energy density υ(E) within a drum filled with concrete computed with OpenMC and SPARC.

depicts the relative deviation of the thermal partial (n,γ) cross sections computed with SPARC in comparison to OpenMC as a histogram. The average deviation of the values is only 0.10%. It jumps into the eye that nearly all elements lie in the same bin. This is plausible since the microscopic capture cross sections for nearly all elements show the same functional behavior in this energy range, namely, one divided by the square root of the energy. In the thermal energy range, the partial (n,γ) cross sections computed with SPARC and OpenMC are therefore almost indistinguishable. The results for the resonant partial (n,γ) cross sections are shown in . Because of the higher model error in the neutron flux computed by SPARC, the average deviation for the resonant partial (n,γ) cross sections is also higher, namely, 7.55%. Since every element has a unique resonance structure and therefore unique cross sections in this energy range, the deviations for the different elements are spread out over a broad energy range, which is a strong contrast when compared with the results for thermal energies.

Fig. 4. Relative deviation of partial (n,γ) cross sections for a drum filled with concrete computed with SPARC from OpenMC.

Fig. 4. Relative deviation of partial (n,γ) cross sections for a drum filled with concrete computed with SPARC from OpenMC.

When the integration in EquationEqs. (4) and (Equation5) are, as mentioned before, adjusted to compute the cross sections for the thermal and resonant energy ranges separately, one sees that the resonant energy density υres(E) computed with SPARC overestimates OpenMC for smaller energies and underestimates OpenMC for higher energies. Since smaller energies have a higher impact on the resulting partial cross sections due to higher cross sections and a higher absolute flux, an overestimation for smaller energies leads directly to an overestimation of the quantity of interest and therefore the asymmetry as seen in .

In order to judge the suitability of SPARC to compute partial (n,γ) cross sections for the application at hand, it is necessary to introduce the reaction rate R, given by

(7) R=σEγ,thermΦtherm+σEγ,resΦres,(7)

where

σEγ,therm, σEγ,res ==

thermal and resonant partial (n,γ) cross sections, respectively

Φtherm, Φres ==

thermal and resonant neutron fluxes, respectively.

The fluxes for this study are simulated; they will be determined by measurement devices in the final evaluation process.

The relative deviation of the reaction rates computed with the partial (n,γ) cross sections from SPARC compared to the values computed with OpenMC for a filled drum is depicted in . The average deviation is only 0.62%, with all outliers limited to the low single-digit percentage range. This result shows that despite the model error of SPARC when describing nonthermal neutrons, the (n,γ) reaction rates can be computed with a comparatively small error. Since the deviation between SPARC and OpenMC is much smaller than the expected measurement error of about 5% to 10%, the computation of the partial (n,γ) cross sections can be carried out with the SP1 method as implemented in SPARC as an alternative to the costly Monte Carlo computations.

Fig. 5. Relative deviation of reaction rates computed with the partial (n,γ) cross sections from SPARC in comparison to the values computed with OpenMC for a drum filled with concrete.

Fig. 5. Relative deviation of reaction rates computed with the partial (n,γ) cross sections from SPARC in comparison to the values computed with OpenMC for a drum filled with concrete.

In order to ensure that σEγ not only can be averaged over the complete drum but also can be computed locally, the same computations have been carried out for a single sector, as depicted in . The average errors in the (n,γ) cross sections in the thermal energy region were −0.31% and −6.57% in the resonant energy region. The error in the reaction rates was −0.56%. Thus, one sees a slight underestimation of the reaction rate whereas for the whole drum, one sees a slight overestimation.

Since SPARC also enables SP3 computations, it was investigated whether a further improvement of the results using the higher-order approximation is possible. With the SP3 approximation, the error in the thermal partial (n,γ) cross sections is 0.11% while it is 8.59% for resonant energies. The error in the reaction rates is 0.71%. The error is slightly higher than the error obtained using the SP1 approximation. However, since the deviations between both results are much smaller than the model error, one should see both methods as yielding indifferent results. Hereinafter, only the SP1 approximation will be used due to the better performance.

In this study, the group structure with 44 energy groups from SCALE (CitationRef. 27) is used for the energy discretization. In order to test its suitability, a convergence study has been performed where the SPARC results are compared to the OpenMC results for the respective group structure. The results are shown in . The energy group structures with 22 groups and 88 groups are generated via uniting two adjacent energy groups and dividing one group into two in a way that the logarithmized width of both groups is the same, respectively. The errors in the thermal and resonant cross partial (n,γ) cross sections, ΔσEγ,therm and ΔσEγ,res, decrease with a finer discretization, as well as the error in the reaction rates ΔR. However, while the error is significantly reduced when using 44 groups instead of 22, the improvement by using 88 groups instead of 44 is minor. This shows that 44 groups are a good choice, and further refinement does not reduce the error enough to justify the additional computational effort.

TABLE I Influence of the Energy Discretization

The results for the other light material, polyethylene, are qualitatively identical to the results for concrete. The sole difference is the smaller deviation of the SPARC results in comparison to OpenMC, with an average deviation of 0.03% for σEγ,therm and of 6.90% for σEγ,res. The reaction rates differ by 0.39%.

IV.B. Metallic Materials

The partial (n,γ) cross sections for lead are depicted in . For thermal neutron energies, the average error in the cross sections is 0.87%. Nearly all elements are concentrated in a small area around the average error, yielding a picture similar to concrete, although with a larger error. For resonant neutron energies, the error is 14.87%, with the results for the different elements distributed over a broad area. depicts the error in the reaction rates. The average error is 3.74%, with many elements concentrated in the same bin, similar to concrete. In general, the error is larger than for concrete, but this is expected. Lead has a much higher absorption cross section, leading to a higher error of the neutron fluxes computed with the SP1 approximation and therefore the partial (n,γ) cross sections and reaction rates.

Fig. 6. Relative deviation of partial (n,γ) cross sections for a drum filled with lead computed with SPARC in comparison to OpenMC.

Fig. 6. Relative deviation of partial (n,γ) cross sections for a drum filled with lead computed with SPARC in comparison to OpenMC.

Fig. 7. Relative deviation of reaction rates computed with the partial (n,γ) cross sections from SPARC in comparison to the values computed with OpenMC for a drum filled with lead.

Fig. 7. Relative deviation of reaction rates computed with the partial (n,γ) cross sections from SPARC in comparison to the values computed with OpenMC for a drum filled with lead.

For cast iron, one obtains results that are qualitatively similar to the results for lead. However, the errors are larger. The errors in the partial (n,γ) cross sections for the thermal energy range are 1.70% and 15.40% for the resonant energy range, resulting in an error of 8.35% for the reaction rates.

IV.C. Performance Comparison

The computations were executed on a cluster consisting of AMD Ryzen Threadrippers 2950X, each with 32 cores. Since the OpenMC computations were executed on nine nodes to reduce the total run time while the SPARC computations were executed on one node, the run times were converted to core hours. OpenMC as a Monte Carlo code scales nearly linear; hence, the internode communication has only a minor influence on the run time, and the core hours are well comparable. The run times are shown . The minimum speedup factor of SPARC in comparison to OpenMC is 20 for polyethylene, getting as large as 32 for lead. This shows that SPARC features a significant performance improvement over OpenMC as it is needed for this application.

TABLE II Run Time Comparison in Core Hours

V. CONCLUSION

This paper shows that partial (n,γ) cross sections for PGNAA can be computed with neutron fluxes obtained using the SP1 and SP3 equations. SPARC, an internal neutron transport code based on open source libraries, was developed and validated with OpenMC. Despite the relatively high model error for higher neutron energies, the partial (n,γ) cross sections can be computed with sufficient accuracy. This is caused by the fact that the model error in the thermal energy range is minimal and the reaction rates are strongly dominated by the contribution of thermal neutrons. For lighter filling materials, the most relevant application case, the average error in the reaction rates is less than a single percent. For metallic materials the error is larger, e.g., about 4% or 8% in the test cases discussed above; however, it is still acceptable. In practice, the drums will typically be filled with compacted scrap metal instead of massive metal like in these test cases. Because of technical limitations, for example, in-drum compactors, the resulting density is expected to be smaller than in the test cases shown. It is expected that the error will also be smaller. The SPARC run time in the test cases shown was 20 to 30 times better than OpenMC using the same compute cluster. It is expected that the performance advantage of SPARC will be even larger when applied to the more detailed models of the measurement facility that will be used for future computations.

Acknowledgments

The work of K. K. was performed within the project QUANTOM, funded by the German Federal Ministry of Education and Research under grant number 15S9406B. The work of A. J. was partially performed within the European Regional Development Fund project ZEBRA, funded by the European Union and the federal state of North Rhine-Westphalia under the funding number EFRE-0800566. We would also like to acknowledge the contribution of Jannick Wolters and Christopher Helmes in the development of the SPARC code, which was used to produce the numerical results for this paper.

Disclosure Statement

No potential conflict of interest was reported by the authors.

References