Commit 87a9ff49 authored by Salko Jr, Robert's avatar Salko Jr, Robert
Browse files

This adds support for tube flow and several other features

Description:
Summary of changes:
* Adds support for modeling flow on the inside of a tube
* Expose the remainder of the steady-state stopping criteria to the user
* Add ability for the user to set the solver
* Fixed a bug in the edits setter method
* Fixed a bug in one of the internal boundary condition getter procedures
* Fixed a bug in the mass/energy balance units
* Fixed a bug in the procedure that allows the user to enable/disable the conduction model
* Added ability for user to enable/disable entrainment/de-entrainment model
* Added ability to set HTC under-relaxation coefficient through SubKit and added test
* Added a new method for setting initial conditions, which can take mdot/massflux and temp/enthalpy. Having one method for setting any combination of initial conditions will be easier for the user to remember and keeps the code cleaner.

Gitlab Ticket # - 4582
parent 9a25fd05
Loading
Loading
Loading
Loading
+244 −45
Original line number Diff line number Diff line
@@ -95,6 +95,9 @@ class Model:
        me.editInterval = []
        me.tEnd = []

        # Modeling options
        me.enableEntrainment = False

        # Edit options
        me.editOptions = collections.OrderedDict()
        me.editOptions["rod_edits"] = False
