various small improvements in documentation - ceiled the MCC3 atom density -...

various small improvements in documentation - ceiled the MCC3 atom density - returned RTFLUX for SE2ANL

CONTRIBUTING.md

@@ -364,7 +364,7 @@ Below is an extended example of the `postrun` block. The main role of this block
 
         # TODO: implement post-processing logic
 
-        self.user_object.results_mycode = 0 # provide some results that we think users will always be interested. This will be used for unit-cell (check results do not change) and for
+        self.user_object.results['mycode'] = 0 # provide some results that we think users will always be interested. This will be used for unit-tests (check results do not change) and can be extracted from Python execution
 
         out.writelines("\n")
         out.close()

PyARCUserObject.py

@@ -101,11 +101,11 @@ class UserObject(object):
         self.results["gamsor_results"] = {} # summarizes all the GAMSOR power results for every core region
         # results["gamsor_results"][timestep][region_id] = [neutron power, gamma ppwer, total power]
         self.results["dassh"] = {}
-        # self.results["dassh"][timestep]["Max_Fuel"] = max_fuel_temp (K)
-        # self.results["dassh"][timestep]["Max_Clad"] = max_clad_temp (K)
-        # self.results["dassh"][timestep]["Peak_Outlet"] = peak_outlet_temp (K)
-        # self.results["dassh"][timestep]["Bulk_Outlet"] = bulk_outlet_temp (K)
-        # self.results["dassh"][timestep]["Max_Total_PD"] = max_total_pressure_drop (K)
+        # results["dassh"][timestep]["Max_Fuel"] = max_fuel_temp (K)
+        # results["dassh"][timestep]["Max_Clad"] = max_clad_temp (K)
+        # results["dassh"][timestep]["Peak_Outlet"] = peak_outlet_temp (K)
+        # results["dassh"][timestep]["Bulk_Outlet"] = bulk_outlet_temp (K)
+        # results["dassh"][timestep]["Max_Total_PD"] = max_total_pressure_drop (K)
         self.results["mycode"] = None
 
         self.dif3d_input_blocks = {}

PyGAMSOR.py

@@ -236,6 +236,8 @@ class GAMSOR(object):
 #                    self.user_object.print_timelines(filename)
                     if "NHFLX0" == filename:
                         shutil.copy("NHFLX0", "../../"+output_files_name+"/GHFLX0")
+                    if "RTFLUX" == filename:
+                        shutil.copy("RTFLUX", "../../"+output_files_name+"/RTFLUX")
 
         if gamsor_crash is False:
             os.chdir("../step3/")

PyMCC3.py

@@ -568,7 +568,7 @@ class MCC3(object):
             for (isotope_name, isotope_type, cell_id, pert_id_new, time_step_new) in list_isotopes_subassembly:
                 isotope_value = self.user_object.isotope_hom_mixture_data[("R", 0)][cell_id][(isotope_name, isotope_type)]
                 if isotope_type != "LUMPED ELEMENT":
-                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], float(isotope_value[0]), float(isotope_value[1]))
+                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], max(float(isotope_value[0]), mcc3_min), float(isotope_value[1]))  # atom density must be larger than mcc3_min
             del list_isotopes_subassembly
             mcc3_inp.writelines(t_composition+"\n")
 
@@ -584,6 +584,7 @@ class MCC3(object):
         lumped = self.user_object.lumped
         lumped_mixing_table = self.user_object.lumped_mixing_table
         sub_assembly = self.user_object.cell_definition[mcc3_id][0]
+        mcc3_min = 1e-23 # minimum atom density allowed in mcc3
         gamma_calc = False
         if calc_id == "GAMMA":
             gamma_calc = True
@@ -849,7 +850,7 @@ class MCC3(object):
             for (isotope_name, isotope_type, cell_id, calc_id_new, time_step_new) in list_isotopes_subassembly:
                 isotope_value = self.user_object.isotope_hom_mixture_data[(calc_id, time_step)][cell_id][(isotope_name, isotope_type)]
                 if isotope_type != "LUMPED ELEMENT":
-                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], float(isotope_value[0]), float(isotope_value[1]))
+                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], max(float(isotope_value[0]), mcc3_min), float(isotope_value[1])) # atom density must be larger than mcc3_min
 
         elif mcc3_cell.automated_hete_treatment is not None:
             t_composition = ""
