add regression test for ce physics

.gitignore

@@ -55,3 +55,5 @@ docs/build
 docs/source/pythonapi/generated/
 
 *.csv
+
+dummy_nuclide.h5

mcdc/input_.py

@@ -322,6 +322,9 @@ def material(
         # Set ID
         card.ID = len(global_.input_deck.materials)
 
+        # Default values
+        card.J = 6
+
         # Set the nuclides
         for i in range(N_nuclide):
             nuc_name = nuclides[i][0]
@@ -335,6 +338,9 @@ def material(
                 # Set ID
                 nuc_card.ID = len(global_.input_deck.nuclides)
 
+                # Default values
+                nuc_card.J = 6
+
                 # Check if the nuclide is available in the nuclear data library
                 dir_name = os.getenv("MCDC_XSLIB")
                 if dir_name == None:

test/regression/slab_ce/input.py

@@ -0,0 +1,72 @@
+import numpy as np
+import os, h5py
+
+import mcdc
+
+
+# Set the XS library directory
+os.environ["MCDC_XSLIB"] = os.getcwd()
+
+# Create the dummy nuclide
+with h5py.File("dummy_nuclide.h5", 'w') as f:
+    f['A'] = 1.0
+
+    f['E_xs'] = np.array([0.0, 1.0 - 1E-6, 1.0 + 1E-6, 2E7])
+    f['capture'] = np.array([0.01344, 0.01344, 0.00384, 0.00384])
+    f['fission'] = np.array([0.06912, 0.06912, 0.00619, 0.00619])
+    f['scatter'] = np.array([0.26304, 0.26304, 0.15024, 0.15024])
+
+    f['E_nu_p'] = np.array([0.0, 1.0 - 1E-6, 1.0 + 1E-6, 2E7])
+    f['nu_p'] = np.array([2.5, 2.5, 2.7, 2.7])
+
+    f['E_chi_p'] = np.array([0.0, 1E5, 2E7])
+    f['chi_p'] = np.array([0.0, 0.0, 1.0])
+
+    f['decay_rate'] = np.zeros(6)
+
+    f['E_nu_d'] = np.array([0.0, 2E7])
+    f['nu_d'] = np.zeros((6,2))
+
+    f['E_chi_d1'] = np.zeros(0)
+    f['E_chi_d2'] = np.zeros(0)
+    f['E_chi_d3'] = np.zeros(0)
+    f['E_chi_d4'] = np.zeros(0)
+    f['E_chi_d5'] = np.zeros(0)
+    f['E_chi_d6'] = np.zeros(0)
+    f['chi_d1'] = np.zeros(0)
+    f['chi_d2'] = np.zeros(0)
+    f['chi_d3'] = np.zeros(0)
+    f['chi_d4'] = np.zeros(0)
+    f['chi_d5'] = np.zeros(0)
+    f['chi_d6'] = np.zeros(0)
+
+# Create the material
+dummy_material = mcdc.material(
+    [
+        ["dummy_nuclide", 1.0],
+    ]
+)
+
+# Set surfaces
+s1 = mcdc.surface("plane-x", x=0.0, bc="reflective")
+s2 = mcdc.surface("plane-x", x=2.0, bc="vacuum")
+
+# Set cells
+mcdc.cell(+s1 & -s2, dummy_material)
+
+# =============================================================================
+# Set source
+# =============================================================================
+
+mcdc.source(
+    x=[0.95, 1.05],
+    energy=np.array([[0.9, 1.1], [1.0, 1.0]]),
+)
+
+# =============================================================================
+# Set tally, setting, and run mcdc
+# =============================================================================
+
+mcdc.tally(scores=["flux"], x=np.linspace(0.0, 2.0, 21), E=np.array([0.0, 1.0, 20e6]))
+mcdc.setting(N_particle=1e3)
+mcdc.run()