@@ -112,9 +115,17 @@ class Model:
        me.chanMap = None
        me.symOpt = None

        # Initial conditions
        me.mdotInit = None
        me.mdotInitialWasSpecified = False
        me.massFluxInit = None
        me.tempInit = None
        me.enthalpyInit = None
        me.pressureInit = None

        # Steady state convergence criteria
        me.steadyConvergenceCriteria = {'massBalance': 0.1,
                                        'energyBalance': 0.1,
        me.steadyConvergenceCriteria = {'massBalance': 1e-3,
                                        'energyBalance': 1e-3,
                                        'maxIterations': 10000,
                                        'p_l2': 1e-4,
                                        'p_linf': 1e-3,
@@ -127,7 +138,21 @@ class Model:
                                        'Tcool_l2': 1e-4,
                                        'Tcool_linf': 1e-3,
                                        'Tsolid_l2': 1e-4,
                                        'Tsolid_linf': 1e-3
                                        'Tsolid_linf': 1e-3,
                                        'void_l2_abs': 1e-4,
                                        'void_linf_abs': 1e-4,
                                        'p_l2_abs': 1e-4,
                                        'p_linf_abs': 1e-4,
                                        'vl_l2_abs': 1e-4,
                                        'vl_linf_abs': 1e-4,
                                        'vv_l2_abs': 1e-4,
                                        'vv_linf_abs': 1e-4,
                                        'vd_l2_abs': 1e-4,
                                        'vd_linf_abs': 1e-4,
                                        'Tcool_l2_abs': 1e-4,
                                        'Tcool_linf_abs': 1e-4,
                                        'Tsolid_l2_abs': 1e-4,
                                        'Tsolid_linf_abs': 1e-4
                                        }

        me.conductorOptions = {'nc': 1,  # 1-D conduction
@@ -147,6 +172,7 @@ class Model:

        me.rtwfp = 100.0

        me.setInitialCalled = False
        me.setInitialConditionsCalled = False
        me.setInitialConditionsMassFluxCalled = False

@@ -162,6 +188,9 @@ class Model:
        # Single-phase turbulent mixing coefficient
        me.beta = 0.037

        # Under-relaxation coefficients
        me.betaHTC = None

    @property
    def is_nodal(self):
        """ Access private is_nodal attr """
@@ -185,15 +214,23 @@ class Model:
        """
        me.title = str(title)

    def setModelOptions(me, beta=None):
    def setModelOptions(me, beta=None, enableEntrainment=None, betaHTC=None):
        """ Set various modeling options

        Any model that can be changed via this method has a default.
        This only needs to be called if the user wishes to change a deafult modeling option.

        Args:
           beta (float): Single-phase turbulent mixing coefficietn
           enableEntrainment (bool): Set to True to enable entrainment/de-entrainment of droplets
           betaHTC (float): Under-relaxation coefficient for fluid/solid heat transfer coefficient
        """

        if beta is not None:
            me.beta = beta
        if enableEntrainment is not None:
            me.enableEntrainment = enableEntrainment
        me.betaHTC = betaHTC

    def setFrictionModel(me, model):
        """ Sets the friction model option
@@ -261,18 +298,24 @@ class Model:

        Must be called to set the initial conditions

        NOTE: This is a deprecated procedure and will be removed in a future release.  Use
        `setInitial` instead.

        Args:
           mdotInit (float): Initial mass flow rate in the system
           tempInit (float): Initial temperature in the model
           pressureInit (float): Initial (and reference) pressure in the system

        """
        me.mdotInit = mdotInit
        # Also set massFlux in case we need it
        me.massFlux = me.mdotInit / me.sections[1].getFlowArea()
        me.pressureInit = pressureInit
        me.tempInit = tempInit
        me.setInitialConditionsMassFluxCalled = False
        if me.setInitialCalled:
            raise RuntimeError("`setInitial` already called")
        if me.setInitialConditionsCalled:
            raise RuntimeError("`setInitialConditions` already called")
        if me.setInitialConditionsMassFluxCalled:
            raise RuntimeError("`setInitialConditionsMassFlux` already called")
        warnings.warn(
            "`setInitialConditions` has been deprecated.  Switch to using `setInitial` instead.", Warning)
        me.setInitial(mdot=mdotInit, temp=tempInit, pressure=pressureInit)
        me.setInitialConditionsCalled = True

    def setInitialConditionsMassFlux(me, massFlux, tempInit, pressureInit):
@@ -283,20 +326,102 @@ class Model:
        1. Will not specify boundary conditions for individual channels
        2. Will have nonzero initial flow

        NOTE: This is a deprecated feature and will be removed in the future.  Use `setInitial` instead.

        Args:
           massFlux (float): Initial mass flux in the system
           tempInit (float): Initial temperature in the model
           pressureInit (float): Initial (and reference) pressure in the system

        """
        me.massFlux = massFlux
        # Also set mdot in case we need it
        me.mdotInit = me.massFlux * me.sections[1].getFlowArea()
        me.pressureInit = pressureInit
        me.tempInit = tempInit
        me.setInitialConditionsCalled = False
        if me.setInitialCalled:
            raise RuntimeError("`setInitial` already called")
        if me.setInitialConditionsCalled:
            raise RuntimeError("`setInitialConditions` already called")
        if me.setInitialConditionsMassFluxCalled:
            raise RuntimeError("`setInitialConditionsMassFlux` already called")
        warnings.warn(
            "`setInitialConditionsMassFlux` has been deprecated.  Switch to using `setInitial` instead..", Warning)
        me.setInitial(massflux=massFluxInit,
                      temp=tempInit, pressure=pressureInit)
        me.setInitialConditionsMassFluxCalled = True

    def setInitial(me, **kwargs):
        """Set model initial conditions.

        Must be called to set the initial conditions
        Three types of parameters must be set:
           1. Flow
           2. Energy
           3. Pressure

        **Flow**

        For flow, you can pass either `mdot` to specify it as absolute mass flow rate in units
        of kg/s or `massfluxInit` to specify it in units of kg/m**2/s.

        If creating a model for execution in parallel, `massfluxInit` must be used if the flow is not
        being specified individually for the channels and if flow is nonzero.

        **Energy**

        For energy, consider that there are both fluid and solid regions, and both need to have their initial
        energy specified.

        For the fluid, you can pass either `temp` to specify it as a temperature in units of Celsius
        or `enthalpy` to specify it as an enthalpy in units of kJ/kg.

        The solid does not necessarily have to have a temperature specified.  If `temp` was specified for
        the fluid, and nothing is specified for the solid, the solid will just take the temperature specified
        by `temp` as its initial temperature.  If `enthalpy` was specified for the fluid, then the solid needs
        to have its initial temperature specified by passing `temp_solid`.  Units are also in Celsius.  It is
        allowable to pass both `temp` and `temp_solid`, in which case the fluid and solid regions will have
        different initial temperatures.

        **Pressure**

        For pressure, simply pass pressure in bar.

        Args:
           mdot (float): Initial mass flow rate in the system
           massflux (float): Initial mass flux in the system
           temp (float): Initial temperature in the model
           enthalpy (float): Initial enthalpy in the model
           pressure (float): Initial (and reference) pressure in the system

        """
        if me.setInitialCalled:
            raise RuntimeError("`setInitial` already called")
        assert not ('mdot' in kwargs and 'massflux' in kwargs)
        assert 'mdot' in kwargs or 'massflux' in kwargs
        assert not ('temp' in kwargs and 'enthalpy' in kwargs)
        assert 'temp' in kwargs or 'enthalpy' in kwargs
        if 'enthalpy' in kwargs:
            assert 'temp_solid' in kwargs
        assert 'pressure' in kwargs

        if 'mdot' in kwargs:
            me.mdotInit = kwargs['mdot']
            me.massFluxInit = me.mdotInit / me.sections[1].getFlowArea()
            me.mdotInitialWasSpecified = True
        else:
            me.massFluxInit = kwargs['massflux']
            me.mdotInit = me.massFluxInit * me.sections[1].getFlowArea()

        if 'temp' in kwargs:
            me.tempInit = kwargs['temp']
        else:
            me.enthalpyInit = kwargs['enthalpy']

        if 'temp_solid' in kwargs:
            me.solidTempInit = kwargs['temp_solid']
        else:
            me.solidTempInit = kwargs['temp']

        me.pressureInit = kwargs['pressure']

        me.setInitialCalled = True

    def setAveragePower(me, qprime, dhfrac=0.0):
        """ Set the average linear heat rate applied to the model.

@@ -570,6 +695,13 @@ class Model:
            me.rodsConnectedToChannel[chID] = []
        me.rodsConnectedToChannel[chID].append(solidID)

    def hasOuterConnections(me, solidID):
        """ Returns True if this solid has channels connected on outside surface """
        if len(me.channelsOuter[solidID]) > 0:
            return True
        else:
            return False

    def getOuterConnections(me, solidID, section=None):
        """ Use to get connections to a solid object.

@@ -598,6 +730,13 @@ class Model:
                percents.append(me.percentsOuter[solidID][i])
        return chans, percents

    def hasInnerConnections(me, solidID):
        """ Returns True if this solid has channels connected on inside surface """
        if len(me.channelsInner[solidID]) > 0:
            return True
        else:
            return False

    def getInnerConnections(me, solidID):
        """ Get inner conenctions to a solid.

@@ -608,14 +747,12 @@ class Model:
           (list of channel IDs connected to the inside surface of the solid, list of percentages of each connection)

        """
        # CTF deck needs 8 connections, with zeros representing no
        # connection, so fill out the connections list to size 8.'''
        fullList = [0] * 8
        fractList = [0.0] * 8
        chans = []
        percents = []
        for i, ch in enumerate(me.channelsInner[solidID]):
            fullList[i] = ch
            fractList[i] = me.percentsInner[solidID][i]
        return (fullList, fractList)
            chans.append(ch)
            percents.append(me.percentsInner[solidID][i])
        return (chans, percents)

    def getRodsConnectedToChannel(me, chID):
        """ Returns a list of all rod IDs (user defined) connected to the passed channel ID.
@@ -813,10 +950,29 @@ class Model:
        me.dtMax.append(dtMax)
        me.rtwfp = 1.0

    def setSolver(me, solver):
        """ Sets the solver to use for solving the pressure matrix.

        Options are:
           direct - Super-LU is used if PETSc is available.  If PETSc is not available, solve time will be slow
           bicgstab - Serial simulations only.  Will work with or without PETSc.
           bicgstab_par - Serial or parallel simulations.  Requires PETSc.
        """
        assert solver in ['direct', 'bicgstab', 'bicgstab_par']
        if solver == 'direct':
            me.solver = 0
        elif solver == 'bicgstab':
            me.solver = 3
        else:
            me.solver = 5

    def setSteadyState(me, energyBalance=None, massBalance=None, maxIterations=None, rtwfp=None,
                       p_l2=None, p_linf=None, Tcool_l2=None, Tcool_linf=None,
                       Tsolid_l2=None, Tsolid_linf=None, vl_l2=None, vl_linf=None, vv_l2=None, vv_linf=None,
                       vd_l2=None, vd_linf=None, dtmin=None, dtmax=None):
                       vd_l2=None, vd_linf=None, void_l2_abs=None, void_linf_abs=None, p_l2_abs=None, p_linf_abs=None, Tcool_l2_abs=None,
                       Tcool_linf_abs=None, Tsolid_l2_abs=None, Tsolid_linf_abs=None, vl_l2_abs=None,
                       vl_linf_abs=None, vv_l2_abs=None, vv_linf_abs=None, vd_l2_abs=None, vd_linf_abs=None,
                       dtmin=None, dtmax=None):
        """ Call to specify the model as steady state.

        Throws exception if transient groups specified already.
@@ -844,6 +1000,20 @@ class Model:
           vv_linf (float): l-infinity-norm criteria on vapor velocity
           vd_l2 (float): l2-norm criteria on droplet velocity
           vd_linf (float): l-infinity-norm criteria on droplet velocity
           void_l2_abs (float): absolute l2-norm criteria on void [-]
           void_linf_abs (float): absolute l-infinity-norm criteria on void [-]
           p_l2_abs (float): absolute l2-norm criteria on pressure [bar]
           p_linf_abs (float): absolute l-infinity-norm criteria on pressure [bar]
           Tcool_l2_abs (float): absolute l2-norm criteria on coolant temperature [C]
           Tcool_linf_abs (float): absolute l-infinity-norm criteria on coolant temperature [C]
           Tsolid_l2_abs (float): absolute l2-norm criteria on solid temperature [C]
           Tsolid_linf_abs (float): absolute l-infinity-norm criteria on solid temperature [C]
           vl_l2_abs (float): absolute l2-norm criteria on liquid velocity [m/s]
           vl_linf_abs (float): absolute l-infinity-norm criteria on liquid velocity [m/s]
           vv_l2_abs (float): absolute l2-norm criteria on vapor velocity [m/s]
           vv_linf_abs (float): absolute l-infinity-norm criteria on vapor velocity [m/s]
           vd_l2_abs (float): absolute l2-norm criteria on droplet velocity [m/s]
           vd_linf_abs (float): absolute l-infinity-norm criteria on droplet velocity [m/s]
           dtmin (float): Minimum allowable timestep size (steady-state is modeled as a pseudo-transient)
           dtmax (float): Maximum allowable timestep size (steady-state is modeled as a pseudo-transient)

@@ -889,6 +1059,34 @@ class Model:
            me.steadyConvergenceCriteria['vd_l2'] = vd_l2
        if vd_linf:
            me.steadyConvergenceCriteria['vd_linf'] = vd_linf
        if void_l2_abs:
            me.steadyConvergenceCriteria['void_l2_abs'] = void_l2_abs
        if void_linf_abs:
            me.steadyConvergenceCriteria['void_linf_abs'] = void_linf_abs
        if p_l2_abs:
            me.steadyConvergenceCriteria['p_l2_abs'] = p_l2_abs
        if p_linf_abs:
            me.steadyConvergenceCriteria['p_linf_abs'] = p_linf_abs
        if Tcool_l2_abs:
            me.steadyConvergenceCriteria['Tcool_l2_abs'] = Tcool_l2_abs
        if Tcool_linf_abs:
            me.steadyConvergenceCriteria['Tcool_linf_abs'] = Tcool_linf_abs
        if Tsolid_l2_abs:
            me.steadyConvergenceCriteria['Tsolid_l2_abs'] = Tsolid_l2_abs
        if Tsolid_linf_abs:
            me.steadyConvergenceCriteria['Tsolid_linf_abs'] = Tsolid_linf_abs
        if vl_l2_abs:
            me.steadyConvergenceCriteria['vl_l2_abs'] = vl_l2_abs
        if vl_linf_abs:
            me.steadyConvergenceCriteria['vl_linf_abs'] = vl_linf_abs
        if vv_l2_abs:
            me.steadyConvergenceCriteria['vv_l2_abs'] = vv_l2_abs
        if vv_linf_abs:
            me.steadyConvergenceCriteria['vv_linf_abs'] = vv_linf_abs
        if vd_l2_abs:
            me.steadyConvergenceCriteria['vd_l2_abs'] = vd_l2_abs
        if vd_linf_abs:
            me.steadyConvergenceCriteria['vd_linf_abs'] = vd_linf_abs
        if dtmin:
            me.dtMin_steady = dtmin
        if dtmax:
@@ -913,20 +1111,20 @@ class Model:
              end of the simulation.

        """
        if nc:
        if nc is not None:
            me.conductorOptions['nc'] = nc
            me.conductorOptionsSet['nc'] = True

        valid = ['none', 'w3', 'biasi', 'groeneveld']

        if chfModel:
        if chfModel is not None:
            if chfModel not in valid:
                raise ValueError("The selected CHF model: " +
                                 str(chfModel) + " is not valid")
            me.conductorOptions['chfModel'] = chfModel
            me.conductorOptionsSet['chfModel'] = True

        if postChfCheck:
        if postChfCheck is not None:
            if postChfCheck not in valid:
                raise ValueError("The selected CHF model: " +
                                 str(postChfCheck) + " is not valid")
@@ -947,19 +1145,19 @@ class Model:
            dnbEdits (bool): Set  to True to produce the CTF .dnb file, which contains information about DNB (default False).

        """
        if chanVTK:
        if chanVTK is not None:
            me.editOptions['fluid_vtk'] = chanVTK
        if rodVTK:
        if rodVTK is not None:
            me.editOptions['rod_vtk'] = rodVTK
        if ctfHDF5:
        if ctfHDF5 is not None:
            me.editOptions['native_hdf5'] = ctfHDF5
        if veracsHDF5:
        if veracsHDF5 is not None:
            me.editOptions['hdf5'] = veracsHDF5
        if chanASCII:
        if chanASCII is not None:
            me.editOptions['chan_edits'] = chanASCII
        if rodEdits:
        if rodEdits is not None:
            me.editOptions['rod_edits'] = rodEdits
        if dnbEdits:
        if dnbEdits is not None:
            me.editOptions['dnb_edits'] = dnbEdits
        return None

@@ -1322,12 +1520,8 @@ class Model:
                    else:
                        pinObj.gapConductance = 0.0

        # Cannot call them both
        assert(
            not (me.setInitialConditionsCalled and me.setInitialConditionsMassFluxCalled))

        # If neither were called, try and infer them from boundary conditions set on the channels
        if not (me.setInitialConditionsCalled or me.setInitialConditionsMassFluxCalled):
        if not me.setInitialCalled:
            mdot, temp, pressure = me._inferInitialConditionsFromBC()
            if mdot and temp and pressure:
                if me.modelParallel:
@@ -1336,13 +1530,13 @@ class Model:
                        raise RuntimeError("A parallel model with more than 1 section is being modeled.  setInitialConditionsMassFlux \
                     was not called and it cannot be infered from the channel boundary conditions that were set.")
                    else:
                        me.setInitialConditionsMassFlux(
                            mdot / me.sections.values()[0].getFlowArea(), temp, pressure)
                        me.setInitial(
                            massflux=mdot / me.sections.values()[0].getFlowArea(), temp=temp, pressure=pressure)
                else:
                    me.setInitialConditions(mdot, temp, pressure)
                    me.setInitial(mdot=mdot, temp=temp, pressure=pressure)
            else:
                raise RuntimeError(
                    "setInitialConditions was not called and the channels do not all have inlet mass/temp and outlet pressure BCs with uniform values, so cannot infer the correct initial conditions")
                    "`setInitial` was not called and the channels do not all have inlet mass/temp and outlet pressure BCs with uniform values, so cannot infer the correct initial conditions")

        # If the inlet mass flow rate distribution was non-uniform, the initial absolute mass flow rate
        # must be defined in CTF, mass flux set to zero, and inlet mass flow rate set to zero
@@ -1469,6 +1663,7 @@ class Model:
            for i, bc in enumerate(ch.bcLevels):
                if ch.bcLevels[i] == 1 and isinstance(ch.bcs[i], BoundaryCondition.MassTemperature):
                    return (ch.bcs[i].mdot, ch.bcs[i].T)
            return None, None

        def getPressureBCOutlet(ch, topLev):
            """ Returns the pressure of the outlet BC if the channel has one"""
@@ -1476,6 +1671,7 @@ class Model:
                if ch.bcLevels[i] == topLev and (isinstance(ch.bcs[i], BoundaryCondition.PressureTemperature) or
                                                 isinstance(ch.bcs[i], BoundaryCondition.PressureEnthalpy)):
                    return ch.bcs[i].pressure
            return None, None

        secID = me.sections.keys()[0]
        mdotInlet = 0.0
@@ -1641,6 +1837,9 @@ class Model:
            for ch in allChan:
                sections.append(me.getSectionOfChannel(ch))
            unique = list(set(sections))
            if len(innerChan) > 0 and len(unique) > 1:
                raise RuntimeError("For solidID: " + str(solidID) +
                                   " has internal connections and exists in multiple sections, which is not currently allowed.")
            # Loop through each section connected to this rod
            for u in unique:
                # In each section, ensure that all percentages of connections add up to 1.0
+55 −20

File changed.

Preview size limit exceeded, changes collapsed.

+4 −4
Original line number Diff line number Diff line
add_test(subkit_unit_tests test_unit_tests.sh ${CMAKE_CURRENT_BINARY_DIR}/../conda/bin/python)
add_test(subkit_input_tests test_input_mode.sh ${CMAKE_CURRENT_BINARY_DIR}/../conda/bin/python)
add_test(subkit_script_mode_tests test_script_mode.sh ${CMAKE_CURRENT_BINARY_DIR}/../conda/bin/python)
add_test(subkit_script_tests test_scripts.sh ${CMAKE_CURRENT_BINARY_DIR}/../conda/bin/python)
add_test(subkit_unit_tests test_unit_tests.sh)
add_test(subkit_input_tests test_input_mode.sh ../SubKit/build/InpBuild.py)
add_test(subkit_script_mode_tests test_script_mode.sh)
add_test(subkit_script_tests test_scripts.sh)

tests/ctfTraining/ge3x3/3x3.inp

deleted100644 → 0
+0 −582

File deleted.

Preview size limit exceeded, changes collapsed.

+3 −3
Original line number Diff line number Diff line
@@ -569,7 +569,7 @@ end 14
*Card 19.2 - absolute stopping criteria [pressure in psia|bar, velocity in ft/s|m/s,
** and temperature in F|C]
** LIAVOID  LIAPRESS  LIATCOOL LIATSOLID     LIAVL     LIAVV     LIAVD
 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04
 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04
*Card 19.2   [%]
** ENERGYBAL     MASSBAL
   1.000000e-01   1.000000e-01
@@ -579,4 +579,4 @@ end 14
 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04
** and temperature in F|C]
** L2AVOID  L2APRESS  L2ATCOOL L2ATSOLID     L2AVL     L2AVV     L2AVD
 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04
 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04 1.000e-04
Loading