Nuclear Inst. and Methods in Physics Research, A

AspartoftheexpansionoftheDirectAcceleratedGeometryMonteCarlo(DAGMC)toolkittosupportotherMonte Carlo codes, FluDAG (FLUKA integrated with DAGMC) was developed. There has been increasing demand from the high energy physics community regarding Computer Aided Design (CAD) geometry support in Monte Carlo codes. In this paper, the development and validation of FluDAG is discussed and its application to a number of high energy physics experiments is demonstrated, along with its validity relative to native FLUKA calculations.


Introduction
Several modern nuclear, particle physics or other high energy applications such as JET, ITER, LHC or others impose significant challenges on the geometric representation of models used in the Monte Carlo (MC) radiation transport process. Specifically, there usually exists a detailed engineering model typically created with manufacturing, assembly, or other engineering analysis in mind. In tandem there are also usually several physics analysis models created to ensure optimal experimental performance, for example signal-to-noise ratio in the case of particle detectors or nuclear heating performance in the case of tokamaks. It is the case that usually the physics geometric model lags behind the engineering design by several months, for several reasons (1) CAD model preparation is typically done by several CAD analysts working in parallel to make modifications, (2) it takes significant human effort to clean and simplify (defeature) these complex CAD models, (3) conversion of even simplified CAD is error prone, tedious and slow. There is also increasing demand from nuclear analysts to use CAD in their analysis usually driven by the availability of these engineering CAD models and the desire to include as many relevant details in their analysis model so as to better reflect reality.
There have been several approaches taken to solve the problem of the use of CAD models in MC simulations, these can be broken down into two branches of solution; (1) translation approaches and (2) direct use of CAD. A non exhaustive list of translation approaches includes method, and the one employed by DAGMC, is to use high resolution tessellated (faceted) representations of the CAD surfaces. Doing this eliminates the need to find higher order roots during navigation and a number of accelerations are used to offset the high number of facets (triangles) required for high resolution representations. The penalty here is memory usage, but this is easily overcome with the large amounts of memory available on the workstations of today.
DAGMC has been widely used within the fusion neutronics community on several analyses for the international experiment ITER [11][12][13] and on a number of other fusion relevant calculations [14,15]. One of the major benefits of DAGMC is the ability to use it as a common geometric base for multi-physics analysis, for example [16]. DAGMC has now been integrated with the FLUKA radiation transport package [17,18] to form FluDAG, a fully featured version of FLUKA that is able to transport particles in CAD geometries.

DAGMC
DAGMC is an open source collection of C++ based software libraries that allows for efficient ray tracing of CAD based geometries. It is built on the back of a mesh based library MOAB [19] which is used to store the DAGMC geometry. As discussed previously, DAGMC geometry is composed entirely of a large number of triangular facets. It would be computationally prohibitive to perform navigation with a large number of triangles without relying upon algorithms that eliminate the need to test each facet. DAGMC has employed acceleration techniques to make it computationally competitive with native geometry implementations. These acceleration techniques are as follows: 1. Imprinting is an operation performed within the CGM CAD engine where surfaces that are coincident have their curves imprinted into the opposing surface, merging then takes the two topologically equivalent surfaces and unifies them into one definition. This provides an acceleration because determining the next volume entered upon crossing a surface is always a (1) operation. 2. By representing the surfaces of the CAD model as a collection of triangles, each ray query is reduced to a planar ray intersection, which individually is quick, and is ultimately handled using the Plücker [20] ray triangle intersection test. 3. The DAGMC geometry is composed of collections of volumes, with surfaces stored as children sets of volumes, where the facets (triangles) are members of the surface set. For each surface in the problem a tree of oriented bounding boxes is built, these trees are then used to determine which triangle a specific ray hits in ( ( )) time.
Using the accelerations above, DAGMC is able to be competitive with a native Combinatorial Solid Geometry (CSG) MC calculation, typically within a factor of 2-3 for a moderately complex geometry. It should be noted that this difference in timing decreases with increasing geometric complexity. It is argued that despite the slower calculation time, an overall reduction in analysis time is possible due to the speed at which a CAD model can be made ready for analysis. This has also been anecdotally realised in several analyses.

