Added regression test for photon production

src/input_xml.F90

@@ -3641,7 +3641,9 @@ contains
     allocate(nuclides(n_nuclides))
     allocate(elements(n_elements))
     allocate(sab_tables(n_sab_tables))
-    if (electron_treatment == ELECTRON_TTB) allocate(ttb(n_materials))
+    if (photon_transport .and. electron_treatment == ELECTRON_TTB) then
+      allocate(ttb(n_materials))
+    end if
 
     ! Read cross sections
     do i = 1, size(materials)

src/physics.F90

@@ -1724,7 +1724,7 @@ contains
     integer :: i
 
     ! Sample the number of photons produced
-    nu_t =  p % wgt / keff * micro_xs(i_nuclide) % photon_prod / &
+    nu_t =  p % wgt * micro_xs(i_nuclide) % photon_prod / &
          micro_xs(i_nuclide) % total
     if (prn() > nu_t - int(nu_t)) then
       nu = int(nu_t)

tests/regression_tests/photon_production/inputs_true.dat

@@ -0,0 +1,50 @@
+<?xml version='1.0' encoding='utf-8'?>
+<geometry>
+  <cell id="13" material="void" region="-9 10 -11" universe="9" />
+  <cell id="14" material="13" region="-9 11 -12" universe="9" />
+  <cell id="15" material="void" region="~(-9 10 -12)" universe="9" />
+  <surface boundary="vacuum" coeffs="0.0 0.0 1e-06" id="9" type="x-cylinder" />
+  <surface boundary="vacuum" coeffs="-1.0" id="10" type="x-plane" />
+  <surface coeffs="1.0" id="11" type="x-plane" />
+  <surface boundary="vacuum" coeffs="11.0" id="12" type="x-plane" />
+</geometry>
+<?xml version='1.0' encoding='utf-8'?>
+<materials>
+  <material id="13">
+    <density units="g/cm3" value="2.6989" />
+    <nuclide ao="1.0" name="Al27" />
+  </material>
+</materials>
+<?xml version='1.0' encoding='utf-8'?>
+<settings>
+  <run_mode>fixed source</run_mode>
+  <particles>10000</particles>
+  <batches>1</batches>
+  <source strength="1.0">
+    <space type="point">
+      <parameters>0 0 0</parameters>
+    </space>
+    <angle reference_uvw="1.0 0.0 0.0" type="monodirectional" />
+    <energy type="discrete">
+      <parameters>14.0 1.0</parameters>
+    </energy>
+  </source>
+  <electron_treatment>ttb</electron_treatment>
+  <photon_transport>true</photon_transport>
+  <cutoff>
+    <energy_photon>1000.0</energy_photon>
+  </cutoff>
+</settings>
+<?xml version='1.0' encoding='utf-8'?>
+<tallies>
+  <filter id="1" type="cell">
+    <bins>14</bins>
+  </filter>
+  <filter id="2" type="particle">
+    <bins>2</bins>
+  </filter>
+  <tally id="1">
+    <filters>1 2</filters>
+    <scores>flux</scores>
+  </tally>
+</tallies>

tests/regression_tests/photon_production/results_true.dat

@@ -0,0 +1,3 @@
+tally 1:
+sum = 1.371553E-08
+sum_sq = 1.881158E-16

tests/regression_tests/photon_production/test.py

@@ -0,0 +1,70 @@
+from math import pi
+
+import numpy as np
+import openmc
+
+from tests.testing_harness import PyAPITestHarness
+
+
+class SourceTestHarness(PyAPITestHarness):
+    def _build_inputs(self):
+        mat = openmc.Material()
+        mat.set_density('g/cm3', 2.6989)
+        mat.add_nuclide('Al27', 1.0)
+        materials = openmc.Materials([mat])
+        materials.export_to_xml()
+
+        cyl = openmc.XCylinder(boundary_type='vacuum', R=1.0e-6)
+        x_plane_left = openmc.XPlane(boundary_type='vacuum', x0=-1.0)
+        x_plane_center = openmc.XPlane(boundary_type='transmission', x0=1.0)
+        x_plane_right = openmc.XPlane(boundary_type='vacuum', x0=11.0)
+
+        inner_cyl_left = openmc.Cell()
+        inner_cyl_right = openmc.Cell()
+        outer_cyl = openmc.Cell()
+
+        inner_cyl_left.region = -cyl & +x_plane_left & -x_plane_center
+        inner_cyl_right.region = -cyl & +x_plane_center & -x_plane_right
+        outer_cyl.region = ~(-cyl & +x_plane_left & -x_plane_right)
+        inner_cyl_right.fill = mat
+        geometry = openmc.Geometry([inner_cyl_left, inner_cyl_right, outer_cyl])
+        geometry.export_to_xml()
+
+        source = openmc.Source()
+        source.space = openmc.stats.Point((0,0,0))
+        source.angle = openmc.stats.Monodirectional()
+        source.energy = openmc.stats.Discrete([14.0], [1.0])
+        source.particle = 'neutron'
+
+        settings = openmc.Settings()
+        settings.particles = 10000
+        settings.run_mode = 'fixed source'
+        settings.batches = 1
+        settings.photon_transport = True
+        settings.electron_treatment = 'ttb'
+        settings.cutoff = {'energy_photon' : 1000.0}
+        settings.source = source
+        settings.export_to_xml()
+
+        cell_filter = openmc.CellFilter(inner_cyl_right)
+        particle_filter = openmc.ParticleFilter('photon')
+        tally = openmc.Tally()
+        tally.filters = [cell_filter, particle_filter]
+        tally.scores = ['flux']
+        tallies = openmc.Tallies([tally])
+        tallies.export_to_xml()
+
+    def _get_results(self):
+        with openmc.StatePoint(self._sp_name) as sp:
+            outstr = ''
+            t = sp.get_tally()
+            outstr += 'tally {}:\n'.format(t.id)
+            outstr += 'sum = {:12.6E}\n'.format(t.sum[0, 0, 0])
+            outstr += 'sum_sq = {:12.6E}\n'.format(t.sum_sq[0, 0, 0])
+
+            return outstr
+
+
+def test_source():
+    harness = SourceTestHarness('statepoint.1.h5')
+    harness.main()