@@ -872,9 +873,9 @@ class MCC3(object):
                             isotope_value = self.user_object.isotope_het_mixture_data[(calc_id, time_step)][cell_id][mat_name][(isotope_name, isotope_type)]
                             if isotope_type != "LUMPED ELEMENT":
                                 if total_density < min_density_cpm:
-                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], min_density_cpm*float(isotope_value[0])/total_density, float(isotope_value[1]))
+                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], min_density_cpm* max(float(isotope_value[0]), mcc3_min)/total_density, float(isotope_value[1])) # atom density must be larger than mcc3_min
                                 else:
-                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], float(isotope_value[0]), float(isotope_value[1]))
+                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3],  max(float(isotope_value[0]), mcc3_min), float(isotope_value[1])) # atom density must be larger than mcc3_min
 
                 del total_density, mat_name
             if mcc3_cell.external_composition is not None:
@@ -884,7 +885,7 @@ class MCC3(object):
                 for (isotope_name, isotope_type, cell_id, calc_id_new, time_step_new) in list_isotopes_subassembly_ext:
                     isotope_value = self.user_object.isotope_hom_mixture_data[(calc_id, time_step)][cell_id][(isotope_name, isotope_type)]
                     if isotope_type != "LUMPED ELEMENT":
-                        t_composition += "\n%6s %6s %6.5e %6.1f "%(isotope_value[4], isotope_value[3], float(isotope_value[0]), float(isotope_value[1]))
+                        t_composition += "\n%6s %6s %6.5e %6.1f "%(isotope_value[4], isotope_value[3], max(float(isotope_value[0]), mcc3_min), float(isotope_value[1]))  # atom density must be larger than mcc3_min
 
         else:
             t_composition = ""
@@ -910,9 +911,9 @@ class MCC3(object):
                             isotope_value = self.user_object.isotope_het_mixture_data[(calc_id, time_step)][cell_id][mat_name][(isotope_name, isotope_type)]
                             if isotope_type != "LUMPED ELEMENT":
                                 if total_density < min_density_cpm:
-                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], min_density_cpm*float(isotope_value[0])/total_density, float(isotope_value[1]))
+                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], min_density_cpm*max(float(isotope_value[0]), mcc3_min)/total_density, float(isotope_value[1]))
                                 else:
-                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], float(isotope_value[0]), float(isotope_value[1]))
+                                    t_composition += "\n%6s %6s %6.5e %6.1f"%(isotope_value[4], isotope_value[3], max(float(isotope_value[0]), mcc3_min), float(isotope_value[1])) # atom density must be larger than mcc3_min
                 del total_density, mat_name
 
             if mcc3_cell.hete_cell_options.external_composition is not None:
@@ -922,7 +923,7 @@ class MCC3(object):
                 for (isotope_name, isotope_type, cell_id, calc_id_new, time_step_new) in list_isotopes_subassembly_ext:
                     isotope_value = self.user_object.isotope_hom_mixture_data[(calc_id, time_step)][cell_id][(isotope_name, isotope_type)]
                     if isotope_type != "LUMPED ELEMENT":
-                        t_composition += "\n%6s %6s %6.5e %6.1f "%(isotope_value[4], isotope_value[3], float(isotope_value[0]), float(isotope_value[1]))
+                        t_composition += "\n%6s %6s %6.5e %6.1f "%(isotope_value[4], isotope_value[3], max(float(isotope_value[0]), mcc3_min), float(isotope_value[1]))  # atom density must be larger than mcc3_min
 
         mcc3_inp.writelines("\n$material"+i_source_region+t_composition_name+i_externalspectrum+t_composition+"\n/")
 
@@ -955,7 +956,7 @@ class MCC3(object):
                                 if isotope_name.lower() == constituent_id and isotope_type == "CONSTITUENT":
                                     isotope_it_value = self.user_object.isotope_hom_mixture_data[(calc_id, time_step)][cell_id][(isotope_name, isotope_type)]
                                     dif3d_isotope_id = isotope_it_value[3]
