allowing no cross_section in volume calculation

src/cross_sections.cpp

@@ -115,7 +115,7 @@ void read_cross_sections_xml(pugi::xml_node root)
   if (!check_for_node(root, "cross_sections")) {
     // No cross_sections.xml file specified in settings.xml, check
     // environment variable
-    if (settings::run_CE) {
+    if (settings::run_CE and settings::run_mode != RunMode::VOLUME) {
       char* envvar = std::getenv("OPENMC_CROSS_SECTIONS");
       if (!envvar) {
         fatal_error(
@@ -127,6 +127,20 @@ void read_cross_sections_xml(pugi::xml_node root)
           "information on how to set up data libraries.");
       }
       settings::path_cross_sections = envvar;
+    } else if (settings::run_mode == RunMode::VOLUME) {
+      char* envvar = std::getenv("OPENMC_CROSS_SECTIONS");
+      if (!envvar) {
+        warning("No cross_sections.xml file was specified in "
+                "materials.xml or in the OPENMC_CROSS_SECTIONS"
+                " environment variable. For the full extend of"
+                " volume calculation OpenMC needs such a file to identify"
+                " where to find data libraries. Some volume calculation"
+                " capability will be disable.");
+        settings::path_cross_sections = "";
+        return;
+      } else {
+        settings::path_cross_sections = envvar;
+      }
     } else {
       char* envvar = std::getenv("OPENMC_MG_CROSS_SECTIONS");
       if (!envvar) {
@@ -330,7 +344,8 @@ void read_ce_cross_sections_xml()
 
 void finalize_cross_sections()
 {
-  if (settings::run_mode != RunMode::PLOTTING) {
+  if (settings::run_mode != RunMode::PLOTTING and
+      settings::path_cross_sections != "") {
     simulation::time_read_xs.start();
     if (settings::run_CE) {
       // Determine desired temperatures for each nuclide and S(a,b) table

src/material.cpp

@@ -219,17 +219,22 @@ Material::Material(pugi::xml_node node)
   if (settings::photon_transport)
     element_.reserve(n);
 
+  // Loading data is not necessary in plotting mode or in volume mode when no
+  // cross section source is specified.
+  bool require_data = !(settings::run_mode == RunMode::PLOTTING ||
+                        (settings::run_mode == RunMode::VOLUME &&
+                          settings::path_cross_sections.empty()));
+
   for (int i = 0; i < n; ++i) {
     const auto& name {names[i]};
 
     // Check that this nuclide is listed in the nuclear data library
     // (cross_sections.xml for CE and the MGXS HDF5 for MG)
-    if (settings::run_mode != RunMode::PLOTTING) {
+    if (require_data) {
       LibraryKey key {Library::Type::neutron, name};
       if (data::library_map.find(key) == data::library_map.end()) {
-        fatal_error("Could not find nuclide " + name +
-                    " in the "
-                    "nuclear data library.");
+        fatal_error(
+          "Could not find nuclide " + name + " in the nuclear data library.");
       }
     }
 
@@ -249,7 +254,7 @@ Material::Material(pugi::xml_node node)
       std::string element = to_element(name);
 
       // Make sure photon cross section data is available
-      if (settings::run_mode != RunMode::PLOTTING) {
+      if (require_data) {
         LibraryKey key {Library::Type::photon, element};
         if (data::library_map.find(key) == data::library_map.end()) {
           fatal_error(
@@ -328,7 +333,7 @@ Material::Material(pugi::xml_node node)
 
       // Check that the thermal scattering table is listed in the
       // cross_sections.xml file
-      if (settings::run_mode != RunMode::PLOTTING) {
+      if (require_data) {
         LibraryKey key {Library::Type::thermal, name};
         if (data::library_map.find(key) == data::library_map.end()) {
           fatal_error("Could not find thermal scattering data " + name +
@@ -1082,20 +1087,29 @@ void Material::to_hdf5(hid_t group) const
   vector<std::string> nuc_names;
   vector<std::string> macro_names;
   vector<double> nuc_densities;
-  if (settings::run_CE) {
-    for (int i = 0; i < nuclide_.size(); ++i) {
-      int i_nuc = nuclide_[i];
-      nuc_names.push_back(data::nuclides[i_nuc]->name_);
-      nuc_densities.push_back(atom_density_(i));
-    }
-  } else {
-    for (int i = 0; i < nuclide_.size(); ++i) {
-      int i_nuc = nuclide_[i];
-      if (data::mg.nuclides_[i_nuc].awr != MACROSCOPIC_AWR) {
-        nuc_names.push_back(data::mg.nuclides_[i_nuc].name);
+
+  // Loading data is not necessary in plotting mode or in volume mode when no
+  // cross section source is specified.
+  bool require_data = !(settings::run_mode == RunMode::PLOTTING ||
+                        (settings::run_mode == RunMode::VOLUME &&
+                          settings::path_cross_sections.empty()));
+
+  if (require_data) {
+    if (settings::run_CE) {
+      for (int i = 0; i < nuclide_.size(); ++i) {
+        int i_nuc = nuclide_[i];
+        nuc_names.push_back(data::nuclides[i_nuc]->name_);
         nuc_densities.push_back(atom_density_(i));
-      } else {
-        macro_names.push_back(data::mg.nuclides_[i_nuc].name);
+      }
+    } else {
+      for (int i = 0; i < nuclide_.size(); ++i) {
+        int i_nuc = nuclide_[i];
+        if (data::mg.nuclides_[i_nuc].awr != MACROSCOPIC_AWR) {
+          nuc_names.push_back(data::mg.nuclides_[i_nuc].name);
+          nuc_densities.push_back(atom_density_(i));
+        } else {
+          macro_names.push_back(data::mg.nuclides_[i_nuc].name);
+        }
       }
     }
   }
@@ -1111,7 +1125,7 @@ void Material::to_hdf5(hid_t group) const
     write_dataset(material_group, "macroscopics", macro_names);
   }
 
-  if (!thermal_tables_.empty()) {
+  if (!thermal_tables_.empty() and require_data) {
     vector<std::string> sab_names;
     for (const auto& table : thermal_tables_) {
       sab_names.push_back(data::thermal_scatt[table.index_table]->name_);

src/tallies/tally.cpp

@@ -855,7 +855,9 @@ void read_tallies_xml(pugi::xml_node root)
   read_meshes(root);
 
   // We only need the mesh info for plotting
-  if (settings::run_mode == RunMode::PLOTTING)
+  if (settings::run_mode == RunMode::PLOTTING or
+      (settings::run_mode == RunMode::VOLUME and
+        settings::path_cross_sections == ""))
     return;
 
   // Read data for tally derivatives

src/volume_calc.cpp

@@ -325,12 +325,16 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
           if (i_material == MATERIAL_VOID)
             continue;
 
-          const auto& mat = model::materials[i_material];
-          for (int k = 0; k < mat->nuclide_.size(); ++k) {
-            // Accumulate nuclide density
-            int i_nuclide = mat->nuclide_[k];
-            atoms(i_nuclide, 0) += mat->atom_density_[k] * f;
-            atoms(i_nuclide, 1) += std::pow(mat->atom_density_[k], 2) * var_f;
+          // Guards against atom density calculation when no cross_section is
+          // presents
+          if (!settings::path_cross_sections.empty()) {
+            const auto& mat = model::materials[i_material];
+            for (int k = 0; k < mat->nuclide_.size(); ++k) {
+              // Accumulate nuclide density
+              int i_nuclide = mat->nuclide_[k];
+              atoms(i_nuclide, 0) += mat->atom_density_[k] * f;
+              atoms(i_nuclide, 1) += std::pow(mat->atom_density_[k], 2) * var_f;
+            }
           }
         }
 
@@ -365,17 +369,21 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
           }
         }
 
-        for (int j = 0; j < n_nuc; ++j) {
-          // Determine total number of atoms. At this point, we have values in
-          // atoms/b-cm. To get to atoms we multiply by 10^24 V.
-          double mean = 1.0e24 * volume_sample * atoms(j, 0);
-          double stdev = 1.0e24 * volume_sample * std::sqrt(atoms(j, 1));
-
-          // Convert full arrays to vectors
-          if (mean > 0.0) {
-            result.nuclides.push_back(j);
-            result.atoms.push_back(mean);
-            result.uncertainty.push_back(stdev);
+        // Guards against atom density calculation when no cross_section is
+        // presents
+        if (!settings::path_cross_sections.empty()) {
+          for (int j = 0; j < n_nuc; ++j) {
+            // Determine total number of atoms. At this point, we have values in
+            // atoms/b-cm. To get to atoms we multiply by 10^24 V.
+            double mean = 1.0e24 * volume_sample * atoms(j, 0);
+            double stdev = 1.0e24 * volume_sample * std::sqrt(atoms(j, 1));
+
+            // Convert full arrays to vectors
+            if (mean > 0.0) {
+              result.nuclides.push_back(j);
+              result.atoms.push_back(mean);
+              result.uncertainty.push_back(stdev);
+            }
           }
         }
       }
@@ -481,15 +489,17 @@ void VolumeCalculation::to_hdf5(
                                           : data::mg.nuclides_[i_nuc].name);
     }
 
-    // Create array of total # of atoms with uncertainty for each nuclide
-    xt::xtensor<double, 2> atom_data({n_nuc, 2});
-    xt::view(atom_data, xt::all(), 0) = xt::adapt(result.atoms);
-    xt::view(atom_data, xt::all(), 1) = xt::adapt(result.uncertainty);
-
-    // Write results
-    write_dataset(group_id, "nuclides", nucnames);
-    write_dataset(group_id, "atoms", atom_data);
+    // Guards against atom density calculation when no cross_section is presents
+    if (!settings::path_cross_sections.empty()) {
+      // Create array of total # of atoms with uncertainty for each nuclide
+      xt::xtensor<double, 2> atom_data({n_nuc, 2});
+      xt::view(atom_data, xt::all(), 0) = xt::adapt(result.atoms);
+      xt::view(atom_data, xt::all(), 1) = xt::adapt(result.uncertainty);
 
+      // Write results
+      write_dataset(group_id, "nuclides", nucnames);
+      write_dataset(group_id, "atoms", atom_data);
+    }
     close_group(group_id);
   }
 

tests/unit_tests/test_volume.py

@@ -43,3 +43,79 @@ def test_no_bcs(run_in_tmpdir):
 
     model.settings.volume_calculations = [vc]
     model.calculate_volumes()
+
+
+def test_volume_no_cross_section(run_in_tmpdir):
+    # store config to allow reset
+    openmc_config = openmc.config['cross_sections']
+
+    # removing path to cross section
+    del openmc.config['cross_sections']
+
+    # setting the simulation
+    uo2 = openmc.Material(1, "uo2")
+    mat = openmc.Material()
+    # Add nuclides to uo2
+    uo2.add_nuclide('U235', 0.03)
+    uo2.add_nuclide('U238', 0.97)
+    uo2.add_nuclide('O16', 2.0)
+    uo2.set_density('g/cm3', 10.0)
+    zirconium = openmc.Material(name="zirconium")
+    zirconium.add_element('Zr', 1.0)
+    zirconium.set_density('g/cm3', 6.6)
+
+    water = openmc.Material(name="h2o")
+    water.add_nuclide('H1', 2.0)
+    water.add_nuclide('O16', 1.0)
+    water.set_density('g/cm3', 1.0)
+    water.add_s_alpha_beta('c_H_in_H2O')
+
+    materials = openmc.Materials([uo2, zirconium, water])
+    materials.export_to_xml()
+
+    fuel_outer_radius = openmc.ZCylinder(r=0.39)
+    clad_inner_radius = openmc.ZCylinder(r=0.40)
+    clad_outer_radius = openmc.ZCylinder(r=0.46)
+    fuel_region = -fuel_outer_radius
+    gap_region = +fuel_outer_radius & -clad_inner_radius
+    clad_region = +clad_inner_radius & -clad_outer_radius
+
+    fuel = openmc.Cell(name='fuel')
+    fuel.fill = uo2
+    fuel.region = fuel_region
+
+    gap = openmc.Cell(name='air gap')
+    gap.region = gap_region
+
+    clad = openmc.Cell(name='clad')
+    clad.fill = zirconium
+    clad.region = clad_region
+
+    pitch = 1.26
+    left = openmc.XPlane(-pitch/2, boundary_type='reflective')
+    right = openmc.XPlane(pitch/2, boundary_type='reflective')
+    bottom = openmc.YPlane(-pitch/2, boundary_type='reflective')
+    top = openmc.YPlane(pitch/2, boundary_type='reflective')
+
+    lower_left = (-pitch/2, -pitch/2, -pitch/2)
+    upper_right = (pitch/2, pitch/2, pitch/2)
+
+    water_region = +left & -right & +bottom & -top & +clad_outer_radius
+
+    moderator = openmc.Cell(name='moderator')
+    moderator.fill = water
+    moderator.region = water_region
+
+    root_universe = openmc.Universe(cells=(fuel, gap, clad, moderator))
+    geometry = openmc.Geometry(root_universe)
+    geometry.export_to_xml()
+    volumes_calc = openmc.VolumeCalculation([fuel, gap, clad, moderator], 100, lower_left, upper_right)
+    settings = openmc.Settings()
+    settings.volume_calculations = [volumes_calc]
+    settings.export_to_xml()
+
+    # This should not fail
+    openmc.calculate_volumes()
+
+    #restore config:
+    openmc.config['cross_sections'] = openmc_config