Merge pull request #130 from clemekay/main

Merging variance deconvolution feature into main.

mcdc/__init__.py

@@ -20,6 +20,7 @@
     weight_roulette,
     IC_generator,
     dsm,
+    uq,
     print_card,
     reset_cards,
 )

mcdc/card.py

@@ -31,6 +31,7 @@ def reset(self):
             "total": False,
             "current": False,
             "eddington": False,
+            "exit": False,
             "mesh": make_card_mesh(),
         }
 
@@ -116,6 +117,14 @@ def reset(self):
             "IC_cycle_stretch": 1.0,
             "branchless_collision": False,
             "dsm_order": 1,
+            "uq": False,
+        }
+
+        self.uq_deltas = {
+            "tag": "Uq",
+            "nuclides": [],
+            "materials": [],
+            "surfaces": [],
         }
 
 
@@ -152,6 +161,7 @@ def make_card_nuclide(G, J):
     card["sensitivity"] = False
     card["sensitivity_ID"] = 0
     card["dsm_Np"] = 1.0
+    card["uq"] = False
     return card
 
 
@@ -176,6 +186,7 @@ def make_card_material(N_nuclide, G, J):
     card["chi_s"] = np.zeros([G, G])
     card["chi_p"] = np.zeros([G, G])
     card["sensitivity"] = False
+    card["uq"] = False
     return card
 
 
@@ -284,3 +295,17 @@ def make_card_mesh():
         "azi": np.array([-PI, PI]),
         "g": np.array([-INF, INF]),
     }
+
+
+def make_card_uq():
+    return {
+        "tag": "t",
+        "ID": -1,
+        "key": "k",
+        "mean": 0.0,
+        "delta": 0.0,
+        "distribution": "d",
+        "rng_seed": 0,
+        "group": False,
+        "group_group": False,
+    }

mcdc/constant.py

@@ -48,3 +48,4 @@
 SEED_SPLIT_SOURCE_PRECURSOR = nb.uint64(0x546F6464)
 SEED_SPLIT_BANK = nb.uint64(0x5279616E)
 SEED_SPLIT_PARTICLE = nb.uint64(0)
+SEED_SPLIT_UQ = nb.uint(0x5368656261)

mcdc/input_.py

@@ -17,6 +17,7 @@
     make_card_universe,
     make_card_lattice,
     make_card_source,
+    make_card_uq,
 )
 from mcdc.constant import (
     GYRATION_RADIUS_ALL,
@@ -1279,6 +1280,83 @@ def dsm(order=1):
     card["dsm_order"] = order
 
 
+def uq(**kw):
+    def append_card(delta_card, global_tag):
+        delta_card["distribution"] = dist
+        delta_card["flags"] = []
+        for key in kw.keys():
+            check_support(parameter["tag"] + " parameter", key, parameter_list, False)
+            delta_card["flags"].append(key)
+            delta_card[key] = kw[key]
+        mcdc.input_deck.uq_deltas[global_tag].append(delta_card)
+
+    mcdc.input_deck.technique["uq"] = True
+    # Make sure N_batch > 1
+    if mcdc.input_deck.setting["N_batch"] <= 1:
+        print_error("Must set N_batch>1 with mcdc.setting() prior to mcdc.uq() call.")
+
+    # Check uq parameter
+    parameter_ = check_support(
+        "uq parameter",
+        list(kw)[0],
+        ["nuclide", "material", "surface", "source"],
+        False,
+    )
+    parameter = kw[parameter_]
+    del kw[parameter_]
+    parameter["uq"] = True
+
+    # Confirm supplied distribution
+    check_requirement("uq", kw, ["distribution"])
+    dist = check_support("distribution", kw["distribution"], ["uniform"], False)
+    del kw["distribution"]
+
+    # Only remaining keywords should be the parameter delta(s)
+
+    if parameter["tag"] == "Material":
+        parameter_list = [
+            "capture",
+            "scatter",
+            "fission",
+            "nu_s",
+            "nu_p",
+            "nu_d",
+            "chi_p",
+            "chi_d",
+            "speed",
+            "decay",
+        ]
+        global_tag = "materials"
+        if parameter["N_nuclide"] == 1:
+            nuc_card = make_card_nuclide(parameter["G"], parameter["J"])
+            nuc_card["ID"] = parameter["nuclide_IDs"][0]
+            append_card(nuc_card, "nuclides")
+        delta_card = make_card_material(
+            parameter["N_nuclide"], parameter["G"], parameter["J"]
+        )
+        for name in ["ID", "nuclide_IDs", "nuclide_densities"]:
+            delta_card[name] = parameter[name]
+    elif parameter["tag"] == "Nuclide":
+        parameter_list = [
+            "capture",
+            "scatter",
+            "fission",
+            "nu_s",
+            "nu_p",
+            "nu_d",
+            "chi_p",
+            "chi_d",
+            "speed",
+            "decay",
+        ]
+        global_tag = "nuclides"
+        delta_card = make_card_nuclide(parameter["G"], parameter["J"])
+        delta_card["ID"] = parameter["ID"]
+    append_card(delta_card, global_tag)
+    # elif parameter['tag'] is 'Surface':
+    # elif parameter['tag'] is 'Source':
+
+
 # ==============================================================================
 # Util
 # ==============================================================================