initial commit of dassh_power pathway to write pin total power distribution with no temperature sweep

Author: MilosAtz

dassh/__main__.py

@@ -14,17 +14,18 @@
 # permissions and limitations under the License.
 ########################################################################
 """
-date: 2022-01-20
+date: 2022-12-19
 author: matz
-Main DASSH calculation procedure
+Main DASSH calculation procedures
 """
 ########################################################################
 import os
 import sys
-import dassh
+import numpy as np
 import argparse
 import cProfile
 import logging
+import dassh
 _log_info = 20  # logging levels must be int
 
 
@@ -279,5 +280,66 @@ def plot():
     dassh_logger.log(_log_info, 'DASSH_PLOT execution complete')
 
 
+def integrate_pin_power(args=None):
+    """Set up DASSH Reactor object, integrate pin power, and write to CSV"""
+    # Get input file from command line arguments
+    parser = argparse.ArgumentParser(description='Process DASSH cmd')
+    parser.add_argument('inputfile',
+                        metavar='inputfile',
+                        help='The input file to run with DASSH')
+    parser.add_argument('--save_reactor',
+                        action='store_true',
+                        help='Save DASSH Reactor object after sweep')
+    args = parser.parse_args()
+
+    # Initiate logger
+    print(dassh._ascii._ascii_title)
+    in_path = os.path.split(args.inputfile)[0]
+    dassh_logger = dassh.logged_class.init_root_logger(in_path, 'dassh_power')
+
+    # Pre-processing
+    # Read input file and set up DASSH input object
+    dassh_logger.log(_log_info, f'Reading input: {args.inputfile}')
+    dassh_input = dassh.DASSHPower_Input(args.inputfile)
+
+    # Initialize the Reactor object
+    reactor = dassh.Reactor(dassh_input, write_output=False)
+
+    # Generate pin power distributions
+    dassh_logger.log(_log_info, 'Generating pin power distributions...')
+    asm_ids = []
+    n_pins = []
+    integrated_pin_powers = []
+    for a in reactor.assemblies:
+        if a.has_rodded:
+            asm_ids.append(a.id)
+            n_pins.append(a.rodded.n_pin)
+            integrated_pin_powers.append(
+                dassh.power._integrate_pin_power(a.power))
+
+    # Save reactor if desired
+    dassh_logger.log(_log_info, 'Saving data')
+    if args.save_reactor:
+        if sys.version_info < (3, 7):
+            handlers = dassh_logger.handlers[:]
+            for handler in handlers:
+                handler.close()
+                dassh_logger.removeHandler(handler)
+        reactor.save()
+        if sys.version_info < (3, 7):
+            dassh_logger = dassh.logged_class.init_root_logger(
+                os.path.split(dassh_logger._root_logfile_path)[0],
+                'dassh', 'a+')
+
+    # Write distributions to CSV
+    arr_to_write = np.zeros((max(n_pins) + 1, len(asm_ids)))
+    arr_to_write[0] = asm_ids
+    for col in range(arr_to_write.shape[1]):
+        arr_to_write[1:(n_pins[col] + 1), col] = integrated_pin_powers[col]
+    outpath = os.path.join(dassh_input.path, 'total_pin_power.csv')
+    np.savetxt(outpath, arr_to_write, delimiter=',')
+    dassh_logger.log(_log_info, 'DASSH_POWER execution complete')
+
+
 if __name__ == '__main__':
     main()

dassh/power.py

@@ -14,7 +14,7 @@
 # permissions and limitations under the License.
 ########################################################################
 """
-date: 2022-03-16
+date: 2022-12-19
 author: matz
 Generate power distributions in assembly components based on neutron
 flux; object to assign to individual assemblies
@@ -1508,82 +1508,33 @@ def _integrate(pp_pins, pp_duct, pp_cool, n_terms):
     return avg_power[:, 1] - avg_power[:, 0]
 
 
-########################################################################
+def _integrate_pin_power(asm_power):
+    """Integrate the pin power distributions
 