DAGMC workflow
The DAGMC workflow has evolved over the course of its development, however the common factor is the generation of a CAD based solid model. The most common route for DAGMC based analysis is to begin by importing the solid model SpaceClaim where most CAD repairs, cleaning, and defeaturing takes place. The model is then exported to ACIS [21] (*.sat) format (although the workflow does support STEP [22] as well as several other CAD formats) and imported into Cubit [23] or Trelis [24]. A DAGMC-based analysis allows a number of attributes of the geometry to be defined within the geometry file. These characteristics generally relate to the physical properties of the volume, for example its material definition or boundary conditions. The overall workflow is shown in Fig. 1.
The Trelis/Cubit tool is intrinsic in the preparation of DAGMC geometries, as it is used to mark up and produce the faceted geometry. All of the infrastructure required to produce the DAGMC geometries are distributed as plugin objects for Trelis/Cubit. There are a number of standalone command line tools that are run sequentially on a model following faceting. We run make_watertight [25] to seal models to ensure no topological weaknesses exist.

University of wisconsin unified workflow
The DAGMC workflow is supported in several MC codes and as the number of supported codes grew it became difficult to translate already existent MC code metadata into the newly supported code form, for example MCNP material descriptions to Geant4 materials. Thus the University of Wisconsin Unified Workflow (UW) 2 was developed. The workflow encodes material descriptions in code agnostic form storing densities, nuclide mass or atom fractions, comments and other extensible metadata. These material objects are translated either a priori or at runtime depending upon the use case. Particular care is given to nuclides concentration given their particular impact to the neutron based simulations. It is also possible, but in a much more limited fashion to define tallies which can be stored in code agnostic form, and again translated to a native format when possible. At the current time, support for tallies/scoring is much less developed than material handling. This lack of support is due in part to the complexity of possible scoring mechanisms and geometric fidelity of the specific MC code, e.g. MCNP allows surface based scoring methods whereas FLUKA & Geant4 offer only boundary crossing scores.

FluDAG development
Development of FluDAG was initiated by the National Aeronautics and Space Administration (NASA) having interest in a CAD based workflow for space radiation analysis. Given the complex nature of space craft and the ready availability of CAD models there is an obvious need for the ability to use CAD models in MC calculations. A schematic of DAGMC interaction with FLUKA is shown in Fig. 2. The FluDAG layer is used to store all additional state that is required for full functioning of the DAGMC layer, including the direction and position of the particle during the last call and if the last step was directly onto a boundary between regions. This state was introduced to handle some specific issues that were encountered during the development of FluDAG that was not expected based on experience of the integration of DAGMC with MCNP5 [26] or Tripoli4 [27].
One of the robustness features of DAGMC is the RayHistory object. It is an object used to store the facets crossed in the current ray direction, and is up to the code author to reset the state of RayHistory appropriately. When DAGMC functions are called and are passed a RayHistory object the routines ignore hits from facets already in the object. Thus the RayHistory should be reset at any direction change, since logically we should be allowed to hit the same facet again.
Electron (along with some other charged particles) transport is different from neutral particle transport due to the underlying physics. One of the differences of electron transport (when using condensed history approaches) is that electrons are allowed to change direction at boundary crossings, in DAGMC terms it means that one must be careful about resetting the expected state of particles and checking for this condition. Neutral particles cannot change direction on a boundary (unless undergoing optical reflection) and therefore do not encounter this problem. FLUKA also takes a 'sensing' step to determine how far away potential boundaries are, which must also be handled by the FluDAG layer and be understood to not be a true geometric step and thus no ray state should be retained.
DAGMC had never been previously used for codes that supported magnetic field tracking. FLUKA has specific geometry routines for point inclusion and tracking when magnetic field tracking is on. What this again means from the DAGMC perspective is that the DAGMC wrappers must check to ensure that the RayHistory state is not perturbed incorrectly, since underlying these calls are DAGMC rayfire calls. In conventional CSG geometry repeated calls to the ray-object intersection call will always give the same result, as will DAGMC. However, in the case of DAGMC the RayHistory object can be used to store previous surface hits, potentially being used to preclude the same facet being intersected more than once. Consider the case of an energetic electron traversing a region of uniform magnetic field. As the electron traverses the orbit around its guiding centre, computationally this is performed by breaking the orbit into a number of angular sub-steps. At the beginning of each sub-step a rayfire is preformed and the distance to the boundary queried. It is clear to imagine that if there are very many sub-steps in a single orbit, or if there are a few facets in the surfaces of a given volume, then it is very likely the same facet will be intersected. Hence it is critical that the history object is reset if the particle does not take the full physics limited step.
Integration of DAGMC with FLUKA requires that the FluDAG wrapper code honours the FluGG [28] (FLUKA with Geant Geometry) API. By utilising the FluGG API there a number of benefits, (1) consistency of external navigation tools is maintained across different codebases, (2) integration of DAGMC into FLUKA bears significant similarity to other packages reducing maintenance burden, (3) Benchmarking of the API has already been performed by the FLUKA team to have confidence that the API is sufficiently verbose

