This adds support for tube flow and several other features (cd6dcfdf) · Commits · Gurecky, William / SubKit

SubKit/build/Model.py

+114 −16

Original line number	Diff line number	Diff line
		@@ -127,7 +127,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
		@@ -570,6 +584,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.

		@@ -592,12 +613,22 @@ class Model:
		sections = [section]
		chans = []
		percents = []
		if not me.hasOuterConnections(solidID):
		raise RuntimeError("Solid ID: " + str(solidID) +
		" does not have any channels connected to outside surface.")
		for i, ch in enumerate(me.channelsOuter[solidID]):
		if me.getSectionOfChannel(ch) in sections:
		chans.append(ch)
		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 +639,15 @@ 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
		if not me.hasInnerConnections(solidID):
		raise RuntimeError("The solidID: " + str(solidID) +
		" does not have inner connections")
		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 +845,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 +895,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 +954,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:
		@@ -947,19 +1040,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

		@@ -1469,6 +1562,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 +1570,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 +1736,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

SubKit/build/genDeck.py

+31 −7

Original line number	Diff line number	Diff line
		@@ -381,8 +381,11 @@ def writeDeck(model, filename):
		group8Data.append("*Card 8.3\n")
		group8Data.append(
		"**NSCH PIE NSCH PIE NSCH PIE NSCH PIE NSCH PIE NSCH PIE NSCH PIE NSCH PIE\n")
		card84Written = False
		for pinID in heated.keys():
		intID = heatedIntIDs[pinID]
		if model.hasInnerConnections(pinID):
		intID = intID * -1
		pin = heated[pinID]
		if (pin.symtype == 4):
		sym = 0.25
		@@ -397,7 +400,11 @@ def writeDeck(model, filename):
		model.solidTypeIDs[pinID]), iaxp, 0, 0, pin.mult, pin.gapConductance,
		len(sections), 0.0, 0.0, sym, 0, 0, 1, 1))
		for section in sections:
		if model.hasOuterConnections(pinID):
		chans, percents = model.getOuterConnections(pinID, section)
		else:
		# Rods with internal connection can only exist in one section, which will have been asserted already
		chans, percents = model.getInnerConnections(pinID)
		if len(chans) > 8:
		raise RuntimeError(
		"CTF does not currently support more than 8 channels connected to one solid, see solid " + str(pinID))
		@@ -407,8 +414,25 @@ def writeDeck(model, filename):
		percents.extend([0.0])
		else:
		chans[i] = chMap.getIndex(chans[i])
		# If there are inside connections, these need to be entered on Card 8.4 and zeros enterd for the channel
		# index on Card 8.3
		if model.hasInnerConnections(pinID):
		chansOuter = [0] * len(chans)
		chansInner = chans
		else:
		chansOuter = chans

		group8Data.append("{0:6d}{1:6.3f}{2:6d}{3:6.3f}{4:6d}{5:6.3f}{6:6d}{7:6.3f}{8:6d}{9:6.3f}{10:6d}{11:6.3f}{12:6d}{13:6.3f}{14:6d}{15:6.3f}\n"
		.format(chans[0], percents[0], chans[1], percents[1], chans[2], percents[2], chans[3], percents[3], chans[4], percents[4], chans[5], percents[5], chans[6], percents[6], chans[7], percents[7]))
		.format(chansOuter[0], percents[0], chansOuter[1], percents[1], chansOuter[2], percents[2], chansOuter[3], percents[3], chansOuter[4], percents[4], chansOuter[5], percents[5], chansOuter[6], percents[6], chansOuter[7], percents[7]))

		if model.hasInnerConnections(pinID):

		if not card84Written:
		group8Data.append("*Card 8.4\n")
		group8Data.append(
		"**NISCH NISCH NISCH NISCH NISCH NISCH NISCH NISCH\n")
		group8Data.append("{:7d}{:7d}{:7d}{:7d}{:7d}{:7d}{:7d}{:7d}\n".format(
		chansInner[0], chansInner[1], chansInner[2], chansInner[3], chansInner[4], chansInner[5], chansInner[6], chansInner[7]))

		group8Data.append("*Card 8.5\n")
		group8Data.append(
		@@ -841,12 +865,12 @@ def writeDeck(model, filename):
		group19Data.append("** and temperature in F\|C]\n")
		group19Data.append(
		"** LIAVOID LIAPRESS LIATCOOL LIATSOLID LIAVL LIAVV LIAVD\n")
		group19Data.append(
		" 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04\n")
		group19Data.append("{:10.3e}{:10.3e}{:10.3e}{:10.3e}{:10.3e}{:10.3e}{:10.3e}\n".format(
		dic['void_linf_abs'], dic['p_linf_abs'], dic['Tcool_linf_abs'], dic['Tsolid_linf_abs'], dic['vl_linf_abs'], dic['vv_linf_abs'], dic['vd_linf_abs']))
		group19Data.append("*Card 19.2 [%]\n")
		group19Data.append("** ENERGYBAL MASSBAL\n")
		group19Data.append("{:15.6e}{:15.6e}\n".format(
		dic['energyBalance'], dic['massBalance']))
		dic['energyBalance'] * 100, dic['massBalance'] * 100))
		group19Data.append("*Card 19.4 [%]\n")
		group19Data.append(
		"** L2PRESS L2TCOOL L2TSOLID L2VL L2VV L2VD\n")
		@@ -857,8 +881,8 @@ def writeDeck(model, filename):
		group19Data.append("** and temperature in F\|C]\n")
		group19Data.append(
		"** L2AVOID L2APRESS L2ATCOOL L2ATSOLID L2AVL L2AVV L2AVD\n")
		group19Data.append(
		" 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04 1.000E-04\n")
		group19Data.append("{:10.3e}{:10.3e}{:10.3e}{:10.3e}{:10.3e}{:10.3e}{:10.3e}\n".format(
		dic['void_l2_abs'], dic['p_l2_abs'], dic['Tcool_l2_abs'], dic['Tsolid_l2_abs'], dic['vl_l2_abs'], dic['vv_l2_abs'], dic['vd_l2_abs']))

		def makeGroupBanner(title):
		assert(isinstance(title, str))

