Merge branch 'develop' into photon-documentation

.travis.yml

@@ -13,6 +13,8 @@ addons:
       - libmpich-dev
       - libhdf5-serial-dev
       - libhdf5-mpich-dev
+      - libblas-dev
+      - liblapack-dev
 cache:
   directories:
   - $HOME/nndc_hdf5
@@ -27,6 +29,7 @@ env:
     - OPENMC_CROSS_SECTIONS=$HOME/nndc_hdf5/cross_sections.xml
     - OPENMC_ENDF_DATA=$HOME/endf-b-vii.1
     - OPENMC_MULTIPOLE_LIBRARY=$HOME/WMP_Library
+    - LD_LIBRARY_PATH=$HOME/MOAB/lib:$HOME/DAGMC/lib
     - PATH=$PATH:$HOME/NJOY2016/build
     - DISPLAY=:99.0
     - COVERALLS_PARALLEL=true
@@ -35,6 +38,8 @@ env:
     - OMP=y MPI=n PHDF5=n
     - OMP=n MPI=y PHDF5=n
     - OMP=n MPI=y PHDF5=y
+    - OMP=n MPI=y PHDF5=y DAGMC=y
+    - OMP=y MPI=y PHDF5=y DAGMC=y        
 notifications:
   webhooks: https://coveralls.io/webhook?repo_token=$COVERALLS_REPO_TOKEN
 install:

CMakeLists.txt

@@ -19,6 +19,9 @@ option(debug    "Compile with debug flags"                       OFF)
 option(optimize "Turn on all compiler optimization flags"        OFF)
 option(coverage "Compile with coverage analysis flags"           OFF)
 option(mpif08   "Use Fortran 2008 MPI interface"                 OFF)
+option(dagmc    "Enable support for DAGMC (CAD) geometry"        OFF)
+
+# Maximum number of nested coordinates levels
 set(maxcoord 10 CACHE STRING "Maximum number of nested coordinate levels")
 
 #===============================================================================
@@ -36,10 +39,26 @@ if(MPI_ENABLED AND mpif08)
   message("-- Using Fortran 2008 MPI bindings")
 endif()
 
+#===============================================================================
+# DAGMC Geometry Support - need DAGMC/MOAB
+#===============================================================================
+if(dagmc)
+  find_package(DAGMC REQUIRED)
+  if(NOT DAGMC_FOUND)
+    message(FATAL_ERROR "Could not find DAGMC installation")
+  endif()
+  link_directories(${DAGMC_LIBRARY_DIRS})
+endif()
+
 #===============================================================================
 # HDF5 for binary output
 #===============================================================================
 
+# Allow user to specify HDF5_ROOT
+if (NOT (CMAKE_VERSION VERSION_LESS 3.12))
+  cmake_policy(SET CMP0074 NEW)
+endif()
+
 # Unfortunately FindHDF5.cmake will always prefer a serial HDF5 installation
 # over a parallel installation if both appear on the user's PATH. To get around
 # this, we check for the environment variable HDF5_ROOT and if it exists, use it