Validation
There are a number of key tests to stress the geometry responses from FluDAG. In terms of high energy physics problems, electrons and neutrons represent the end state of majority of histories, and represents the majority of the CPU time. Indeed for most problems energy is mostly transferred to electrons. A number of the tests performed originated from the development of FluGG and these are reproduced since these tests proved pathological for the case of electron transport. The key comparison in all of the following tests is that when comparing native FLUKA calculations to the exact FluDAG geometry the results are statistically equivalent, i.e. agreement to within ±3 .

Aluminium and gold
The AlAuAl (aluminium-gold-aluminium) benchmark originated from the development of FluGG. This benchmark is specifically designed to stress electron transport in very thin layers, especially for the previously mentioned electron behaviour at geometric boundaries. The thin layers with very different densities and atomic numbers strongly affect the underlying electron transport, and thereby photon production by bremßtrahlung, pair production etc. The geometry of the setup is shown in Fig. 3.
The geometry is composed of 3 material regions, a thin aluminium layer, a thin gold layer, followed by a aluminium layer. The source is a 1 MeV pencil electron beam pointed in the positive z direction, with particles starting at 10.0 cm downstream from the first layer. The CAD model for the FluDAG was created by exporting the native FLUKA geometry to MCNP format, then using mcnp2cad [29] the CAD model was created. There were 20 batches of calculation, where each batch contained 5.0 × 10 5 primary histories. In each material zone the track length estimate of electron flux is determined, as shown in Fig. 4. The non-physical undulating artifacts present in Fig. 5 are due to the choice of electron physics used, the key take away point is that the features are exactly reproduced in each code across the whole energy range at each stage. This benchmark provided guidance and showed shortcomings in the original development of FluDAG regarding particle direction changes on boundaries, highlighting that despite the simplicity of this benchmark it exposed aberrant electron behaviour, which was subsequently remedied to produce the results shown here.

Magnetic fields, spheres and cylinders
During the development of FluGG a specific test case gave particularly pathological behaviour for electrons crossing boundaries between cells. The test is subsequently known as ''MagnSph'' and the geometry is shown in Fig. 6. This test is particularly pathological for electron transport due to some of the peculiarities of electron transport as discussed previously regarding direction changes on boundaries. The native geometry is composed of cylinders and spheres which are numerically coincident at their extrema, in this case the cylinders have a radius of 0.5 cm and are offset such that adjacent cylinders centres are 1.0 cm from one another.
This benchmark highlights one of the major differences in the geometric capabilities between traditional combinatorial solid geometry and CAD geometry. It is not possible for CAD engines to represent the exact mathematical relationships to the precision possible in CSG. The CSG engines within MC codes have been tuned for accuracy using double precision arithmetic specifically to avoid problems with roundoff. CAD engines are tolerant to a given fixed precision which is many orders of magnitude larger than a double precision number, for example the ACIS CAD kernel treats entities closer than 1.0 × 10 −6 cm as being coincident. It is not possible to resolve numerically touching entities where the numerical precision required is higher than that of the CAD engine, for example ACIS can only distinguish vertices as being distinct entities when they are greater than 1.0×10 −6 cm apart. In this instance we found that the problem as originally defined resulted in several imprint and merge issues. Reducing the radii of the cylinders from 0.5 to 0.49999 and the radii of the spheres from 0.3 to 0.29999 resulted in no issues regarding imprinting or merging. However, changing the radii to 0.499999 and 0.299999 respectively had the same issues as the unmodified model. Thus, the final model used for the benchmark calculations were has radii of cylinders and spheres of 0.49999 and 0.29999 respectively. This of course means that the geometries are not quite identical, but for the purposes of this benchmark are treated as being equivalent. The particle source is a square cross sectioned beam of side 1 cm, with the primary source particle being positrons at 1 GeV, starting at 2, −4, −1 cm directed into the positive z direction as indicated by the face of the brown volume (for colour see the online version) in Fig. 6. There is a uniform magnetic field of 60 T directed along the x direction.
The results of electron transport are shown in Fig. 7, no statistically significant artifacts can be seen between the FLUKA and FluDAG simulations, thus showing that the FluDAG magnetic tracking routines are performing as expected. The resultant secondary photon distribution is shown in Fig. 8, the photon distribution should somewhat follow the electron distribution since the primary particles are positrons. The energy deposition, shown in Fig. 9, shows excellent agreement between the FLUKA and FluDAG implementations.