-    # def renormalize(self, z, dz):
-    #     """Need to normalize the power based on the discretization that
-    #     DASSH will make during the sweep in order to accurately deliver
-    #     the correct power
-    #
-    #     Parameters
-    #     ----------
-    #     z : numpy.ndarray
-    #         Absolute axial points (m)
-    #     dz : numpy.ndarray
-    #         Mesh step sizes (m)
-    #
-    #     Returns
-    #     -------
-    #     numpy.ndarray
-    #         Array of correction factors for each fine mesh interval
-    #
-    #     """
-    #     # Preprocess the axial mesh points
-    #     z = z[1:] * 100  # m --> cm
-    #     z = z - 0.5 * dz  # get the axial mesh midpoints
-    #     # Get the fine mesh intervals and the transformed z values for
-    #     # each of the axial mesh points where power is evaluated.
-    #     # kf = get_kfint2(self, z)
-    #     # zm = transform_z2(self, kf, z)
-    #     kf = self.get_kfint2(z)
-    #     zm = self.transform_z2(kf, z)
-    #     zm = np.expand_dims(zm, 1)
-    #
-    #     # Raise those transformed z values to the monomial exponent
-    #     z_exp = np.power(zm, np.arange(self.n_terms))
-    #     z_exp = np.expand_dims(z_exp, 1)
-    #
-    #     # Multiply the coefficients through and sum to get the total
-    #     # power in each axial step: do this for pins, duct, and coolant
-    #     mult_terms = self.pin_power[kf] * z_exp
-    #     sum_terms = np.sum(mult_terms, axis=(1, 2))
-    #     mult_terms = self.duct_power[kf] * z_exp
-    #     sum_terms += np.sum(mult_terms, axis=(1, 2))
-    #     mult_terms = self.coolant_power[kf] * z_exp
-    #     sum_terms += np.sum(mult_terms, axis=(1, 2))
-    #
-    #     # Dimension is equal to n_step x 1; equals total in each step
-    #     sum_terms = sum_terms * 100  # W/cm --> W/m
-    #
-    #     # Normalize each step linear power by step size
-    #     dz_fint = self.z_finemesh[1:] - self.z_finemesh[:-1]
-    #     dz_fint = dz_fint[kf]  # cm
-    #     sum_terms = sum_terms * dz / dz_fint  # W/m * m / cm = W/cm
-    #
-    #     # Determine the corrective factor in each kfint
-    #     x = np.zeros(len(self.z_finemesh) - 1)
-    #     for i in range(len(x)):
-    #         # self.avg_power has units of W/cm; so does sum_terms
-    #         x[i] = self.avg_power[i] / np.sum(sum_terms[kf == i])
-    #         # x[i] = self.avg_power[i] * 100 / np.average(sum_terms[kf == i])
-    #     return x
+    Parameters
+    ----------
+    asm_power : DASSH AssemblyPower object
+
+    Returns
+    -------
+    numpy.ndarray
+        Total power in each pin
 
+    """
+    z_bnds = np.array([-0.5, 0.5])
+    z_bnds = z_bnds.reshape(2, 1)
+    int_exponents = np.arange(1, asm_power.n_terms + 1)
+    z_int = np.power(z_bnds, int_exponents)
+    # Integrate in each region, evaluate at lower/upper bound
+    # shape is n_region x n_pin x 2 (upper/lower bound)
+    integrated = np.dot(asm_power.pin_power / int_exponents, z_int.T)
+    # Take the difference across the region
+    # shape is n_region x n_pin
+    diff_across_region = integrated[:, :, 1] - integrated[:, :, 0]
+    # multiply linear power by z-bounds and sum to get total power
+    # in each pin
+    dz_finemesh = asm_power.z_finemesh[1:] - asm_power.z_finemesh[:-1]
+    power_per_pin = np.dot(dz_finemesh, diff_across_region)
+    return power_per_pin
 
 ########################################################################
-# Old
-# def calculate_total_power(self):
-#     """Calculate the total power (W) produced by the assembly using
-#     the linear power (W/m) shape functions."""
-#
-#     dz = [self.z_finemesh[k] - self.z_finemesh[k - 1]
-#           for k in range(1, len(self.z_finemesh))]
-#     p_total = 0.0
-#     for k in range(len(dz)):
-#         zc = self.z_finemesh[k] + dz[k] / 2  # center of mesh cell
-#         if k < self.k_bnds[0] or k >= self.k_bnds[1]:
-#             p_total += dz[k] * self.get_refl_power(zc)
-#         else:
-#             p_total += dz[k] * (np.sum(self.get_pin_power(zc))
-#                                 + np.sum(self.get_duct_power(zc))
-#                                 + np.sum(self.get_coolant_power(zc)))
-#     return p_total

setup.py

@@ -14,7 +14,7 @@
 # permissions and limitations under the License.
 ########################################################################
 """
-date: 2022-01-20
+date: 2022-12-16
 author: Milos Atz
 """
 ########################################################################
@@ -130,7 +130,8 @@ def check_dependencies():
                        'Operating System :: MacOS'],
           entry_points={'console_scripts': [
               'dassh = dassh.__main__:main',
-              'dassh_plot = dassh.__main__:plot']}
+              'dassh_plot = dassh.__main__:plot',
+              'dassh_power = dassh.__main__:integrate_pin_power']}
           )