@@ -290,6 +309,7 @@ add_library(libopenmc SHARED
   src/algorithm.F90
   src/bank_header.F90
   src/api.F90
+  src/dagmc_header.F90
   src/cmfd_data.F90
   src/cmfd_execute.F90
   src/cmfd_header.F90
@@ -378,6 +398,7 @@ add_library(libopenmc SHARED
   src/tallies/tally_header.F90
   src/tallies/trigger.F90
   src/tallies/trigger_header.F90
+  src/dagmc.cpp
   src/cell.cpp
   src/cmfd_execute.cpp
   src/distribution.cpp
@@ -423,6 +444,7 @@ add_library(libopenmc SHARED
   src/thermal.cpp
   src/xml_interface.cpp
   src/xsdata.cpp)
+
 set_target_properties(libopenmc PROPERTIES
   OUTPUT_NAME openmc
   LINKER_LANGUAGE Fortran)
@@ -468,10 +490,15 @@ endif()
 target_link_libraries(libopenmc ${ldflags} ${HDF5_LIBRARIES} pugixml
                       faddeeva xtensor)
 
+if(dagmc)
+  target_compile_definitions(libopenmc PRIVATE DAGMC)
+  target_link_libraries(libopenmc ${DAGMC_LIBRARIES})
+  target_include_directories(libopenmc PRIVATE ${DAGMC_INCLUDE_DIRS})
+endif()
+
 #===============================================================================
 # openmc executable
 #===============================================================================
-
 add_executable(openmc src/main.cpp)
 target_compile_options(openmc PRIVATE ${cxxflags})
 target_link_libraries(openmc libopenmc)

Dockerfile

@@ -0,0 +1,39 @@
+FROM ubuntu:latest
+
+# Setup environment variables for Docker image
+ENV FC=/usr/bin/mpif90 CC=/usr/bin/mpicc CXX=/usr/bin/mpicxx \
+    PATH=/opt/openmc/bin:/opt/NJOY2016/build:$PATH \
+    LD_LIBRARY_PATH=/opt/openmc/lib:$LD_LIBRARY_PATH \
+    OPENMC_CROSS_SECTIONS=/root/nndc_hdf5/cross_sections.xml \
+    OPENMC_MULTIPOLE_LIBRARY=/root/WMP_Library \
+    OPENMC_ENDF_DATA=/root/endf-b-vii.1
+
+# Install dependencies from Debian package manager
+RUN apt-get update -y && \
+    apt-get upgrade -y && \
+    apt-get install -y python3-pip && \
+    apt-get install -y wget git emacs && \
+    apt-get install -y gfortran g++ cmake && \
+    apt-get install -y mpich libmpich-dev && \
+    apt-get install -y libhdf5-serial-dev libhdf5-mpich-dev && \
+    apt-get install -y imagemagick && \
+    apt-get autoremove
+
+# Update system-provided pip
+RUN pip3 install --upgrade pip
+
+# Clone and install NJOY2016
+RUN git clone https://github.com/njoy/NJOY2016 /opt/NJOY2016 && \
+    cd /opt/NJOY2016 && \
+    mkdir build && cd build && \
+    cmake -Dstatic=on .. && make 2>/dev/null && make install
+
+# Clone and install OpenMC
+RUN git clone https://github.com/openmc-dev/openmc.git /opt/openmc && \
+    cd /opt/openmc && mkdir -p build && cd build && \
+    cmake  -Doptimize=on -DHDF5_PREFER_PARALLEL=on .. && \
+    make && make install && \
+    cd .. && pip install -e .[test]
+
+# Download cross sections (NNDC and WMP) and ENDF data needed by test suite
+RUN ./opt/openmc/tools/ci/download-xs.sh

cmake/Modules/FindDAGMC.cmake

@@ -0,0 +1,18 @@
+# Try to find DAGMC
+#
+# Once done this will define
+#
+#  DAGMC_FOUND - system has DAGMC
+#  DAGMC_INCLUDE_DIRS - the DAGMC include directory
+#  DAGMC_LIBRARIES - Link these to use DAGMC
+#  DAGMC_DEFINITIONS - Compiler switches required for using DAGMC
+
+find_path(DAGMC_CMAKE_CONFIG NAMES DAGMCConfig.cmake
+          HINTS ${DAGMC_ROOT} $ENV{DAGMC_ROOT}
+          PATHS ENV LD_LIBRARY_PATH
+          PATH_SUFFIXES lib Lib cmake lib/cmake
+          NO_DEFAULT_PATH)
+
+message(STATUS "Found DAGMC in ${DAGMC_CMAKE_CONFIG}")
+
+include(${DAGMC_CMAKE_CONFIG}/DAGMCConfig.cmake)

docs/source/devguide/docker.rst

@@ -0,0 +1,57 @@
+.. _devguide_docker:
+
+======================
+Deployment with Docker
+======================
+
+OpenMC can be easily deployed using `Docker <https://www.docker.com/>`_ on any
+Windows, Mac or Linux system. With Docker running, execute the following
+command in the shell to build a `Docker image`_ called ``debian/openmc:latest``:
+
+.. code-block:: sh
+
+    docker build -t debian/openmc:latest https://github.com/openmc-dev/openmc.git#develop
+
+.. note:: This may take 5 -- 10 minutes to run to completion.
+
+This command will execute the instructions in OpenMC's ``Dockerfile`` to
+build a Docker image with OpenMC installed. The image includes OpenMC with
+MPICH and parallel HDF5 in the ``/opt/openmc`` directory, and
+`Miniconda3 <https://conda.io/miniconda.html>`_ with all of the Python
+pre-requisites (NumPy, SciPy, Pandas, etc.) installed. The
+`NJOY2016 <http://www.njoy21.io/NJOY2016/>`_ codebase is installed in
+``/opt/NJOY2016`` to support full functionality and testing of the
+``openmc.data`` Python module. The publicly available nuclear data libraries
+necessary to run OpenMC's test suite -- including NNDC and WMP cross sections
+and ENDF data --  are in the ``/opt/openmc/data directory``, and the
+corresponding :envvar:`OPENMC_CROSS_SECTIONS`,
+:envvar:`OPENMC_MULTIPOLE_LIBRARY`, and :envvar:`OPENMC_ENDF_DATA`
+environment variables are initialized.
+
+After building the Docker image, you can run the following to see the names of
+all images on your machine, including ``debian/openmc:latest``:
+
+.. code-block:: sh
+
+    docker image ls
+
+Now you can run the following to create a `Docker container`_ called
+``my_openmc`` based on the ``debian/openmc:latest`` image:
+
+.. code-block:: sh
+
+    docker run -it --name=my_openmc debian/openmc:latest
+
+This command will open an interactive shell running from within the
+Docker container where you have access to use OpenMC.
+
+.. note:: The ``docker run`` command supports many
+          `options <https://docs.docker.com/engine/reference/commandline/run/>`_
+          for spawning containers -- including `mounting volumes`_ from the
+          host filesystem -- which many users will find useful.
+
+.. _Docker image: https://docs.docker.com/engine/reference/commandline/images/
+.. _Docker container: https://www.docker.com/resources/what-container
+.. _options: https://docs.docker.com/engine/reference/commandline/run/
+.. _mounting volumes: https://docs.docker.com/storage/volumes/
+

docs/source/devguide/index.rst

@@ -18,3 +18,4 @@ other related topics.
     tests
     user-input
     docbuild
+    docker

docs/source/io_formats/fission_energy.rst

@@ -1,53 +0,0 @@
-.. _usersguide_fission_energy:
-
-==================================
-Fission Energy Release File Format
-==================================
-
-This file is a compact HDF5 representation of the ENDF MT=1, MF=458 data (see
-ENDF-102_ for details).  It gives the information needed to compute the energy
-carried away from fission reactions by each reaction product (e.g. fragment
-nuclei, neutrons) which depends on the incident neutron energy.  OpenMC is
-distributed with one of these files under
-data/fission_Q_data_endfb71.h5.  More files of this format can be created from
-ENDF files with the
-``openmc.data.write_compact_458_library`` function.  They can be read with the
-``openmc.data.FissionEnergyRelease.from_compact_hdf5`` class method.
-
-:Attributes: - **comment** (*char[]*) -- An optional text comment
-             - **component order** (*char[][]*) -- An array of strings
-               specifying the order each reaction product occurs in the data
-               arrays.  The components use the 2-3 letter abbreviations
-               specified in ENDF-102 e.g. EFR for fission fragments and ENP for
-               prompt neutrons.
-
-**/<nuclide name>/**
-    Nuclides are named by concatenating their atomic symbol and mass number. For
-    example, 'U235' or 'Pu239'.  Metastable nuclides are appended with an
-    '_m' and their metastable number.  For example, 'Am242_m1'
-
-:Datasets:
-          - **data** (*double[][][]*) -- The energy release coefficients.  The
-             first axis indexes the component type.  The second axis specifies
-             values or uncertainties.  The third axis indexes the polynomial
-             order.  If the data uses the Sher-Beck format, then the last axis
-             will have a length of one and ENDF-102 should be consulted for
-             energy dependence.  Otherwise, the data uses the Madland format
-             which is a polynomial of incident energy.
-
-             For example, if 'EFR' is given first in the **component order**
-             attribute and the data uses the Madland format, then the energy
-             released in the form of fission fragments at an incident energy
-             :math:`E` is given by
-
-             .. math::
-                \text{data}[0, 0, 0] + \text{data}[0, 0, 1] \cdot E
-                + \text{data}[0, 0, 2] \cdot E^2 + \ldots
-              
-             And its uncertainty is
-
-             .. math::
-                \text{data}[0, 1, 0] + \text{data}[0, 1, 1] \cdot E
-                + \text{data}[0, 1, 2] \cdot E^2 + \ldots
-
-.. _ENDF-102: http://www.nndc.bnl.gov/endfdocs/ENDF-102-2012.pdf

docs/source/io_formats/index.rst

@@ -34,7 +34,6 @@ Data Files
    nuclear_data
    mgxs_library
    data_wmp
-   fission_energy
 
 ------------
 Output Files

docs/source/io_formats/summary.rst

@@ -24,6 +24,8 @@ The current version of the summary file format is 6.0.
              - **n_universes** (*int*) -- Number of unique universes in the
                problem.
              - **n_lattices** (*int*) -- Number of lattices in the problem.
+             - **dagmc** (*int*) -- Indicates that a DAGMC geometry was used
+               if present.
 
 **/geometry/cells/cell <uid>/**
 

docs/source/pythonapi/data.rst

@@ -39,7 +39,6 @@ Core Functions
     openmc.data.linearize
     openmc.data.thin
     openmc.data.water_density
-    openmc.data.write_compact_458_library
     openmc.data.zam
 
 Angle-Energy Distributions

docs/source/quickinstall.rst

@@ -58,6 +58,29 @@ are no longer supported.
 .. _Personal Package Archive: https://launchpad.net/~paulromano/+archive/staging
 .. _APT package manager: https://help.ubuntu.com/community/AptGet/Howto
 
+-------------------------------------------
+Installing on Linux/Mac/Windows with Docker
+-------------------------------------------
+
+OpenMC can be easily deployed using `Docker <https://www.docker.com/>`_ on any
+Windows, Mac or Linux system. With Docker running, execute the following
+command in the shell to download and run a `Docker image`_ with the most recent release of OpenMC from `DockerHub <https://hub.docker.com/>`_ called ``openmc/openmc:v0.10.0``:
+
+.. code-block:: sh
+
+    docker run openmc/openmc:v0.10.0
+
+This will take several minutes to run depending on your internet download speed. The command will place you in an interactive shell running in a `Docker container`_ with OpenMC installed.
+
+.. note:: The ``docker run`` command supports many `options`_ for spawning
+          containers -- including `mounting volumes`_ from the host
+	  filesystem -- which many users will find useful.
+
+.. _Docker image: https://docs.docker.com/engine/reference/commandline/images/
+.. _Docker container: https://www.docker.com/resources/what-container
+.. _options: https://docs.docker.com/engine/reference/commandline/run/
+.. _mounting volumes: https://docs.docker.com/storage/volumes/
+
 ---------------------------------------
 Installing from Source on Ubuntu 15.04+
 ---------------------------------------

include/openmc/cell.h

@@ -13,6 +13,9 @@
 #include "openmc/constants.h"
 #include "openmc/position.h"
 
+#ifdef DAGMC
+#include "DagMC.hpp"
+#endif
 
 namespace openmc {
 
@@ -107,9 +110,8 @@ public:
 
   std::vector<int32_t> offset_;  //!< Distribcell offset table
 
-  Cell() {};
-
   explicit Cell(pugi::xml_node cell_node);
+  Cell() {};
 
   //! \brief Determine if a cell contains the particle at a given location.
   //!
@@ -130,21 +132,57 @@ public:
   //! \param on_surface The signed index of a surface that the coordinate is
   //!   known to be on.  This index takes precedence over surface sense
   //!   calculations.
+  virtual bool
+  contains(Position r, Direction u, int32_t on_surface) const = 0;
+
+  //! Find the oncoming boundary of this cell.
+  virtual std::pair<double, int32_t>
+  distance(Position r, Direction u, int32_t on_surface) const = 0;
+
+  //! Write all information needed to reconstruct the cell to an HDF5 group.
+  //! @param group_id An HDF5 group id.
+  virtual void to_hdf5(hid_t group_id) const = 0;
+
+  virtual ~Cell() {}
+};
+ 
+class CSGCell : public Cell
+{
+public:
+
+  CSGCell();
+  
+  explicit CSGCell(pugi::xml_node cell_node);
+
   bool
   contains(Position r, Direction u, int32_t on_surface) const;
 
-  //! Find the oncoming boundary of this cell.
   std::pair<double, int32_t>
   distance(Position r, Direction u, int32_t on_surface) const;
 
-  //! \brief Write cell information to an HDF5 group.
-  //! \param group_id An HDF5 group id.
   void to_hdf5(hid_t group_id) const;
 
+
+  
 protected:
   bool contains_simple(Position r, Direction u, int32_t on_surface) const;
   bool contains_complex(Position r, Direction u, int32_t on_surface) const;
 };
 
+#ifdef DAGMC
+class DAGCell : public Cell
+{
+public:
+  moab::DagMC* dagmc_ptr_;
+  DAGCell();
+
+  std::pair<double, int32_t> distance(Position r, Direction u, int32_t on_surface) const;
+  bool contains(Position r, Direction u, int32_t on_surface) const;
+
+  void to_hdf5(hid_t group_id) const;
+
+};
+#endif
+ 
 } // namespace openmc
 #endif // OPENMC_CELL_H

include/openmc/dagmc.h

@@ -0,0 +1,28 @@
+
+#ifndef OPENMC_DAGMC_H
+#define OPENMC_DAGMC_H
+
+#ifdef DAGMC
+
+#include "DagMC.hpp"
+#include "openmc/cell.h"
+#include "openmc/surface.h"
+
+namespace openmc {
+
+extern moab::DagMC* DAG;
+
+extern "C" void load_dagmc_geometry();
+extern "C" void free_memory_dagmc();
+
+}
+
+#endif // DAGMC
+
+#endif // OPENMC_DAGMC_H
+
+#ifdef DAGMC
+extern "C" constexpr bool dagmc_enabled = true;
+#else
+extern "C" constexpr bool dagmc_enabled = false;
+#endif

include/openmc/endf.h

@@ -29,6 +29,7 @@ bool is_fission(int MT);
 class Function1D {
 public:
   virtual double operator()(double x) const = 0;
+  virtual ~Function1D() = default;
 };
 
 //==============================================================================

include/openmc/position.h

@@ -1,7 +1,8 @@
 #ifndef OPENMC_POSITION_H
 #define OPENMC_POSITION_H
 
-#include <cmath>
+#include <cmath> // for sqrt
+#include <stdexcept> // for out_of_range
 #include <vector>
 
 namespace openmc {
@@ -24,11 +25,14 @@ struct Position {
   Position& operator-=(double);
   Position& operator*=(Position);
   Position& operator*=(double);
+
   const double& operator[](int i) const {
     switch (i) {
       case 0: return x;
       case 1: return y;
       case 2: return z;
+      default:
+        throw std::out_of_range{"Index in Position must be between 0 and 2."};
     }
   }
   double& operator[](int i) {
@@ -36,6 +40,8 @@ struct Position {
       case 0: return x;
       case 1: return y;
       case 2: return z;
+      default:
+        throw std::out_of_range{"Index in Position must be between 0 and 2."};
     }
   }
 

include/openmc/settings.h

@@ -45,6 +45,7 @@ extern "C" bool ufs_on;                  //!< uniform fission site method on?
 extern "C" bool urr_ptables_on;          //!< use unresolved resonance prob. tables?
 extern "C" bool write_all_tracks;        //!< write track files for every particle?
 extern "C" bool write_initial_source;    //!< write out initial source file?
+extern "C" bool dagmc;                   //!< indicator of DAGMC geometry
 
 // Paths to various files
 extern std::string path_cross_sections;   //!< path to cross_sections.xml
@@ -86,7 +87,6 @@ extern "C" int trigger_batch_interval;   //!< Batch interval for triggers
 extern "C" int verbosity;                //!< How verbose to make output
 extern "C" double weight_cutoff;         //!< Weight cutoff for Russian roulette
 extern "C" double weight_survive;        //!< Survival weight after Russian roulette
-
 } // namespace settings
 
 //! Read settings from XML file

include/openmc/surface.h

@@ -12,6 +12,9 @@
 #include "openmc/constants.h"
 #include "openmc/position.h"
 
+#ifdef DAGMC
+#include "DagMC.hpp"
+#endif
 
 namespace openmc {
 
@@ -65,6 +68,7 @@ public:
   std::vector<int> neighbor_neg_; //!< List of cells on negative side
 
   explicit Surface(pugi::xml_node surf_node);
+  Surface();
 
   virtual ~Surface() {}
 
@@ -104,12 +108,40 @@ public:
   //! Write all information needed to reconstruct the surface to an HDF5 group.
   //! \param group_id An HDF5 group id.
   //TODO: this probably needs to include i_periodic for PeriodicSurface
+  virtual void to_hdf5(hid_t group_id) const = 0;
+
+};
+
+class CSGSurface : public Surface
+{
+public:
+  explicit CSGSurface(pugi::xml_node surf_node);
+  CSGSurface();
+
   void to_hdf5(hid_t group_id) const;
 
 protected:
   virtual void to_hdf5_inner(hid_t group_id) const = 0;
 };
 
+//==============================================================================
+//! A `Surface` representing a DAGMC-based surface in DAGMC.
+//==============================================================================
+#ifdef DAGMC
+class DAGSurface : public Surface
+{
+public:
+  moab::DagMC* dagmc_ptr_;
+  DAGSurface();
+  double evaluate(Position r) const;
+  double distance(Position r, Direction u, bool coincident) const;
+  Direction normal(Position r) const;
+  //! Get the bounding box of this surface.
+  BoundingBox bounding_box() const;
+
+  void to_hdf5(hid_t group_id) const;
+};
+#endif
 //==============================================================================
 //! A `Surface` that supports periodic boundary conditions.
 //!
@@ -118,7 +150,7 @@ protected:
 //! `XPlane`-`YPlane` pairs.
 //==============================================================================
 
-class PeriodicSurface : public Surface
+class PeriodicSurface : public CSGSurface
 {
 public:
   int i_periodic_{C_NONE};    //!< Index of corresponding periodic surface
@@ -227,7 +259,7 @@ public:
 //! \f$(y - y_0)^2 + (z - z_0)^2 - R^2 = 0\f$
 //==============================================================================
 
-class SurfaceXCylinder : public Surface
+class SurfaceXCylinder : public CSGSurface
 {
   double y0_, z0_, radius_;
 public:
@@ -245,7 +277,7 @@ public:
 //! \f$(x - x_0)^2 + (z - z_0)^2 - R^2 = 0\f$
 //==============================================================================
 
-class SurfaceYCylinder : public Surface
+class SurfaceYCylinder : public CSGSurface
 {
   double x0_, z0_, radius_;
 public:
@@ -263,7 +295,7 @@ public:
 //! \f$(x - x_0)^2 + (y - y_0)^2 - R^2 = 0\f$
 //==============================================================================
 
-class SurfaceZCylinder : public Surface
+class SurfaceZCylinder : public CSGSurface
 {
   double x0_, y0_, radius_;
 public:
@@ -281,7 +313,7 @@ public:
 //! \f$(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2 - R^2 = 0\f$
 //==============================================================================
 
-class SurfaceSphere : public Surface
+class SurfaceSphere : public CSGSurface
 {
   double x0_, y0_, z0_, radius_;
 public:
@@ -299,7 +331,7 @@ public:
 //! \f$(y - y_0)^2 + (z - z_0)^2 - R^2 (x - x_0)^2 = 0\f$
 //==============================================================================
 
-class SurfaceXCone : public Surface
+class SurfaceXCone : public CSGSurface
 {
   double x0_, y0_, z0_, radius_sq_;
 public:
@@ -317,7 +349,7 @@ public:
 //! \f$(x - x_0)^2 + (z - z_0)^2 - R^2 (y - y_0)^2 = 0\f$
 //==============================================================================
 
-class SurfaceYCone : public Surface
+class SurfaceYCone : public CSGSurface
 {
   double x0_, y0_, z0_, radius_sq_;
 public:
@@ -335,7 +367,7 @@ public:
 //! \f$(x - x_0)^2 + (y - y_0)^2 - R^2 (z - z_0)^2 = 0\f$
 //==============================================================================
 
-class SurfaceZCone : public Surface
+class SurfaceZCone : public CSGSurface
 {
   double x0_, y0_, z0_, radius_sq_;
 public:
@@ -352,7 +384,7 @@ public:
 //! \f$A x^2 + B y^2 + C z^2 + D x y + E y z + F x z + G x + H y + J z + K = 0\f$
 //==============================================================================
 
-class SurfaceQuadric : public Surface
+class SurfaceQuadric : public CSGSurface
 {
   // Ax^2 + By^2 + Cz^2 + Dxy + Eyz + Fxz + Gx + Hy + Jz + K = 0
   double A_, B_, C_, D_, E_, F_, G_, H_, J_, K_;

openmc/capi/__init__.py

@@ -12,7 +12,7 @@ objects in the :mod:`openmc.capi` subpackage, for example:
 
 """
 
-from ctypes import CDLL
+from ctypes import CDLL, c_bool
 import os
 import sys
 
@@ -38,6 +38,8 @@ else:
     from unittest.mock import Mock
     _dll = Mock()
 
+dagmc_enabled = bool(c_bool.in_dll(_dll, "dagmc_enabled"))
+    
 from .error import *
 from .core import *
 from .nuclide import *

openmc/capi/core.py

@@ -112,7 +112,7 @@ def find_material(xyz):
     if isinstance(mats, openmc.capi.Material):
         return mats
     else:
-        return mats[instance]
+        return mats[instance.value]
 
 def hard_reset():
     """Reset tallies, timers, and pseudo-random number generator state."""

openmc/data/fission_energy.py

@@ -43,10 +43,6 @@ class FissionEnergyRelease(EqualityMixin):
     dataset.  The :meth:`FissionEnergyRelease.from_hdf5` method builds this
     class from the usual OpenMC HDF5 data files.
     :meth:`FissionEnergyRelease.from_endf` uses ENDF-formatted data.
-    :meth:`FissionEnergyRelease.from_compact_hdf5` uses a different HDF5 format
-    that is meant to be compact and store the exact same data as the ENDF
-    format.  Files with this format can be generated with the
-    :func:`openmc.data.write_compact_458_library` function.
 
     References
     ----------

openmc/geometry.py

@@ -266,7 +266,10 @@ class Geometry(object):
             Dictionary mapping cell IDs to :class:`openmc.Cell` instances
 
         """
-        return self.root_universe.get_all_cells()
+        if self.root_universe is not None:
+            return self.root_universe.get_all_cells()
+        else:
+            return []
 
     def get_all_universes(self):
         """Return all universes in the geometry.

openmc/settings.py

@@ -37,6 +37,8 @@ class Settings(object):
         weight assigned to particles that are not killed after Russian
         roulette. Value of energy should be a float indicating energy in eV
         below which particle type will be killed.
+    dagmc : bool
+        Indicate that a CAD-based DAGMC geometry will be used.
     electron_treatment : {'led', 'ttb'}
         Whether to deposit all energy from electrons locally ('led') or create
         secondary bremsstrahlung photons ('ttb').
@@ -224,6 +226,8 @@ class Settings(object):
         self._create_fission_neutrons = None
         self._log_grid_bins = None
 
+        self._dagmc = None
+        
     @property
     def run_mode(self):
         return self._run_mode
@@ -368,6 +372,10 @@ class Settings(object):
     def log_grid_bins(self):
         return self._log_grid_bins
 
+    @property
+    def dagmc(self):
+        return self._dagmc
+    
     @run_mode.setter
     def run_mode(self, run_mode):
         cv.check_value('run mode', run_mode, _RUN_MODES)
@@ -511,6 +519,11 @@ class Settings(object):
         cv.check_type('photon transport', photon_transport, bool)
         self._photon_transport = photon_transport
 
+    @dagmc.setter
+    def dagmc(self, dagmc):
+        cv.check_type('dagmc geometry', dagmc, bool)
+        self._dagmc = dagmc
+        
     @ptables.setter
     def ptables(self, ptables):
         cv.check_type('probability tables', ptables, bool)
@@ -945,6 +958,11 @@ class Settings(object):
             elem = ET.SubElement(root, "log_grid_bins")
             elem.text = str(self._log_grid_bins)
 
+    def _create_dagmc_subelement(self, root):
+        if self._dagmc is not None:
+            elem = ET.SubElement(root, "dagmc")
+            element.text = str(self._dagmc).lower()
+        
     def export_to_xml(self, path='settings.xml'):
         """Export simulation settings to an XML file.
 
@@ -992,7 +1010,8 @@ class Settings(object):
         self._create_volume_calcs_subelement(root_element)
         self._create_create_fission_neutrons_subelement(root_element)
         self._create_log_grid_bins_subelement(root_element)
-
+        self._create_dagmc_subelement(root_element)
+        
         # Clean the indentation in the file to be user-readable
         clean_indentation(root_element)
 

openmc/summary.py

@@ -97,6 +97,9 @@ class Summary(object):
                 self._macroscopics = name.decode()
 
     def _read_geometry(self):
+        if "dagmc" in self._f['geometry'].attrs.keys():
+            return
+
         # Read in and initialize the Materials and Geometry
         self._read_materials()
         self._read_surfaces()

scripts/openmc-convert-mcnp71-data

@@ -25,8 +25,6 @@ parser = argparse.ArgumentParser(
 )
 parser.add_argument('-d', '--destination', default='mcnp_endfb71',
                     help='Directory to create new library in')
-parser.add_argument('-f', '--fission_energy_release',
-                    help='HDF5 file containing fission energy release data')
 parser.add_argument('--libver', choices=['earliest', 'latest'],
                     default='earliest', help="Output HDF5 versioning. Use "
                     "'earliest' for backwards compatibility or 'latest' for "
@@ -64,13 +62,6 @@ for basename, xs_list in sorted(suffixes.items()):
     else:
         data = openmc.data.IncidentNeutron.from_ace(filename, 'mcnp')
 
-        # Add fission energy release data, if available
-        if args.fission_energy_release is not None:
-            fer = openmc.data.FissionEnergyRelease.from_compact_hdf5(
-                args.fission_energy_release, data)
-            if fer is not None:
-                data.fission_energy = fer
-
     # For each higher temperature, add cross sections to the existing table
     for xs in xs_list[1:]:
         filename = '.'.join((basename, xs))

src/api.F90

@@ -29,6 +29,10 @@ module openmc_api
   use timer_header
   use volume_calc,     only: openmc_calculate_volumes
 
+#ifdef DAGMC
+  use dagmc_header,      only: free_memory_dagmc
+#endif
+
   implicit none
 
   private
@@ -149,6 +153,7 @@ contains
     root_universe = -1
     run_CE = .true.
     run_mode = -1
+    dagmc = .false.
     satisfy_triggers = .false.
     call openmc_set_seed(DEFAULT_SEED)
     source_latest = .false.
@@ -325,6 +330,9 @@ contains
     call free_memory_tally_filter()
     call free_memory_tally_derivative()
     call free_memory_bank()
+#ifdef DAGMC
+    call free_memory_dagmc()
+#endif
 
     ! Deallocate CMFD
     call deallocate_cmfd(cmfd)

src/cell.cpp

@@ -15,7 +15,6 @@
 #include "openmc/surface.h"
 #include "openmc/xml_interface.h"
 
-
 namespace openmc {
 
 //==============================================================================
@@ -209,7 +208,9 @@ Universe::to_hdf5(hid_t universes_group) const
 // Cell implementation
 //==============================================================================
 
-Cell::Cell(pugi::xml_node cell_node)
+CSGCell::CSGCell() {} // empty constructor
+  
+CSGCell::CSGCell(pugi::xml_node cell_node)
 {
   if (check_for_node(cell_node, "id")) {
     id_ = std::stoi(get_node_value(cell_node, "id"));
@@ -393,7 +394,7 @@ Cell::Cell(pugi::xml_node cell_node)
 //==============================================================================
 
 bool
-Cell::contains(Position r, Direction u, int32_t on_surface) const
+CSGCell::contains(Position r, Direction u, int32_t on_surface) const
 {
   if (simple_) {
     return contains_simple(r, u, on_surface);
@@ -405,7 +406,7 @@ Cell::contains(Position r, Direction u, int32_t on_surface) const
 //==============================================================================
 
 std::pair<double, int32_t>
-Cell::distance(Position r, Direction u, int32_t on_surface) const
+CSGCell::distance(Position r, Direction u, int32_t on_surface) const
 {
   double min_dist {INFTY};
   int32_t i_surf {std::numeric_limits<int32_t>::max()};
@@ -434,12 +435,12 @@ Cell::distance(Position r, Direction u, int32_t on_surface) const
 //==============================================================================
 
 void
-Cell::to_hdf5(hid_t cells_group) const
+CSGCell::to_hdf5(hid_t cell_group) const
 {
   // Create a group for this cell.
   std::stringstream group_name;
   group_name << "cell " << id_;
-  auto group = create_group(cells_group, group_name);
+  auto group = create_group(cell_group, group_name);
 
   if (!name_.empty()) {
     write_string(group, "name", name_, false);
@@ -513,7 +514,7 @@ Cell::to_hdf5(hid_t cells_group) const
 //==============================================================================
 
 bool
-Cell::contains_simple(Position r, Direction u, int32_t on_surface) const
+CSGCell::contains_simple(Position r, Direction u, int32_t on_surface) const
 {
   for (int32_t token : rpn_) {
     if (token < OP_UNION) {
@@ -537,7 +538,7 @@ Cell::contains_simple(Position r, Direction u, int32_t on_surface) const
 //==============================================================================
 
 bool
-Cell::contains_complex(Position r, Direction u, int32_t on_surface) const
+CSGCell::contains_complex(Position r, Direction u, int32_t on_surface) const
 {
   // Make a stack of booleans.  We don't know how big it needs to be, but we do
   // know that rpn.size() is an upper-bound.
@@ -585,6 +586,50 @@ Cell::contains_complex(Position r, Direction u, int32_t on_surface) const
   }
 }
 
+//==============================================================================
+// DAGMC Cell implementation
+//==============================================================================
+#ifdef DAGMC
+DAGCell::DAGCell() : Cell{} {};
+
+std::pair<double, int32_t>
+DAGCell::distance(Position r, Direction u, int32_t on_surface) const
+{
+  moab::ErrorCode rval;
+  moab::EntityHandle vol = dagmc_ptr_->entity_by_id(3, id_);
+  moab::EntityHandle hit_surf;
+  double dist;
+  double pnt[3] = {r.x, r.y, r.z};
+  double dir[3] = {u.x, u.y, u.z};
+  rval = dagmc_ptr_->ray_fire(vol, pnt, dir, hit_surf, dist);
+  MB_CHK_ERR_CONT(rval);
+  int surf_idx;
+  if (hit_surf != 0) {
+    surf_idx = dagmc_ptr_->index_by_handle(hit_surf);
+  } else {  // indicate that particle is lost
+    surf_idx = -1;
+  }
+  
+  return {dist, surf_idx};
+}
+  
+bool DAGCell::contains(Position r, Direction u, int32_t on_surface) const
+{
+  moab::ErrorCode rval;
+  moab::EntityHandle vol = dagmc_ptr_->entity_by_id(3, id_);
+
+  int result = 0;
+  double pnt[3] = {r.x, r.y, r.z};
+  double dir[3] = {u.x, u.y, u.z};
+  rval = dagmc_ptr_->point_in_volume(vol, pnt, result, dir);
+  MB_CHK_ERR_CONT(rval);  
+  return result;
+}
+
+void DAGCell::to_hdf5(hid_t group_id) const { return; }
+
+#endif
+
 //==============================================================================
 // Non-method functions
 //==============================================================================
@@ -601,7 +646,7 @@ read_cells(pugi::xml_node* node)
   // Loop over XML cell elements and populate the array.
   cells.reserve(n_cells);
   for (pugi::xml_node cell_node: node->children("cell")) {
-    cells.push_back(new Cell(cell_node));
+    cells.push_back(new CSGCell(cell_node));
   }
 
   // Populate the Universe vector and map.
@@ -725,6 +770,21 @@ extern "C" {
 
   int cell_type(Cell* c) {return c->type_;}
 
+#ifdef DAGMC
+
+  int32_t next_cell(DAGCell* cur_cell, DAGSurface* surf_xed )
+  {
+    moab::EntityHandle surf = surf_xed->dagmc_ptr_->entity_by_id(2,surf_xed->id_);
+    moab::EntityHandle vol = cur_cell->dagmc_ptr_->entity_by_id(3,cur_cell->id_);
+
+    moab::EntityHandle new_vol;
+    cur_cell->dagmc_ptr_->next_vol(surf, vol, new_vol);
+
+    return cur_cell->dagmc_ptr_->index_by_handle(new_vol);
+  }
+
+#endif
+
   int32_t cell_universe(Cell* c) {return c->universe_;}
 
   int32_t cell_fill(Cell* c) {return c->fill_;}
@@ -753,7 +813,7 @@ extern "C" {
   {
     cells.reserve(cells.size() + n);
     for (int32_t i = 0; i < n; i++) {
-      cells.push_back(new Cell());
+      cells.push_back(new CSGCell());
     }
     n_cells = cells.size();
   }

src/dagmc.cpp

@@ -0,0 +1,144 @@
+
+#include "openmc/dagmc.h"
+#include "openmc/error.h"
+#include "openmc/string_functions.h"
+#include "openmc/settings.h"
+#include "openmc/geometry.h"
+
+#include <string>
+#include <sstream>
+#include <algorithm>
+
+#ifdef DAGMC
+
+namespace openmc {
+
+moab::DagMC* DAG;
+
+void load_dagmc_geometry()
+{
+  if (!DAG) {
+    DAG = new moab::DagMC();
+  }
+
+  int32_t dagmc_univ_id = 0; // universe is always 0 for DAGMC
+
+  moab::ErrorCode rval = DAG->load_file("dagmc.h5m");
+  MB_CHK_ERR_CONT(rval);
+
+  rval = DAG->init_OBBTree();
+  MB_CHK_ERR_CONT(rval);
+
+  std::vector<std::string> prop_keywords;
+  prop_keywords.push_back("mat");
+  prop_keywords.push_back("boundary");
+
+  std::map<std::string, std::string> ph;
+  DAG->parse_properties(prop_keywords, ph, ":");
+  MB_CHK_ERR_CONT(rval);
+
+  // initialize cell objects
+  n_cells = DAG->num_entities(3);
+
+  // Allocate the cell overlap count if necessary.
+  if (settings::check_overlaps) overlap_check_count.resize(n_cells, 0);
+
+  for (int i = 0; i < n_cells; i++) {
+    moab::EntityHandle vol_handle = DAG->entity_by_index(3, i+1);
+
+    // set cell ids using global IDs
+    DAGCell* c = new DAGCell();
+    c->id_ = DAG->id_by_index(3, i+1);
+    c->dagmc_ptr_ = DAG;
+    c->universe_ = dagmc_univ_id; // set to zero for now
+    c->fill_ = C_NONE; // no fill, single universe
+
+    cells.push_back(c);
+    cell_map[c->id_] = c->id_;
+
+    // Populate the Universe vector and dict
+    auto it = universe_map.find(dagmc_univ_id);
+    if (it == universe_map.end()) {
+      universes.push_back(new Universe());
+      universes.back()-> id_ = dagmc_univ_id;
+      universes.back()->cells_.push_back(i);
+      universe_map[dagmc_univ_id] = universes.size() - 1;
+    } else {
+      universes[it->second]->cells_.push_back(i);
+    }
+
+    if (DAG->is_implicit_complement(vol_handle)) {
+      // assuming implicit complement is void for now
+      c->material_.push_back(MATERIAL_VOID);
+      continue;
+    }
+
+    if (DAG->has_prop(vol_handle, "mat")){
+      std::string mat_value;
+      rval = DAG->prop_value(vol_handle, "mat", mat_value);
+      MB_CHK_ERR_CONT(rval);
+      to_lower(mat_value);
+
+      if (mat_value == "void" || mat_value == "vacuum") {
+        c->material_.push_back(MATERIAL_VOID);
+      } else {
+        c->material_.push_back(std::stoi(mat_value));
+      }
+    } else {
+      std::stringstream err_msg;
+      err_msg << "Volume " << c->id_ << " has no material assignment.";
+      fatal_error(err_msg.str());
+    }
+  }
+
+  // initialize surface objects
+  n_surfaces = DAG->num_entities(2);
+  surfaces.resize(n_surfaces);
+
+  for (int i = 0; i < n_surfaces; i++) {
+    moab::EntityHandle surf_handle = DAG->entity_by_index(2, i+1);
+
+    // set cell ids using global IDs
+    DAGSurface* s = new DAGSurface();
+    s->id_ = DAG->id_by_index(2, i+1);
+    s->dagmc_ptr_ = DAG;
+
+    if (DAG->has_prop(surf_handle, "boundary")) {
+      std::string bc_value;
+      rval = DAG->prop_value(surf_handle, "boundary", bc_value);
+      MB_CHK_ERR_CONT(rval);
+      to_lower(bc_value);
+
+      if (bc_value == "transmit" || bc_value == "transmission") {
+        s->bc_ = BC_TRANSMIT;
+      } else if (bc_value == "vacuum") {
+        s->bc_ = BC_VACUUM;
+      } else if (bc_value == "reflective" || bc_value == "reflect" || bc_value == "reflecting") {
+        s->bc_ = BC_REFLECT;
+      } else if (bc_value == "periodic") {
+        fatal_error("Periodic boundary condition not supported in DAGMC.");
+      } else {
+        std::stringstream err_msg;
+        err_msg << "Unknown boundary condition \"" << s->bc_
+                << "\" specified on surface " << s->id_;
+        fatal_error(err_msg);
+      }
+    } else {   // if no BC property is found, set to transmit
+      s->bc_ = BC_TRANSMIT;
+    }
+
+    // add to global array and map
+    surfaces[i] = s;
+    surface_map[s->id_] = s->id_;
+  }
+
+  return;
+}
+
+void free_memory_dagmc()
+{
+  delete DAG;
+}
+
+}
+#endif

src/dagmc_header.F90

@@ -0,0 +1,20 @@
+#ifdef DAGMC
+
+module dagmc_header
+
+  use, intrinsic :: ISO_C_BINDING
+
+  implicit none
+
+  interface
+    subroutine load_dagmc_geometry() bind(C)
+    end subroutine load_dagmc_geometry
+
+    subroutine free_memory_dagmc() bind(C)
+    end subroutine free_memory_dagmc
+
+  end interface
+
+end module dagmc_header
+
+#endif

src/distribution.cpp

@@ -4,6 +4,7 @@
 #include <cmath>     // for sqrt, floor, max
 #include <iterator>  // for back_inserter
 #include <numeric>   // for accumulate
+#include <stdexcept> // for runtime_error
 #include <string>    // for string, stod
 
 #include "openmc/error.h"
@@ -44,7 +45,7 @@ double Discrete::sample() const
       c += p_[i];
       if (xi < c) return x_[i];
     }
-    // throw exception?
+    throw std::runtime_error{"Error when sampling probability mass function."};
   } else {
     return x_[0];
   }
@@ -248,19 +249,21 @@ UPtrDist distribution_from_xml(pugi::xml_node node)
   std::string type = get_node_value(node, "type", true, true);
 
   // Allocate extension of Distribution
+  UPtrDist dist;
   if (type == "uniform") {
-    return UPtrDist{new Uniform(node)};
+    dist = UPtrDist{new Uniform(node)};
   } else if (type == "maxwell") {
-    return UPtrDist{new Maxwell(node)};
+    dist = UPtrDist{new Maxwell(node)};
   } else if (type == "watt") {
-    return UPtrDist{new Watt(node)};
+    dist = UPtrDist{new Watt(node)};
   } else if (type == "discrete") {
-    return UPtrDist{new Discrete(node)};
+    dist = UPtrDist{new Discrete(node)};
   } else if (type == "tabular") {
-    return UPtrDist{new Tabular(node)};
+    dist = UPtrDist{new Tabular(node)};
   } else {
     openmc::fatal_error("Invalid distribution type: " + type);
   }
+  return dist;
 }
 
 } // namespace openmc

src/eigenvalue.cpp

@@ -146,7 +146,7 @@ extern "C" double entropy_c(int i)
   return entropy.at(i - 1);
 }
 
-extern "C" double entropy_clear()
+extern "C" void entropy_clear()
 {
   entropy.clear();
 }

src/endf.cpp

@@ -3,6 +3,7 @@
 #include <algorithm> // for copy
 #include <cmath>     // for log, exp
 #include <iterator>  // for back_inserter
+#include <stdexcept> // for runtime_error
 
 #include "xtensor/xarray.hpp"
 #include "xtensor/xview.hpp"
@@ -19,17 +20,23 @@ namespace openmc {
 
 Interpolation int2interp(int i)
 {
+  // TODO: We are ignoring specification of two-dimensional interpolation
+  // schemes (method of corresponding points and unit base interpolation). Those
+  // should be accounted for in the distribution classes somehow.
+
   switch (i) {
-  case 1:
+  case 1: case 11: case 21:
     return Interpolation::histogram;
-  case 2:
+  case 2: case 12: case 22:
     return Interpolation::lin_lin;
-  case 3:
+  case 3: case 13: case 23:
     return Interpolation::lin_log;
-  case 4:
+  case 4: case 14: case 24:
     return Interpolation::log_lin;
-  case 5:
+  case 5: case 15: case 25:
     return Interpolation::log_log;
+  default:
+    throw std::runtime_error{"Invalid interpolation code."};
   }
 }
 
@@ -142,6 +149,8 @@ double Tabulated1D::operator()(double x) const
   case Interpolation::log_log:
     r = log(x/x0)/log(x1/x0);
     return y0*exp(r*log(y1/y0));
+  default:
+    throw std::runtime_error{"Invalid interpolation scheme."};
   }
 }
 

src/geometry.F90

@@ -55,10 +55,34 @@ module geometry
 
     subroutine neighbor_lists() bind(C)
     end subroutine neighbor_lists
+
+#ifdef DAGMC
+
+    function next_cell_c(current_cell, surface_crossed) &
+         bind(C, name="next_cell") result(new_cell)
+      import C_PTR, C_INT32_T
+      type(C_PTR), intent(in), value :: current_cell
+      type(C_PTR), intent(in), value :: surface_crossed
+      integer(C_INT32_T)             :: new_cell
+    end function next_cell_c
+
+#endif
+
   end interface
 
 contains
 
+#ifdef DAGMC
+
+  function next_cell(c, s) result(new_cell)
+    type(Cell), intent(in) :: c
+    type(Surface), intent(in) :: s
+    integer :: new_cell
+    new_cell = next_cell_c(c%ptr, s%ptr)
+  end function next_cell
+
+#endif
+
 !===============================================================================
 ! FIND_CELL determines what cell a source particle is in within a particular
 ! universe. If the base universe is passed, the particle should be found as long

src/geometry_aux.cpp

@@ -389,6 +389,7 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
         return path.str();
       }
     }
+    throw std::runtime_error{"Error determining distribcell path."};
   }
 }
 

src/input_xml.F90

@@ -11,6 +11,9 @@ module input_xml
   use error,            only: fatal_error, warning, write_message, openmc_err_msg
   use geometry,         only: neighbor_lists
   use geometry_header
+#ifdef DAGMC
+  use dagmc_header
+#endif
   use hdf5_interface
   use list_header,      only: ListChar, ListInt, ListReal
   use material_header
@@ -176,7 +179,9 @@ contains
 
     ! After reading input and basic geometry setup is complete, build lists of
     ! neighboring cells for efficient tracking
-    call neighbor_lists()
+    if (.not. dagmc) then
+      call neighbor_lists()
+    end if
 
     ! Assign temperatures to cells that don't have temperatures already assigned
     call assign_temperatures()
@@ -347,6 +352,63 @@ contains
 
   end subroutine read_settings_xml_f
 
+
+#ifdef DAGMC
+
+!===============================================================================
+! READ_GEOMETRY_DAGMC reads data from a DAGMC .h5m file, checking
+! for material properties and surface boundary conditions
+! some universe information is spoofed for now
+!===============================================================================
+
+  subroutine read_geometry_dagmc()
+
+    integer :: i, j
+    integer :: univ_id
+    integer :: n_cells_in_univ
+    logical :: file_exists
+    character(MAX_LINE_LEN) :: filename
+    type(Cell),     pointer :: c
+    type(VectorInt) :: univ_ids      ! List of all universe IDs
+    type(DictIntInt) :: cells_in_univ_dict ! Used to count how many cells each
+                                           ! universe contains
+
+    ! Check if dagmc.h5m exists
+    filename = trim(path_input) // "dagmc.h5m"
+    inquire(FILE=filename, EXIST=file_exists)
+    if (.not. file_exists) then
+      call fatal_error("Geometry DAGMC file '" // trim(filename) // "' does not &
+           &exist!")
+    end if
+
+    call write_message("Reading DAGMC geometry...", 5)
+    call load_dagmc_geometry()
+    call allocate_surfaces()
+    call allocate_cells()
+
+    ! setup universe data structs
+    do i = 1, n_cells
+      c => cells(i)
+      ! additional metadata spoofing
+      univ_id = c % universe()
+
+      if (.not. cells_in_univ_dict % has(univ_id)) then
+        n_universes = n_universes + 1
+        n_cells_in_univ = 1
+        call universe_dict % set(univ_id, n_universes)
+        call univ_ids % push_back(univ_id)
+      else
+        n_cells_in_univ = 1 + cells_in_univ_dict % get(univ_id)
+      end if
+      call cells_in_univ_dict % set(univ_id, n_cells_in_univ)
+    end do
+
+    root_universe = find_root_universe()
+
+  end subroutine read_geometry_dagmc
+
+#endif
+
 !===============================================================================
 ! READ_GEOMETRY_XML reads data from a geometry.xml file and parses it, checking
 ! for errors and placing properly-formatted data in the right data structures
@@ -374,6 +436,12 @@ contains
     type(VectorInt) :: univ_ids      ! List of all universe IDs
     type(DictIntInt) :: cells_in_univ_dict ! Used to count how many cells each
                                            ! universe contains
+#ifdef DAGMC
+    if (dagmc) then
+      call read_geometry_dagmc()
+      return
+    end if
+#endif
 
     ! Display output message
     call write_message("Reading geometry XML file...", 5)
@@ -401,7 +469,6 @@ contains
 
     ! Allocate surfaces array
     allocate(surfaces(n_surfaces))
-
     do i = 1, n_surfaces
       surfaces(i) % ptr = surface_pointer(i - 1);
 
@@ -559,6 +626,40 @@ contains
 
   end subroutine read_geometry_xml
 
+  subroutine allocate_surfaces()
+    integer :: i
+
+    ! Allocate surfaces array
+    allocate(surfaces(n_surfaces))
+
+    do i = 1, n_surfaces
+      surfaces(i) % ptr = surface_pointer(i - 1);
+      ! Add surface to dictionary
+      call surface_dict % set(surfaces(i) % id(), i)
+    end do
+
+  end subroutine allocate_surfaces
+
+  subroutine allocate_cells()
+    integer :: i
+    type(Cell), pointer :: c
+
+    ! Allocate cells array
+    allocate(cells(n_cells))
+
+    do i = 1, n_cells
+      c => cells(i)
+      c % ptr = cell_pointer(i - 1)
+      ! Check to make sure 'id' hasn't been used
+      if (cell_dict % has(c % id())) then
+        call fatal_error("Two or more cells use the same unique ID: " &
+               // to_str(c % id()))
+      end if
+      ! Add cell to dictionary
+      call cell_dict % set(c % id(), i)
+    end do
+  end subroutine allocate_cells
+
 !===============================================================================
 ! READ_MATERIAL_XML reads data from a materials.xml file and parses it, checking
 ! for errors and placing properly-formatted data in the right data structures

src/mesh.cpp

@@ -162,12 +162,15 @@ int RegularMesh::get_bin(Position r) const
 
 int RegularMesh::get_bin_from_indices(const int* ijk) const
 {
-  if (n_dimension_ == 1) {
-    return ijk[0];
-  } else if (n_dimension_ == 2) {
-    return (ijk[1] - 1)*shape_[0] + ijk[0];
-  } else if (n_dimension_ == 3) {
-    return ((ijk[2] - 1)*shape_[1] + (ijk[1] - 1))*shape_[0] + ijk[0];
+  switch (n_dimension_) {
+    case 1:
+      return ijk[0];
+    case 2:
+      return (ijk[1] - 1)*shape_[0] + ijk[0];
+    case 3:
+      return ((ijk[2] - 1)*shape_[1] + (ijk[1] - 1))*shape_[0] + ijk[0];
+    default:
+      throw std::runtime_error{"Invalid number of mesh dimensions"};
   }
 }
 
@@ -206,6 +209,8 @@ bool RegularMesh::intersects(Position r0, Position r1) const
       return intersects_2d(r0, r1);
     case 3:
       return intersects_3d(r0, r1);
+    default:
+      throw std::runtime_error{"Invalid number of mesh dimensions."};
   }
 }
 

src/relaxng/settings.rnc

@@ -57,6 +57,8 @@ element settings {
 
   element ptables { xsd:boolean }? &
 
+  element dagmc { xsd:boolean }? &
+
   element run_cmfd { xsd:boolean }? &
 
   element run_mode { xsd:string }? &

src/relaxng/settings.rng

@@ -254,6 +254,11 @@
         <data type="boolean"/>
       </element>
     </optional>
+    <optional>
+      <element name="dagmc">
+        <data type="boolean"/>
+      </element>
+    </optional>
     <optional>
       <element name="run_cmfd">
         <data type="boolean"/>

src/settings.F90

@@ -80,6 +80,9 @@ module settings
   ! Mode to run in (fixed source, eigenvalue, plotting, etc)
   integer(C_INT), bind(C) :: run_mode
 
+  ! flag for use of DAGMC geometry
+  logical(C_BOOL), bind(C) :: dagmc
+
   ! Restart run
   logical(C_BOOL), bind(C) :: restart_run
 

src/settings.cpp

@@ -56,7 +56,8 @@ bool ufs_on                  {false};
 bool urr_ptables_on          {true};
 bool write_all_tracks        {false};
 bool write_initial_source    {false};
-
+bool dagmc                   {false};
+  
 std::string path_cross_sections;
 std::string path_input;
 std::string path_multipole;
@@ -207,6 +208,17 @@ void read_settings_xml()
     verbosity = std::stoi(get_node_value(root, "verbosity"));
   }
 
+  // DAGMC geometry check
+  if (check_for_node(root, "dagmc")) {
+    dagmc = get_node_value_bool(root, "dagmc");
+  }
+
+#ifndef DAGMC
+  if (dagmc) {
+    fatal_error("DAGMC mode unsupported for this build of OpenMC");
+  }
+#endif
+  
   // To this point, we haven't displayed any output since we didn't know what
   // the verbosity is. Now that we checked for it, show the title if necessary
   if (openmc_master) {
@@ -396,6 +408,13 @@ void read_settings_xml()
 #endif
   }
 
+#ifdef _OPENMP
+  if (dagmc && omp_get_max_threads() > 1) {
+    warning("Forcing number of threads to 1 for DAGMC simulation.");
+    omp_set_num_threads(1);
+  }
+#endif
+  
   // ==========================================================================
   // EXTERNAL SOURCE
 

src/summary.F90

@@ -82,9 +82,9 @@ contains
     integer(HID_T) :: nuclide_group
     integer(HID_T) :: macro_group
     integer :: i
-    character(12), allocatable :: nuc_names(:)
-    character(12), allocatable :: macro_names(:)
-    real(8), allocatable :: awrs(:)
+    character(kind=C_CHAR, len=20), allocatable :: nuc_names(:)
+    character(kind=C_CHAR, len=20), allocatable :: macro_names(:)
+    real(C_DOUBLE), allocatable :: awrs(:)
     integer :: num_nuclides
     integer :: num_macros
     integer :: j
@@ -172,8 +172,8 @@ contains
     integer :: k
     integer :: n
     integer :: err
-    character(20), allocatable :: nuc_names(:)
-    character(20), allocatable :: macro_names(:)
+    character(kind=C_CHAR, len=20), allocatable :: nuc_names(:)
+    character(kind=C_CHAR, len=20), allocatable :: macro_names(:)
     real(8) :: volume
     real(8), allocatable :: nuc_densities(:)
     integer :: num_nuclides

src/summary.cpp

@@ -2,13 +2,22 @@
 #include "openmc/hdf5_interface.h"
 #include "openmc/lattice.h"
 #include "openmc/surface.h"
-
+#include "openmc/settings.h"
 
 namespace openmc {
 
 extern "C" void
 write_geometry(hid_t file_id) {
+
   auto geom_group = create_group(file_id, "geometry");
+  
+#ifdef DAGMC
+  if (settings::dagmc) {
+    write_attribute(geom_group, "dagmc", 1);
+    return;
+  }
+#endif
+
   write_attribute(geom_group, "n_cells", cells.size());
   write_attribute(geom_group, "n_surfaces", surfaces.size());
   write_attribute(geom_group, "n_universes", universes.size());

src/surface.cpp

@@ -138,6 +138,8 @@ void read_coeffs(pugi::xml_node surf_node, int surf_id, double &c1, double &c2,
 // Surface implementation
 //==============================================================================
 
+Surface::Surface() {} // empty constructor
+
 Surface::Surface(pugi::xml_node surf_node)
 {
   if (check_for_node(surf_node, "id")) {
@@ -207,8 +209,11 @@ Surface::reflect(Position r, Direction u) const
   return u -= (2.0 * projection / magnitude) * n;
 }
 
+CSGSurface::CSGSurface() : Surface{} {};
+CSGSurface::CSGSurface(pugi::xml_node surf_node) : Surface{surf_node} {};
+
 void
-Surface::to_hdf5(hid_t group_id) const
+CSGSurface::to_hdf5(hid_t group_id) const
 {
   std::string group_name {"surface "};
   group_name += std::to_string(id_);
@@ -239,12 +244,62 @@ Surface::to_hdf5(hid_t group_id) const
   close_group(surf_group);
 }
 
+//==============================================================================
+// DAGSurface implementation
+//==============================================================================
+#ifdef DAGMC
+DAGSurface::DAGSurface() : Surface{} {} // empty constructor
+
+double DAGSurface::evaluate(Position r) const
+{
+  return 0.0;
+}
+
+double
+DAGSurface::distance(Position r, Direction u, bool coincident) const
+{
+  moab::ErrorCode rval;
+  moab::EntityHandle surf = dagmc_ptr_->entity_by_id(2, id_);
+  moab::EntityHandle hit_surf;
+  double dist;
+  double pnt[3] = {r.x, r.y, r.z};
+  double dir[3] = {u.x, u.y, u.z};
+  rval = dagmc_ptr_->ray_fire(surf, pnt, dir, hit_surf, dist, NULL, 0, 0);
+  MB_CHK_ERR_CONT(rval);
+  if (dist < 0.0) dist = INFTY;
+  return dist;
+}
+
+Direction DAGSurface::normal(Position r) const
+{
+  moab::ErrorCode rval;
+  Direction u;
+  moab::EntityHandle surf = dagmc_ptr_->entity_by_id(2, id_);
+  double pnt[3] = {r.x, r.y, r.z};
+  double dir[3] = {u.x, u.y, u.z};
+  rval = dagmc_ptr_->get_angle(surf, pnt, dir);
+  MB_CHK_ERR_CONT(rval);
+  return u;
+}
+
+BoundingBox DAGSurface::bounding_box() const
+{
+  moab::ErrorCode rval;
+  moab::EntityHandle surf = dagmc_ptr_->entity_by_id(2, id_);
+  double min[3], max[3];
+  rval = dagmc_ptr_->getobb(surf, min, max);
+  MB_CHK_ERR_CONT(rval);
+  return {min[0], max[0], min[1], max[1], min[2], max[2]};
+}
+
+void DAGSurface::to_hdf5(hid_t group_id) const {}
+#endif
 //==============================================================================
 // PeriodicSurface implementation
 //==============================================================================
 
 PeriodicSurface::PeriodicSurface(pugi::xml_node surf_node)
-  : Surface {surf_node}
+  : CSGSurface {surf_node}
 {
   if (check_for_node(surf_node, "periodic_surface_id")) {
     i_periodic_ = std::stoi(get_node_value(surf_node, "periodic_surface_id"));
@@ -580,7 +635,7 @@ axis_aligned_cylinder_normal(Position r, double offset1, double offset2)
 //==============================================================================
 
 SurfaceXCylinder::SurfaceXCylinder(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, y0_, z0_, radius_);
 }
@@ -614,7 +669,7 @@ void SurfaceXCylinder::to_hdf5_inner(hid_t group_id) const
 //==============================================================================
 
 SurfaceYCylinder::SurfaceYCylinder(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, x0_, z0_, radius_);
 }
@@ -647,7 +702,7 @@ void SurfaceYCylinder::to_hdf5_inner(hid_t group_id) const
 //==============================================================================
 
 SurfaceZCylinder::SurfaceZCylinder(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, x0_, y0_, radius_);
 }
@@ -680,7 +735,7 @@ void SurfaceZCylinder::to_hdf5_inner(hid_t group_id) const
 //==============================================================================
 
 SurfaceSphere::SurfaceSphere(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, x0_, y0_, z0_, radius_);
 }
@@ -833,7 +888,7 @@ axis_aligned_cone_normal(Position r, double offset1, double offset2,
 //==============================================================================
 
 SurfaceXCone::SurfaceXCone(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, x0_, y0_, z0_, radius_sq_);
 }
@@ -866,7 +921,7 @@ void SurfaceXCone::to_hdf5_inner(hid_t group_id) const
 //==============================================================================
 
 SurfaceYCone::SurfaceYCone(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, x0_, y0_, z0_, radius_sq_);
 }
@@ -899,7 +954,7 @@ void SurfaceYCone::to_hdf5_inner(hid_t group_id) const
 //==============================================================================
 
 SurfaceZCone::SurfaceZCone(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, x0_, y0_, z0_, radius_sq_);
 }
@@ -932,7 +987,7 @@ void SurfaceZCone::to_hdf5_inner(hid_t group_id) const
 //==============================================================================
 
 SurfaceQuadric::SurfaceQuadric(pugi::xml_node surf_node)
-  : Surface(surf_node)
+  : CSGSurface(surf_node)
 {
   read_coeffs(surf_node, id_, A_, B_, C_, D_, E_, F_, G_, H_, J_, K_);
 }

src/tallies/tally_filter_universe.F90

@@ -78,7 +78,7 @@ contains
 
     allocate(universe_ids(size(this % universes)))
     do i = 1, size(this % universes)
-      universe_ids(i) = universe_id(this % universes(i)-1)
+      universe_ids(i) = universe_id(this % universes(i))
     end do
     call write_dataset(filter_group, "bins", universe_ids)
   end subroutine to_statepoint_universe

src/thermal.cpp

@@ -1,7 +1,7 @@
 #include "openmc/thermal.h"
 
 #include <algorithm> // for sort, move, min, max, find
-#include <cmath>     // for round, sqrt, fabs
+#include <cmath>     // for round, sqrt, abs
 #include <sstream>   // for stringstream
 
 #include "xtensor/xarray.hpp"
@@ -80,7 +80,7 @@ ThermalScattering::ThermalScattering(hid_t group, const std::vector<double>& tem
 
       auto i_closest = xt::argmin(xt::abs(temps_available - T))[0];
       auto temp_actual = temps_available[i_closest];
-      if (std::fabs(temp_actual - T) < tolerance) {
+      if (std::abs(temp_actual - T) < tolerance) {
         if (std::find(temps_to_read.begin(), temps_to_read.end(), std::round(temp_actual))
             == temps_to_read.end()) {
           temps_to_read.push_back(std::round(temp_actual));
@@ -156,7 +156,7 @@ ThermalScattering::calculate_xs(double E, double sqrtkT, int* i_temp,
   if (settings::temperature_method == TEMPERATURE_NEAREST) {
     // If using nearest temperature, do linear search on temperature
     for (i = 0; i < kTs_.size(); ++i) {
-      if (abs(kTs_[i] - kT) < K_BOLTZMANN*settings::temperature_tolerance) {
+      if (std::abs(kTs_[i] - kT) < K_BOLTZMANN*settings::temperature_tolerance) {
         break;
       }
     }
@@ -574,7 +574,7 @@ ThermalData::sample(const NuclideMicroXS* micro_xs, double E,
   // Because of floating-point roundoff, it may be possible for mu to be
   // outside of the range [-1,1). In these cases, we just set mu to exactly
   // -1 or 1
-  if (std::fabs(*mu) > 1.0) *mu = std::copysign(1.0, *mu);
+  if (std::abs(*mu) > 1.0) *mu = std::copysign(1.0, *mu);
 
 }
 

src/tracking.F90

@@ -7,6 +7,10 @@ module tracking
   use geometry_header,    only: cells
   use geometry,           only: find_cell, distance_to_boundary, cross_lattice,&
                                 check_cell_overlap
+#ifdef DAGMC
+  use geometry,           only: next_cell
+#endif
+
   use material_header,    only: materials, Material
   use message_passing
   use mgxs_interface
@@ -309,6 +313,7 @@ contains
     real(8) :: norm       ! "norm" of surface normal
     real(8) :: xyz(3)     ! Saved global coordinate
     integer :: i_surface  ! index in surfaces
+    integer :: i_cell     ! index of new cell
     logical :: rotational ! if rotational periodic BC applied
     logical :: found      ! particle found in universe?
     class(Surface), pointer :: surf
@@ -467,6 +472,21 @@ contains
     ! ==========================================================================
     ! SEARCH NEIGHBOR LISTS FOR NEXT CELL
 
+#ifdef DAGMC
+    if (dagmc) then
+      i_cell = next_cell(cells(p % last_cell(1) + 1), surfaces(abs(p % surface)))
+      ! save material and temp
+      p % last_material = p % material
+      p % last_sqrtkT = p % sqrtKT
+      ! set new cell value
+      p % coord(1) % cell = i_cell-1 ! decrement for C++ indexing
+      p % cell_instance = 1
+      p % material = cells(i_cell) % material(1)
+      p % sqrtKT = cells(i_cell) % sqrtKT(1)
+      return
+    end if
+#endif
+
     call find_cell(p, found, p % surface)
     if (found) return
 

src/xml_interface.F90

@@ -12,6 +12,7 @@ module xml_interface
   public :: check_for_node
   public :: get_node_list
   public :: get_node_value
+  public :: get_node_value_bool
   public :: get_node_array
   public :: node_value_string
   public :: node_word_count

src/xml_interface.cpp

@@ -53,6 +53,7 @@ get_node_value_bool(pugi::xml_node node, const char* name)
             << node.name() << "\" XML node";
     fatal_error(err_msg);
   }
+  return false;
 }
 
 } // namespace openmc

tests/regression_tests/dagmc/materials.xml

@@ -0,0 +1,16 @@
+<?xml version="1.0"?>
+<materials>
+
+  <material id="40">
+    <density value="11" units="g/cc" />
+    <nuclide name="U235" ao="1.0" />
+  </material>
+
+  <material id="41">
+    <density value="1.0" units="g/cc" />
+    <nuclide name="H1" ao="2.0" />
+    <nuclide name="O16" ao="1.0" />
+    <sab name="c_H_in_H2O"/>
+  </material>
+  
+</materials>

tests/regression_tests/dagmc/results_true.dat

@@ -0,0 +1,5 @@
+k-combined:
+1.115067E+00 5.423808E-02
+tally 1:
+8.543144E+00
+1.530584E+01

tests/regression_tests/dagmc/settings.xml

@@ -0,0 +1,16 @@
+<?xml version="1.0"?>
+<settings>
+
+  <run_mode>eigenvalue</run_mode>
+  <dagmc>true</dagmc>
+  <batches>5</batches>
+  <inactive>0</inactive>
+  <particles>100</particles>
+  <!-- Starting source -->
+  <source>
+    <space type="box">
+      <parameters>-4 -4 -4  4  4  4</parameters>
+    </space>
+  </source>
+
+</settings>

tests/regression_tests/dagmc/tallies.xml

@@ -0,0 +1,13 @@
+<?xml version="1.0"?>
+<tallies>
+
+  <filter id="1" type="cell">
+    <bins>1</bins>
+  </filter>
+
+  <tally id="1">
+    <filters>1</filters>
+    <scores>total </scores>
+  </tally>
+
+</tallies>

tests/regression_tests/dagmc/test.py

@@ -0,0 +1,12 @@
+from tests.testing_harness import TestHarness
+import os
+import pytest
+import openmc
+
+pytestmark = pytest.mark.skipif(
+    not openmc.capi.dagmc_enabled,
+    reason="DAGMC CAD geometry is not enabled.")
+
+def test_dagmc():
+    harness = TestHarness('statepoint.5.h5')
+    harness.main()    

tools/ci/download-xs.sh

@@ -0,0 +1,19 @@
+#!/bin/bash
+set -ex
+
+# Download NNDC HDF5 data
+if [[ ! -e $HOME/nndc_hdf5/cross_sections.xml ]]; then
+    wget https://anl.box.com/shared/static/a0eflty17atnpd0pp7460exagr3nuhm7.xz -O - | tar -C $HOME -xvJ
+fi
+
+# Download ENDF/B-VII.1 distribution
+ENDF=$HOME/endf-b-vii.1/
+if [[ ! -d $ENDF/neutrons || ! -d $ENDF/photoat || ! -d $ENDF/atomic_relax ]]; then
+    wget https://anl.box.com/shared/static/4kd2gxnf4gtk4w1c8eua5fsua22kvgjb.xz -O - | tar -C $HOME -xvJ
+fi
+
+# Download multipole library
+if [[ ! -e $HOME/WMP_Library/092235.h5 ]]; then
+    wget https://github.com/mit-crpg/WMP_Library/releases/download/v1.0/WMP_Library_v1.0.tar.gz
+    tar -C $HOME -xzvf WMP_Library_v1.0.tar.gz
+fi

tools/ci/travis-before-script.sh

@@ -5,19 +5,5 @@ set -ex
 # https://docs.travis-ci.com/user/gui-and-headless-browsers/#Using-xvfb-to-Run-Tests-That-Require-a-GUI
 sh -e /etc/init.d/xvfb start
 
-# Download NNDC HDF5 data
-if [[ ! -e $HOME/nndc_hdf5/cross_sections.xml ]]; then
-    wget https://anl.box.com/shared/static/a0eflty17atnpd0pp7460exagr3nuhm7.xz -O - | tar -C $HOME -xvJ
-fi
-
-# Download ENDF/B-VII.1 distribution
-ENDF=$HOME/endf-b-vii.1/
-if [[ ! -d $ENDF/neutrons || ! -d $ENDF/photoat || ! -d $ENDF/atomic_relax ]]; then
-    wget https://anl.box.com/shared/static/4kd2gxnf4gtk4w1c8eua5fsua22kvgjb.xz -O - | tar -C $HOME -xvJ
-fi
-
-# Download multipole library
-if [[ ! -e $HOME/WMP_Library/092235.h5 ]]; then
-    wget https://github.com/mit-crpg/WMP_Library/releases/download/v1.0/WMP_Library_v1.0.tar.gz
-    tar -C $HOME -xzvf WMP_Library_v1.0.tar.gz
-fi
+# Download NNDC HDF5 data, ENDF/B-VII.1 distribution, multipole library
+sh ./tools/ci/download-xs.sh

tools/ci/travis-install-dagmc.sh

@@ -0,0 +1,36 @@
+
+#!/bin/bash
+set -ex
+
+# MOAB Variables
+MOAB_BRANCH='Version5.0'
+MOAB_REPO='https://bitbucket.org/fathomteam/moab/'
+MOAB_INSTALL_DIR=$HOME/MOAB/
+
+# DAGMC Variables
+DAGMC_BRANCH='develop'
+DAGMC_REPO='https://github.com/svalinn/dagmc'
+DAGMC_INSTALL_DIR=$HOME/DAGMC/
+
+CURRENT_DIR=$(pwd)
+
+# MOAB Install
+cd $HOME
+mkdir MOAB && cd MOAB
+git clone -b $MOAB_BRANCH $MOAB_REPO
+mkdir build && cd build
+cmake ../moab -DENABLE_HDF5=ON -DCMAKE_INSTALL_PREFIX=$MOAB_INSTALL_DIR
+make -j && make -j test install
+rm -rf $HOME/MOAB/moab
+export LD_LIBRARY_PATH=$MOAB_INSTALL_DIR/lib:$LD_LIBRARY_PATH
+
+# DAGMC Install
+mkdir DAGMC && cd DAGMC
+git clone -b $DAGMC_BRANCH $DAGMC_REPO
+mkdir build && cd build
+cmake ../dagmc -DBUILD_TALLY=ON -DCMAKE_INSTALL_PREFIX=$DAGMC_INSTALL_DIR
+make -j install
+rm -rf $HOME/DAGMC/dagmc
+export LD_LIBRARY_PATH=$DAGMC_INSTALL_DIR/lib:$LD_LIBRARY_PATH
+
+cd $CURRENT_DIR

tools/ci/travis-install.py

@@ -2,7 +2,6 @@ import os
 import shutil
 import subprocess
 
-
 def which(program):
     def is_exe(fpath):
         return os.path.isfile(fpath) and os.access(fpath, os.X_OK)
@@ -20,7 +19,7 @@ def which(program):
     return None
 
 
-def install(omp=False, mpi=False, phdf5=False):
+def install(omp=False, mpi=False, phdf5=False, dagmc=False):
     # Create build directory and change to it
     shutil.rmtree('build', ignore_errors=True)
     os.mkdir('build')
@@ -48,23 +47,26 @@ def install(omp=False, mpi=False, phdf5=False):
     else:
         cmake_cmd.append('-DHDF5_PREFER_PARALLEL=OFF')
 
+    if dagmc:
+        cmake_cmd.append('-Ddagmc=ON')
+        
     # Build and install
     cmake_cmd.append('..')
     print(' '.join(cmake_cmd))
     subprocess.check_call(cmake_cmd)
-    subprocess.check_call(['make', '-j'])
+    subprocess.check_call(['make', '-j4'])
     subprocess.check_call(['sudo', 'make', 'install'])
 
-
 def main():
     # Convert Travis matrix environment variables into arguments for install()
     omp = (os.environ.get('OMP') == 'y')
     mpi = (os.environ.get('MPI') == 'y')
     phdf5 = (os.environ.get('PHDF5') == 'y')
 
-    # Build and install
-    install(omp, mpi, phdf5)
+    dagmc = (os.environ.get('DAGMC') == 'y')
 
+    # Build and install
+    install(omp, mpi, phdf5, dagmc)
 
 if __name__ == '__main__':
     main()

tools/ci/travis-install.sh

@@ -4,6 +4,11 @@ set -ex
 # Install NJOY 2016
 ./tools/ci/travis-install-njoy.sh
 
+# Install DAGMC if needed
+if [[ $DAGMC = 'y' ]]; then
+    ./tools/ci/travis-install-dagmc.sh
+fi
+
 # Upgrade pip before doing anything else
 pip install --upgrade pip