-                            mcc3_inp.writelines("  %s   %6.5e\n"%(dif3d_isotope_id, afrac))
+                            mcc3_inp.writelines("  %s   %6.5e\n"%(dif3d_isotope_id, max(afrac, mcc3_min)))  # atom density must be larger than mcc3_min
         mcc3_inp.writelines("/")
 
         del list_isotopes_subassembly

README.md

@@ -41,6 +41,7 @@ This module is developed under the leadership of ANL with the assistance and con
     * [Contact](#contact) - contact information for questions
 8. [Tutorial](/tutorial/README.md#tutorial) - Links to tutorial, sample inputs and benchmark problems
 9. [Unit Tests](/test/README.md#unit-tests) - Testing pre/post processing and execution logic of the ARC codes
+9. [Results](#results) - Structure of the results being returned by PyARC as python object
 
 ## License
 
@@ -286,3 +287,47 @@ The PyARC module relies on the following sub-modules:
 * The PyARC module also relies on the [__PySCL__](/PySCL/README.md) module that is developed by ORNL to provide the standard composition library (SCL) in a consistent way between different codes
 * The [__plot2d__](/utilities_plot2d/README.md) provides 2D plots to display geometry and some results based on summary file
 * The [__covmat__](/utilities_covmat/README.md) enable covariance matrix generation from nuclear data uncertainties on reactivity coefficients - obtained through a series of PERSENT simulations
+
+
+# Results
+
+When executing PyARC withing a python scrip or through [__WATTS__](https://github.com/watts-dev/watts), the main calculation results can be obtained within a dictionary called `user_object.results`. The following of this dictionary is the following:
+
+```javascript
+# DIF3D results
+results["keff_dif3d"][timestep] = keff
+results["dif3d_results"][timestep][region_ID] = [Volume [cm3], Neutron Power [W], Max power density [Wm3], Int Flux [n-cm/s],Peak Flux[n/cm2-s],Int Fast flux [n-cm/s],Peak Fast Flux [n/cm2-s]] 
+results["keff_rzmflx"][timestep] = keff # this is TWODANT keff result
+results["keff_nodal"][timestep] = keff  # this is PROTEUS-Nodal keff result
+
+# REBUS results
+results["rebus_inventory"][day][(area_ID, isotope)] = mass # [kg]
+results["rebus_results"]["PEAK BURNUP"] = burnup # [GWd/T]
+results["rebus_results"]["PEAK FAST FLUENCE"] = fast_fluence # [n/cm2]
+results["rebus_results"]["ENRICHMENT"] = rebus_enrichment_result # [%]
+results["rebus_results"]["AVG BURNUP"] = average_discharge_burnup # [GWd/T]
+results["rebus_results"]["EXCESS REACTIVITY"] = excess_reactivity # [pcm]
+results["rebus_results"]["CORE LIFE"] = core_lifetime # [EFPD]
+
+# MCC3 results
+results["keff_mcc3"][(R, num_timestep, MCC3 step number, cell_ID, region))] = kinf
+
+# PERSENT results
+results["kpert_persent"][pert_ID] = total_perturbation
+results["kuq_persent"][sens_id] = total_uncertainty
+results["kuq_persent"][sens_id, sensitivity_name] = total_sensitivity
+
+# CHANNEL results - summarizes all the results (GAMSOR and PERSENT) summed over the different channels
+results["channel_results"][output_file][CHAN#][Axial_ID] =  [weight perturbed, result (power oreactivity coeff), [list of regions considered]]
+
+# GAMSOR results
+results["power_gamsor"] = power # provides just the total GAMSOR power (for unit-test purposes)
+results["gamsor_results"][timestep][region_id] = [neutron power, gamma ppwer, total power] # summarizes all the GAMSOR power results for every core region
+
+# DASSH results
+results["dassh"][timestep]["Max_Fuel"] = max_fuel_temp (K)
+results["dassh"][timestep]["Max_Clad"] = max_clad_temp (K)
+results["dassh"][timestep]["Peak_Outlet"] = peak_outlet_temp (K)
+results["dassh"][timestep]["Bulk_Outlet"] = bulk_outlet_temp (K)
+results["dassh"][timestep]["Max_Total_PD"] = max_total_pressure_drop (K)
+```
\ No newline at end of file