tests/test_script_mode.sh

+1 −0

Original line number	Diff line number	Diff line
		@@ -58,6 +58,7 @@ runScriptTest "msrTraining/fullCore/base" "genModel.py" "$rebaseline"
		runScriptTest "msrTraining/pincell" "genModel.py" "$rebaseline"
		runScriptTest "ctfTraining/nodal_lof" "lof.py" "$rebaseline"
		runScriptTest "ctfTraining/pincell" "pincell.py" "$rebaseline"
		runScriptTest "tubeFlow" "make_deck.py" "$rebaseline"

		# Put these at the end because they take forever
		runScriptTest "mcfr" "coreRegion.py" "$rebaseline"

tests/tubeFlow/deck.inp.gold

0 → 100644

+246 −0

File added.

Preview size limit exceeded, changes collapsed.

tests/tubeFlow/make_deck.py

0 → 100755

+67 −0

Original line number	Diff line number	Diff line
		# This test is testing the following new capabilities:
		# - The remainder of the steady-state stopping criteria (absolute units) have been added and are set through
		# this test
		# - Ability to select the pressure matrix solver was exposed and added here
		# - Ability to model flow inside a tube was added and tested here
		# - There was a bug in the setEdits procedure where if the user set something to False which defaults to True,
		# it would ignore it and use the True value
		# - There was a bug in the boundary condition getter in the Model class, which was getting hit by this script
		# and crashing the code
		# - The mass/energy balance units were wrong, which was fixed here.

		import sys
		sys.path.insert(0, '../../../')
		import SubKit.build as skb
		import SubKit.utils.UnitConversions as units
		import math


		def main():
		nLev = 320
		model = skb.Model()
		section = skb.Section(nLev=nLev, height=9.0)
		section.addChannel(chID=1, area=7.8540e-5, pw=3.1416e-2)
		model.addSection(secID=1, section=section)

		axial = [0.00000e+00, 4.87999e+00, 4.88000e+00,
		8.90000e+00, 8.90001e+00, 9.00000e+00]
		power = [0.0, 0.0, 1.0, 1.0, 0.0, 0.0]
		model.addAxialPowerShape(powID=1, axial=axial, power=power)

		# Define the outer tube that the fluid is flowing through

		# Tube thickness does not really matter, so just pick something
		tubeThick = 0.01 # m
		tubeInner = 0.01 # m

		outerTube = skb.Tube(d_outer=tubeInner + tubeThick, d_inner=tubeInner)
		model.addSolidType(solidTypeID=1, geoObj=outerTube)
		heatedObj = skb.HeatedSolid(outerTube, 1)
		model.addSolid(solidID=1, solidTypeID=1, solidObject=heatedObj)
		model.addSolidInnerChan(solidID=1, chID=1, percent=1.0)

		# Apply the boundary condition
		model.addBoundaryCondition(bc=skb.MassEnthalpy(
		mdot=7.85398e-02, h=3.05035e+02), chID=1, lev=1)
		model.addBoundaryCondition(bc=skb.PressureEnthalpy(
		pressure=30.0, h=3.05035e+02), chID=1, lev=nLev + 2)

		model.setAveragePower(qprime=31.41593)
		model.setInitialConditionsMassFlux(1000.0, 233.85845, 30.0)
		model.setSolver('direct')
		model.setFluidProperties('asme_1968')
		model.setEditOptions(chanVTK=False, rodVTK=False, ctfHDF5=True, veracsHDF5=False, chanASCII=False,
		rodEdits=False, dnbEdits=False)
		model.setSteadyState(energyBalance=1e-4, massBalance=1e-4, maxIterations=200000, rtwfp=1.0,
		p_l2=1e2, p_linf=1e-3, Tcool_l2=1e2, Tcool_linf=1e-3, Tsolid_l2=1e2, Tsolid_linf=1e99,
		vl_l2=1e2, vl_linf=1e-3, vv_l2=1e2, vv_linf=1e-3, vd_l2=1e2, vd_linf=1e-3,
		void_l2_abs=1e2, void_linf_abs=1e-3, p_l2_abs=1e2, p_linf_abs=1e-3, Tcool_l2_abs=1e2,
		Tcool_linf_abs=1e-3, Tsolid_l2_abs=1e2, Tsolid_linf_abs=1e99, vl_l2_abs=1e2, vl_linf_abs=1e-3,
		vv_l2_abs=1e2, vv_linf_abs=1e-3, vd_l2_abs=1e2, vd_linf_abs=1e-3,
		dtmin=1e-7, dtmax=1e-4)

		model.generateModel()


		if __name__ == "__main__":
		main()