Normalize fixed source calculations by total external source strength

src/tallies/tally_header.F90

@@ -9,7 +9,8 @@ module tally_header
   use dict_header,         only: DictIntInt
   use message_passing,     only: n_procs
   use nuclide_header,      only: nuclide_dict
-  use settings,            only: reduce_tallies
+  use settings,            only: reduce_tallies, run_mode
+  use source_header,       only: external_source
   use stl_vector,          only: VectorInt
   use string,              only: to_lower, to_f_string, str_to_int, to_str
   use tally_filter_header, only: TallyFilterContainer, filters, n_filters
@@ -161,6 +162,7 @@ contains
 
     integer :: i, j
     real(C_DOUBLE) :: val
+    real(C_DOUBLE) :: total_source
 
     ! Increment number of realizations
     if (reduce_tallies) then
@@ -169,10 +171,17 @@ contains
       this % n_realizations = this % n_realizations + n_procs
     end if
 
+    ! Calculate total source strength for normalization
+    if (run_mode == MODE_FIXEDSOURCE) then
+      total_source = sum(external_source(:) % strength)
+    else
+      total_source = ONE
+    end if
+
     ! Accumulate each result
     do j = 1, size(this % results, 3)
       do i = 1, size(this % results, 2)
-        val = this % results(RESULT_VALUE, i, j)/total_weight
+        val = this % results(RESULT_VALUE, i, j)/total_weight * total_source
         this % results(RESULT_VALUE, i, j) = ZERO
 
         this % results(RESULT_SUM, i, j) = &

tests/test_fixed_source/geometry.xml

@@ -1,8 +0,0 @@
-<?xml version="1.0"?>
-<geometry>
-
-  <!-- Sphere with radius 10 -->
-  <surface id="1" type="sphere" coeffs="0 0 0 10" boundary="vacuum"/>
-  <cell id="1" material="1" region="-1" />
-    
-</geometry>

tests/test_fixed_source/inputs_true.dat

@@ -0,0 +1,31 @@
+<?xml version='1.0' encoding='utf-8'?>
+<geometry>
+  <cell id="1" material="1" region="-1" universe="1" />
+  <surface boundary="vacuum" coeffs="0.0 0.0 0.0 10.0" id="1" type="sphere" />
+</geometry>
+<?xml version='1.0' encoding='utf-8'?>
+<materials>
+  <material id="1">
+    <density units="g/cc" value="7.5" />
+    <nuclide ao="1.0" name="O16" />
+    <nuclide ao="0.0001" name="U238" />
+  </material>
+</materials>
+<?xml version='1.0' encoding='utf-8'?>
+<settings>
+  <run_mode>fixed source</run_mode>
+  <particles>100</particles>
+  <batches>10</batches>
+  <source strength="10.0">
+    <space type="point">
+      <parameters>0.0 0.0 0.0</parameters>
+    </space>
+  </source>
+  <temperature_default>294</temperature_default>
+</settings>
+<?xml version='1.0' encoding='utf-8'?>
+<tallies>
+  <tally id="1">
+    <scores>flux</scores>
+  </tally>
+</tallies>

tests/test_fixed_source/materials.xml

@@ -1,10 +0,0 @@
-<?xml version="1.0"?>
-<materials>
-
-  <material id="1">
-    <density value="7.5" units="g/cc" />
-    <nuclide name="O16" ao="1.0" />
-    <nuclide name="U238" ao="0.0001" />
-  </material>
-
-</materials>

tests/test_fixed_source/results_true.dat

@@ -1,6 +1,6 @@
 tally 1:
-4.448476E+02
-1.984373E+04
+4.448476E+03
+1.984373E+06
 leakage:
 9.790000E+00
 9.588300E+00

tests/test_fixed_source/settings.xml

@@ -1,14 +0,0 @@
-<?xml version="1.0"?>
-<settings>
-
-  <run_mode>fixed source</run_mode>
-  <batches>10</batches>
-  <particles>100</particles>
-
-  <temperature_default>294</temperature_default>
-
-  <source>
-    <space type="point" parameters="0 0 0" />
-  </source>
-
-</settings>

tests/test_fixed_source/tallies.xml

@@ -1,8 +0,0 @@
-<?xml version="1.0"?>
-<tallies>
-
-  <tally id="1">
-    <scores>flux</scores>
-  </tally>
-
-</tallies>

tests/test_fixed_source/test_fixed_source.py

@@ -5,17 +5,18 @@ import os
 import sys
 import numpy as np
 sys.path.insert(0, os.pardir)
-from testing_harness import TestHarness
-from openmc import StatePoint
+from testing_harness import PyAPITestHarness
+import openmc
+import openmc.stats
 
 
-class FixedSourceTestHarness(TestHarness):
+class FixedSourceTestHarness(PyAPITestHarness):
     def _get_results(self):
         """Digest info in the statepoint and return as a string."""
         # Read the statepoint file.
         statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name))[0]
         outstr = ''
-        with StatePoint(statepoint) as sp:
+        with openmc.StatePoint(statepoint) as sp:
             # Write out tally data.
             for i, tally_ind in enumerate(sp.tallies):
                 tally = sp.tallies[tally_ind]
@@ -36,5 +37,28 @@ class FixedSourceTestHarness(TestHarness):
 
 
 if __name__ == '__main__':
-    harness = FixedSourceTestHarness('statepoint.10.h5')
+    mat = openmc.Material()
+    mat.add_nuclide('O16', 1.0)
+    mat.add_nuclide('U238', 0.0001)
+    mat.set_density('g/cc', 7.5)
+
+    surf = openmc.Sphere(R=10.0, boundary_type='vacuum')
+    cell = openmc.Cell(fill=mat, region=-surf)
+
+    model = openmc.model.Model()
+    model.geometry.root_universe = openmc.Universe(cells=[cell])
+    model.materials.append(mat)
+
+    model.settings.run_mode = 'fixed source'
+    model.settings.batches = 10
+    model.settings.particles = 100
+    model.settings.temperature = {'default': 294}
+    model.settings.source = openmc.Source(space=openmc.stats.Point(),
+                                          strength=10.0)
+
+    tally = openmc.Tally()
+    tally.scores = ['flux']
+    model.tallies.append(tally)
+
+    harness = FixedSourceTestHarness('statepoint.10.h5', model)
     harness.main()