ATIC
The Advanced Thin Ionisation Calorimeter (ATIC) is a high altitude balloon based cosmic ray detector. The detector is composed of several layers of Bismuth Germanate (BiGeO) scintillator, interspersed between them are silicon matrices in order to determine the charge of particles traversing the layer. During the development of the detector several models were created, including native FLUKA models and Geant4 models. The FLUKA model was used as the basis of conversion for this analysis, using the ''export to MCNP'' feature in Flair, then using the mcnp2cad tool to produce an ACIS file.
This benchmark was used as evidence that the development work done for FluDAG was correct, it was decided by NASA that this geometry, shown in Fig. 10, would be used to validate the FluDAG code.
The primary incident particle here was a pencil beam of 1 GeV protons, the proton flux resultant from this source is shown in Fig. 11 and extracting the flux along the centre of Fig. 11 is shown in Fig. 12. There are no regions of disagreement, specifically in the lineout there is random undulation around unity but always within the error bars.
Similarly the resulting energy deposition of all primary and secondary particles is shown in Fig. 13, again with a lineout of data along the midplane is shown in Fig. 14. Similarly to the proton flux shown in Fig. 11 there is no disagreement between FLUKA & FLUDAG, and the lineout shows good agreement, even in regions of greater error.
The results were shown to give excellent agreement between FLUKA & FluDAG. For all the particles examined there was agreement within the statistical errors with no outliers beyond 3 . This along with the previous benchmarks showed that there is a direct agreement between FLUKA & FluDAG when equivalent geometries are used.

nTOF
The n_TOF (Neutron Time Of Flight) facility at CERN is a pulsed neutron source designed to study neutron-nucleus interactions for a wide range of neutron energies up to several GeV. The study of neutron interactions is of critical importance to several fields including nuclear technology, nuclear fusion, and accelerator driven systems. Bunches of ∼10 12 protons impinge on a largely lead target. Combinations of direct (p,xn) and (n,xn) interactions ultimately produce some 300 neutrons per incident proton. These neutrons are moderated using water at the rear of the target and then directed along a number of beamlines. These  Flair ''export to MCNP'', then using the mcnp2cad tool to convert the MCNP geometry to CAD. All the calculations run for this section were run on at the Centre for High Throughput Computing (CHTC) system, for 5000 primary histories per calculation, 5 calculations per batch, and 1000 groups of simulation, for a total of 25 million primaries.
A detailed comparison of the neutron spectrum leaving the rear of the device was performed. In the three cases a one way angular neutron flux was performed considering the neutrons crossing within 15 • of the surface normal and was binned into 1000 logarithmically spaced bins between 1.0 × 10 −14 GeV and 20 GeV. The purpose in this instance is to a) perform a like-for-like comparison between the native and direct translated geometries and b) determine the effect on downstream neutron production when a geometry representative of the true geometry is used.
The results for the comparison between FLUKA and FluDAG-Simplified are shown in Fig. 16. The comparison between the native FLUKA and detailed FluDAG comparison is shown in Fig. 17. Good agreement is seen at the higher energies above 400 MeV, however at lower energies there is a statistically significant difference. This difference is due to the heterogeneity of the CAD geometry in comparison to the native geometry, this is confirmed by the statistically equivalent results of the FluDAG-Simplified geometry results shown in Fig. 16.
The comparisons shown for the nTOF geometry reconfirms that for directly equivalent geometries produce statistically equivalent results for FLUKA and FluDAG. However, this example also demonstrates the importance of including sufficient details in your model for the analysis at hand.

Conclusions
In this paper the validation performed on FluDAG was discussed. The development of FluDAG was discussed, highlighting important lessons that were learned during the process, specifically regarding electron transport and charged particle transport in magnetic fields. Two benchmarks were reported on as part of FluDAG development, one covering very thin metallic layers (with thicknesses of less than 20 μm) named 'AlAuAl', the other was more geometrically complex with magnetic field enabled, and both were shown to give agreement to within statistical error with native FLUKA geometry and navigation. The complex geometry of the Advanced Thin Ionisation Calorimeter was used to compare FLUKA and FluDAG and were shown to give excellent agreement throughout the model to within statistical error. As an example of high energy physics use, the neutron time of flight (nTOF) facility was analysed and when using identical geometries were shown to give excellent agreement, however when using a detailed CAD geometry there were some differences beyond statistics, due to the heterogeneity. FluDAG and DAGMC have been shown to be a viable method for performing analyses allowing the use of complex CAD geometries and eliminating the need to approximate as-built geometries to confirm to MC second order primitive limitations. FluDAG development was funded by NASA under the Bioastronautics and the Human Health and Performance contracts with The University of Wisconsin-Madison.