https://mooseframework.inl.gov
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
QuadInterWrapper1PhaseProblem Class Referenceabstract

Quadrilateral interwrapper solver. More...

#include <QuadInterWrapper1PhaseProblem.h>

Inheritance diagram for QuadInterWrapper1PhaseProblem:
[legend]

Public Types

enum  Direction { Direction::TO_EXTERNAL_APP, Direction::FROM_EXTERNAL_APP }
 
enum  CoverageCheckMode {
  CoverageCheckMode::FALSE, CoverageCheckMode::TRUE, CoverageCheckMode::OFF, CoverageCheckMode::ON,
  CoverageCheckMode::SKIP_LIST, CoverageCheckMode::ONLY_LIST
}
 
typedef DataFileName DataFileParameterType
 

Public Member Functions

 QuadInterWrapper1PhaseProblem (const InputParameters &params)
 
virtual void externalSolve () override
 
virtual void syncSolutions (Direction direction) override
 
virtual bool solverSystemConverged (const unsigned int) override
 
virtual void initialSetup () override
 
virtual void solve (unsigned int nl_sys_num=0) override final
 
virtual void addExternalVariables ()
 
virtual libMesh::EquationSystemses () override
 
virtual MooseMeshmesh () override
 
virtual const MooseMeshmesh () const override
 
const MooseMeshmesh (bool use_displaced) const override
 
void setCoordSystem (const std::vector< SubdomainName > &blocks, const MultiMooseEnum &coord_sys)
 
void setAxisymmetricCoordAxis (const MooseEnum &rz_coord_axis)
 
void setCoupling (Moose::CouplingType type)
 
Moose::CouplingType coupling () const
 
void setCouplingMatrix (std::unique_ptr< libMesh::CouplingMatrix > cm, const unsigned int nl_sys_num)
 
void setCouplingMatrix (libMesh::CouplingMatrix *cm, const unsigned int nl_sys_num)
 
const libMesh::CouplingMatrixcouplingMatrix (const unsigned int nl_sys_num) const override
 
void setNonlocalCouplingMatrix ()
 
bool areCoupled (const unsigned int ivar, const unsigned int jvar, const unsigned int nl_sys_num) const
 
bool hasUOAuxStateCheck () const
 
bool checkingUOAuxState () const
 
void trustUserCouplingMatrix ()
 
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase *> > & couplingEntries (const THREAD_ID tid, const unsigned int nl_sys_num)
 
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase *> > & nonlocalCouplingEntries (const THREAD_ID tid, const unsigned int nl_sys_num)
 
virtual bool hasVariable (const std::string &var_name) const override
 
bool hasSolverVariable (const std::string &var_name) const
 
virtual const MooseVariableFieldBasegetVariable (const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const override
 
virtual const MooseVariableFieldBasegetVariable (const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const =0
 
virtual MooseVariableFieldBasegetVariable (const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY)
 
virtual MooseVariableFieldBasegetVariable (const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY)
 
MooseVariableFieldBasegetActualFieldVariable (const THREAD_ID tid, const std::string &var_name) override
 
virtual MooseVariablegetStandardVariable (const THREAD_ID tid, const std::string &var_name) override
 
virtual VectorMooseVariablegetVectorVariable (const THREAD_ID tid, const std::string &var_name) override
 
virtual ArrayMooseVariablegetArrayVariable (const THREAD_ID tid, const std::string &var_name) override
 
virtual bool hasScalarVariable (const std::string &var_name) const override
 
virtual MooseVariableScalargetScalarVariable (const THREAD_ID tid, const std::string &var_name) override
 
virtual libMesh::SystemgetSystem (const std::string &var_name) override
 
virtual void setActiveElementalMooseVariables (const std::set< MooseVariableFEBase * > &moose_vars, const THREAD_ID tid) override
 
virtual void clearActiveElementalMooseVariables (const THREAD_ID tid) override
 
virtual void clearActiveFEVariableCoupleableMatrixTags (const THREAD_ID tid) override
 
virtual void clearActiveFEVariableCoupleableVectorTags (const THREAD_ID tid) override
 
virtual void setActiveFEVariableCoupleableVectorTags (std::set< TagID > &vtags, const THREAD_ID tid) override
 
virtual void setActiveFEVariableCoupleableMatrixTags (std::set< TagID > &mtags, const THREAD_ID tid) override
 
virtual void clearActiveScalarVariableCoupleableMatrixTags (const THREAD_ID tid) override
 
virtual void clearActiveScalarVariableCoupleableVectorTags (const THREAD_ID tid) override
 
virtual void setActiveScalarVariableCoupleableVectorTags (std::set< TagID > &vtags, const THREAD_ID tid) override
 
virtual void setActiveScalarVariableCoupleableMatrixTags (std::set< TagID > &mtags, const THREAD_ID tid) override
 
virtual void createQRules (libMesh::QuadratureType type, libMesh::Order order, libMesh::Order volume_order=libMesh::INVALID_ORDER, libMesh::Order face_order=libMesh::INVALID_ORDER, SubdomainID block=Moose::ANY_BLOCK_ID, bool allow_negative_qweights=true)
 
void bumpVolumeQRuleOrder (libMesh::Order order, SubdomainID block)
 
void bumpAllQRuleOrder (libMesh::Order order, SubdomainID block)
 
unsigned int getMaxQps () const
 
libMesh::Order getMaxScalarOrder () const
 
void checkNonlocalCoupling ()
 
void checkUserObjectJacobianRequirement (THREAD_ID tid)
 
void setVariableAllDoFMap (const std::vector< const MooseVariableFEBase * > &moose_vars)
 
const std::vector< const MooseVariableFEBase *> & getUserObjectJacobianVariables (const THREAD_ID tid) const
 
virtual Assemblyassembly (const THREAD_ID tid, const unsigned int sys_num) override
 
virtual const Assemblyassembly (const THREAD_ID tid, const unsigned int sys_num) const override
 
virtual std::vector< VariableName > getVariableNames ()
 
void checkDuplicatePostprocessorVariableNames ()
 
void timestepSetup () override
 
void customSetup (const ExecFlagType &exec_type) override
 
void residualSetup () override
 
void jacobianSetup () override
 
virtual void prepare (const Elem *elem, const THREAD_ID tid) override
 
virtual void prepare (const Elem *elem, unsigned int ivar, unsigned int jvar, const std::vector< dof_id_type > &dof_indices, const THREAD_ID tid) override
 
virtual void prepareFace (const Elem *elem, const THREAD_ID tid) override
 
virtual void setCurrentSubdomainID (const Elem *elem, const THREAD_ID tid) override
 
virtual void setNeighborSubdomainID (const Elem *elem, unsigned int side, const THREAD_ID tid) override
 
virtual void setNeighborSubdomainID (const Elem *elem, const THREAD_ID tid)
 
virtual void prepareAssembly (const THREAD_ID tid) override
 
virtual void addGhostedElem (dof_id_type elem_id) override
 
virtual void addGhostedBoundary (BoundaryID boundary_id) override
 
virtual void ghostGhostedBoundaries () override
 
virtual void sizeZeroes (unsigned int size, const THREAD_ID tid)
 
virtual bool reinitDirac (const Elem *elem, const THREAD_ID tid) override
 
virtual void reinitElem (const Elem *elem, const THREAD_ID tid) override
 
virtual void reinitElemPhys (const Elem *elem, const std::vector< Point > &phys_points_in_elem, const THREAD_ID tid) override
 
void reinitElemFace (const Elem *elem, unsigned int side, BoundaryID, const THREAD_ID tid)
 
virtual void reinitElemFace (const Elem *elem, unsigned int side, const THREAD_ID tid) override
 
virtual void reinitLowerDElem (const Elem *lower_d_elem, const THREAD_ID tid, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr) override
 
virtual void reinitNode (const Node *node, const THREAD_ID tid) override
 
virtual void reinitNodeFace (const Node *node, BoundaryID bnd_id, const THREAD_ID tid) override
 
virtual void reinitNodes (const std::vector< dof_id_type > &nodes, const THREAD_ID tid) override
 
virtual void reinitNodesNeighbor (const std::vector< dof_id_type > &nodes, const THREAD_ID tid) override
 
virtual void reinitNeighbor (const Elem *elem, unsigned int side, const THREAD_ID tid) override
 
virtual void reinitNeighborPhys (const Elem *neighbor, unsigned int neighbor_side, const std::vector< Point > &physical_points, const THREAD_ID tid) override
 
virtual void reinitNeighborPhys (const Elem *neighbor, const std::vector< Point > &physical_points, const THREAD_ID tid) override
 
virtual void reinitElemNeighborAndLowerD (const Elem *elem, unsigned int side, const THREAD_ID tid) override
 
virtual void reinitScalars (const THREAD_ID tid, bool reinit_for_derivative_reordering=false) override
 
virtual void reinitOffDiagScalars (const THREAD_ID tid) override
 
virtual void getDiracElements (std::set< const Elem * > &elems) override
 
virtual void clearDiracInfo () override
 
virtual void subdomainSetup (SubdomainID subdomain, const THREAD_ID tid)
 
virtual void neighborSubdomainSetup (SubdomainID subdomain, const THREAD_ID tid)
 
virtual void newAssemblyArray (std::vector< std::shared_ptr< SolverSystem >> &solver_systems)
 
virtual void initNullSpaceVectors (const InputParameters &parameters, std::vector< std::shared_ptr< NonlinearSystemBase >> &nl)
 
virtual void init () override
 
virtual void solveLinearSystem (const unsigned int linear_sys_num, const Moose::PetscSupport::PetscOptions *po=nullptr)
 
virtual void setException (const std::string &message)
 
virtual bool hasException ()
 
virtual void checkExceptionAndStopSolve (bool print_message=true)
 
virtual unsigned int nNonlinearIterations (const unsigned int nl_sys_num) const override
 
virtual unsigned int nLinearIterations (const unsigned int nl_sys_num) const override
 
virtual Real finalNonlinearResidual (const unsigned int nl_sys_num) const override
 
virtual bool computingPreSMOResidual (const unsigned int nl_sys_num) const override
 
virtual std::string solverTypeString (unsigned int solver_sys_num=0)
 
virtual bool startedInitialSetup ()
 
virtual void onTimestepBegin () override
 
virtual void onTimestepEnd () override
 
virtual Realtime () const
 
virtual RealtimeOld () const
 
virtual inttimeStep () const
 
virtual Realdt () const
 
virtual RealdtOld () const
 
Real getTimeFromStateArg (const Moose::StateArg &state) const
 
virtual void transient (bool trans)
 
virtual bool isTransient () const override
 
virtual void addTimeIntegrator (const std::string &type, const std::string &name, InputParameters &parameters)
 
virtual void addPredictor (const std::string &type, const std::string &name, InputParameters &parameters)
 
virtual void copySolutionsBackwards ()
 
virtual void advanceState ()
 
virtual void restoreSolutions ()
 
virtual void saveOldSolutions ()
 
virtual void restoreOldSolutions ()
 
void needSolutionState (unsigned int oldest_needed, Moose::SolutionIterationType iteration_type)
 
virtual void outputStep (ExecFlagType type)
 
virtual void postExecute ()
 
void forceOutput ()
 
virtual void initPetscOutputAndSomeSolverSettings ()
 
Moose::PetscSupport::PetscOptionsgetPetscOptions ()
 
void logAdd (const std::string &system, const std::string &name, const std::string &type, const InputParameters &params) const
 
virtual void addFunction (const std::string &type, const std::string &name, InputParameters &parameters)
 
virtual bool hasFunction (const std::string &name, const THREAD_ID tid=0)
 
virtual FunctiongetFunction (const std::string &name, const THREAD_ID tid=0)
 
virtual void addMeshDivision (const std::string &type, const std::string &name, InputParameters &params)
 
MeshDivisiongetMeshDivision (const std::string &name, const THREAD_ID tid=0) const
 
virtual void addConvergence (const std::string &type, const std::string &name, InputParameters &parameters)
 
virtual ConvergencegetConvergence (const std::string &name, const THREAD_ID tid=0) const
 
virtual const std::vector< std::shared_ptr< Convergence > > & getConvergenceObjects (const THREAD_ID tid=0) const
 
virtual bool hasConvergence (const std::string &name, const THREAD_ID tid=0) const
 
bool needToAddDefaultNonlinearConvergence () const
 
bool needToAddDefaultMultiAppFixedPointConvergence () const
 
void setNeedToAddDefaultNonlinearConvergence ()
 
void setNeedToAddDefaultMultiAppFixedPointConvergence ()
 
bool hasSetMultiAppFixedPointConvergenceName () const
 
virtual void addDefaultNonlinearConvergence (const InputParameters &params)
 
virtual bool onlyAllowDefaultNonlinearConvergence () const
 
void addDefaultMultiAppFixedPointConvergence (const InputParameters &params)
 
virtual void addLineSearch (const InputParameters &)
 
virtual void lineSearch ()
 
LineSearchgetLineSearch () override
 
virtual void addDistribution (const std::string &type, const std::string &name, InputParameters &parameters)
 
virtual DistributiongetDistribution (const std::string &name)
 
virtual void addSampler (const std::string &type, const std::string &name, InputParameters &parameters)
 
virtual SamplergetSampler (const std::string &name, const THREAD_ID tid=0)
 
NonlinearSystemBasegetNonlinearSystemBase (const unsigned int sys_num)
 
const NonlinearSystemBasegetNonlinearSystemBase (const unsigned int sys_num) const
 
void setCurrentNonlinearSystem (const unsigned int nl_sys_num)
 
NonlinearSystemBasecurrentNonlinearSystem ()
 
const NonlinearSystemBasecurrentNonlinearSystem () const
 
virtual const SystemBasesystemBaseNonlinear (const unsigned int sys_num) const override
 
virtual SystemBasesystemBaseNonlinear (const unsigned int sys_num) override
 
virtual const SystemBasesystemBaseSolver (const unsigned int sys_num) const override
 
virtual SystemBasesystemBaseSolver (const unsigned int sys_num) override
 
virtual const SystemBasesystemBaseAuxiliary () const override
 
virtual SystemBasesystemBaseAuxiliary () override
 
virtual NonlinearSystemgetNonlinearSystem (const unsigned int sys_num)
 
virtual const SystemBasegetSystemBase (const unsigned int sys_num) const
 
virtual SystemBasegetSystemBase (const unsigned int sys_num)
 
LinearSystemgetLinearSystem (unsigned int sys_num)
 
const LinearSystemgetLinearSystem (unsigned int sys_num) const
 
SolverSystemgetSolverSystem (unsigned int sys_num)
 
const SolverSystemgetSolverSystem (unsigned int sys_num) const
 
void setCurrentLinearSystem (unsigned int sys_num)
 
LinearSystemcurrentLinearSystem ()
 
const LinearSystemcurrentLinearSystem () const
 
virtual const SystemBasesystemBaseLinear (unsigned int sys_num) const override
 
virtual SystemBasesystemBaseLinear (unsigned int sys_num) override
 
virtual void addVariable (const std::string &var_type, const std::string &var_name, InputParameters &params)
 
virtual void addKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addHDGKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addNodalKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addScalarKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addBoundaryCondition (const std::string &bc_name, const std::string &name, InputParameters &parameters)
 
virtual void addConstraint (const std::string &c_name, const std::string &name, InputParameters &parameters)
 
virtual void setInputParametersFEProblem (InputParameters &parameters)
 
virtual void addAuxVariable (const std::string &var_type, const std::string &var_name, InputParameters &params)
 
virtual void addAuxVariable (const std::string &var_name, const libMesh::FEType &type, const std::set< SubdomainID > *const active_subdomains=NULL)
 
virtual void addAuxArrayVariable (const std::string &var_name, const libMesh::FEType &type, unsigned int components, const std::set< SubdomainID > *const active_subdomains=NULL)
 
virtual void addAuxScalarVariable (const std::string &var_name, libMesh::Order order, Real scale_factor=1., const std::set< SubdomainID > *const active_subdomains=NULL)
 
virtual void addAuxKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addAuxScalarKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
AuxiliarySystemgetAuxiliarySystem ()
 
virtual void addDiracKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addDGKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addFVKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addLinearFVKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addFVBC (const std::string &fv_bc_name, const std::string &name, InputParameters &parameters)
 
virtual void addLinearFVBC (const std::string &fv_bc_name, const std::string &name, InputParameters &parameters)
 
virtual void addFVInterfaceKernel (const std::string &fv_ik_name, const std::string &name, InputParameters &parameters)
 
virtual void addInterfaceKernel (const std::string &kernel_name, const std::string &name, InputParameters &parameters)
 
virtual void addInitialCondition (const std::string &ic_name, const std::string &name, InputParameters &parameters)
 
virtual void addFVInitialCondition (const std::string &ic_name, const std::string &name, InputParameters &parameters)
 
void projectSolution ()
 
unsigned short getCurrentICState ()
 
void projectInitialConditionOnCustomRange (libMesh::ConstElemRange &elem_range, ConstBndNodeRange &bnd_node_range)
 
virtual void addMaterial (const std::string &material_name, const std::string &name, InputParameters &parameters)
 
virtual void addMaterialHelper (std::vector< MaterialWarehouse * > warehouse, const std::string &material_name, const std::string &name, InputParameters &parameters)
 
virtual void addInterfaceMaterial (const std::string &material_name, const std::string &name, InputParameters &parameters)
 
virtual void addFunctorMaterial (const std::string &functor_material_name, const std::string &name, InputParameters &parameters)
 
void prepareMaterials (const std::unordered_set< unsigned int > &consumer_needed_mat_props, const SubdomainID blk_id, const THREAD_ID tid)
 
void reinitMaterials (SubdomainID blk_id, const THREAD_ID tid, bool swap_stateful=true)
 
void reinitMaterialsFace (SubdomainID blk_id, const THREAD_ID tid, bool swap_stateful=true, const std::deque< MaterialBase * > *reinit_mats=nullptr)
 
void reinitMaterialsNeighbor (SubdomainID blk_id, const THREAD_ID tid, bool swap_stateful=true, const std::deque< MaterialBase * > *reinit_mats=nullptr)
 
void reinitMaterialsBoundary (BoundaryID boundary_id, const THREAD_ID tid, bool swap_stateful=true, const std::deque< MaterialBase * > *reinit_mats=nullptr)
 
void reinitMaterialsInterface (BoundaryID boundary_id, const THREAD_ID tid, bool swap_stateful=true)
 
virtual void swapBackMaterials (const THREAD_ID tid)
 
virtual void swapBackMaterialsFace (const THREAD_ID tid)
 
virtual void swapBackMaterialsNeighbor (const THREAD_ID tid)
 
void setActiveMaterialProperties (const std::unordered_set< unsigned int > &mat_prop_ids, const THREAD_ID tid)
 
bool hasActiveMaterialProperties (const THREAD_ID tid) const
 
void clearActiveMaterialProperties (const THREAD_ID tid)
 
std::vector< std::shared_ptr< T > > addObject (const std::string &type, const std::string &name, InputParameters &parameters, const bool threaded=true, const std::string &var_param_name="variable")
 
virtual void addPostprocessor (const std::string &pp_name, const std::string &name, InputParameters &parameters)
 
virtual void addVectorPostprocessor (const std::string &pp_name, const std::string &name, InputParameters &parameters)
 
virtual void addReporter (const std::string &type, const std::string &name, InputParameters &parameters)
 
const ReporterDatagetReporterData () const
 
ReporterDatagetReporterData (ReporterData::WriteKey)
 
virtual std::vector< std::shared_ptr< UserObject > > addUserObject (const std::string &user_object_name, const std::string &name, InputParameters &parameters)
 
const ExecuteMooseObjectWarehouse< UserObject > & getUserObjects () const
 
T & getUserObject (const std::string &name, unsigned int tid=0) const
 
const UserObjectgetUserObjectBase (const std::string &name, const THREAD_ID tid=0) const
 
const PositionsgetPositionsObject (const std::string &name) const
 
bool hasUserObject (const std::string &name) const
 
bool hasPostprocessorValueByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueByName (const PostprocessorName &name, std::size_t t_index=0) const
 
void setPostprocessorValueByName (const PostprocessorName &name, const PostprocessorValue &value, std::size_t t_index=0)
 
bool hasPostprocessor (const std::string &name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueByName (const std::string &object_name, const std::string &vector_name, std::size_t t_index=0) const
 
void setVectorPostprocessorValueByName (const std::string &object_name, const std::string &vector_name, const VectorPostprocessorValue &value, std::size_t t_index=0)
 
const VectorPostprocessorgetVectorPostprocessorObjectByName (const std::string &object_name, const THREAD_ID tid=0) const
 
virtual void addDamper (const std::string &damper_name, const std::string &name, InputParameters &parameters)
 
void setupDampers ()
 
bool hasDampers ()
 
virtual void addIndicator (const std::string &indicator_name, const std::string &name, InputParameters &parameters)
 
virtual void addMarker (const std::string &marker_name, const std::string &name, InputParameters &parameters)
 
virtual void addMultiApp (const std::string &multi_app_name, const std::string &name, InputParameters &parameters)
 
std::shared_ptr< MultiAppgetMultiApp (const std::string &multi_app_name) const
 
std::vector< std::shared_ptr< Transfer > > getTransfers (ExecFlagType type, Transfer::DIRECTION direction) const
 
std::vector< std::shared_ptr< Transfer > > getTransfers (Transfer::DIRECTION direction) const
 
const ExecuteMooseObjectWarehouse< Transfer > & getMultiAppTransferWarehouse (Transfer::DIRECTION direction) const
 
void execMultiAppTransfers (ExecFlagType type, Transfer::DIRECTION direction)
 
bool execMultiApps (ExecFlagType type, bool auto_advance=true)
 
void finalizeMultiApps ()
 
void incrementMultiAppTStep (ExecFlagType type)
 
void advanceMultiApps (ExecFlagType type)
 
void finishMultiAppStep (ExecFlagType type, bool recurse_through_multiapp_levels=false)
 
void backupMultiApps (ExecFlagType type)
 
void restoreMultiApps (ExecFlagType type, bool force=false)
 
Real computeMultiAppsDT (ExecFlagType type)
 
virtual void addTransfer (const std::string &transfer_name, const std::string &name, InputParameters &parameters)
 
void execTransfers (ExecFlagType type)
 
Real computeResidualL2Norm (NonlinearSystemBase &sys)
 
Real computeResidualL2Norm (LinearSystem &sys)
 
virtual Real computeResidualL2Norm ()
 
virtual void computeResidualSys (libMesh::NonlinearImplicitSystem &sys, const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual)
 
void computeResidual (libMesh::NonlinearImplicitSystem &sys, const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual)
 
virtual void computeResidual (const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual, const unsigned int nl_sys_num)
 
void computeResidualAndJacobian (const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual, libMesh::SparseMatrix< libMesh::Number > &jacobian)
 
virtual void computeResidualTag (const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual, TagID tag)
 
virtual void computeResidualType (const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual, TagID tag)
 
virtual void computeResidualInternal (const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual, const std::set< TagID > &tags)
 
virtual void computeResidualTags (const std::set< TagID > &tags)
 
virtual void computeJacobianSys (libMesh::NonlinearImplicitSystem &sys, const NumericVector< libMesh::Number > &soln, libMesh::SparseMatrix< libMesh::Number > &jacobian)
 
virtual void computeJacobian (const NumericVector< libMesh::Number > &soln, libMesh::SparseMatrix< libMesh::Number > &jacobian, const unsigned int nl_sys_num)
 
virtual void computeJacobianTag (const NumericVector< libMesh::Number > &soln, libMesh::SparseMatrix< libMesh::Number > &jacobian, TagID tag)
 
virtual void computeJacobianInternal (const NumericVector< libMesh::Number > &soln, libMesh::SparseMatrix< libMesh::Number > &jacobian, const std::set< TagID > &tags)
 
virtual void computeJacobianTags (const std::set< TagID > &tags)
 
virtual void computeJacobianBlocks (std::vector< JacobianBlock * > &blocks, const unsigned int nl_sys_num)
 
virtual void computeJacobianBlock (libMesh::SparseMatrix< libMesh::Number > &jacobian, libMesh::System &precond_system, unsigned int ivar, unsigned int jvar)
 
virtual void computeLinearSystemSys (libMesh::LinearImplicitSystem &sys, libMesh::SparseMatrix< libMesh::Number > &system_matrix, NumericVector< libMesh::Number > &rhs, const bool compute_gradients=true)
 
void computeLinearSystemTags (const NumericVector< libMesh::Number > &soln, const std::set< TagID > &vector_tags, const std::set< TagID > &matrix_tags, const bool compute_gradients=true)
 
virtual Real computeDamping (const NumericVector< libMesh::Number > &soln, const NumericVector< libMesh::Number > &update)
 
virtual bool shouldUpdateSolution ()
 
virtual bool updateSolution (NumericVector< libMesh::Number > &vec_solution, NumericVector< libMesh::Number > &ghosted_solution)
 
virtual void predictorCleanup (NumericVector< libMesh::Number > &ghosted_solution)
 
virtual void computeBounds (libMesh::NonlinearImplicitSystem &sys, NumericVector< libMesh::Number > &lower, NumericVector< libMesh::Number > &upper)
 
virtual void computeNearNullSpace (libMesh::NonlinearImplicitSystem &sys, std::vector< NumericVector< libMesh::Number > * > &sp)
 
virtual void computeNullSpace (libMesh::NonlinearImplicitSystem &sys, std::vector< NumericVector< libMesh::Number > * > &sp)
 
virtual void computeTransposeNullSpace (libMesh::NonlinearImplicitSystem &sys, std::vector< NumericVector< libMesh::Number > * > &sp)
 
virtual void computePostCheck (libMesh::NonlinearImplicitSystem &sys, const NumericVector< libMesh::Number > &old_soln, NumericVector< libMesh::Number > &search_direction, NumericVector< libMesh::Number > &new_soln, bool &changed_search_direction, bool &changed_new_soln)
 
virtual void computeIndicatorsAndMarkers ()
 
virtual void computeIndicators ()
 
virtual void computeMarkers ()
 
virtual void addResidual (const THREAD_ID tid) override
 
virtual void addResidualNeighbor (const THREAD_ID tid) override
 
virtual void addResidualLower (const THREAD_ID tid) override
 
virtual void addResidualScalar (const THREAD_ID tid=0)
 
virtual void cacheResidual (const THREAD_ID tid) override
 
virtual void cacheResidualNeighbor (const THREAD_ID tid) override
 
virtual void addCachedResidual (const THREAD_ID tid) override
 
virtual void addCachedResidualDirectly (NumericVector< libMesh::Number > &residual, const THREAD_ID tid)
 
virtual void setResidual (NumericVector< libMesh::Number > &residual, const THREAD_ID tid) override
 
virtual void setResidual (libMesh::NumericVector< libMesh::Number > &residual, const THREAD_ID tid)=0
 
virtual void setResidualNeighbor (NumericVector< libMesh::Number > &residual, const THREAD_ID tid) override
 
virtual void setResidualNeighbor (libMesh::NumericVector< libMesh::Number > &residual, const THREAD_ID tid)=0
 
virtual void addJacobian (const THREAD_ID tid) override
 
virtual void addJacobianNeighbor (const THREAD_ID tid) override
 
virtual void addJacobianNeighbor (libMesh::SparseMatrix< libMesh::Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, std::vector< dof_id_type > &dof_indices, std::vector< dof_id_type > &neighbor_dof_indices, const std::set< TagID > &tags, const THREAD_ID tid) override
 
virtual void addJacobianNeighbor (libMesh::SparseMatrix< libMesh::Number > &jacobian, unsigned int ivar, unsigned int jvar, const libMesh::DofMap &dof_map, std::vector< dof_id_type > &dof_indices, std::vector< dof_id_type > &neighbor_dof_indices, const std::set< TagID > &tags, const THREAD_ID tid)=0
 
virtual void addJacobianNeighborLowerD (const THREAD_ID tid) override
 
virtual void addJacobianLowerD (const THREAD_ID tid) override
 
virtual void addJacobianBlockTags (libMesh::SparseMatrix< libMesh::Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, std::vector< dof_id_type > &dof_indices, const std::set< TagID > &tags, const THREAD_ID tid)
 
virtual void addJacobianScalar (const THREAD_ID tid=0)
 
virtual void addJacobianOffDiagScalar (unsigned int ivar, const THREAD_ID tid=0)
 
virtual void cacheJacobian (const THREAD_ID tid) override
 
virtual void cacheJacobianNeighbor (const THREAD_ID tid) override
 
virtual void addCachedJacobian (const THREAD_ID tid) override
 
virtual void prepareShapes (unsigned int var, const THREAD_ID tid) override
 
virtual void prepareFaceShapes (unsigned int var, const THREAD_ID tid) override
 
virtual void prepareNeighborShapes (unsigned int var, const THREAD_ID tid) override
 
virtual void addDisplacedProblem (std::shared_ptr< DisplacedProblem > displaced_problem)
 
virtual std::shared_ptr< const DisplacedProblemgetDisplacedProblem () const
 
virtual std::shared_ptr< DisplacedProblemgetDisplacedProblem ()
 
virtual void updateGeomSearch (GeometricSearchData::GeometricSearchType type=GeometricSearchData::ALL) override
 
virtual void updateMortarMesh ()
 
void createMortarInterface (const std::pair< BoundaryID, BoundaryID > &primary_secondary_boundary_pair, const std::pair< SubdomainID, SubdomainID > &primary_secondary_subdomain_pair, bool on_displaced, bool periodic, const bool debug, const bool correct_edge_dropping, const Real minimum_projection_angle)
 
const std::unordered_map< std::pair< BoundaryID, BoundaryID >, AutomaticMortarGeneration > & getMortarInterfaces (bool on_displaced) const
 
virtual void possiblyRebuildGeomSearchPatches ()
 
virtual GeometricSearchDatageomSearchData () override
 
void setRestartFile (const std::string &file_name)
 
const MaterialPropertyRegistrygetMaterialPropertyRegistry () const
 
const InitialConditionWarehousegetInitialConditionWarehouse () const
 
const FVInitialConditionWarehousegetFVInitialConditionWarehouse () const
 
SolverParamssolverParams (unsigned int solver_sys_num=0)
 
const SolverParamssolverParams (unsigned int solver_sys_num=0) const
 
Adaptivityadaptivity ()
 
virtual void initialAdaptMesh ()
 
virtual bool adaptMesh ()
 
unsigned int getNumCyclesCompleted ()
 
bool hasInitialAdaptivity () const
 
bool hasInitialAdaptivity () const
 
void initXFEM (std::shared_ptr< XFEMInterface > xfem)
 
std::shared_ptr< XFEMInterfacegetXFEM ()
 
bool haveXFEM ()
 
virtual bool updateMeshXFEM ()
 
virtual void meshChanged (bool intermediate_change, bool contract_mesh, bool clean_refinement_flags)
 
void notifyWhenMeshChanges (MeshChangedInterface *mci)
 
void notifyWhenMeshDisplaces (MeshDisplacedInterface *mdi)
 
void initElementStatefulProps (const libMesh::ConstElemRange &elem_range, const bool threaded)
 
virtual void checkProblemIntegrity ()
 
void registerRandomInterface (RandomInterface &random_interface, const std::string &name)
 
void setConstJacobian (bool state)
 
void setKernelCoverageCheck (CoverageCheckMode mode)
 
void setKernelCoverageCheck (bool flag)
 
void setKernelCoverageCheck (CoverageCheckMode mode)
 
void setMaterialCoverageCheck (CoverageCheckMode mode)
 
void setMaterialCoverageCheck (bool flag)
 
void setMaterialCoverageCheck (CoverageCheckMode mode)
 
void setParallelBarrierMessaging (bool flag)
 
void setVerboseProblem (bool verbose)
 
bool verboseMultiApps () const
 
void parentOutputPositionChanged ()
 
unsigned int subspaceDim (const std::string &prefix) const
 
const MaterialWarehousegetMaterialWarehouse () const
 
const MaterialWarehousegetRegularMaterialsWarehouse () const
 
const MaterialWarehousegetDiscreteMaterialWarehouse () const
 
const MaterialWarehousegetInterfaceMaterialsWarehouse () const
 
std::shared_ptr< MaterialBasegetMaterial (std::string name, Moose::MaterialDataType type, const THREAD_ID tid=0, bool no_warn=false)
 
MaterialDatagetMaterialData (Moose::MaterialDataType type, const THREAD_ID tid=0) const
 
bool restoreOriginalNonzeroPattern () const
 
bool errorOnJacobianNonzeroReallocation () const
 
void setErrorOnJacobianNonzeroReallocation (bool state)
 
bool preserveMatrixSparsityPattern () const
 
void setPreserveMatrixSparsityPattern (bool preserve)
 
bool ignoreZerosInJacobian () const
 
void setIgnoreZerosInJacobian (bool state)
 
bool acceptInvalidSolution () const
 
bool allowInvalidSolution () const
 
bool showInvalidSolutionConsole () const
 
bool immediatelyPrintInvalidSolution () const
 
bool hasTimeIntegrator () const
 
virtual void execute (const ExecFlagType &exec_type)
 
virtual void executeAllObjects (const ExecFlagType &exec_type)
 
virtual ExecutorgetExecutor (const std::string &name)
 
virtual void computeUserObjects (const ExecFlagType &type, const Moose::AuxGroup &group)
 
virtual void computeUserObjectByName (const ExecFlagType &type, const Moose::AuxGroup &group, const std::string &name)
 
void needsPreviousNewtonIteration (bool state)
 
bool needsPreviousNewtonIteration () const
 
ExecuteMooseObjectWarehouse< Control > & getControlWarehouse ()
 
void executeControls (const ExecFlagType &exec_type)
 
void executeSamplers (const ExecFlagType &exec_type)
 
virtual void updateActiveObjects ()
 
void reportMooseObjectDependency (MooseObject *a, MooseObject *b)
 
ExecuteMooseObjectWarehouse< MultiApp > & getMultiAppWarehouse ()
 
bool hasJacobian () const
 
bool constJacobian () const
 
void addOutput (const std::string &, const std::string &, InputParameters &)
 
TheWarehousetheWarehouse () const
 
void setSNESMFReuseBase (bool reuse, bool set_by_user)
 
bool useSNESMFReuseBase ()
 
void skipExceptionCheck (bool skip_exception_check)
 
bool isSNESMFReuseBaseSetbyUser ()
 
bool & petscOptionsInserted ()
 
PetscOptions & petscOptionsDatabase ()
 
virtual void setUDotRequested (const bool u_dot_requested)
 
virtual void setUDotDotRequested (const bool u_dotdot_requested)
 
virtual void setUDotOldRequested (const bool u_dot_old_requested)
 
virtual void setUDotDotOldRequested (const bool u_dotdot_old_requested)
 
virtual bool uDotRequested ()
 
virtual bool uDotDotRequested ()
 
virtual bool uDotOldRequested ()
 
virtual bool uDotDotOldRequested ()
 
void haveADObjects (bool have_ad_objects) override
 
virtual void haveADObjects (bool have_ad_objects)
 
bool haveADObjects () const
 
bool haveADObjects () const
 
bool shouldSolve () const
 
const MortarDatamortarData () const
 
MortarDatamortarData ()
 
virtual bool hasNeighborCoupling () const
 
virtual bool hasMortarCoupling () const
 
void computingNonlinearResid (bool computing_nonlinear_residual) final
 
bool computingNonlinearResid () const
 
virtual void computingNonlinearResid (const bool computing_nonlinear_residual)
 
bool computingNonlinearResid () const
 
void setCurrentlyComputingResidual (bool currently_computing_residual) final
 
void numGridSteps (unsigned int num_grid_steps)
 
void uniformRefine ()
 
void automaticScaling (bool automatic_scaling) override
 
virtual void automaticScaling (bool automatic_scaling)
 
bool automaticScaling () const
 
bool automaticScaling () const
 
virtual void reinitElemFaceRef (const Elem *elem, unsigned int side, Real tolerance, const std::vector< Point > *const pts, const std::vector< Real > *const weights=nullptr, const THREAD_ID tid=0) override
 
virtual void reinitNeighborFaceRef (const Elem *neighbor_elem, unsigned int neighbor_side, Real tolerance, const std::vector< Point > *const pts, const std::vector< Real > *const weights=nullptr, const THREAD_ID tid=0) override
 
bool fvBCsIntegrityCheck () const
 
void fvBCsIntegrityCheck (bool fv_bcs_integrity_check)
 
void getFVMatsAndDependencies (SubdomainID block_id, std::vector< std::shared_ptr< MaterialBase >> &face_materials, std::vector< std::shared_ptr< MaterialBase >> &neighbor_materials, std::set< MooseVariableFieldBase * > &variables, const THREAD_ID tid)
 
void resizeMaterialData (Moose::MaterialDataType data_type, unsigned int nqp, const THREAD_ID tid)
 
bool haveDisplaced () const override final
 
bool hasLinearConvergenceObjects () const
 
void setNonlinearConvergenceNames (const std::vector< ConvergenceName > &convergence_names)
 
void setLinearConvergenceNames (const std::vector< ConvergenceName > &convergence_names)
 
void setMultiAppFixedPointConvergenceName (const ConvergenceName &convergence_name)
 
const std::vector< ConvergenceName > & getNonlinearConvergenceNames () const
 
const std::vector< ConvergenceName > & getLinearConvergenceNames () const
 
const ConvergenceName & getMultiAppFixedPointConvergenceName () const
 
void computingScalingJacobian (bool computing_scaling_jacobian)
 
bool computingScalingJacobian () const override final
 
void computingScalingResidual (bool computing_scaling_residual)
 
bool computingScalingResidual () const override final
 
MooseAppCoordTransformcoordTransform ()
 
virtual std::size_t numNonlinearSystems () const override
 
virtual std::size_t numLinearSystems () const override
 
virtual std::size_t numSolverSystems () const override
 
bool isSolverSystemNonlinear (const unsigned int sys_num)
 
virtual unsigned int currentNlSysNum () const override
 
virtual unsigned int currentLinearSysNum () const override
 
virtual unsigned int nlSysNum (const NonlinearSystemName &nl_sys_name) const override
 
unsigned int linearSysNum (const LinearSystemName &linear_sys_name) const override
 
unsigned int solverSysNum (const SolverSystemName &solver_sys_name) const override
 
unsigned int systemNumForVariable (const VariableName &variable_name) const
 
bool getFailNextNonlinearConvergenceCheck () const
 
bool getFailNextSystemConvergenceCheck () const
 
void setFailNextNonlinearConvergenceCheck ()
 
void setFailNextSystemConvergenceCheck ()
 
void resetFailNextNonlinearConvergenceCheck ()
 
void resetFailNextSystemConvergenceCheck ()
 
void setExecutionPrinting (const ExecFlagEnum &print_exec)
 
bool shouldPrintExecution (const THREAD_ID tid) const
 
void reinitMortarUserObjects (BoundaryID primary_boundary_id, BoundaryID secondary_boundary_id, bool displaced)
 
virtual const std::vector< VectorTag > & currentResidualVectorTags () const override
 
void setCurrentResidualVectorTags (const std::set< TagID > &vector_tags)
 
void clearCurrentResidualVectorTags ()
 
void clearCurrentJacobianMatrixTags ()
 
virtual void needFV () override
 
virtual bool haveFV () const override
 
virtual bool hasNonlocalCoupling () const override
 
bool identifyVariableGroupsInNL () const
 
virtual void setCurrentLowerDElem (const Elem *const lower_d_elem, const THREAD_ID tid) override
 
virtual void setCurrentBoundaryID (BoundaryID bid, const THREAD_ID tid) override
 
const std::vector< NonlinearSystemName > & getNonlinearSystemNames () const
 
const std::vector< LinearSystemName > & getLinearSystemNames () const
 
const std::vector< SolverSystemName > & getSolverSystemNames () const
 
virtual const libMesh::CouplingMatrixnonlocalCouplingMatrix (const unsigned i) const override
 
virtual bool checkNonlocalCouplingRequirement () const override
 
virtual Moose::FEBackend feBackend () const
 
const bool & currentlyComputingResidual () const
 
const bool & currentlyComputingResidual () const
 
virtual bool nlConverged (const unsigned int nl_sys_num)
 
virtual bool converged (const unsigned int sys_num)
 
bool defaultGhosting ()
 
virtual TagID addVectorTag (const TagName &tag_name, const Moose::VectorTagType type=Moose::VECTOR_TAG_RESIDUAL)
 
void addNotZeroedVectorTag (const TagID tag)
 
bool vectorTagNotZeroed (const TagID tag) const
 
virtual const VectorTaggetVectorTag (const TagID tag_id) const
 
std::vector< VectorTaggetVectorTags (const std::set< TagID > &tag_ids) const
 
virtual const std::vector< VectorTag > & getVectorTags (const Moose::VectorTagType type=Moose::VECTOR_TAG_ANY) const
 
virtual TagID getVectorTagID (const TagName &tag_name) const
 
virtual TagName vectorTagName (const TagID tag) const
 
virtual bool vectorTagExists (const TagID tag_id) const
 
virtual bool vectorTagExists (const TagName &tag_name) const
 
virtual unsigned int numVectorTags (const Moose::VectorTagType type=Moose::VECTOR_TAG_ANY) const
 
virtual Moose::VectorTagType vectorTagType (const TagID tag_id) const
 
virtual TagID addMatrixTag (TagName tag_name)
 
virtual TagID getMatrixTagID (const TagName &tag_name) const
 
virtual TagName matrixTagName (TagID tag)
 
virtual bool matrixTagExists (const TagName &tag_name) const
 
virtual bool matrixTagExists (TagID tag_id) const
 
virtual unsigned int numMatrixTags () const
 
virtual std::map< TagName, TagID > & getMatrixTags ()
 
virtual bool hasLinearVariable (const std::string &var_name) const
 
virtual bool hasAuxiliaryVariable (const std::string &var_name) const
 
virtual const std::set< MooseVariableFieldBase *> & getActiveElementalMooseVariables (const THREAD_ID tid) const
 
virtual bool hasActiveElementalMooseVariables (const THREAD_ID tid) const
 
Moose::CoordinateSystemType getCoordSystem (SubdomainID sid) const
 
unsigned int getAxisymmetricRadialCoord () const
 
virtual DiracKernelInfodiracKernelInfo ()
 
void reinitNeighborLowerDElem (const Elem *elem, const THREAD_ID tid=0)
 
void reinitMortarElem (const Elem *elem, const THREAD_ID tid=0)
 
virtual void storeSubdomainMatPropName (SubdomainID block_id, const std::string &name)
 
virtual void storeBoundaryMatPropName (BoundaryID boundary_id, const std::string &name)
 
virtual void storeSubdomainZeroMatProp (SubdomainID block_id, const MaterialPropertyName &name)
 
virtual void storeBoundaryZeroMatProp (BoundaryID boundary_id, const MaterialPropertyName &name)
 
virtual void storeSubdomainDelayedCheckMatProp (const std::string &requestor, SubdomainID block_id, const std::string &name)
 
virtual void storeBoundaryDelayedCheckMatProp (const std::string &requestor, BoundaryID boundary_id, const std::string &name)
 
virtual void checkBlockMatProps ()
 
virtual void checkBoundaryMatProps ()
 
virtual void markMatPropRequested (const std::string &)
 
virtual bool isMatPropRequested (const std::string &prop_name) const
 
void addConsumedPropertyName (const MooseObjectName &obj_name, const std::string &prop_name)
 
const std::map< MooseObjectName, std::set< std::string > > & getConsumedPropertyMap () const
 
virtual std::set< SubdomainIDgetMaterialPropertyBlocks (const std::string &prop_name)
 
virtual std::vector< SubdomainName > getMaterialPropertyBlockNames (const std::string &prop_name)
 
virtual bool hasBlockMaterialProperty (SubdomainID block_id, const std::string &prop_name)
 
virtual std::set< BoundaryIDgetMaterialPropertyBoundaryIDs (const std::string &prop_name)
 
virtual std::vector< BoundaryName > getMaterialPropertyBoundaryNames (const std::string &prop_name)
 
virtual bool hasBoundaryMaterialProperty (BoundaryID boundary_id, const std::string &prop_name)
 
virtual std::set< dof_id_type > & ghostedElems ()
 
const bool & currentlyComputingJacobian () const
 
void setCurrentlyComputingJacobian (const bool currently_computing_jacobian)
 
const bool & currentlyComputingResidualAndJacobian () const
 
void setCurrentlyComputingResidualAndJacobian (bool currently_computing_residual_and_jacobian)
 
virtual bool safeAccessTaggedMatrices () const
 
virtual bool safeAccessTaggedVectors () const
 
const std::set< TagID > & getActiveScalarVariableCoupleableVectorTags (const THREAD_ID tid) const
 
const std::set< TagID > & getActiveScalarVariableCoupleableMatrixTags (const THREAD_ID tid) const
 
const std::set< TagID > & getActiveFEVariableCoupleableVectorTags (const THREAD_ID tid) const
 
const std::set< TagID > & getActiveFEVariableCoupleableMatrixTags (const THREAD_ID tid) const
 
void addAlgebraicGhostingFunctor (libMesh::GhostingFunctor &algebraic_gf, bool to_mesh=true)
 
void addCouplingGhostingFunctor (libMesh::GhostingFunctor &coupling_gf, bool to_mesh=true)
 
void removeAlgebraicGhostingFunctor (libMesh::GhostingFunctor &algebraic_gf)
 
void removeCouplingGhostingFunctor (libMesh::GhostingFunctor &coupling_gf)
 
void hasScalingVector (const unsigned int nl_sys_num)
 
void clearAllDofIndices ()
 
const Moose::Functor< T > & getFunctor (const std::string &name, const THREAD_ID tid, const std::string &requestor_name, bool requestor_is_ad)
 
bool hasFunctor (const std::string &name, const THREAD_ID tid) const
 
bool hasFunctorWithType (const std::string &name, const THREAD_ID tid) const
 
void addFunctor (const std::string &name, const Moose::FunctorBase< T > &functor, const THREAD_ID tid)
 
const Moose::FunctorBase< T > & addPiecewiseByBlockLambdaFunctor (const std::string &name, PolymorphicLambda my_lammy, const std::set< ExecFlagType > &clearance_schedule, const MooseMesh &mesh, const std::set< SubdomainID > &block_ids, const THREAD_ID tid)
 
void setFunctorOutput (bool set_output)
 
void registerUnfilledFunctorRequest (T *functor_interface, const std::string &functor_name, const THREAD_ID tid)
 
void reinitFVFace (const THREAD_ID tid, const FaceInfo &fi)
 
void preparePRefinement ()
 
bool doingPRefinement () const
 
bool havePRefinement () const
 
MooseVariableFEBasegetVariableHelper (const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type, Moose::VarFieldType expected_var_field_type, const std::vector< T > &systems, const SystemBase &aux) const
 
void _setCLIOption ()
 
virtual void terminateSolve ()
 
virtual bool isSolveTerminationRequested () const
 
const ConsoleStreamconsole () const
 
virtual bool enabled () const
 
std::shared_ptr< MooseObjectgetSharedPtr ()
 
std::shared_ptr< const MooseObjectgetSharedPtr () const
 
MooseAppgetMooseApp () const
 
const std::string & type () const
 
virtual const std::string & name () const
 
std::string typeAndName () const
 
std::string errorPrefix (const std::string &error_type) const
 
void callMooseError (std::string msg, const bool with_prefix) const
 
MooseObjectParameterName uniqueParameterName (const std::string &parameter_name) const
 
const InputParametersparameters () const
 
MooseObjectName uniqueName () const
 
const T & getParam (const std::string &name) const
 
std::vector< std::pair< T1, T2 > > getParam (const std::string &param1, const std::string &param2) const
 
const T * queryParam (const std::string &name) const
 
const T & getRenamedParam (const std::string &old_name, const std::string &new_name) const
 
getCheckedPointerParam (const std::string &name, const std::string &error_string="") const
 
bool isParamValid (const std::string &name) const
 
bool isParamSetByUser (const std::string &nm) const
 
void paramError (const std::string &param, Args... args) const
 
void paramWarning (const std::string &param, Args... args) const
 
void paramInfo (const std::string &param, Args... args) const
 
void connectControllableParams (const std::string &parameter, const std::string &object_type, const std::string &object_name, const std::string &object_parameter) const
 
void mooseError (Args &&... args) const
 
void mooseErrorNonPrefixed (Args &&... args) const
 
void mooseDocumentedError (const std::string &repo_name, const unsigned int issue_num, Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseInfo (Args &&... args) const
 
std::string getDataFileName (const std::string &param) const
 
std::string getDataFileNameByName (const std::string &relative_path) const
 
std::string getDataFilePath (const std::string &relative_path) const
 
PerfGraphperfGraph ()
 
const libMesh::ConstElemRangegetEvaluableElementRange ()
 
const libMesh::ConstElemRangegetEvaluableElementRange ()
 
const libMesh::ConstElemRangegetNonlinearEvaluableElementRange ()
 
const libMesh::ConstElemRangegetNonlinearEvaluableElementRange ()
 
const libMesh::ConstElemRangegetCurrentAlgebraicElementRange ()
 
const libMesh::ConstElemRangegetCurrentAlgebraicElementRange ()
 
const libMesh::ConstNodeRangegetCurrentAlgebraicNodeRange ()
 
const libMesh::ConstNodeRangegetCurrentAlgebraicNodeRange ()
 
const ConstBndNodeRangegetCurrentAlgebraicBndNodeRange ()
 
const ConstBndNodeRangegetCurrentAlgebraicBndNodeRange ()
 
void setCurrentAlgebraicElementRange (libMesh::ConstElemRange *range)
 
void setCurrentAlgebraicElementRange (libMesh::ConstElemRange *range)
 
void setCurrentAlgebraicNodeRange (libMesh::ConstNodeRange *range)
 
void setCurrentAlgebraicNodeRange (libMesh::ConstNodeRange *range)
 
void setCurrentAlgebraicBndNodeRange (ConstBndNodeRange *range)
 
void setCurrentAlgebraicBndNodeRange (ConstBndNodeRange *range)
 
void allowOutput (bool state)
 
void allowOutput (bool state)
 
void allowOutput (bool state)
 
void allowOutput (bool state)
 
bool hasMultiApps () const
 
bool hasMultiApps (ExecFlagType type) const
 
bool hasMultiApps () const
 
bool hasMultiApps (ExecFlagType type) const
 
bool hasMultiApp (const std::string &name) const
 
bool hasMultiApp (const std::string &name) const
 
const AutomaticMortarGenerationgetMortarInterface (const std::pair< BoundaryID, BoundaryID > &primary_secondary_boundary_pair, const std::pair< SubdomainID, SubdomainID > &primary_secondary_subdomain_pair, bool on_displaced) const
 
AutomaticMortarGenerationgetMortarInterface (const std::pair< BoundaryID, BoundaryID > &primary_secondary_boundary_pair, const std::pair< SubdomainID, SubdomainID > &primary_secondary_subdomain_pair, bool on_displaced)
 
const AutomaticMortarGenerationgetMortarInterface (const std::pair< BoundaryID, BoundaryID > &primary_secondary_boundary_pair, const std::pair< SubdomainID, SubdomainID > &primary_secondary_subdomain_pair, bool on_displaced) const
 
AutomaticMortarGenerationgetMortarInterface (const std::pair< BoundaryID, BoundaryID > &primary_secondary_boundary_pair, const std::pair< SubdomainID, SubdomainID > &primary_secondary_subdomain_pair, bool on_displaced)
 
const MaterialPropertyStoragegetMaterialPropertyStorage ()
 
const MaterialPropertyStoragegetMaterialPropertyStorage ()
 
const MaterialPropertyStoragegetBndMaterialPropertyStorage ()
 
const MaterialPropertyStoragegetBndMaterialPropertyStorage ()
 
const MaterialPropertyStoragegetNeighborMaterialPropertyStorage ()
 
const MaterialPropertyStoragegetNeighborMaterialPropertyStorage ()
 
const MooseObjectWarehouse< Indicator > & getIndicatorWarehouse ()
 
const MooseObjectWarehouse< Indicator > & getIndicatorWarehouse ()
 
const MooseObjectWarehouse< InternalSideIndicatorBase > & getInternalSideIndicatorWarehouse ()
 
const MooseObjectWarehouse< InternalSideIndicatorBase > & getInternalSideIndicatorWarehouse ()
 
const MooseObjectWarehouse< Marker > & getMarkerWarehouse ()
 
const MooseObjectWarehouse< Marker > & getMarkerWarehouse ()
 
bool needBoundaryMaterialOnSide (BoundaryID bnd_id, const THREAD_ID tid)
 
bool needBoundaryMaterialOnSide (BoundaryID bnd_id, const THREAD_ID tid)
 
bool needInterfaceMaterialOnSide (BoundaryID bnd_id, const THREAD_ID tid)
 
bool needInterfaceMaterialOnSide (BoundaryID bnd_id, const THREAD_ID tid)
 
bool needSubdomainMaterialOnSide (SubdomainID subdomain_id, const THREAD_ID tid)
 
bool needSubdomainMaterialOnSide (SubdomainID subdomain_id, const THREAD_ID tid)
 
const ExecFlagTypegetCurrentExecuteOnFlag () const
 
const ExecFlagTypegetCurrentExecuteOnFlag () const
 
void setCurrentExecuteOnFlag (const ExecFlagType &)
 
void setCurrentExecuteOnFlag (const ExecFlagType &)
 
const Parallel::Communicator & comm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static InputParameters validParams ()
 
static void selectVectorTagsFromSystem (const SystemBase &system, const std::vector< VectorTag > &input_vector_tags, std::set< TagID > &selected_tags)
 
static void selectMatrixTagsFromSystem (const SystemBase &system, const std::map< TagName, TagID > &input_matrix_tags, std::set< TagID > &selected_tags)
 
static void objectSetupHelper (const std::vector< T * > &objects, const ExecFlagType &exec_flag)
 
static void objectSetupHelper (const std::vector< T * > &objects, const ExecFlagType &exec_flag)
 
static void objectExecuteHelper (const std::vector< T * > &objects)
 
static void objectExecuteHelper (const std::vector< T * > &objects)
 

Public Attributes

std::map< std::string, std::vector< dof_id_type > > _var_dof_map
 
const ConsoleStream _console
 
std::vector< Real_real_zero
 
std::vector< VariableValue_scalar_zero
 
std::vector< VariableValue_zero
 
std::vector< VariablePhiValue_phi_zero
 
std::vector< MooseArray< ADReal > > _ad_zero
 
std::vector< VariableGradient_grad_zero
 
std::vector< MooseArray< ADRealVectorValue > > _ad_grad_zero
 
std::vector< VariablePhiGradient_grad_phi_zero
 
std::vector< VariableSecond_second_zero
 
std::vector< MooseArray< ADRealTensorValue > > _ad_second_zero
 
std::vector< VariablePhiSecond_second_phi_zero
 
std::vector< Point > _point_zero
 
std::vector< VectorVariableValue_vector_zero
 
std::vector< VectorVariableCurl_vector_curl_zero
 

Protected Member Functions

virtual Real computeFrictionFactor (Real Re) override
 Returns friction factor. More...
 
virtual void computeWijFromSolve (int iblock)
 Computes diversion crossflow per gap for block with index iblock Block is a partition of the whole domain. More...
 
virtual void computeSumWij (int iblock)
 Computes net diversion crossflow per channel for block iblock. More...
 
virtual void computeMdot (int iblock)
 Computes mass flow per channel for block iblock. More...
 
virtual void computeWijPrime (int iblock)
 Computes turbulent crossflow per gap for block iblock. More...
 
virtual void computeDP (int iblock)
 Computes Pressure Drop per channel for block iblock. More...
 
virtual void computeP (int iblock)
 Computes Pressure per channel for block iblock. More...
 
virtual void computeh (int iblock)
 Computes Enthalpy per channel for block iblock. More...
 
virtual void computeT (int iblock)
 Computes Temperature per channel for block iblock. More...
 
virtual void computeRho (int iblock)
 Computes Density per channel for block iblock. More...
 
virtual void computeMu (int iblock)
 Computes Viscosity per channel for block iblock. More...
 
virtual void computeWij (int iblock)
 Computes cross fluxes for block iblock. More...
 
virtual Real computeAddedHeat (unsigned int i_ch, unsigned int iz)
 Computes added heat for channel i_ch and cell iz. More...
 
virtual libMesh::DenseVector< RealresidualFunction (int iblock, libMesh::DenseVector< Real > solution)
 Computes Residual per gap for block iblock. More...
 
virtual PetscErrorCode petscSnesSolver (int iblock, const libMesh::DenseVector< Real > &solution, libMesh::DenseVector< Real > &root)
 Computes solution of nonlinear equation using snes. More...
 
virtual PetscErrorCode implicitPetscSolve (int iblock)
 Computes implicit solve using PetSc. More...
 
virtual void initializeSolution ()
 Function to initialize the solution fields. More...
 
virtual PetscScalar computeInterpolationCoefficients (PetscScalar Peclet=0.0)
 Functions that computes the interpolation scheme given the Peclet number. More...
 
virtual PetscScalar computeInterpolatedValue (PetscScalar topValue, PetscScalar botValue, PetscScalar Peclet=0.0)
 
PetscErrorCode cleanUp ()
 
virtual PetscErrorCode createPetscVector (Vec &v, PetscInt n)
 Petsc Functions. More...
 
virtual PetscErrorCode createPetscMatrix (Mat &M, PetscInt n, PetscInt m)
 
template<class T >
PetscErrorCode populateVectorFromDense (Vec &x, const T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
 
template<class T >
PetscErrorCode populateDenseFromVector (const Vec &x, T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
 
template<class T >
PetscErrorCode populateVectorFromHandle (Vec &x, const T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
 
template<class T >
PetscErrorCode populateSolutionChan (const Vec &x, T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
 
template<class T >
PetscErrorCode populateSolutionGap (const Vec &x, T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
 
virtual void meshChanged ()
 
MooseVariableFieldBasegetVariableHelper (const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type, Moose::VarFieldType expected_var_field_type, const std::vector< T > &nls, const SystemBase &aux) const
 
void createTagVectors ()
 
void createTagSolutions ()
 
virtual void meshDisplaced ()
 
void computeSystems (const ExecFlagType &type)
 
bool duplicateVariableCheck (const std::string &var_name, const libMesh::FEType &type, bool is_aux, const std::set< SubdomainID > *const active_subdomains)
 
void computeUserObjectsInternal (const ExecFlagType &type, const Moose::AuxGroup &group, TheWarehouse::Query &query)
 
void checkDisplacementOrders ()
 
void checkUserObjects ()
 
void checkDependMaterialsHelper (const std::map< SubdomainID, std::vector< std::shared_ptr< MaterialBase >>> &materials_map)
 
void checkCoordinateSystems ()
 
void reinitBecauseOfGhostingOrNewGeomObjects (bool mortar_changed=false)
 
void addObjectParamsHelper (InputParameters &params, const std::string &object_name, const std::string &var_param_name="variable")
 
bool verifyVectorTags () const
 
void markFamilyPRefinement (const InputParameters &params)
 
PerfID registerTimedSection (const std::string &section_name, const unsigned int level) const
 
PerfID registerTimedSection (const std::string &section_name, const unsigned int level, const std::string &live_message, const bool print_dots=true) const
 
std::string timedSectionName (const std::string &section_name) const
 
T & declareRestartableData (const std::string &data_name, Args &&... args)
 
ManagedValue< T > declareManagedRestartableDataWithContext (const std::string &data_name, void *context, Args &&... args)
 
const T & getRestartableData (const std::string &data_name) const
 
T & declareRestartableDataWithContext (const std::string &data_name, void *context, Args &&... args)
 
T & declareRecoverableData (const std::string &data_name, Args &&... args)
 
T & declareRestartableDataWithObjectName (const std::string &data_name, const std::string &object_name, Args &&... args)
 
T & declareRestartableDataWithObjectNameWithContext (const std::string &data_name, const std::string &object_name, void *context, Args &&... args)
 
std::string restartableName (const std::string &data_name) const
 

Protected Attributes

const InterWrapperMesh_subchannel_mesh
 
unsigned int _n_blocks
 number of axial blocks More...
 
libMesh::DenseMatrix< Real > & _Wij
 
libMesh::DenseMatrix< Real_Wij_old
 
libMesh::DenseMatrix< Real_WijPrime
 
libMesh::DenseMatrix< Real_Wij_residual_matrix
 
const Real _g_grav
 
const Real_kij
 
const unsigned int _n_cells
 
const unsigned int _n_gaps
 
const unsigned int _n_assemblies
 
const unsigned int _n_channels
 
const unsigned int _block_size
 
std::vector< Real_z_grid
 axial location of nodes More...
 
Real _one
 
Real _TR
 Flag that activates or deactivates the transient parts of the equations we solve by multiplication. More...
 
const bool _compute_density
 Flag that activates or deactivates the calculation of density. More...
 
const bool _compute_viscosity
 Flag that activates or deactivates the calculation of viscosity. More...
 
const bool _compute_power
 Flag that informs if we need to solve the Enthalpy/Temperature equations or not. More...
 
const bool _pin_mesh_exist
 Flag that informs if there is a pin mesh or not. More...
 
bool _converged
 Variable that informs whether we exited external solve with a converged solution or not. More...
 
const Real_dt
 Time step. More...
 
const Real_P_out
 Outlet Pressure. More...
 
const Real_beta
 Thermal diffusion coefficient used in turbulent crossflow. More...
 
const Real_CT
 Turbulent modeling parameter used in axial momentum equation. More...
 
const Real_P_tol
 Convergence tolerance for the pressure loop in external solve. More...
 
const Real_T_tol
 Convergence tolerance for the temperature loop in external solve. More...
 
const int_T_maxit
 Maximum iterations for the inner temperature loop. More...
 
const PetscReal & _rtol
 The relative convergence tolerance, (relative decrease) for the ksp linear solver. More...
 
const PetscReal & _atol
 The absolute convergence tolerance for the ksp linear solver. More...
 
const PetscReal & _dtol
 The divergence tolerance for the ksp linear solver. More...
 
const PetscInt & _maxit
 The maximum number of iterations to use for the ksp linear solver. More...
 
const MooseEnum _interpolation_scheme
 The interpolation method used in constructing the systems. More...
 
const bool _implicit_bool
 Flag to define the usage of a implicit or explicit solution. More...
 
const bool _staggered_pressure_bool
 Flag to define the usage of staggered or collocated pressure. More...
 
const bool _segregated_bool
 Segregated solve. More...
 
const bool _monolithic_thermal_bool
 Thermal monolithic bool. More...
 
const SinglePhaseFluidProperties_fp
 Solutions handles and link to TH tables properties. More...
 
std::unique_ptr< SolutionHandle_mdot_soln
 
std::unique_ptr< SolutionHandle_SumWij_soln
 
std::unique_ptr< SolutionHandle_P_soln
 
std::unique_ptr< SolutionHandle_DP_soln
 
std::unique_ptr< SolutionHandle_h_soln
 
std::unique_ptr< SolutionHandle_T_soln
 
std::unique_ptr< SolutionHandle_Tpin_soln
 
std::unique_ptr< SolutionHandle_rho_soln
 
std::unique_ptr< SolutionHandle_mu_soln
 
std::unique_ptr< SolutionHandle_S_flow_soln
 
std::unique_ptr< SolutionHandle_w_perim_soln
 
std::unique_ptr< SolutionHandle_q_prime_soln
 
Mat _mc_sumWij_mat
 Mass conservation Mass conservation - sum of cross fluxes. More...
 
Vec _Wij_vec
 
Vec _prod
 
Vec _prodp
 
Mat _mc_axial_convection_mat
 Mass conservation - axial convection. More...
 
Vec _mc_axial_convection_rhs
 
Mat _amc_turbulent_cross_flows_mat
 Mass conservation - density time derivative No implicit matrix. More...
 
Vec _amc_turbulent_cross_flows_rhs
 
Mat _amc_time_derivative_mat
 Axial momentum conservation - time derivative. More...
 
Vec _amc_time_derivative_rhs
 
Mat _amc_advective_derivative_mat
 Axial momentum conservation - advective (Eulerian) derivative. More...
 
Vec _amc_advective_derivative_rhs
 
Mat _amc_cross_derivative_mat
 Axial momentum conservation - cross flux derivative. More...
 
Vec _amc_cross_derivative_rhs
 
Mat _amc_friction_force_mat
 Axial momentum conservation - friction force. More...
 
Vec _amc_friction_force_rhs
 
Vec _amc_gravity_rhs
 Axial momentum conservation - buoyancy force No implicit matrix. More...
 
Mat _amc_pressure_force_mat
 Axial momentum conservation - pressure force. More...
 
Vec _amc_pressure_force_rhs
 
Mat _amc_sys_mdot_mat
 Axial momentum system matrix. More...
 
Vec _amc_sys_mdot_rhs
 
Mat _cmc_time_derivative_mat
 Cross momentum Cross momentum conservation - time derivative. More...
 
Vec _cmc_time_derivative_rhs
 
Mat _cmc_advective_derivative_mat
 Cross momentum conservation - advective (Eulerian) derivative. More...
 
Vec _cmc_advective_derivative_rhs
 
Mat _cmc_friction_force_mat
 Cross momentum conservation - friction force. More...
 
Vec _cmc_friction_force_rhs
 
Mat _cmc_pressure_force_mat
 Cross momentum conservation - pressure force. More...
 
Vec _cmc_pressure_force_rhs
 
Mat _cmc_sys_Wij_mat
 Lateral momentum system matrix. More...
 
Vec _cmc_sys_Wij_rhs
 
Vec _cmc_Wij_channel_dummy
 
Mat _hc_time_derivative_mat
 Enthalpy Enthalpy conservation - time derivative. More...
 
Vec _hc_time_derivative_rhs
 
Mat _hc_advective_derivative_mat
 Enthalpy conservation - advective (Eulerian) derivative;. More...
 
Vec _hc_advective_derivative_rhs
 
Mat _hc_cross_derivative_mat
 Enthalpy conservation - cross flux derivative. More...
 
Vec _hc_cross_derivative_rhs
 
Vec _hc_added_heat_rhs
 Enthalpy conservation - source and sink. More...
 
Mat _hc_sys_h_mat
 System matrices. More...
 
Vec _hc_sys_h_rhs
 
PetscInt _global_counter = 0
 No implicit matrix. More...
 
PetscScalar _added_K = 0.0
 Added resistances for monolithic convergence. More...
 
PetscScalar _added_K_old = 1000.0
 
PetscScalar _max_sumWij
 
PetscScalar _max_sumWij_new
 
PetscScalar _correction_factor = 1.0
 
MooseMesh_mesh
 
bool _initialized
 
std::optional< std::vector< ConvergenceName > > _nonlinear_convergence_names
 
std::optional< std::vector< ConvergenceName > > _linear_convergence_names
 
std::optional< ConvergenceName > _multiapp_fixed_point_convergence_name
 
std::set< TagID_fe_vector_tags
 
std::set< TagID_fe_matrix_tags
 
std::set< TagID_linear_vector_tags
 
std::set< TagID_linear_matrix_tags
 
const bool & _solve
 
bool _transient
 
Real_time
 
Real_time_old
 
int_t_step
 
Real_dt_old
 
bool _need_to_add_default_nonlinear_convergence
 
bool _need_to_add_default_multiapp_fixed_point_convergence
 
const std::vector< LinearSystemName > _linear_sys_names
 
const std::size_t _num_linear_sys
 
std::vector< std::shared_ptr< LinearSystem > > _linear_systems
 
std::map< LinearSystemName, unsigned int_linear_sys_name_to_num
 
LinearSystem_current_linear_sys
 
const bool _using_default_nl
 
const std::vector< NonlinearSystemName > _nl_sys_names
 
const std::size_t _num_nl_sys
 
std::vector< std::shared_ptr< NonlinearSystemBase > > _nl
 
std::map< NonlinearSystemName, unsigned int_nl_sys_name_to_num
 
NonlinearSystemBase_current_nl_sys
 
SolverSystem_current_solver_sys
 
std::vector< std::shared_ptr< SolverSystem > > _solver_systems
 
std::map< SolverVariableName, unsigned int_solver_var_to_sys_num
 
std::map< SolverSystemName, unsigned int_solver_sys_name_to_num
 
std::vector< SolverSystemName > _solver_sys_names
 
std::shared_ptr< AuxiliarySystem_aux
 
Moose::CouplingType _coupling
 
std::vector< std::unique_ptr< libMesh::CouplingMatrix > > _cm
 
std::map< std::string, unsigned int_subspace_dim
 
std::vector< std::vector< std::unique_ptr< Assembly > > > _assembly
 
MooseObjectWarehouse< MeshDivision_mesh_divisions
 
MooseObjectWarehouse< Function_functions
 
MooseObjectWarehouse< Convergence_convergences
 
MooseObjectWarehouse< KernelBase_nonlocal_kernels
 
MooseObjectWarehouse< IntegratedBCBase_nonlocal_integrated_bcs
 
MaterialPropertyRegistry _material_prop_registry
 
MaterialPropertyStorage_material_props
 
MaterialPropertyStorage_bnd_material_props
 
MaterialPropertyStorage_neighbor_material_props
 
MooseObjectWarehouse< Marker_markers
 
ReporterData _reporter_data
 
ExecuteMooseObjectWarehouse< UserObject_all_user_objects
 
ExecuteMooseObjectWarehouse< MultiApp_multi_apps
 
ExecuteMooseObjectWarehouse< TransientMultiApp_transient_multi_apps
 
ExecuteMooseObjectWarehouse< Transfer_transfers
 
ExecuteMooseObjectWarehouse< Transfer_to_multi_app_transfers
 
ExecuteMooseObjectWarehouse< Transfer_from_multi_app_transfers
 
ExecuteMooseObjectWarehouse< Transfer_between_multi_app_transfers
 
std::map< std::string, std::unique_ptr< RandomData > > _random_data_objects
 
std::vector< std::unordered_map< SubdomainID, bool > > _block_mat_side_cache
 
std::vector< std::unordered_map< BoundaryID, bool > > _bnd_mat_side_cache
 
std::vector< std::unordered_map< BoundaryID, bool > > _interface_mat_side_cache
 
std::vector< MeshChangedInterface *> _notify_when_mesh_changes
 
std::vector< MeshDisplacedInterface *> _notify_when_mesh_displaces
 
Adaptivity _adaptivity
 
unsigned int _cycles_completed
 
std::shared_ptr< XFEMInterface_xfem
 
MooseMesh_displaced_mesh
 
std::shared_ptr< DisplacedProblem_displaced_problem
 
GeometricSearchData _geometric_search_data
 
MortarData _mortar_data
 
bool _reinit_displaced_elem
 
bool _reinit_displaced_face
 
bool _reinit_displaced_neighbor
 
bool _input_file_saved
 
bool _has_dampers
 
bool _has_constraints
 
bool _snesmf_reuse_base
 
bool _skip_exception_check
 
bool _snesmf_reuse_base_set_by_user
 
bool _has_initialized_stateful
 
bool _const_jacobian
 
bool _has_jacobian
 
bool _needs_old_newton_iter
 
bool _previous_nl_solution_required
 
bool _has_nonlocal_coupling
 
bool _calculate_jacobian_in_uo
 
std::vector< std::vector< const MooseVariableFEBase *> > _uo_jacobian_moose_vars
 
std::vector< unsigned char > _has_active_material_properties
 
std::vector< SolverParams_solver_params
 
CoverageCheckMode _kernel_coverage_check
 
std::vector< SubdomainName > _kernel_coverage_blocks
 
const bool _boundary_restricted_node_integrity_check
 
const bool _boundary_restricted_elem_integrity_check
 
CoverageCheckMode _material_coverage_check
 
std::vector< SubdomainName > _material_coverage_blocks
 
bool _fv_bcs_integrity_check
 
const bool _material_dependency_check
 
const bool _uo_aux_state_check
 
unsigned int _max_qps
 
libMesh::Order _max_scalar_order
 
bool _has_time_integrator
 
bool _has_exception
 
bool _parallel_barrier_messaging
 
MooseEnum _verbose_setup
 
bool _verbose_multiapps
 
bool _verbose_restore
 
std::string _exception_message
 
ExecFlagType _current_execute_on_flag
 
ExecuteMooseObjectWarehouse< Control_control_warehouse
 
Moose::PetscSupport::PetscOptions _petsc_options
 
PetscOptions _petsc_option_data_base
 
bool _is_petsc_options_inserted
 
std::shared_ptr< LineSearch_line_search
 
std::unique_ptr< libMesh::ConstElemRange_evaluable_local_elem_range
 
std::unique_ptr< libMesh::ConstElemRange_nl_evaluable_local_elem_range
 
std::unique_ptr< libMesh::ConstElemRange_aux_evaluable_local_elem_range
 
std::unique_ptr< libMesh::ConstElemRange_current_algebraic_elem_range
 
std::unique_ptr< libMesh::ConstNodeRange_current_algebraic_node_range
 
std::unique_ptr< ConstBndNodeRange_current_algebraic_bnd_node_range
 
bool _using_ad_mat_props
 
unsigned short _current_ic_state
 
const bool _use_hash_table_matrix_assembly
 
std::map< TagName, TagID_matrix_tag_name_to_tag_id
 
std::map< TagID, TagName > _matrix_tag_id_to_tag_name
 
Factory_factory
 
DiracKernelInfo _dirac_kernel_info
 
std::map< SubdomainID, std::set< std::string > > _map_block_material_props
 
std::map< BoundaryID, std::set< std::string > > _map_boundary_material_props
 
std::map< SubdomainID, std::set< MaterialPropertyName > > _zero_block_material_props
 
std::map< BoundaryID, std::set< MaterialPropertyName > > _zero_boundary_material_props
 
std::set< std::string > _material_property_requested
 
std::vector< std::set< MooseVariableFieldBase *> > _active_elemental_moose_variables
 
std::vector< unsigned int_has_active_elemental_moose_variables
 
std::vector< std::set< TagID > > _active_fe_var_coupleable_matrix_tags
 
std::vector< std::set< TagID > > _active_fe_var_coupleable_vector_tags
 
std::vector< std::set< TagID > > _active_sc_var_coupleable_matrix_tags
 
std::vector< std::set< TagID > > _active_sc_var_coupleable_vector_tags
 
bool _default_ghosting
 
std::set< dof_id_type_ghosted_elems
 
bool _currently_computing_jacobian
 
bool _currently_computing_residual_and_jacobian
 
bool _computing_nonlinear_residual
 
bool _currently_computing_residual
 
bool _safe_access_tagged_matrices
 
bool _safe_access_tagged_vectors
 
bool _have_ad_objects
 
std::unordered_set< TagID_not_zeroed_tagged_vectors
 
bool _cli_option_found
 
bool _color_output
 
bool _termination_requested
 
const bool & _enabled
 
MooseApp_app
 
const std::string _type
 
const std::string _name
 
const InputParameters_pars
 
ActionFactory_action_factory
 
MooseApp_pg_moose_app
 
const std::string _prefix
 
MooseApp_restartable_app
 
const std::string _restartable_system_name
 
const THREAD_ID _restartable_tid
 
const bool _restartable_read_only
 
InitialConditionWarehouse _ics
 
FVInitialConditionWarehouse _fv_ics
 
ScalarInitialConditionWarehouse _scalar_ics
 
MaterialWarehouse _materials
 
MaterialWarehouse _interface_materials
 
MaterialWarehouse _discrete_materials
 
MaterialWarehouse _all_materials
 
MooseObjectWarehouse< Indicator_indicators
 
MooseObjectWarehouse< InternalSideIndicatorBase_internal_side_indicators
 
std::map< SubdomainID, std::multimap< std::string, std::string > > _map_block_material_props_check
 
std::map< BoundaryID, std::multimap< std::string, std::string > > _map_boundary_material_props_check
 
const Parallel::Communicator & _communicator
 

Detailed Description

Quadrilateral interwrapper solver.

Definition at line 19 of file QuadInterWrapper1PhaseProblem.h.

Constructor & Destructor Documentation

◆ QuadInterWrapper1PhaseProblem()

QuadInterWrapper1PhaseProblem::QuadInterWrapper1PhaseProblem ( const InputParameters params)

Definition at line 23 of file QuadInterWrapper1PhaseProblem.C.

25 {
26 }
InterWrapper1PhaseProblem(const InputParameters &params)

Member Function Documentation

◆ cleanUp()

PetscErrorCode InterWrapper1PhaseProblem::cleanUp ( )
protectedinherited

Definition at line 253 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::~InterWrapper1PhaseProblem().

254 {
256  // We need to clean up the petsc matrices/vectors
257  // Mass conservation components
258  LibmeshPetscCall(MatDestroy(&_mc_sumWij_mat));
259  LibmeshPetscCall(VecDestroy(&_Wij_vec));
260  LibmeshPetscCall(VecDestroy(&_prod));
261  LibmeshPetscCall(VecDestroy(&_prodp));
262  LibmeshPetscCall(MatDestroy(&_mc_axial_convection_mat));
263  LibmeshPetscCall(VecDestroy(&_mc_axial_convection_rhs));
264 
265  // Axial momentum conservation components
266  LibmeshPetscCall(MatDestroy(&_amc_turbulent_cross_flows_mat));
267  LibmeshPetscCall(VecDestroy(&_amc_turbulent_cross_flows_rhs));
268  LibmeshPetscCall(MatDestroy(&_amc_time_derivative_mat));
269  LibmeshPetscCall(VecDestroy(&_amc_time_derivative_rhs));
270  LibmeshPetscCall(MatDestroy(&_amc_advective_derivative_mat));
271  LibmeshPetscCall(VecDestroy(&_amc_advective_derivative_rhs));
272  LibmeshPetscCall(MatDestroy(&_amc_cross_derivative_mat));
273  LibmeshPetscCall(VecDestroy(&_amc_cross_derivative_rhs));
274  LibmeshPetscCall(MatDestroy(&_amc_friction_force_mat));
275  LibmeshPetscCall(VecDestroy(&_amc_friction_force_rhs));
276  LibmeshPetscCall(VecDestroy(&_amc_gravity_rhs));
277  LibmeshPetscCall(MatDestroy(&_amc_pressure_force_mat));
278  LibmeshPetscCall(VecDestroy(&_amc_pressure_force_rhs));
279  LibmeshPetscCall(MatDestroy(&_amc_sys_mdot_mat));
280  LibmeshPetscCall(VecDestroy(&_amc_sys_mdot_rhs));
281 
282  // Lateral momentum conservation components
283  LibmeshPetscCall(MatDestroy(&_cmc_time_derivative_mat));
284  LibmeshPetscCall(VecDestroy(&_cmc_time_derivative_rhs));
285  LibmeshPetscCall(MatDestroy(&_cmc_advective_derivative_mat));
286  LibmeshPetscCall(VecDestroy(&_cmc_advective_derivative_rhs));
287  LibmeshPetscCall(MatDestroy(&_cmc_friction_force_mat));
288  LibmeshPetscCall(VecDestroy(&_cmc_friction_force_rhs));
289  LibmeshPetscCall(MatDestroy(&_cmc_pressure_force_mat));
290  LibmeshPetscCall(VecDestroy(&_cmc_pressure_force_rhs));
291  LibmeshPetscCall(MatDestroy(&_cmc_sys_Wij_mat));
292  LibmeshPetscCall(VecDestroy(&_cmc_sys_Wij_rhs));
293  LibmeshPetscCall(VecDestroy(&_cmc_Wij_channel_dummy));
294 
295  // Energy conservation components
296  LibmeshPetscCall(MatDestroy(&_hc_time_derivative_mat));
297  LibmeshPetscCall(VecDestroy(&_hc_time_derivative_rhs));
298  LibmeshPetscCall(MatDestroy(&_hc_advective_derivative_mat));
299  LibmeshPetscCall(VecDestroy(&_hc_advective_derivative_rhs));
300  LibmeshPetscCall(MatDestroy(&_hc_cross_derivative_mat));
301  LibmeshPetscCall(VecDestroy(&_hc_cross_derivative_rhs));
302  LibmeshPetscCall(VecDestroy(&_hc_added_heat_rhs));
303  LibmeshPetscCall(MatDestroy(&_hc_sys_h_mat));
304  LibmeshPetscCall(VecDestroy(&_hc_sys_h_rhs));
305 
306  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
307 }
Mat _amc_turbulent_cross_flows_mat
Mass conservation - density time derivative No implicit matrix.
Mat _cmc_friction_force_mat
Cross momentum conservation - friction force.
Mat _amc_advective_derivative_mat
Axial momentum conservation - advective (Eulerian) derivative.
Mat _cmc_pressure_force_mat
Cross momentum conservation - pressure force.
Mat _mc_sumWij_mat
Mass conservation Mass conservation - sum of cross fluxes.
PetscFunctionBegin
Mat _mc_axial_convection_mat
Mass conservation - axial convection.
Mat _hc_time_derivative_mat
Enthalpy Enthalpy conservation - time derivative.
Mat _amc_cross_derivative_mat
Axial momentum conservation - cross flux derivative.
Mat _hc_advective_derivative_mat
Enthalpy conservation - advective (Eulerian) derivative;.
Mat _amc_time_derivative_mat
Axial momentum conservation - time derivative.
Mat _cmc_advective_derivative_mat
Cross momentum conservation - advective (Eulerian) derivative.
Vec _hc_added_heat_rhs
Enthalpy conservation - source and sink.
Mat _cmc_sys_Wij_mat
Lateral momentum system matrix.
Mat _amc_pressure_force_mat
Axial momentum conservation - pressure force.
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)
Mat _hc_cross_derivative_mat
Enthalpy conservation - cross flux derivative.
Vec _amc_gravity_rhs
Axial momentum conservation - buoyancy force No implicit matrix.
Mat _cmc_time_derivative_mat
Cross momentum Cross momentum conservation - time derivative.
Mat _amc_sys_mdot_mat
Axial momentum system matrix.
Mat _amc_friction_force_mat
Axial momentum conservation - friction force.

◆ computeAddedHeat()

Real InterWrapper1PhaseProblem::computeAddedHeat ( unsigned int  i_ch,
unsigned int  iz 
)
protectedvirtualinherited

Computes added heat for channel i_ch and cell iz.

Definition at line 1477 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::computeh().

1478 {
1479  auto dz = _z_grid[iz] - _z_grid[iz - 1];
1480  if (_pin_mesh_exist)
1481  {
1482  auto heat_rate_in = 0.0;
1483  auto heat_rate_out = 0.0;
1484  for (auto i_pin : _subchannel_mesh.getChannelPins(i_ch))
1485  {
1486  auto * node_in = _subchannel_mesh.getPinNode(i_pin, iz - 1);
1487  auto * node_out = _subchannel_mesh.getPinNode(i_pin, iz);
1488  heat_rate_out += 0.25 * (*_q_prime_soln)(node_out);
1489  heat_rate_in += 0.25 * (*_q_prime_soln)(node_in);
1490  }
1491  return (heat_rate_in + heat_rate_out) * dz / 2.0;
1492  }
1493  else
1494  {
1495  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1496  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1497  return ((*_q_prime_soln)(node_out) + (*_q_prime_soln)(node_in)) * dz / 2.0;
1498  }
1499 }
virtual const std::vector< unsigned int > & getChannelPins(unsigned int i_chan) const =0
Return a vector of pin indices for a given channel index.
std::vector< Real > _z_grid
axial location of nodes
virtual Node * getPinNode(unsigned int i_pin, unsigned iz) const =0
Get the pin mesh node for a given pin index and elevation index.
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
std::unique_ptr< SolutionHandle > _q_prime_soln
const bool _pin_mesh_exist
Flag that informs if there is a pin mesh or not.

◆ computeDP()

void InterWrapper1PhaseProblem::computeDP ( int  iblock)
protectedvirtualinherited

Computes Pressure Drop per channel for block iblock.

Time derivative term

Advective derivative term

Cross derivative term

Friction term

Gravity force

Assembling system

Reimplemented in TriInterWrapper1PhaseProblem.

Definition at line 793 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::residualFunction().

794 {
795  unsigned int last_node = (iblock + 1) * _block_size;
796  unsigned int first_node = iblock * _block_size + 1;
797 
798  if (!_implicit_bool)
799  {
800  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
801  {
802  auto k_grid = _subchannel_mesh.getKGrid();
803  auto dz = _z_grid[iz] - _z_grid[iz - 1];
804  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
805  {
806  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
807  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
808  auto rho_in = (*_rho_soln)(node_in);
809  auto rho_out = (*_rho_soln)(node_out);
810  auto mu_in = (*_mu_soln)(node_in);
811  auto S = (*_S_flow_soln)(node_in);
812  auto w_perim = (*_w_perim_soln)(node_in);
813  // hydraulic diameter in the i direction
814  auto Dh_i = 4.0 * S / w_perim;
815  auto time_term = _TR * ((*_mdot_soln)(node_out)-_mdot_soln->old(node_out)) * dz / _dt -
816  dz * 2.0 * (*_mdot_soln)(node_out) * (rho_out - _rho_soln->old(node_out)) /
817  rho_in / _dt;
818  auto mass_term1 =
819  std::pow((*_mdot_soln)(node_out), 2.0) * (1.0 / S / rho_out - 1.0 / S / rho_in);
820  auto mass_term2 = -2.0 * (*_mdot_soln)(node_out) * (*_SumWij_soln)(node_out) / S / rho_in;
821  auto crossflow_term = 0.0;
822  auto turbulent_term = 0.0;
823  unsigned int counter = 0;
824  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
825  {
826  auto chans = _subchannel_mesh.getGapChannels(i_gap);
827  unsigned int ii_ch = chans.first;
828  unsigned int jj_ch = chans.second;
829  auto * node_in_i = _subchannel_mesh.getChannelNode(ii_ch, iz - 1);
830  auto * node_in_j = _subchannel_mesh.getChannelNode(jj_ch, iz - 1);
831  auto * node_out_i = _subchannel_mesh.getChannelNode(ii_ch, iz);
832  auto * node_out_j = _subchannel_mesh.getChannelNode(jj_ch, iz);
833  auto rho_i = (*_rho_soln)(node_in_i);
834  auto rho_j = (*_rho_soln)(node_in_j);
835  auto Si = (*_S_flow_soln)(node_in_i);
836  auto Sj = (*_S_flow_soln)(node_in_j);
837  Real u_star = 0.0;
838  // figure out donor axial velocity
839  if (_Wij(i_gap, iz) > 0.0)
840  u_star = (*_mdot_soln)(node_out_i) / Si / rho_i;
841  else
842  u_star = (*_mdot_soln)(node_out_j) / Sj / rho_j;
843 
844  crossflow_term +=
845  _subchannel_mesh.getCrossflowSign(i_ch, counter) * _Wij(i_gap, iz) * u_star;
846 
847  turbulent_term += _WijPrime(i_gap, iz) * (2 * (*_mdot_soln)(node_out) / rho_in / S -
848  (*_mdot_soln)(node_out_j) / Sj / rho_j -
849  (*_mdot_soln)(node_out_i) / Si / rho_i);
850  counter++;
851  }
852  turbulent_term *= _CT;
853  auto Re = (((*_mdot_soln)(node_in) / S) * Dh_i / mu_in);
854  auto fi = computeFrictionFactor(Re);
855  auto ki = k_grid[i_ch][iz - 1];
856  auto friction_term = (fi * dz / Dh_i + ki) * 0.5 *
857  (std::pow((*_mdot_soln)(node_out), 2.0)) /
858  (S * (*_rho_soln)(node_out));
859  auto gravity_term = _g_grav * (*_rho_soln)(node_out)*dz * S;
860  auto DP = std::pow(S, -1.0) * (time_term + mass_term1 + mass_term2 + crossflow_term +
861  turbulent_term + friction_term + gravity_term); // Pa
862  _DP_soln->set(node_out, DP);
863  }
864  }
865  }
866  else
867  {
868  LibmeshPetscCall(MatZeroEntries(_amc_time_derivative_mat));
869  LibmeshPetscCall(MatZeroEntries(_amc_advective_derivative_mat));
870  LibmeshPetscCall(MatZeroEntries(_amc_cross_derivative_mat));
871  LibmeshPetscCall(MatZeroEntries(_amc_friction_force_mat));
872  LibmeshPetscCall(VecZeroEntries(_amc_time_derivative_rhs));
873  LibmeshPetscCall(VecZeroEntries(_amc_advective_derivative_rhs));
874  LibmeshPetscCall(VecZeroEntries(_amc_cross_derivative_rhs));
875  LibmeshPetscCall(VecZeroEntries(_amc_friction_force_rhs));
876  LibmeshPetscCall(VecZeroEntries(_amc_gravity_rhs));
877  LibmeshPetscCall(MatZeroEntries(_amc_sys_mdot_mat));
878  LibmeshPetscCall(VecZeroEntries(_amc_sys_mdot_rhs));
879  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
880  {
881  auto k_grid = _subchannel_mesh.getKGrid();
882  auto dz = _z_grid[iz] - _z_grid[iz - 1];
883  auto iz_ind = iz - first_node;
884  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
885  {
886  // inlet and outlet nodes
887  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
888  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
889 
890  // interpolation weight coefficient
891  PetscScalar Pe = 0.5;
893 
894  // inlet, outlet, and interpolated density
895  auto rho_in = (*_rho_soln)(node_in);
896  auto rho_out = (*_rho_soln)(node_out);
897  auto rho_interp = computeInterpolatedValue(rho_out, rho_in, Pe);
898 
899  // inlet, outlet, and interpolated viscosity
900  auto mu_in = (*_mu_soln)(node_in);
901  auto mu_out = (*_mu_soln)(node_out);
902  auto mu_interp = computeInterpolatedValue(mu_out, mu_in, Pe);
903 
904  // inlet, outlet, and interpolated axial surface area
905  auto S_in = (*_S_flow_soln)(node_in);
906  auto S_out = (*_S_flow_soln)(node_out);
907  auto S_interp = computeInterpolatedValue(S_out, S_in, Pe);
908 
909  // inlet, outlet, and interpolated wetted perimeter
910  auto w_perim_in = (*_w_perim_soln)(node_in);
911  auto w_perim_out = (*_w_perim_soln)(node_out);
912  auto w_perim_interp = computeInterpolatedValue(w_perim_out, w_perim_in, Pe);
913 
914  // hydraulic diameter in the i direction
915  auto Dh_i = 4.0 * S_interp / w_perim_interp;
916 
918  if (iz == first_node)
919  {
920  PetscScalar value_vec_tt = -1.0 * _TR * alpha * (*_mdot_soln)(node_in)*dz / _dt;
921  PetscInt row_vec_tt = i_ch + _n_channels * iz_ind;
922  LibmeshPetscCall(
923  VecSetValues(_amc_time_derivative_rhs, 1, &row_vec_tt, &value_vec_tt, ADD_VALUES));
924  }
925  else
926  {
927  PetscInt row_tt = i_ch + _n_channels * iz_ind;
928  PetscInt col_tt = i_ch + _n_channels * (iz_ind - 1);
929  PetscScalar value_tt = _TR * alpha * dz / _dt;
930  LibmeshPetscCall(MatSetValues(
931  _amc_time_derivative_mat, 1, &row_tt, 1, &col_tt, &value_tt, INSERT_VALUES));
932  }
933 
934  // Adding diagonal elements
935  PetscInt row_tt = i_ch + _n_channels * iz_ind;
936  PetscInt col_tt = i_ch + _n_channels * iz_ind;
937  PetscScalar value_tt = _TR * (1.0 - alpha) * dz / _dt;
938  LibmeshPetscCall(MatSetValues(
939  _amc_time_derivative_mat, 1, &row_tt, 1, &col_tt, &value_tt, INSERT_VALUES));
940 
941  // Adding RHS elements
942  PetscScalar mdot_old_interp =
943  computeInterpolatedValue(_mdot_soln->old(node_out), _mdot_soln->old(node_in), Pe);
944  PetscScalar value_vec_tt = _TR * mdot_old_interp * dz / _dt;
945  PetscInt row_vec_tt = i_ch + _n_channels * iz_ind;
946  LibmeshPetscCall(
947  VecSetValues(_amc_time_derivative_rhs, 1, &row_vec_tt, &value_vec_tt, ADD_VALUES));
948 
950  if (iz == first_node)
951  {
952  PetscScalar value_vec_at = std::pow((*_mdot_soln)(node_in), 2.0) / (S_in * rho_in);
953  PetscInt row_vec_at = i_ch + _n_channels * iz_ind;
954  LibmeshPetscCall(VecSetValues(
955  _amc_advective_derivative_rhs, 1, &row_vec_at, &value_vec_at, ADD_VALUES));
956  }
957  else
958  {
959  PetscInt row_at = i_ch + _n_channels * iz_ind;
960  PetscInt col_at = i_ch + _n_channels * (iz_ind - 1);
961  PetscScalar value_at = -1.0 * std::abs((*_mdot_soln)(node_in)) / (S_in * rho_in);
962  LibmeshPetscCall(MatSetValues(
963  _amc_advective_derivative_mat, 1, &row_at, 1, &col_at, &value_at, INSERT_VALUES));
964  }
965 
966  // Adding diagonal elements
967  PetscInt row_at = i_ch + _n_channels * iz_ind;
968  PetscInt col_at = i_ch + _n_channels * iz_ind;
969  PetscScalar value_at = std::abs((*_mdot_soln)(node_out)) / (S_out * rho_out);
970  LibmeshPetscCall(MatSetValues(
971  _amc_advective_derivative_mat, 1, &row_at, 1, &col_at, &value_at, INSERT_VALUES));
972 
974  unsigned int counter = 0;
975  unsigned int cross_index = iz; // iz-1;
976  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
977  {
978  auto chans = _subchannel_mesh.getGapChannels(i_gap);
979  unsigned int ii_ch = chans.first;
980  unsigned int jj_ch = chans.second;
981  auto * node_in_i = _subchannel_mesh.getChannelNode(ii_ch, iz - 1);
982  auto * node_in_j = _subchannel_mesh.getChannelNode(jj_ch, iz - 1);
983  auto * node_out_i = _subchannel_mesh.getChannelNode(ii_ch, iz);
984  auto * node_out_j = _subchannel_mesh.getChannelNode(jj_ch, iz);
985  auto rho_i =
986  computeInterpolatedValue((*_rho_soln)(node_out_i), (*_rho_soln)(node_in_i), Pe);
987  auto rho_j =
988  computeInterpolatedValue((*_rho_soln)(node_out_j), (*_rho_soln)(node_in_j), Pe);
989  auto S_i =
990  computeInterpolatedValue((*_S_flow_soln)(node_out_i), (*_S_flow_soln)(node_in_i), Pe);
991  auto S_j =
992  computeInterpolatedValue((*_S_flow_soln)(node_out_j), (*_S_flow_soln)(node_in_j), Pe);
993  auto u_star = 0.0;
994  // figure out donor axial velocity
995  if (_Wij(i_gap, cross_index) > 0.0)
996  {
997  if (iz == first_node)
998  {
999  u_star = (*_mdot_soln)(node_in_i) / S_i / rho_i;
1000  PetscScalar value_vec_ct = -1.0 * alpha *
1001  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1002  _Wij(i_gap, cross_index) * u_star;
1003  PetscInt row_vec_ct = i_ch + _n_channels * iz_ind;
1004  LibmeshPetscCall(VecSetValues(
1005  _amc_cross_derivative_rhs, 1, &row_vec_ct, &value_vec_ct, ADD_VALUES));
1006  }
1007  else
1008  {
1009  PetscScalar value_ct = alpha * _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1010  _Wij(i_gap, cross_index) / S_i / rho_i;
1011  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1012  PetscInt col_ct = ii_ch + _n_channels * (iz_ind - 1);
1013  LibmeshPetscCall(MatSetValues(
1014  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1015  }
1016  PetscScalar value_ct = (1.0 - alpha) *
1017  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1018  _Wij(i_gap, cross_index) / S_i / rho_i;
1019  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1020  PetscInt col_ct = ii_ch + _n_channels * iz_ind;
1021  LibmeshPetscCall(MatSetValues(
1022  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1023  }
1024  else if (_Wij(i_gap, cross_index) < 0.0) // _Wij=0 operations not necessary
1025  {
1026  if (iz == first_node)
1027  {
1028  u_star = (*_mdot_soln)(node_in_j) / S_j / rho_j;
1029  PetscScalar value_vec_ct = -1.0 * alpha *
1030  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1031  _Wij(i_gap, cross_index) * u_star;
1032  PetscInt row_vec_ct = i_ch + _n_channels * iz_ind;
1033  LibmeshPetscCall(VecSetValues(
1034  _amc_cross_derivative_rhs, 1, &row_vec_ct, &value_vec_ct, ADD_VALUES));
1035  }
1036  else
1037  {
1038  PetscScalar value_ct = alpha * _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1039  _Wij(i_gap, cross_index) / S_j / rho_j;
1040  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1041  PetscInt col_ct = jj_ch + _n_channels * (iz_ind - 1);
1042  LibmeshPetscCall(MatSetValues(
1043  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1044  }
1045  PetscScalar value_ct = (1.0 - alpha) *
1046  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1047  _Wij(i_gap, cross_index) / S_j / rho_j;
1048  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1049  PetscInt col_ct = jj_ch + _n_channels * iz_ind;
1050  LibmeshPetscCall(MatSetValues(
1051  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1052  }
1053 
1054  if (iz == first_node)
1055  {
1056  PetscScalar value_vec_ct = -2.0 * alpha *
1057  (*_mdot_soln)(node_in)*_WijPrime(i_gap, cross_index) /
1058  (rho_interp * S_interp);
1059  value_vec_ct +=
1060  alpha * (*_mdot_soln)(node_in_j)*_WijPrime(i_gap, cross_index) / (rho_j * S_j);
1061  value_vec_ct +=
1062  alpha * (*_mdot_soln)(node_in_i)*_WijPrime(i_gap, cross_index) / (rho_i * S_i);
1063  PetscInt row_vec_ct = i_ch + _n_channels * iz_ind;
1064  LibmeshPetscCall(
1065  VecSetValues(_amc_cross_derivative_rhs, 1, &row_vec_ct, &value_vec_ct, ADD_VALUES));
1066  }
1067  else
1068  {
1069  PetscScalar value_center_ct =
1070  2.0 * alpha * _WijPrime(i_gap, cross_index) / (rho_interp * S_interp);
1071  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1072  PetscInt col_ct = i_ch + _n_channels * (iz_ind - 1);
1073  LibmeshPetscCall(MatSetValues(
1074  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_center_ct, ADD_VALUES));
1075 
1076  PetscScalar value_left_ct =
1077  -1.0 * alpha * _WijPrime(i_gap, cross_index) / (rho_j * S_j);
1078  row_ct = i_ch + _n_channels * iz_ind;
1079  col_ct = jj_ch + _n_channels * (iz_ind - 1);
1080  LibmeshPetscCall(MatSetValues(
1081  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_left_ct, ADD_VALUES));
1082 
1083  PetscScalar value_right_ct =
1084  -1.0 * alpha * _WijPrime(i_gap, cross_index) / (rho_i * S_i);
1085  row_ct = i_ch + _n_channels * iz_ind;
1086  col_ct = ii_ch + _n_channels * (iz_ind - 1);
1087  LibmeshPetscCall(MatSetValues(
1088  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_right_ct, ADD_VALUES));
1089  }
1090 
1091  PetscScalar value_center_ct =
1092  2.0 * (1.0 - alpha) * _WijPrime(i_gap, cross_index) / (rho_interp * S_interp);
1093  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1094  PetscInt col_ct = i_ch + _n_channels * iz_ind;
1095  LibmeshPetscCall(MatSetValues(
1096  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_center_ct, ADD_VALUES));
1097 
1098  PetscScalar value_left_ct =
1099  -1.0 * (1.0 - alpha) * _WijPrime(i_gap, cross_index) / (rho_j * S_j);
1100  row_ct = i_ch + _n_channels * iz_ind;
1101  col_ct = jj_ch + _n_channels * iz_ind;
1102  LibmeshPetscCall(MatSetValues(
1103  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_left_ct, ADD_VALUES));
1104 
1105  PetscScalar value_right_ct =
1106  -1.0 * (1.0 - alpha) * _WijPrime(i_gap, cross_index) / (rho_i * S_i);
1107  row_ct = i_ch + _n_channels * iz_ind;
1108  col_ct = ii_ch + _n_channels * iz_ind;
1109  LibmeshPetscCall(MatSetValues(
1110  _amc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_right_ct, ADD_VALUES));
1111 
1112  counter++;
1113  }
1114 
1116  PetscScalar mdot_interp =
1117  computeInterpolatedValue((*_mdot_soln)(node_out), (*_mdot_soln)(node_in), Pe);
1118  auto Re = ((mdot_interp / S_interp) * Dh_i / mu_interp);
1119  auto fi = computeFrictionFactor(Re);
1120  auto ki = computeInterpolatedValue(k_grid[i_ch][iz], k_grid[i_ch][iz - 1], Pe);
1121  auto coef = (fi * dz / Dh_i + ki) * 0.5 * std::abs((*_mdot_soln)(node_out)) /
1122  (S_interp * rho_interp);
1123  if (iz == first_node)
1124  {
1125  PetscScalar value_vec = -1.0 * alpha * coef * (*_mdot_soln)(node_in);
1126  PetscInt row_vec = i_ch + _n_channels * iz_ind;
1127  LibmeshPetscCall(
1128  VecSetValues(_amc_friction_force_rhs, 1, &row_vec, &value_vec, ADD_VALUES));
1129  }
1130  else
1131  {
1132  PetscInt row = i_ch + _n_channels * iz_ind;
1133  PetscInt col = i_ch + _n_channels * (iz_ind - 1);
1134  PetscScalar value = alpha * coef;
1135  LibmeshPetscCall(
1136  MatSetValues(_amc_friction_force_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
1137  }
1138 
1139  // Adding diagonal elements
1140  PetscInt row = i_ch + _n_channels * iz_ind;
1141  PetscInt col = i_ch + _n_channels * iz_ind;
1142  PetscScalar value = (1.0 - alpha) * coef;
1143  LibmeshPetscCall(
1144  MatSetValues(_amc_friction_force_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
1145 
1147  PetscScalar value_vec = -1.0 * _g_grav * rho_interp * dz * S_interp;
1148  PetscInt row_vec = i_ch + _n_channels * iz_ind;
1149  LibmeshPetscCall(VecSetValues(_amc_gravity_rhs, 1, &row_vec, &value_vec, ADD_VALUES));
1150  }
1151  }
1153  LibmeshPetscCall(MatZeroEntries(_amc_sys_mdot_mat));
1154  LibmeshPetscCall(VecZeroEntries(_amc_sys_mdot_rhs));
1155  LibmeshPetscCall(MatAssemblyBegin(_amc_time_derivative_mat, MAT_FINAL_ASSEMBLY));
1156  LibmeshPetscCall(MatAssemblyEnd(_amc_time_derivative_mat, MAT_FINAL_ASSEMBLY));
1157  LibmeshPetscCall(MatAssemblyBegin(_amc_advective_derivative_mat, MAT_FINAL_ASSEMBLY));
1158  LibmeshPetscCall(MatAssemblyEnd(_amc_advective_derivative_mat, MAT_FINAL_ASSEMBLY));
1159  LibmeshPetscCall(MatAssemblyBegin(_amc_cross_derivative_mat, MAT_FINAL_ASSEMBLY));
1160  LibmeshPetscCall(MatAssemblyEnd(_amc_cross_derivative_mat, MAT_FINAL_ASSEMBLY));
1161  LibmeshPetscCall(MatAssemblyBegin(_amc_friction_force_mat, MAT_FINAL_ASSEMBLY));
1162  LibmeshPetscCall(MatAssemblyEnd(_amc_friction_force_mat, MAT_FINAL_ASSEMBLY));
1163  LibmeshPetscCall(MatAssemblyBegin(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1164  LibmeshPetscCall(MatAssemblyEnd(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1165  // Matrix
1166 #if !PETSC_VERSION_LESS_THAN(3, 15, 0)
1167  LibmeshPetscCall(
1168  MatAXPY(_amc_sys_mdot_mat, 1.0, _amc_time_derivative_mat, UNKNOWN_NONZERO_PATTERN));
1169  LibmeshPetscCall(MatAssemblyBegin(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1170  LibmeshPetscCall(MatAssemblyEnd(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1171  LibmeshPetscCall(
1172  MatAXPY(_amc_sys_mdot_mat, 1.0, _amc_advective_derivative_mat, UNKNOWN_NONZERO_PATTERN));
1173  LibmeshPetscCall(MatAssemblyBegin(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1174  LibmeshPetscCall(MatAssemblyEnd(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1175  LibmeshPetscCall(
1176  MatAXPY(_amc_sys_mdot_mat, 1.0, _amc_friction_force_mat, UNKNOWN_NONZERO_PATTERN));
1177 #else
1178  LibmeshPetscCall(
1179  MatAXPY(_amc_sys_mdot_mat, 1.0, _amc_time_derivative_mat, DIFFERENT_NONZERO_PATTERN));
1180  LibmeshPetscCall(MatAssemblyBegin(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1181  LibmeshPetscCall(MatAssemblyEnd(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1182  LibmeshPetscCall(
1183  MatAXPY(_amc_sys_mdot_mat, 1.0, _amc_advective_derivative_mat, DIFFERENT_NONZERO_PATTERN));
1184  LibmeshPetscCall(MatAssemblyBegin(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1185  LibmeshPetscCall(MatAssemblyEnd(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1186  LibmeshPetscCall(
1187  MatAXPY(_amc_sys_mdot_mat, 1.0, _amc_friction_force_mat, DIFFERENT_NONZERO_PATTERN));
1188 #endif
1189  LibmeshPetscCall(MatAssemblyBegin(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1190  LibmeshPetscCall(MatAssemblyEnd(_amc_sys_mdot_mat, MAT_FINAL_ASSEMBLY));
1191  _console << "Block: " << iblock << " - Linear momentum conservation matrix assembled"
1192  << std::endl;
1193  // RHS
1194  LibmeshPetscCall(VecAXPY(_amc_sys_mdot_rhs, 1.0, _amc_time_derivative_rhs));
1195  LibmeshPetscCall(VecAXPY(_amc_sys_mdot_rhs, 1.0, _amc_advective_derivative_rhs));
1196  LibmeshPetscCall(VecAXPY(_amc_sys_mdot_rhs, 1.0, _amc_friction_force_rhs));
1197  LibmeshPetscCall(VecAXPY(_amc_sys_mdot_rhs, 1.0, _amc_gravity_rhs));
1198 
1199  if (_segregated_bool)
1200  {
1201  // Assembly the matrix system
1202  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
1203  _prod, *_mdot_soln, first_node, last_node, _n_channels));
1204  Vec ls;
1205  LibmeshPetscCall(VecDuplicate(_amc_sys_mdot_rhs, &ls));
1206  LibmeshPetscCall(MatMult(_amc_sys_mdot_mat, _prod, ls));
1207  LibmeshPetscCall(VecAXPY(ls, -1.0, _amc_sys_mdot_rhs));
1208  PetscScalar * xx;
1209  LibmeshPetscCall(VecGetArray(ls, &xx));
1210  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1211  {
1212  auto iz_ind = iz - first_node;
1213  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1214  {
1215  // Setting nodes
1216  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1217  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1218 
1219  // inlet, outlet, and interpolated axial surface area
1220  auto S_in = (*_S_flow_soln)(node_in);
1221  auto S_out = (*_S_flow_soln)(node_out);
1222  auto S_interp = computeInterpolatedValue(S_out, S_in);
1223 
1224  // Setting solutions
1225  if (S_interp != 0)
1226  {
1227  auto DP = std::pow(S_interp, -1.0) * xx[iz_ind * _n_channels + i_ch];
1228  _DP_soln->set(node_out, DP);
1229  }
1230  else
1231  {
1232  auto DP = 0.0;
1233  _DP_soln->set(node_out, DP);
1234  }
1235  }
1236  }
1237  LibmeshPetscCall(VecDestroy(&ls));
1238  }
1239  }
1240 }
virtual PetscScalar computeInterpolatedValue(PetscScalar topValue, PetscScalar botValue, PetscScalar Peclet=0.0)
virtual const std::vector< std::vector< Real > > & getKGrid() const
Get axial cell location and value of loss coefficient.
const bool _segregated_bool
Segregated solve.
Real _TR
Flag that activates or deactivates the transient parts of the equations we solve by multiplication...
std::unique_ptr< SolutionHandle > _S_flow_soln
Mat _amc_advective_derivative_mat
Axial momentum conservation - advective (Eulerian) derivative.
virtual const std::pair< unsigned int, unsigned int > & getGapChannels(unsigned int i_gap) const =0
Return a pair of inter-wrapper indices for a given gap index.
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
std::unique_ptr< SolutionHandle > _SumWij_soln
std::unique_ptr< SolutionHandle > _rho_soln
std::unique_ptr< SolutionHandle > _mdot_soln
std::vector< Real > _z_grid
axial location of nodes
Mat _amc_cross_derivative_mat
Axial momentum conservation - cross flux derivative.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
static const std::string S
Definition: NS.h:163
Mat _amc_time_derivative_mat
Axial momentum conservation - time derivative.
const InterWrapperMesh & _subchannel_mesh
virtual PetscScalar computeInterpolationCoefficients(PetscScalar Peclet=0.0)
Functions that computes the interpolation scheme given the Peclet number.
std::unique_ptr< SolutionHandle > _DP_soln
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string alpha
Definition: NS.h:134
virtual Real computeFrictionFactor(Real Re)=0
Returns friction factor.
const ConsoleStream _console
libMesh::DenseMatrix< Real > & _Wij
libMesh::DenseMatrix< Real > _WijPrime
MooseUnits pow(const MooseUnits &, int)
virtual const std::vector< unsigned int > & getChannelGaps(unsigned int i_chan) const =0
Return a vector of gap indices for a given channel index.
Vec _amc_gravity_rhs
Axial momentum conservation - buoyancy force No implicit matrix.
const Real & _CT
Turbulent modeling parameter used in axial momentum equation.
Mat _amc_sys_mdot_mat
Axial momentum system matrix.
Mat _amc_friction_force_mat
Axial momentum conservation - friction force.
virtual const Real & getCrossflowSign(unsigned int i_chan, unsigned int i_local) const =0
Return a signs for the cross flow given a inter-wrapper index and local neighbor index.

◆ computeFrictionFactor()

double QuadInterWrapper1PhaseProblem::computeFrictionFactor ( Real  Re)
overrideprotectedvirtual

Returns friction factor.

Implements InterWrapper1PhaseProblem.

Definition at line 29 of file QuadInterWrapper1PhaseProblem.C.

30 {
31  Real a, b;
32  if (Re < 1)
33  {
34  return 64.0;
35  }
36  else if (Re >= 1 and Re < 5000)
37  {
38  a = 64.0;
39  b = -1.0;
40  }
41  else if (Re >= 5000 and Re < 30000)
42  {
43  a = 0.316;
44  b = -0.25;
45  }
46  else
47  {
48  a = 0.184;
49  b = -0.20;
50  }
51  return a * std::pow(Re, b);
52 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MooseUnits pow(const MooseUnits &, int)

◆ computeh()

void InterWrapper1PhaseProblem::computeh ( int  iblock)
protectedvirtualinherited

Computes Enthalpy per channel for block iblock.

Time derivative term

Advective derivative term

Cross derivative term

Added heat enthalpy

Assembling system

Reimplemented in TriInterWrapper1PhaseProblem.

Definition at line 1502 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and InterWrapper1PhaseProblem::implicitPetscSolve().

1503 {
1504  unsigned int last_node = (iblock + 1) * _block_size;
1505  unsigned int first_node = iblock * _block_size + 1;
1506 
1507  if (iblock == 0)
1508  {
1509  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1510  {
1511  auto * node = _subchannel_mesh.getChannelNode(i_ch, 0);
1512  auto h_out = _fp->h_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node));
1513  if (h_out < 0)
1514  {
1515  mooseError(
1516  name(), " : Calculation of negative Enthalpy h_out = : ", h_out, " Axial Level= : ", 0);
1517  }
1518  _h_soln->set(node, h_out);
1519  }
1520  }
1521 
1522  if (!_implicit_bool)
1523  {
1524  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1525  {
1526  auto dz = _z_grid[iz] - _z_grid[iz - 1];
1527  auto heated_length = _subchannel_mesh.getHeatedLength();
1528  auto unheated_length_entry = _subchannel_mesh.getHeatedLengthEntry();
1529  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1530  {
1531  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1532  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1533  auto mdot_in = (*_mdot_soln)(node_in);
1534  auto h_in = (*_h_soln)(node_in); // J/kg
1535  auto volume = dz * (*_S_flow_soln)(node_in);
1536  auto mdot_out = (*_mdot_soln)(node_out);
1537  auto h_out = 0.0;
1538  Real sumWijh = 0.0;
1539  Real sumWijPrimeDhij = 0.0;
1540  Real added_enthalpy;
1541  if (_z_grid[iz] > unheated_length_entry &&
1542  _z_grid[iz] <= unheated_length_entry + heated_length)
1543  added_enthalpy = computeAddedHeat(i_ch, iz);
1544  else
1545  added_enthalpy = 0.0;
1546 
1547  // Calculate sum of crossflow into channel i from channels j around i
1548  unsigned int counter = 0;
1549  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
1550  {
1551  auto chans = _subchannel_mesh.getGapChannels(i_gap);
1552  unsigned int ii_ch = chans.first;
1553  // i is always the smallest and first index in the mapping
1554  unsigned int jj_ch = chans.second;
1555  auto * node_in_i = _subchannel_mesh.getChannelNode(ii_ch, iz - 1);
1556  auto * node_in_j = _subchannel_mesh.getChannelNode(jj_ch, iz - 1);
1557  // Define donor enthalpy
1558  auto h_star = 0.0;
1559  if (_Wij(i_gap, iz) > 0.0)
1560  h_star = (*_h_soln)(node_in_i);
1561  else if (_Wij(i_gap, iz) < 0.0)
1562  h_star = (*_h_soln)(node_in_j);
1563  // take care of the sign by applying the map, use donor cell
1564  sumWijh += _subchannel_mesh.getCrossflowSign(i_ch, counter) * _Wij(i_gap, iz) * h_star;
1565  sumWijPrimeDhij += _WijPrime(i_gap, iz) * (2 * (*_h_soln)(node_in) -
1566  (*_h_soln)(node_in_j) - (*_h_soln)(node_in_i));
1567  counter++;
1568  }
1569 
1570  h_out = (mdot_in * h_in - sumWijh - sumWijPrimeDhij + added_enthalpy +
1571  _TR * _rho_soln->old(node_out) * _h_soln->old(node_out) * volume / _dt) /
1572  (mdot_out + _TR * (*_rho_soln)(node_out)*volume / _dt);
1573 
1574  if (h_out < 0)
1575  {
1576  mooseError(name(),
1577  " : Calculation of negative Enthalpy h_out = : ",
1578  h_out,
1579  " Axial Level= : ",
1580  iz);
1581  }
1582  _h_soln->set(node_out, h_out); // J/kg
1583  }
1584  }
1585  }
1586  else
1587  {
1588  _console << "Setting matrices" << std::endl;
1589  LibmeshPetscCall(MatZeroEntries(_hc_time_derivative_mat));
1590  LibmeshPetscCall(MatZeroEntries(_hc_advective_derivative_mat));
1591  LibmeshPetscCall(MatZeroEntries(_hc_cross_derivative_mat));
1592  LibmeshPetscCall(VecZeroEntries(_hc_time_derivative_rhs));
1593  LibmeshPetscCall(VecZeroEntries(_hc_advective_derivative_rhs));
1594  LibmeshPetscCall(VecZeroEntries(_hc_cross_derivative_rhs));
1595  LibmeshPetscCall(VecZeroEntries(_hc_added_heat_rhs));
1596  LibmeshPetscCall(MatZeroEntries(_hc_sys_h_mat));
1597  LibmeshPetscCall(VecZeroEntries(_hc_sys_h_rhs));
1598 
1599  _console << "Starting enthalpy assembly" << std::endl;
1600  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1601  {
1602  auto dz = _z_grid[iz] - _z_grid[iz - 1];
1603  auto heated_length = _subchannel_mesh.getHeatedLength();
1604  auto unheated_length_entry = _subchannel_mesh.getHeatedLengthEntry();
1605  auto iz_ind = iz - first_node;
1606  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1607  {
1608  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1609  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1610  auto volume = dz * (*_S_flow_soln)(node_in);
1611 
1612  // interpolation weight coefficient
1613  PetscScalar Pe = 0.5;
1615 
1617  if (iz == first_node)
1618  {
1619  PetscScalar value_vec_tt =
1620  -1.0 * _TR * alpha * (*_rho_soln)(node_in) * (*_h_soln)(node_in)*volume / _dt;
1621  PetscInt row_vec_tt = i_ch + _n_channels * iz_ind;
1622  LibmeshPetscCall(
1623  VecSetValues(_hc_time_derivative_rhs, 1, &row_vec_tt, &value_vec_tt, ADD_VALUES));
1624  }
1625  else
1626  {
1627  PetscInt row_tt = i_ch + _n_channels * iz_ind;
1628  PetscInt col_tt = i_ch + _n_channels * (iz_ind - 1);
1629  PetscScalar value_tt = _TR * alpha * (*_rho_soln)(node_in)*volume / _dt;
1630  LibmeshPetscCall(MatSetValues(
1631  _hc_time_derivative_mat, 1, &row_tt, 1, &col_tt, &value_tt, INSERT_VALUES));
1632  }
1633 
1634  // Adding diagonal elements
1635  PetscInt row_tt = i_ch + _n_channels * iz_ind;
1636  PetscInt col_tt = i_ch + _n_channels * iz_ind;
1637  PetscScalar value_tt = _TR * (1.0 - alpha) * (*_rho_soln)(node_out)*volume / _dt;
1638  LibmeshPetscCall(MatSetValues(
1639  _hc_time_derivative_mat, 1, &row_tt, 1, &col_tt, &value_tt, INSERT_VALUES));
1640 
1641  // Adding RHS elements
1642  PetscScalar rho_old_interp =
1643  computeInterpolatedValue(_rho_soln->old(node_out), _rho_soln->old(node_in), Pe);
1644  PetscScalar h_old_interp =
1645  computeInterpolatedValue(_h_soln->old(node_out), _h_soln->old(node_in), Pe);
1646 
1647  PetscScalar value_vec_tt = _TR * rho_old_interp * h_old_interp * volume / _dt;
1648  PetscInt row_vec_tt = i_ch + _n_channels * iz_ind;
1649  LibmeshPetscCall(
1650  VecSetValues(_hc_time_derivative_rhs, 1, &row_vec_tt, &value_vec_tt, ADD_VALUES));
1651 
1653  if (iz == first_node)
1654  {
1655  PetscScalar value_vec_at = (*_mdot_soln)(node_in) * (*_h_soln)(node_in);
1656  PetscInt row_vec_at = i_ch + _n_channels * iz_ind;
1657  LibmeshPetscCall(VecSetValues(
1658  _hc_advective_derivative_rhs, 1, &row_vec_at, &value_vec_at, ADD_VALUES));
1659  }
1660  else
1661  {
1662  PetscInt row_at = i_ch + _n_channels * iz_ind;
1663  PetscInt col_at = i_ch + _n_channels * (iz_ind - 1);
1664  PetscScalar value_at = -1.0 * (*_mdot_soln)(node_in);
1665  LibmeshPetscCall(MatSetValues(
1666  _hc_advective_derivative_mat, 1, &row_at, 1, &col_at, &value_at, INSERT_VALUES));
1667  }
1668 
1669  // Adding diagonal elements
1670  PetscInt row_at = i_ch + _n_channels * iz_ind;
1671  PetscInt col_at = i_ch + _n_channels * iz_ind;
1672  PetscScalar value_at = (*_mdot_soln)(node_out);
1673  LibmeshPetscCall(MatSetValues(
1674  _hc_advective_derivative_mat, 1, &row_at, 1, &col_at, &value_at, INSERT_VALUES));
1675 
1677  unsigned int counter = 0;
1678  unsigned int cross_index = iz; // iz-1;
1679  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
1680  {
1681  auto chans = _subchannel_mesh.getGapChannels(i_gap);
1682  unsigned int ii_ch = chans.first;
1683  unsigned int jj_ch = chans.second;
1684  auto * node_in_i = _subchannel_mesh.getChannelNode(ii_ch, iz - 1);
1685  auto * node_in_j = _subchannel_mesh.getChannelNode(jj_ch, iz - 1);
1686 
1687  PetscScalar h_star;
1688  // figure out donor axial velocity
1689  if (_Wij(i_gap, cross_index) > 0.0)
1690  {
1691  if (iz == first_node)
1692  {
1693  h_star = (*_h_soln)(node_in_i);
1694  PetscScalar value_vec_ct = -1.0 * alpha *
1695  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1696  _Wij(i_gap, cross_index) * h_star;
1697  PetscInt row_vec_ct = i_ch + _n_channels * iz_ind;
1698  LibmeshPetscCall(VecSetValues(
1699  _hc_cross_derivative_rhs, 1, &row_vec_ct, &value_vec_ct, ADD_VALUES));
1700  }
1701  else
1702  {
1703  PetscScalar value_ct = alpha * _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1704  _Wij(i_gap, cross_index);
1705  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1706  PetscInt col_ct = ii_ch + _n_channels * (iz_ind - 1);
1707  LibmeshPetscCall(MatSetValues(
1708  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1709  }
1710  PetscScalar value_ct = (1.0 - alpha) *
1711  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1712  _Wij(i_gap, cross_index);
1713  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1714  PetscInt col_ct = ii_ch + _n_channels * iz_ind;
1715  LibmeshPetscCall(MatSetValues(
1716  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1717  }
1718  else if (_Wij(i_gap, cross_index) < 0.0) // _Wij=0 operations not necessary
1719  {
1720  if (iz == first_node)
1721  {
1722  h_star = (*_h_soln)(node_in_j);
1723  PetscScalar value_vec_ct = -1.0 * alpha *
1724  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1725  _Wij(i_gap, cross_index) * h_star;
1726  PetscInt row_vec_ct = i_ch + _n_channels * iz_ind;
1727  LibmeshPetscCall(VecSetValues(
1728  _hc_cross_derivative_rhs, 1, &row_vec_ct, &value_vec_ct, ADD_VALUES));
1729  }
1730  else
1731  {
1732  PetscScalar value_ct = alpha * _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1733  _Wij(i_gap, cross_index);
1734  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1735  PetscInt col_ct = jj_ch + _n_channels * (iz_ind - 1);
1736  LibmeshPetscCall(MatSetValues(
1737  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1738  }
1739  PetscScalar value_ct = (1.0 - alpha) *
1740  _subchannel_mesh.getCrossflowSign(i_ch, counter) *
1741  _Wij(i_gap, cross_index);
1742  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1743  PetscInt col_ct = jj_ch + _n_channels * iz_ind;
1744  LibmeshPetscCall(MatSetValues(
1745  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_ct, ADD_VALUES));
1746  }
1747 
1748  // Turbulent cross flows
1749  if (iz == first_node)
1750  {
1751  PetscScalar value_vec_ct =
1752  -2.0 * alpha * (*_mdot_soln)(node_in)*_WijPrime(i_gap, cross_index);
1753  value_vec_ct = alpha * (*_mdot_soln)(node_in_j)*_WijPrime(i_gap, cross_index);
1754  value_vec_ct += alpha * (*_mdot_soln)(node_in_i)*_WijPrime(i_gap, cross_index);
1755  PetscInt row_vec_ct = i_ch + _n_channels * iz_ind;
1756  LibmeshPetscCall(
1757  VecSetValues(_hc_cross_derivative_rhs, 1, &row_vec_ct, &value_vec_ct, ADD_VALUES));
1758  }
1759  else
1760  {
1761  PetscScalar value_center_ct = 2.0 * alpha * _WijPrime(i_gap, cross_index);
1762  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1763  PetscInt col_ct = i_ch + _n_channels * (iz_ind - 1);
1764  LibmeshPetscCall(MatSetValues(
1765  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_center_ct, ADD_VALUES));
1766 
1767  PetscScalar value_left_ct = -1.0 * alpha * _WijPrime(i_gap, cross_index);
1768  row_ct = i_ch + _n_channels * iz_ind;
1769  col_ct = jj_ch + _n_channels * (iz_ind - 1);
1770  LibmeshPetscCall(MatSetValues(
1771  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_left_ct, ADD_VALUES));
1772 
1773  PetscScalar value_right_ct = -1.0 * alpha * _WijPrime(i_gap, cross_index);
1774  row_ct = i_ch + _n_channels * iz_ind;
1775  col_ct = ii_ch + _n_channels * (iz_ind - 1);
1776  LibmeshPetscCall(MatSetValues(
1777  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_right_ct, ADD_VALUES));
1778  }
1779 
1780  PetscScalar value_center_ct = 2.0 * (1.0 - alpha) * _WijPrime(i_gap, cross_index);
1781  PetscInt row_ct = i_ch + _n_channels * iz_ind;
1782  PetscInt col_ct = i_ch + _n_channels * iz_ind;
1783  LibmeshPetscCall(MatSetValues(
1784  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_center_ct, ADD_VALUES));
1785 
1786  PetscScalar value_left_ct = -1.0 * (1.0 - alpha) * _WijPrime(i_gap, cross_index);
1787  row_ct = i_ch + _n_channels * iz_ind;
1788  col_ct = jj_ch + _n_channels * iz_ind;
1789  LibmeshPetscCall(MatSetValues(
1790  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_left_ct, ADD_VALUES));
1791 
1792  PetscScalar value_right_ct = -1.0 * (1.0 - alpha) * _WijPrime(i_gap, cross_index);
1793  row_ct = i_ch + _n_channels * iz_ind;
1794  col_ct = ii_ch + _n_channels * iz_ind;
1795  LibmeshPetscCall(MatSetValues(
1796  _hc_cross_derivative_mat, 1, &row_ct, 1, &col_ct, &value_right_ct, ADD_VALUES));
1797 
1798  counter++;
1799  }
1800 
1802  PetscScalar added_enthalpy;
1803  if (_z_grid[iz] > unheated_length_entry &&
1804  _z_grid[iz] <= unheated_length_entry + heated_length)
1805  added_enthalpy = computeAddedHeat(i_ch, iz);
1806  else
1807  added_enthalpy = 0.0;
1808 
1809  PetscInt row_vec_ht = i_ch + _n_channels * iz_ind;
1810  LibmeshPetscCall(
1811  VecSetValues(_hc_added_heat_rhs, 1, &row_vec_ht, &added_enthalpy, ADD_VALUES));
1812  }
1813  }
1814  _console << "Done with enthalpy assembly" << std::endl;
1815 
1817  LibmeshPetscCall(MatAssemblyBegin(_hc_time_derivative_mat, MAT_FINAL_ASSEMBLY));
1818  LibmeshPetscCall(MatAssemblyEnd(_hc_time_derivative_mat, MAT_FINAL_ASSEMBLY));
1819  LibmeshPetscCall(MatAssemblyBegin(_hc_advective_derivative_mat, MAT_FINAL_ASSEMBLY));
1820  LibmeshPetscCall(MatAssemblyEnd(_hc_advective_derivative_mat, MAT_FINAL_ASSEMBLY));
1821  LibmeshPetscCall(MatAssemblyBegin(_hc_cross_derivative_mat, MAT_FINAL_ASSEMBLY));
1822  LibmeshPetscCall(MatAssemblyEnd(_hc_cross_derivative_mat, MAT_FINAL_ASSEMBLY));
1823  LibmeshPetscCall(MatAssemblyBegin(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1824  LibmeshPetscCall(MatAssemblyEnd(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1825  // Matrix
1826 #if !PETSC_VERSION_LESS_THAN(3, 15, 0)
1827  LibmeshPetscCall(MatAXPY(_hc_sys_h_mat, 1.0, _hc_time_derivative_mat, UNKNOWN_NONZERO_PATTERN));
1828  LibmeshPetscCall(MatAssemblyBegin(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1829  LibmeshPetscCall(MatAssemblyEnd(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1830  LibmeshPetscCall(
1831  MatAXPY(_hc_sys_h_mat, 1.0, _hc_advective_derivative_mat, UNKNOWN_NONZERO_PATTERN));
1832  LibmeshPetscCall(MatAssemblyBegin(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1833  LibmeshPetscCall(MatAssemblyEnd(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1834  LibmeshPetscCall(
1835  MatAXPY(_hc_sys_h_mat, 1.0, _hc_cross_derivative_mat, UNKNOWN_NONZERO_PATTERN));
1836 #else
1837  LibmeshPetscCall(
1838  MatAXPY(_hc_sys_h_mat, 1.0, _hc_time_derivative_mat, DIFFERENT_NONZERO_PATTERN));
1839  LibmeshPetscCall(MatAssemblyBegin(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1840  LibmeshPetscCall(MatAssemblyEnd(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1841  LibmeshPetscCall(
1842  MatAXPY(_hc_sys_h_mat, 1.0, _hc_advective_derivative_mat, DIFFERENT_NONZERO_PATTERN));
1843  LibmeshPetscCall(MatAssemblyBegin(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1844  LibmeshPetscCall(MatAssemblyEnd(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1845  LibmeshPetscCall(
1846  MatAXPY(_hc_sys_h_mat, 1.0, _hc_cross_derivative_mat, DIFFERENT_NONZERO_PATTERN));
1847 #endif
1848  LibmeshPetscCall(MatAssemblyBegin(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1849  LibmeshPetscCall(MatAssemblyEnd(_hc_sys_h_mat, MAT_FINAL_ASSEMBLY));
1850  _console << "Block: " << iblock << " - Enthalpy conservation matrix assembled" << std::endl;
1851  // RHS
1852  LibmeshPetscCall(VecAXPY(_hc_sys_h_rhs, 1.0, _hc_time_derivative_rhs));
1853  LibmeshPetscCall(VecAXPY(_hc_sys_h_rhs, 1.0, _hc_advective_derivative_rhs));
1854  LibmeshPetscCall(VecAXPY(_hc_sys_h_rhs, 1.0, _hc_cross_derivative_rhs));
1855  LibmeshPetscCall(VecAXPY(_hc_sys_h_rhs, 1.0, _hc_added_heat_rhs));
1856 
1857  if (_segregated_bool)
1858  {
1859  // Assembly the matrix system
1860  KSP ksploc;
1861  PC pc;
1862  Vec sol;
1863  LibmeshPetscCall(VecDuplicate(_hc_sys_h_rhs, &sol));
1864  LibmeshPetscCall(KSPCreate(PETSC_COMM_WORLD, &ksploc));
1865  LibmeshPetscCall(KSPSetOperators(ksploc, _hc_sys_h_mat, _hc_sys_h_mat));
1866  LibmeshPetscCall(KSPGetPC(ksploc, &pc));
1867  LibmeshPetscCall(PCSetType(pc, PCJACOBI));
1868  LibmeshPetscCall(KSPSetTolerances(ksploc, _rtol, _atol, _dtol, _maxit));
1869  LibmeshPetscCall(KSPSetOptionsPrefix(ksploc, "h_cons_sys_"));
1870  LibmeshPetscCall(KSPSetFromOptions(ksploc));
1871  LibmeshPetscCall(KSPSolve(ksploc, _hc_sys_h_rhs, sol));
1872  PetscScalar * xx;
1873  LibmeshPetscCall(VecGetArray(sol, &xx));
1874 
1875  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1876  {
1877  auto iz_ind = iz - first_node;
1878  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1879  {
1880  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1881  auto h_out = xx[iz_ind * _n_channels + i_ch];
1882  if (h_out < 0)
1883  {
1884  mooseError(name(),
1885  " : Calculation of negative Enthalpy h_out = : ",
1886  h_out,
1887  " Axial Level= : ",
1888  iz);
1889  }
1890  _h_soln->set(node_out, h_out);
1891  }
1892  }
1893  LibmeshPetscCall(KSPDestroy(&ksploc));
1894  LibmeshPetscCall(VecDestroy(&sol));
1895  }
1896  }
1897 }
virtual PetscScalar computeInterpolatedValue(PetscScalar topValue, PetscScalar botValue, PetscScalar Peclet=0.0)
const bool _segregated_bool
Segregated solve.
Real _TR
Flag that activates or deactivates the transient parts of the equations we solve by multiplication...
std::unique_ptr< SolutionHandle > _P_soln
virtual const Real & getHeatedLength() const
Return heated length.
virtual const std::pair< unsigned int, unsigned int > & getGapChannels(unsigned int i_gap) const =0
Return a pair of inter-wrapper indices for a given gap index.
virtual const Real & getHeatedLengthEntry() const
Return unheated length at entry.
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
const Real & _P_out
Outlet Pressure.
const PetscReal & _dtol
The divergence tolerance for the ksp linear solver.
std::unique_ptr< SolutionHandle > _rho_soln
const PetscReal & _rtol
The relative convergence tolerance, (relative decrease) for the ksp linear solver.
const SinglePhaseFluidProperties * _fp
Solutions handles and link to TH tables properties.
virtual const std::string & name() const
Mat _hc_time_derivative_mat
Enthalpy Enthalpy conservation - time derivative.
const PetscInt & _maxit
The maximum number of iterations to use for the ksp linear solver.
std::vector< Real > _z_grid
axial location of nodes
std::unique_ptr< SolutionHandle > _T_soln
Mat _hc_advective_derivative_mat
Enthalpy conservation - advective (Eulerian) derivative;.
const InterWrapperMesh & _subchannel_mesh
virtual PetscScalar computeInterpolationCoefficients(PetscScalar Peclet=0.0)
Functions that computes the interpolation scheme given the Peclet number.
Real volume(const MeshBase &mesh, unsigned int dim=libMesh::invalid_uint)
Vec _hc_added_heat_rhs
Enthalpy conservation - source and sink.
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string alpha
Definition: NS.h:134
void mooseError(Args &&... args) const
std::unique_ptr< SolutionHandle > _h_soln
virtual Real computeAddedHeat(unsigned int i_ch, unsigned int iz)
Computes added heat for channel i_ch and cell iz.
const ConsoleStream _console
libMesh::DenseMatrix< Real > & _Wij
libMesh::DenseMatrix< Real > _WijPrime
Mat _hc_cross_derivative_mat
Enthalpy conservation - cross flux derivative.
virtual const std::vector< unsigned int > & getChannelGaps(unsigned int i_chan) const =0
Return a vector of gap indices for a given channel index.
const PetscReal & _atol
The absolute convergence tolerance for the ksp linear solver.
virtual const Real & getCrossflowSign(unsigned int i_chan, unsigned int i_local) const =0
Return a signs for the cross flow given a inter-wrapper index and local neighbor index.

◆ computeInterpolatedValue()

PetscScalar InterWrapper1PhaseProblem::computeInterpolatedValue ( PetscScalar  topValue,
PetscScalar  botValue,
PetscScalar  Peclet = 0.0 
)
protectedvirtualinherited

Definition at line 343 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::computeDP(), InterWrapper1PhaseProblem::computeh(), InterWrapper1PhaseProblem::computeP(), and InterWrapper1PhaseProblem::computeWij().

346 {
348  return alpha * botValue + (1.0 - alpha) * topValue;
349 }
auto Peclet(const T1 &volume_fraction, const T2 &cp, const T3 &rho, const T4 &vel, const T5 &D_h, const T6 &k)
Compute Peclet number.
Definition: Numerics.h:153
virtual PetscScalar computeInterpolationCoefficients(PetscScalar Peclet=0.0)
Functions that computes the interpolation scheme given the Peclet number.
static const std::string alpha
Definition: NS.h:134

◆ computeInterpolationCoefficients()

PetscScalar InterWrapper1PhaseProblem::computeInterpolationCoefficients ( PetscScalar  Peclet = 0.0)
protectedvirtualinherited

Functions that computes the interpolation scheme given the Peclet number.

Definition at line 316 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::computeDP(), InterWrapper1PhaseProblem::computeh(), InterWrapper1PhaseProblem::computeInterpolatedValue(), InterWrapper1PhaseProblem::computeP(), and InterWrapper1PhaseProblem::computeWij().

317 {
318  if (_interpolation_scheme == "upwind")
319  {
320  return 1.0;
321  }
322  else if (_interpolation_scheme == "downwind")
323  {
324  return 0.0;
325  }
326  else if (_interpolation_scheme == "central_difference")
327  {
328  return 0.5;
329  }
330  else if (_interpolation_scheme == "exponential")
331  {
332  return ((Peclet - 1.0) * std::exp(Peclet) + 1) / (Peclet * (std::exp(Peclet) - 1.) + 1e-10);
333  }
334  else
335  {
336  mooseError(name(),
337  ": Interpolation scheme should be a string: upwind, downwind, central_difference, "
338  "exponential");
339  }
340 }
virtual const std::string & name() const
const MooseEnum _interpolation_scheme
The interpolation method used in constructing the systems.
auto Peclet(const T1 &volume_fraction, const T2 &cp, const T3 &rho, const T4 &vel, const T5 &D_h, const T6 &k)
Compute Peclet number.
Definition: Numerics.h:153
void mooseError(Args &&... args) const

◆ computeMdot()

void InterWrapper1PhaseProblem::computeMdot ( int  iblock)
protectedvirtualinherited

Computes mass flow per channel for block iblock.

Definition at line 566 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::residualFunction().

567 {
568  unsigned int last_node = (iblock + 1) * _block_size;
569  unsigned int first_node = iblock * _block_size + 1;
570  if (!_implicit_bool)
571  {
572  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
573  {
574  auto dz = _z_grid[iz] - _z_grid[iz - 1];
575  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
576  {
577  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
578  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
579  auto volume = dz * (*_S_flow_soln)(node_in);
580  auto time_term = _TR * ((*_rho_soln)(node_out)-_rho_soln->old(node_out)) * volume / _dt;
581  // Wij positive out of i into j;
582  auto mdot_out = (*_mdot_soln)(node_in) - (*_SumWij_soln)(node_out)-time_term;
583  if (mdot_out < 0)
584  {
585  _console << "Wij = : " << _Wij << "\n";
586  mooseError(name(),
587  " : Calculation of negative mass flow mdot_out = : ",
588  mdot_out,
589  " Axial Level= : ",
590  iz,
591  " - Implicit solves are required for recirculating flow.");
592  }
593  _mdot_soln->set(node_out, mdot_out); // kg/sec
594  }
595  }
596  }
597  else
598  {
599  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
600  {
601  auto dz = _z_grid[iz] - _z_grid[iz - 1];
602  auto iz_ind = iz - first_node;
603  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
604  {
605  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
606  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
607  auto volume = dz * (*_S_flow_soln)(node_in);
608 
609  // Adding time derivative to the RHS
610  auto time_term = _TR * ((*_rho_soln)(node_out)-_rho_soln->old(node_out)) * volume / _dt;
611  PetscInt row_vec = i_ch + _n_channels * iz_ind;
612  PetscScalar value_vec = -1.0 * time_term;
613  LibmeshPetscCall(
614  VecSetValues(_mc_axial_convection_rhs, 1, &row_vec, &value_vec, INSERT_VALUES));
615 
616  // Imposing bottom boundary condition or adding of diagonal elements
617  if (iz == first_node)
618  {
619  PetscScalar value_vec = (*_mdot_soln)(node_in);
620  PetscInt row_vec = i_ch + _n_channels * iz_ind;
621  LibmeshPetscCall(
622  VecSetValues(_mc_axial_convection_rhs, 1, &row_vec, &value_vec, ADD_VALUES));
623  }
624  else
625  {
626  PetscInt row = i_ch + _n_channels * iz_ind;
627  PetscInt col = i_ch + _n_channels * (iz_ind - 1);
628  PetscScalar value = -1.0;
629  LibmeshPetscCall(
630  MatSetValues(_mc_axial_convection_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
631  }
632 
633  // Adding diagonal elements
634  PetscInt row = i_ch + _n_channels * iz_ind;
635  PetscInt col = i_ch + _n_channels * iz_ind;
636  PetscScalar value = 1.0;
637  LibmeshPetscCall(
638  MatSetValues(_mc_axial_convection_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
639 
640  // Adding cross flows RHS
641  if (_segregated_bool)
642  {
643  PetscScalar value_vec_2 = -1.0 * (*_SumWij_soln)(node_out);
644  PetscInt row_vec_2 = i_ch + _n_channels * iz_ind;
645  LibmeshPetscCall(
646  VecSetValues(_mc_axial_convection_rhs, 1, &row_vec_2, &value_vec_2, ADD_VALUES));
647  }
648  }
649  }
650  LibmeshPetscCall(MatAssemblyBegin(_mc_axial_convection_mat, MAT_FINAL_ASSEMBLY));
651  LibmeshPetscCall(MatAssemblyEnd(_mc_axial_convection_mat, MAT_FINAL_ASSEMBLY));
652  _console << "Block: " << iblock << " - Mass conservation matrix assembled" << std::endl;
653 
654  if (_segregated_bool)
655  {
656  KSP ksploc;
657  PC pc;
658  Vec sol;
659  LibmeshPetscCall(VecDuplicate(_mc_axial_convection_rhs, &sol));
660  LibmeshPetscCall(KSPCreate(PETSC_COMM_WORLD, &ksploc));
661  LibmeshPetscCall(KSPSetOperators(ksploc, _mc_axial_convection_mat, _mc_axial_convection_mat));
662  LibmeshPetscCall(KSPGetPC(ksploc, &pc));
663  LibmeshPetscCall(PCSetType(pc, PCJACOBI));
664  LibmeshPetscCall(KSPSetTolerances(ksploc, _rtol, _atol, _dtol, _maxit));
665  LibmeshPetscCall(KSPSetOptionsPrefix(ksploc, "mass_sys_"));
666  LibmeshPetscCall(KSPSetFromOptions(ksploc));
667  LibmeshPetscCall(KSPSolve(ksploc, _mc_axial_convection_rhs, sol));
668  LibmeshPetscCall(populateSolutionChan<SolutionHandle>(
669  sol, *_mdot_soln, first_node, last_node, _n_channels));
670  LibmeshPetscCall(VecZeroEntries(_mc_axial_convection_rhs));
671  LibmeshPetscCall(KSPDestroy(&ksploc));
672  LibmeshPetscCall(VecDestroy(&sol));
673  }
674  }
675 }
const bool _segregated_bool
Segregated solve.
Real _TR
Flag that activates or deactivates the transient parts of the equations we solve by multiplication...
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
std::unique_ptr< SolutionHandle > _SumWij_soln
const PetscReal & _dtol
The divergence tolerance for the ksp linear solver.
std::unique_ptr< SolutionHandle > _rho_soln
const PetscReal & _rtol
The relative convergence tolerance, (relative decrease) for the ksp linear solver.
Mat _mc_axial_convection_mat
Mass conservation - axial convection.
std::unique_ptr< SolutionHandle > _mdot_soln
virtual const std::string & name() const
const PetscInt & _maxit
The maximum number of iterations to use for the ksp linear solver.
std::vector< Real > _z_grid
axial location of nodes
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
const InterWrapperMesh & _subchannel_mesh
Real volume(const MeshBase &mesh, unsigned int dim=libMesh::invalid_uint)
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
void mooseError(Args &&... args) const
const ConsoleStream _console
libMesh::DenseMatrix< Real > & _Wij
const PetscReal & _atol
The absolute convergence tolerance for the ksp linear solver.

◆ computeMu()

void InterWrapper1PhaseProblem::computeMu ( int  iblock)
protectedvirtualinherited

Computes Viscosity per channel for block iblock.

Definition at line 1941 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

1942 {
1943  unsigned int last_node = (iblock + 1) * _block_size;
1944  unsigned int first_node = iblock * _block_size + 1;
1945 
1946  if (iblock == 0)
1947  {
1948  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1949  {
1950  auto * node = _subchannel_mesh.getChannelNode(i_ch, 0);
1951  _mu_soln->set(node, _fp->mu_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node)));
1952  }
1953  }
1954 
1955  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1956  {
1957  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1958  {
1959  auto * node = _subchannel_mesh.getChannelNode(i_ch, iz);
1960  _mu_soln->set(node, _fp->mu_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node)));
1961  }
1962  }
1963 }
std::unique_ptr< SolutionHandle > _P_soln
const Real & _P_out
Outlet Pressure.
const SinglePhaseFluidProperties * _fp
Solutions handles and link to TH tables properties.
std::unique_ptr< SolutionHandle > _T_soln
const InterWrapperMesh & _subchannel_mesh
std::unique_ptr< SolutionHandle > _mu_soln
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.

◆ computeP()

void InterWrapper1PhaseProblem::computeP ( int  iblock)
protectedvirtualinherited

Computes Pressure per channel for block iblock.

Definition at line 1243 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::residualFunction().

1244 {
1245  unsigned int last_node = (iblock + 1) * _block_size;
1246  unsigned int first_node = iblock * _block_size + 1;
1247  if (!_implicit_bool)
1248  {
1250  {
1251  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
1252  {
1253  // Calculate pressure in the inlet of the cell assuming known outlet
1254  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1255  {
1256  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1257  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1258  // update Pressure solution
1259  _P_soln->set(node_in, (*_P_soln)(node_out) + (*_DP_soln)(node_out));
1260  }
1261  }
1262  }
1263  else
1264  {
1265  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
1266  {
1267  // Calculate pressure in the inlet of the cell assuming known outlet
1268  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1269  {
1270  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1271  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1272  // update Pressure solution
1273  // Note: assuming uniform axial discretization in the curren code
1274  // We will need to update this later if we allow non-uniform refinements in the axial
1275  // direction
1276  PetscScalar Pe = 0.5;
1278  if (iz == last_node)
1279  {
1280  _P_soln->set(node_in, (*_P_soln)(node_out) + (*_DP_soln)(node_out) / 2.0);
1281  }
1282  else
1283  {
1284  _P_soln->set(node_in,
1285  (*_P_soln)(node_out) + (1.0 - alpha) * (*_DP_soln)(node_out) +
1286  alpha * (*_DP_soln)(node_in));
1287  }
1288  }
1289  }
1290  }
1291  }
1292  else
1293  {
1295  {
1296  LibmeshPetscCall(VecZeroEntries(_amc_pressure_force_rhs));
1297  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
1298  {
1299  auto iz_ind = iz - first_node;
1300  // Calculate pressure in the inlet of the cell assuming known outlet
1301  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1302  {
1303  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1304  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1305 
1306  // inlet, outlet, and interpolated axial surface area
1307  auto S_in = (*_S_flow_soln)(node_in);
1308  auto S_out = (*_S_flow_soln)(node_out);
1309  auto S_interp = computeInterpolatedValue(S_out, S_in);
1310 
1311  // Creating matrix of coefficients and RHS vector
1312  PetscInt row = i_ch + _n_channels * iz_ind;
1313  PetscInt col = i_ch + _n_channels * iz_ind;
1314  PetscScalar value = -1.0 * S_interp;
1315  LibmeshPetscCall(
1316  MatSetValues(_amc_pressure_force_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
1317 
1318  if (iz == last_node)
1319  {
1320  PetscScalar value = -1.0 * (*_P_soln)(node_out)*S_interp;
1321  PetscInt row = i_ch + _n_channels * iz_ind;
1322  LibmeshPetscCall(VecSetValues(_amc_pressure_force_rhs, 1, &row, &value, ADD_VALUES));
1323  }
1324  else
1325  {
1326  PetscInt row = i_ch + _n_channels * iz_ind;
1327  PetscInt col = i_ch + _n_channels * (iz_ind + 1);
1328  PetscScalar value = 1.0 * S_interp;
1329  LibmeshPetscCall(
1330  MatSetValues(_amc_pressure_force_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
1331  }
1332 
1333  if (_segregated_bool)
1334  {
1335  auto dp_out = (*_DP_soln)(node_out);
1336  PetscScalar value_v = -1.0 * dp_out * S_interp;
1337  PetscInt row_v = i_ch + _n_channels * iz_ind;
1338  LibmeshPetscCall(
1339  VecSetValues(_amc_pressure_force_rhs, 1, &row_v, &value_v, ADD_VALUES));
1340  }
1341  }
1342  }
1343  // Solving pressure problem
1344  LibmeshPetscCall(MatAssemblyBegin(_amc_pressure_force_mat, MAT_FINAL_ASSEMBLY));
1345  LibmeshPetscCall(MatAssemblyEnd(_amc_pressure_force_mat, MAT_FINAL_ASSEMBLY));
1346  if (_segregated_bool)
1347  {
1348  KSP ksploc;
1349  PC pc;
1350  Vec sol;
1351  LibmeshPetscCall(VecDuplicate(_amc_pressure_force_rhs, &sol));
1352  LibmeshPetscCall(KSPCreate(PETSC_COMM_WORLD, &ksploc));
1353  LibmeshPetscCall(KSPSetOperators(ksploc, _amc_pressure_force_mat, _amc_pressure_force_mat));
1354  LibmeshPetscCall(KSPGetPC(ksploc, &pc));
1355  LibmeshPetscCall(PCSetType(pc, PCJACOBI));
1356  LibmeshPetscCall(KSPSetTolerances(ksploc, _rtol, _atol, _dtol, _maxit));
1357  LibmeshPetscCall(KSPSetOptionsPrefix(ksploc, "pressure_sys_"));
1358  LibmeshPetscCall(KSPSetFromOptions(ksploc));
1359  LibmeshPetscCall(KSPSolve(ksploc, _amc_pressure_force_rhs, sol));
1360  PetscScalar * xx;
1361  LibmeshPetscCall(VecGetArray(sol, &xx));
1362  // update Pressure solution
1363  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
1364  {
1365  auto iz_ind = iz - first_node;
1366  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1367  {
1368  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1369  PetscScalar value = xx[iz_ind * _n_channels + i_ch];
1370  _P_soln->set(node_in, value);
1371  }
1372  }
1373  LibmeshPetscCall(VecZeroEntries(_amc_pressure_force_rhs));
1374  LibmeshPetscCall(KSPDestroy(&ksploc));
1375  LibmeshPetscCall(VecDestroy(&sol));
1376  }
1377  }
1378  else
1379  {
1380  LibmeshPetscCall(VecZeroEntries(_amc_pressure_force_rhs));
1381  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
1382  {
1383  auto iz_ind = iz - first_node;
1384  // Calculate pressure in the inlet of the cell assuming known outlet
1385  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1386  {
1387  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
1388  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1389 
1390  // inlet, outlet, and interpolated axial surface area
1391  auto S_in = (*_S_flow_soln)(node_in);
1392  auto S_out = (*_S_flow_soln)(node_out);
1393  auto S_interp = computeInterpolatedValue(S_out, S_in);
1394 
1395  // Creating matrix of coefficients
1396  PetscInt row = i_ch + _n_channels * iz_ind;
1397  PetscInt col = i_ch + _n_channels * iz_ind;
1398  PetscScalar value = -1.0 * S_interp;
1399  LibmeshPetscCall(
1400  MatSetValues(_amc_pressure_force_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
1401 
1402  if (iz == last_node)
1403  {
1404  PetscScalar value = -1.0 * (*_P_soln)(node_out)*S_interp;
1405  PetscInt row = i_ch + _n_channels * iz_ind;
1406  LibmeshPetscCall(VecSetValues(_amc_pressure_force_rhs, 1, &row, &value, ADD_VALUES));
1407 
1408  auto dp_out = (*_DP_soln)(node_out);
1409  PetscScalar value_v = -1.0 * dp_out / 2.0 * S_interp;
1410  PetscInt row_v = i_ch + _n_channels * iz_ind;
1411  LibmeshPetscCall(
1412  VecSetValues(_amc_pressure_force_rhs, 1, &row_v, &value_v, ADD_VALUES));
1413  }
1414  else
1415  {
1416  PetscInt row = i_ch + _n_channels * iz_ind;
1417  PetscInt col = i_ch + _n_channels * (iz_ind + 1);
1418  PetscScalar value = 1.0 * S_interp;
1419  LibmeshPetscCall(
1420  MatSetValues(_amc_pressure_force_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
1421 
1422  if (_segregated_bool)
1423  {
1424  auto dp_in = (*_DP_soln)(node_in);
1425  auto dp_out = (*_DP_soln)(node_out);
1426  auto dp_interp = computeInterpolatedValue(dp_out, dp_in);
1427  PetscScalar value_v = -1.0 * dp_interp * S_interp;
1428  PetscInt row_v = i_ch + _n_channels * iz_ind;
1429  LibmeshPetscCall(
1430  VecSetValues(_amc_pressure_force_rhs, 1, &row_v, &value_v, ADD_VALUES));
1431  }
1432  }
1433  }
1434  }
1435  // Solving pressure problem
1436  LibmeshPetscCall(MatAssemblyBegin(_amc_pressure_force_mat, MAT_FINAL_ASSEMBLY));
1437  LibmeshPetscCall(MatAssemblyEnd(_amc_pressure_force_mat, MAT_FINAL_ASSEMBLY));
1438  _console << "Block: " << iblock << " - Axial momentum pressure force matrix assembled"
1439  << std::endl;
1440 
1441  if (_segregated_bool)
1442  {
1443  KSP ksploc;
1444  PC pc;
1445  Vec sol;
1446  LibmeshPetscCall(VecDuplicate(_amc_pressure_force_rhs, &sol));
1447  LibmeshPetscCall(KSPCreate(PETSC_COMM_WORLD, &ksploc));
1448  LibmeshPetscCall(KSPSetOperators(ksploc, _amc_pressure_force_mat, _amc_pressure_force_mat));
1449  LibmeshPetscCall(KSPGetPC(ksploc, &pc));
1450  LibmeshPetscCall(PCSetType(pc, PCJACOBI));
1451  LibmeshPetscCall(KSPSetTolerances(ksploc, _rtol, _atol, _dtol, _maxit));
1452  LibmeshPetscCall(KSPSetOptionsPrefix(ksploc, "axial_mom_sys_"));
1453  LibmeshPetscCall(KSPSetFromOptions(ksploc));
1454  LibmeshPetscCall(KSPSolve(ksploc, _amc_pressure_force_rhs, sol));
1455  PetscScalar * xx;
1456  LibmeshPetscCall(VecGetArray(sol, &xx));
1457  // update Pressure solution
1458  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
1459  {
1460  auto iz_ind = iz - first_node;
1461  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1462  {
1463  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1464  PetscScalar value = xx[iz_ind * _n_channels + i_ch];
1465  _P_soln->set(node_in, value);
1466  }
1467  }
1468  LibmeshPetscCall(VecZeroEntries(_amc_pressure_force_rhs));
1469  LibmeshPetscCall(KSPDestroy(&ksploc));
1470  LibmeshPetscCall(VecDestroy(&sol));
1471  }
1472  }
1473  }
1474 }
virtual PetscScalar computeInterpolatedValue(PetscScalar topValue, PetscScalar botValue, PetscScalar Peclet=0.0)
const bool _segregated_bool
Segregated solve.
std::unique_ptr< SolutionHandle > _P_soln
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
const PetscReal & _dtol
The divergence tolerance for the ksp linear solver.
const PetscReal & _rtol
The relative convergence tolerance, (relative decrease) for the ksp linear solver.
const PetscInt & _maxit
The maximum number of iterations to use for the ksp linear solver.
const bool _staggered_pressure_bool
Flag to define the usage of staggered or collocated pressure.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
const InterWrapperMesh & _subchannel_mesh
virtual PetscScalar computeInterpolationCoefficients(PetscScalar Peclet=0.0)
Functions that computes the interpolation scheme given the Peclet number.
std::unique_ptr< SolutionHandle > _DP_soln
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
static const std::string alpha
Definition: NS.h:134
Mat _amc_pressure_force_mat
Axial momentum conservation - pressure force.
const ConsoleStream _console
const PetscReal & _atol
The absolute convergence tolerance for the ksp linear solver.

◆ computeRho()

void InterWrapper1PhaseProblem::computeRho ( int  iblock)
protectedvirtualinherited

Computes Density per channel for block iblock.

Definition at line 1916 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

1917 {
1918  unsigned int last_node = (iblock + 1) * _block_size;
1919  unsigned int first_node = iblock * _block_size + 1;
1920 
1921  if (iblock == 0)
1922  {
1923  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1924  {
1925  auto * node = _subchannel_mesh.getChannelNode(i_ch, 0);
1926  _rho_soln->set(node, _fp->rho_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node)));
1927  }
1928  }
1929 
1930  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1931  {
1932  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1933  {
1934  auto * node = _subchannel_mesh.getChannelNode(i_ch, iz);
1935  _rho_soln->set(node, _fp->rho_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node)));
1936  }
1937  }
1938 }
std::unique_ptr< SolutionHandle > _P_soln
const Real & _P_out
Outlet Pressure.
std::unique_ptr< SolutionHandle > _rho_soln
const SinglePhaseFluidProperties * _fp
Solutions handles and link to TH tables properties.
std::unique_ptr< SolutionHandle > _T_soln
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.

◆ computeSumWij()

void InterWrapper1PhaseProblem::computeSumWij ( int  iblock)
protectedvirtualinherited

Computes net diversion crossflow per channel for block iblock.

Definition at line 501 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::residualFunction().

502 {
503  unsigned int last_node = (iblock + 1) * _block_size;
504  unsigned int first_node = iblock * _block_size + 1;
505  // Add to solution vector if explicit
506  if (!_implicit_bool)
507  {
508  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
509  {
510  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
511  {
512  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
513  Real sumWij = 0.0;
514  // Calculate sum of crossflow into channel i from channels j around i
515  unsigned int counter = 0;
516  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
517  {
518  sumWij += _subchannel_mesh.getCrossflowSign(i_ch, counter) * _Wij(i_gap, iz);
519  counter++;
520  }
521  // The net crossflow coming out of cell i [kg/sec]
522  _SumWij_soln->set(node_out, sumWij);
523  }
524  }
525  }
526  // Add to matrix if implicit
527  else
528  {
529  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
530  {
531  unsigned int iz_ind = iz - first_node;
532  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
533  {
534  // Calculate sum of crossflow into channel i from channels j around i
535  unsigned int counter = 0;
536  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
537  {
538  PetscInt row = i_ch + _n_channels * iz_ind;
539  PetscInt col = i_gap + _n_gaps * iz_ind;
540  PetscScalar value = _subchannel_mesh.getCrossflowSign(i_ch, counter);
541  LibmeshPetscCall(MatSetValues(_mc_sumWij_mat, 1, &row, 1, &col, &value, INSERT_VALUES));
542  counter++;
543  }
544  }
545  }
546  LibmeshPetscCall(MatAssemblyBegin(_mc_sumWij_mat, MAT_FINAL_ASSEMBLY));
547  LibmeshPetscCall(MatAssemblyEnd(_mc_sumWij_mat, MAT_FINAL_ASSEMBLY));
548  if (_segregated_bool)
549  {
550  Vec loc_prod;
551  Vec loc_Wij;
552  LibmeshPetscCall(VecDuplicate(_amc_sys_mdot_rhs, &loc_prod));
553  LibmeshPetscCall(VecDuplicate(_Wij_vec, &loc_Wij));
555  loc_Wij, _Wij, first_node, last_node + 1, _n_gaps));
556  LibmeshPetscCall(MatMult(_mc_sumWij_mat, loc_Wij, loc_prod));
557  LibmeshPetscCall(populateSolutionChan<SolutionHandle>(
558  loc_prod, *_SumWij_soln, first_node, last_node, _n_channels));
559  LibmeshPetscCall(VecDestroy(&loc_prod));
560  LibmeshPetscCall(VecDestroy(&loc_Wij));
561  }
562  }
563 }
const bool _segregated_bool
Segregated solve.
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
std::unique_ptr< SolutionHandle > _SumWij_soln
Mat _mc_sumWij_mat
Mass conservation Mass conservation - sum of cross fluxes.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
PetscErrorCode populateVectorFromDense(Vec &x, const T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
libMesh::DenseMatrix< Real > & _Wij
virtual const std::vector< unsigned int > & getChannelGaps(unsigned int i_chan) const =0
Return a vector of gap indices for a given channel index.
virtual const Real & getCrossflowSign(unsigned int i_chan, unsigned int i_local) const =0
Return a signs for the cross flow given a inter-wrapper index and local neighbor index.

◆ computeT()

void InterWrapper1PhaseProblem::computeT ( int  iblock)
protectedvirtualinherited

Computes Temperature per channel for block iblock.

Definition at line 1900 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

1901 {
1902  unsigned int last_node = (iblock + 1) * _block_size;
1903  unsigned int first_node = iblock * _block_size + 1;
1904 
1905  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
1906  {
1907  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
1908  {
1909  auto * node = _subchannel_mesh.getChannelNode(i_ch, iz);
1910  _T_soln->set(node, _fp->T_from_p_h((*_P_soln)(node) + _P_out, (*_h_soln)(node)));
1911  }
1912  }
1913 }
std::unique_ptr< SolutionHandle > _P_soln
const Real & _P_out
Outlet Pressure.
const SinglePhaseFluidProperties * _fp
Solutions handles and link to TH tables properties.
std::unique_ptr< SolutionHandle > _T_soln
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
std::unique_ptr< SolutionHandle > _h_soln

◆ computeWij()

void InterWrapper1PhaseProblem::computeWij ( int  iblock)
protectedvirtualinherited

Computes cross fluxes for block iblock.

Assembling system

Definition at line 1966 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::residualFunction().

1967 {
1968  // Cross flow residual
1969  if (!_implicit_bool)
1970  {
1971  unsigned int last_node = (iblock + 1) * _block_size;
1972  unsigned int first_node = iblock * _block_size;
1973 
1974  const Real & pitch = _subchannel_mesh.getPitch();
1975  for (unsigned int iz = first_node + 1; iz < last_node + 1; iz++)
1976  {
1977  auto dz = _z_grid[iz] - _z_grid[iz - 1];
1978  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
1979  {
1980  auto chans = _subchannel_mesh.getGapChannels(i_gap);
1981  unsigned int i_ch = chans.first;
1982  unsigned int j_ch = chans.second;
1983  auto * node_in_i = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
1984  auto * node_out_i = _subchannel_mesh.getChannelNode(i_ch, iz);
1985  auto * node_in_j = _subchannel_mesh.getChannelNode(j_ch, iz - 1);
1986  auto * node_out_j = _subchannel_mesh.getChannelNode(j_ch, iz);
1987  auto rho_i = (*_rho_soln)(node_in_i);
1988  auto rho_j = (*_rho_soln)(node_in_j);
1989  auto Si = (*_S_flow_soln)(node_in_i);
1990  auto Sj = (*_S_flow_soln)(node_in_j);
1991  auto Sij = dz * _subchannel_mesh.getGapWidth(i_gap);
1992  auto Lij = pitch;
1993  // total local form loss in the ij direction
1994  auto friction_term = _kij * _Wij(i_gap, iz) * std::abs(_Wij(i_gap, iz));
1995  auto DPij = (*_P_soln)(node_in_i) - (*_P_soln)(node_in_j);
1996  // Figure out donor cell density
1997  auto rho_star = 0.0;
1998  if (_Wij(i_gap, iz) > 0.0)
1999  rho_star = rho_i;
2000  else if (_Wij(i_gap, iz) < 0.0)
2001  rho_star = rho_j;
2002  else
2003  rho_star = (rho_i + rho_j) / 2.0;
2004 
2005  auto mass_term_out =
2006  (*_mdot_soln)(node_out_i) / Si / rho_i + (*_mdot_soln)(node_out_j) / Sj / rho_j;
2007  auto mass_term_in =
2008  (*_mdot_soln)(node_in_i) / Si / rho_i + (*_mdot_soln)(node_in_j) / Sj / rho_j;
2009  auto term_out = Sij * rho_star * (Lij / dz) * mass_term_out * _Wij(i_gap, iz);
2010  auto term_in = Sij * rho_star * (Lij / dz) * mass_term_in * _Wij(i_gap, iz - 1);
2011  auto inertia_term = term_out - term_in;
2012  auto pressure_term = 2 * std::pow(Sij, 2.0) * DPij * rho_star;
2013  auto time_term =
2014  _TR * 2.0 * (_Wij(i_gap, iz) - _Wij_old(i_gap, iz)) * Lij * Sij * rho_star / _dt;
2015 
2016  _Wij_residual_matrix(i_gap, iz - 1 - iblock * _block_size) =
2017  time_term + friction_term + inertia_term - pressure_term;
2018  }
2019  }
2020  }
2021  else
2022  {
2023 
2024  // Initializing to zero the elements of the lateral momentum assembly
2025  LibmeshPetscCall(MatZeroEntries(_cmc_time_derivative_mat));
2026  LibmeshPetscCall(MatZeroEntries(_cmc_advective_derivative_mat));
2027  LibmeshPetscCall(MatZeroEntries(_cmc_friction_force_mat));
2028  LibmeshPetscCall(MatZeroEntries(_cmc_pressure_force_mat));
2029  LibmeshPetscCall(VecZeroEntries(_cmc_time_derivative_rhs));
2030  LibmeshPetscCall(VecZeroEntries(_cmc_advective_derivative_rhs));
2031  LibmeshPetscCall(VecZeroEntries(_cmc_friction_force_rhs));
2032  LibmeshPetscCall(VecZeroEntries(_cmc_pressure_force_rhs));
2033  LibmeshPetscCall(MatZeroEntries(_cmc_sys_Wij_mat));
2034  LibmeshPetscCall(VecZeroEntries(_cmc_sys_Wij_rhs));
2035 
2036  unsigned int last_node = (iblock + 1) * _block_size;
2037  unsigned int first_node = iblock * _block_size;
2038 
2039  const Real & pitch = _subchannel_mesh.getPitch();
2040  for (unsigned int iz = first_node + 1; iz < last_node + 1; iz++)
2041  {
2042  auto dz = _z_grid[iz] - _z_grid[iz - 1];
2043  auto iz_ind = iz - first_node - 1;
2044  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
2045  {
2046  auto chans = _subchannel_mesh.getGapChannels(i_gap);
2047  unsigned int i_ch = chans.first;
2048  unsigned int j_ch = chans.second;
2049  auto * node_in_i = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
2050  auto * node_out_i = _subchannel_mesh.getChannelNode(i_ch, iz);
2051  auto * node_in_j = _subchannel_mesh.getChannelNode(j_ch, iz - 1);
2052  auto * node_out_j = _subchannel_mesh.getChannelNode(j_ch, iz);
2053 
2054  // inlet, outlet, and interpolated densities
2055  auto rho_i_in = (*_rho_soln)(node_in_i);
2056  auto rho_i_out = (*_rho_soln)(node_out_i);
2057  auto rho_i_interp = computeInterpolatedValue(rho_i_out, rho_i_in);
2058  auto rho_j_in = (*_rho_soln)(node_in_j);
2059  auto rho_j_out = (*_rho_soln)(node_out_j);
2060  auto rho_j_interp = computeInterpolatedValue(rho_j_out, rho_j_in);
2061 
2062  // inlet, outlet, and interpolated areas
2063  auto S_i_in = (*_S_flow_soln)(node_in_i);
2064  auto S_i_out = (*_S_flow_soln)(node_out_i);
2065  auto S_j_in = (*_S_flow_soln)(node_in_j);
2066  auto S_j_out = (*_S_flow_soln)(node_out_j);
2067 
2068  // Cross-sectional gap area
2069  auto Sij = dz * _subchannel_mesh.getGapWidth(i_gap);
2070  auto Lij = pitch;
2071 
2072  // Figure out donor cell density
2073  auto rho_star = 0.0;
2074  if (_Wij(i_gap, iz) > 0.0)
2075  rho_star = rho_i_interp;
2076  else if (_Wij(i_gap, iz) < 0.0)
2077  rho_star = rho_j_interp;
2078  else
2079  rho_star = (rho_i_interp + rho_j_interp) / 2.0;
2080 
2081  // Assembling time derivative
2082  PetscScalar time_factor = _TR * Lij * Sij * rho_star / _dt;
2083  PetscInt row_td = i_gap + _n_gaps * iz_ind;
2084  PetscInt col_td = i_gap + _n_gaps * iz_ind;
2085  PetscScalar value_td = time_factor;
2086  LibmeshPetscCall(MatSetValues(
2087  _cmc_time_derivative_mat, 1, &row_td, 1, &col_td, &value_td, INSERT_VALUES));
2088  PetscScalar value_td_rhs = time_factor * _Wij_old(i_gap, iz);
2089  LibmeshPetscCall(
2090  VecSetValues(_cmc_time_derivative_rhs, 1, &row_td, &value_td_rhs, INSERT_VALUES));
2091 
2092  // Assembling inertial term
2093  auto alpha = 0.0; // hard code downwind;
2094  auto mass_term_out = (*_mdot_soln)(node_out_i) / S_i_out / rho_i_out +
2095  (*_mdot_soln)(node_out_j) / S_j_out / rho_j_out;
2096  auto mass_term_in = (*_mdot_soln)(node_in_i) / S_i_in / rho_i_in +
2097  (*_mdot_soln)(node_in_j) / S_j_in / rho_j_in;
2098  auto term_out = Sij * rho_star * (Lij / dz) * mass_term_out / 2.0;
2099  auto term_in = Sij * rho_star * (Lij / dz) * mass_term_in / 2.0;
2100  if (iz == first_node + 1)
2101  {
2102  PetscInt row_ad = i_gap + _n_gaps * iz_ind;
2103  PetscScalar value_ad = term_in * alpha * _Wij(i_gap, iz - 1);
2104  LibmeshPetscCall(
2105  VecSetValues(_cmc_advective_derivative_rhs, 1, &row_ad, &value_ad, ADD_VALUES));
2106 
2107  PetscInt col_ad = i_gap + _n_gaps * iz_ind;
2108  value_ad = -1.0 * term_in * (1.0 - alpha) + term_out * alpha;
2109  LibmeshPetscCall(MatSetValues(
2110  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2111 
2112  col_ad = i_gap + _n_gaps * (iz_ind + 1);
2113  value_ad = term_out * (1.0 - alpha);
2114  LibmeshPetscCall(MatSetValues(
2115  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2116  }
2117  else if (iz == last_node)
2118  {
2119  PetscInt row_ad = i_gap + _n_gaps * iz_ind;
2120  PetscInt col_ad = i_gap + _n_gaps * (iz_ind - 1);
2121  PetscScalar value_ad = -1.0 * term_in * alpha;
2122  LibmeshPetscCall(MatSetValues(
2123  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2124 
2125  col_ad = i_gap + _n_gaps * iz_ind;
2126  value_ad = -1.0 * term_in * (1.0 - alpha) + term_out * alpha;
2127  LibmeshPetscCall(MatSetValues(
2128  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2129 
2130  value_ad = -1.0 * term_out * (1.0 - alpha) * _Wij(i_gap, iz);
2131  LibmeshPetscCall(
2132  VecSetValues(_cmc_advective_derivative_rhs, 1, &row_ad, &value_ad, ADD_VALUES));
2133  }
2134  else
2135  {
2136  PetscInt row_ad = i_gap + _n_gaps * iz_ind;
2137  PetscInt col_ad = i_gap + _n_gaps * (iz_ind - 1);
2138  PetscScalar value_ad = -1.0 * term_in * alpha;
2139  LibmeshPetscCall(MatSetValues(
2140  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2141 
2142  col_ad = i_gap + _n_gaps * iz_ind;
2143  value_ad = -1.0 * term_in * (1.0 - alpha) + term_out * alpha;
2144  LibmeshPetscCall(MatSetValues(
2145  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2146 
2147  col_ad = i_gap + _n_gaps * (iz_ind + 1);
2148  value_ad = term_out * (1.0 - alpha);
2149  LibmeshPetscCall(MatSetValues(
2150  _cmc_advective_derivative_mat, 1, &row_ad, 1, &col_ad, &value_ad, INSERT_VALUES));
2151  }
2152  // Assembling friction force
2153  PetscInt row_ff = i_gap + _n_gaps * iz_ind;
2154  PetscInt col_ff = i_gap + _n_gaps * iz_ind;
2155  PetscScalar value_ff = _kij * std::abs(_Wij(i_gap, iz)) / 2.0;
2156  LibmeshPetscCall(MatSetValues(
2157  _cmc_friction_force_mat, 1, &row_ff, 1, &col_ff, &value_ff, INSERT_VALUES));
2158 
2159  // Assembling pressure force
2161 
2163  {
2164  PetscScalar pressure_factor = std::pow(Sij, 2.0) * rho_star;
2165  PetscInt row_pf = i_gap + _n_gaps * iz_ind;
2166  PetscInt col_pf = i_ch + _n_channels * iz_ind;
2167  PetscScalar value_pf = -1.0 * alpha * pressure_factor;
2168  LibmeshPetscCall(
2169  MatSetValues(_cmc_pressure_force_mat, 1, &row_pf, 1, &col_pf, &value_pf, ADD_VALUES));
2170  col_pf = j_ch + _n_channels * iz_ind;
2171  value_pf = alpha * pressure_factor;
2172  LibmeshPetscCall(
2173  MatSetValues(_cmc_pressure_force_mat, 1, &row_pf, 1, &col_pf, &value_pf, ADD_VALUES));
2174 
2175  if (iz == last_node)
2176  {
2177  PetscInt row_pf = i_gap + _n_gaps * iz_ind;
2178  PetscScalar value_pf = (1.0 - alpha) * pressure_factor * (*_P_soln)(node_out_i);
2179  LibmeshPetscCall(
2180  VecSetValues(_cmc_pressure_force_rhs, 1, &row_pf, &value_pf, ADD_VALUES));
2181  value_pf = -1.0 * (1.0 - alpha) * pressure_factor * (*_P_soln)(node_out_j);
2182  LibmeshPetscCall(
2183  VecSetValues(_cmc_pressure_force_rhs, 1, &row_pf, &value_pf, ADD_VALUES));
2184  }
2185  else
2186  {
2187  row_pf = i_gap + _n_gaps * iz_ind;
2188  col_pf = i_ch + _n_channels * (iz_ind + 1);
2189  value_pf = -1.0 * (1.0 - alpha) * pressure_factor;
2190  LibmeshPetscCall(MatSetValues(
2191  _cmc_pressure_force_mat, 1, &row_pf, 1, &col_pf, &value_pf, ADD_VALUES));
2192  col_pf = j_ch + _n_channels * (iz_ind + 1);
2193  value_pf = (1.0 - alpha) * pressure_factor;
2194  LibmeshPetscCall(MatSetValues(
2195  _cmc_pressure_force_mat, 1, &row_pf, 1, &col_pf, &value_pf, ADD_VALUES));
2196  }
2197  }
2198  else
2199  {
2200  PetscScalar pressure_factor = std::pow(Sij, 2.0) * rho_star;
2201  PetscInt row_pf = i_gap + _n_gaps * iz_ind;
2202  PetscInt col_pf = i_ch + _n_channels * iz_ind;
2203  PetscScalar value_pf = -1.0 * pressure_factor;
2204  LibmeshPetscCall(
2205  MatSetValues(_cmc_pressure_force_mat, 1, &row_pf, 1, &col_pf, &value_pf, ADD_VALUES));
2206  col_pf = j_ch + _n_channels * iz_ind;
2207  value_pf = pressure_factor;
2208  LibmeshPetscCall(
2209  MatSetValues(_cmc_pressure_force_mat, 1, &row_pf, 1, &col_pf, &value_pf, ADD_VALUES));
2210  }
2211  }
2212  }
2214  LibmeshPetscCall(MatZeroEntries(_cmc_sys_Wij_mat));
2215  LibmeshPetscCall(VecZeroEntries(_cmc_sys_Wij_rhs));
2216  LibmeshPetscCall(MatAssemblyBegin(_cmc_time_derivative_mat, MAT_FINAL_ASSEMBLY));
2217  LibmeshPetscCall(MatAssemblyEnd(_cmc_time_derivative_mat, MAT_FINAL_ASSEMBLY));
2218  LibmeshPetscCall(MatAssemblyBegin(_cmc_advective_derivative_mat, MAT_FINAL_ASSEMBLY));
2219  LibmeshPetscCall(MatAssemblyEnd(_cmc_advective_derivative_mat, MAT_FINAL_ASSEMBLY));
2220  LibmeshPetscCall(MatAssemblyBegin(_cmc_friction_force_mat, MAT_FINAL_ASSEMBLY));
2221  LibmeshPetscCall(MatAssemblyEnd(_cmc_friction_force_mat, MAT_FINAL_ASSEMBLY));
2222  LibmeshPetscCall(MatAssemblyBegin(_cmc_pressure_force_mat, MAT_FINAL_ASSEMBLY));
2223  LibmeshPetscCall(MatAssemblyEnd(_cmc_pressure_force_mat, MAT_FINAL_ASSEMBLY));
2224  LibmeshPetscCall(MatAssemblyBegin(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2225  LibmeshPetscCall(MatAssemblyEnd(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2226  // Matrix
2227 #if !PETSC_VERSION_LESS_THAN(3, 15, 0)
2228  LibmeshPetscCall(
2229  MatAXPY(_cmc_sys_Wij_mat, 1.0, _cmc_time_derivative_mat, UNKNOWN_NONZERO_PATTERN));
2230  LibmeshPetscCall(MatAssemblyBegin(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2231  LibmeshPetscCall(MatAssemblyEnd(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2232  LibmeshPetscCall(
2233  MatAXPY(_cmc_sys_Wij_mat, 1.0, _cmc_advective_derivative_mat, UNKNOWN_NONZERO_PATTERN));
2234  LibmeshPetscCall(MatAssemblyBegin(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2235  LibmeshPetscCall(MatAssemblyEnd(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2236  LibmeshPetscCall(
2237  MatAXPY(_cmc_sys_Wij_mat, 1.0, _cmc_friction_force_mat, UNKNOWN_NONZERO_PATTERN));
2238 #else
2239  LibmeshPetscCall(
2240  MatAXPY(_cmc_sys_Wij_mat, 1.0, _cmc_time_derivative_mat, DIFFERENT_NONZERO_PATTERN));
2241  LibmeshPetscCall(MatAssemblyBegin(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2242  LibmeshPetscCall(MatAssemblyEnd(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2243  LibmeshPetscCall(
2244  MatAXPY(_cmc_sys_Wij_mat, 1.0, _cmc_advective_derivative_mat, DIFFERENT_NONZERO_PATTERN));
2245  LibmeshPetscCall(MatAssemblyBegin(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2246  LibmeshPetscCall(MatAssemblyEnd(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2247  LibmeshPetscCall(
2248  MatAXPY(_cmc_sys_Wij_mat, 1.0, _cmc_friction_force_mat, DIFFERENT_NONZERO_PATTERN));
2249 #endif
2250  LibmeshPetscCall(MatAssemblyBegin(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2251  LibmeshPetscCall(MatAssemblyEnd(_cmc_sys_Wij_mat, MAT_FINAL_ASSEMBLY));
2252  _console << "Block: " << iblock << " - Cross flow system matrix assembled" << std::endl;
2253  _console << "Block: " << iblock << " - Cross flow pressure force matrix assembled" << std::endl;
2254  // RHS
2255  LibmeshPetscCall(VecAXPY(_cmc_sys_Wij_rhs, 1.0, _cmc_time_derivative_rhs));
2256  LibmeshPetscCall(VecAXPY(_cmc_sys_Wij_rhs, 1.0, _cmc_advective_derivative_rhs));
2257  LibmeshPetscCall(VecAXPY(_cmc_sys_Wij_rhs, 1.0, _cmc_friction_force_rhs));
2258 
2259  if (_segregated_bool)
2260  {
2261  // Assembly the matrix system
2262  Vec sol_holder_P;
2263  LibmeshPetscCall(createPetscVector(sol_holder_P, _block_size * _n_gaps));
2264  Vec sol_holder_W;
2265  LibmeshPetscCall(createPetscVector(sol_holder_W, _block_size * _n_gaps));
2266  Vec loc_holder_Wij;
2267  LibmeshPetscCall(createPetscVector(loc_holder_Wij, _block_size * _n_gaps));
2268  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2269  _prodp, *_P_soln, iblock * _block_size, (iblock + 1) * _block_size - 1, _n_channels));
2271  loc_holder_Wij, _Wij, first_node, last_node, _n_gaps));
2272  LibmeshPetscCall(MatMult(_cmc_sys_Wij_mat, _Wij_vec, sol_holder_W));
2273  LibmeshPetscCall(VecAXPY(sol_holder_W, -1.0, _cmc_sys_Wij_rhs));
2274  LibmeshPetscCall(MatMult(_cmc_pressure_force_mat, _prodp, sol_holder_P));
2275  LibmeshPetscCall(VecAXPY(sol_holder_P, -1.0, _cmc_pressure_force_rhs));
2276  LibmeshPetscCall(VecAXPY(sol_holder_W, 1.0, sol_holder_P));
2277  PetscScalar * xx;
2278  LibmeshPetscCall(VecGetArray(sol_holder_W, &xx));
2279  for (unsigned int iz = first_node + 1; iz < last_node + 1; iz++)
2280  {
2281  auto iz_ind = iz - first_node - 1;
2282  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
2283  {
2284  _Wij_residual_matrix(i_gap, iz - 1 - iblock * _block_size) = xx[iz_ind * _n_gaps + i_gap];
2285  }
2286  }
2287  LibmeshPetscCall(VecDestroy(&sol_holder_P));
2288  LibmeshPetscCall(VecDestroy(&sol_holder_W));
2289  LibmeshPetscCall(VecDestroy(&loc_holder_Wij));
2290  }
2291  }
2292 }
virtual PetscScalar computeInterpolatedValue(PetscScalar topValue, PetscScalar botValue, PetscScalar Peclet=0.0)
Mat _cmc_friction_force_mat
Cross momentum conservation - friction force.
const bool _segregated_bool
Segregated solve.
Real _TR
Flag that activates or deactivates the transient parts of the equations we solve by multiplication...
std::unique_ptr< SolutionHandle > _P_soln
virtual const std::pair< unsigned int, unsigned int > & getGapChannels(unsigned int i_gap) const =0
Return a pair of inter-wrapper indices for a given gap index.
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
Mat _cmc_pressure_force_mat
Cross momentum conservation - pressure force.
virtual PetscErrorCode createPetscVector(Vec &v, PetscInt n)
Petsc Functions.
std::unique_ptr< SolutionHandle > _mdot_soln
std::vector< Real > _z_grid
axial location of nodes
const bool _staggered_pressure_bool
Flag to define the usage of staggered or collocated pressure.
static const std::string pitch
const InterWrapperMesh & _subchannel_mesh
virtual PetscScalar computeInterpolationCoefficients(PetscScalar Peclet=0.0)
Functions that computes the interpolation scheme given the Peclet number.
Mat _cmc_advective_derivative_mat
Cross momentum conservation - advective (Eulerian) derivative.
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
virtual Real getGapWidth(unsigned int gap_index) const =0
Return gap width for a given gap index.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string alpha
Definition: NS.h:134
Mat _cmc_sys_Wij_mat
Lateral momentum system matrix.
PetscErrorCode populateVectorFromDense(Vec &x, const T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
libMesh::DenseMatrix< Real > _Wij_old
const ConsoleStream _console
libMesh::DenseMatrix< Real > & _Wij
MooseUnits pow(const MooseUnits &, int)
libMesh::DenseMatrix< Real > _Wij_residual_matrix
virtual const Real & getPitch() const
Return the pitch between 2 inter-wrappers.
Mat _cmc_time_derivative_mat
Cross momentum Cross momentum conservation - time derivative.

◆ computeWijFromSolve()

void InterWrapper1PhaseProblem::computeWijFromSolve ( int  iblock)
protectedvirtualinherited

Computes diversion crossflow per gap for block with index iblock Block is a partition of the whole domain.

Definition at line 468 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

469 {
470  unsigned int last_node = (iblock + 1) * _block_size;
471  unsigned int first_node = iblock * _block_size + 1;
472  // Initial guess, port crossflow of block (iblock) into a vector that will act as my initial guess
473  libMesh::DenseVector<Real> solution_seed(_n_gaps * _block_size, 0.0);
474  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
475  {
476  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
477  {
478  int i = _n_gaps * (iz - first_node) + i_gap; // column wise transfer
479  solution_seed(i) = _Wij(i_gap, iz);
480  }
481  }
482 
483  // Solving the combined lateral momentum equation for Wij using a PETSc solver and update vector
484  // root
486  LibmeshPetscCall(petscSnesSolver(iblock, solution_seed, root));
487 
488  // Assign the solution to the cross-flow matrix
489  int i = 0;
490  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
491  {
492  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
493  {
494  _Wij(i_gap, iz) = root(i);
495  i++;
496  }
497  }
498 }
virtual PetscErrorCode petscSnesSolver(int iblock, const libMesh::DenseVector< Real > &solution, libMesh::DenseVector< Real > &root)
Computes solution of nonlinear equation using snes.
Real root(std::function< Real(Real)> const &f, Real x1, Real x2, Real tol=1.0e-12)
Finds the root of a function using Brent&#39;s method.
Definition: BrentsMethod.C:66
libMesh::DenseMatrix< Real > & _Wij

◆ computeWijPrime()

void InterWrapper1PhaseProblem::computeWijPrime ( int  iblock)
protectedvirtualinherited

Computes turbulent crossflow per gap for block iblock.

Update turbulent crossflow

Definition at line 678 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::residualFunction().

679 {
680  unsigned int last_node = (iblock + 1) * _block_size;
681  unsigned int first_node = iblock * _block_size + 1;
682 
683  if (!_implicit_bool)
684  {
685  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
686  {
687  auto dz = _z_grid[iz] - _z_grid[iz - 1];
688  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
689  {
690  auto chans = _subchannel_mesh.getGapChannels(i_gap);
691  unsigned int i_ch = chans.first;
692  unsigned int j_ch = chans.second;
693  auto * node_in_i = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
694  auto * node_out_i = _subchannel_mesh.getChannelNode(i_ch, iz);
695  auto * node_in_j = _subchannel_mesh.getChannelNode(j_ch, iz - 1);
696  auto * node_out_j = _subchannel_mesh.getChannelNode(j_ch, iz);
697  auto Si_in = (*_S_flow_soln)(node_in_i);
698  auto Sj_in = (*_S_flow_soln)(node_in_j);
699  auto Si_out = (*_S_flow_soln)(node_out_i);
700  auto Sj_out = (*_S_flow_soln)(node_out_j);
701  // crossflow area between channels i,j (dz*gap_width)
702  auto Sij = dz * _subchannel_mesh.getGapWidth(i_gap);
703  // Calculation of Turbulent Crossflow
704  _WijPrime(i_gap, iz) =
705  _beta * 0.5 *
706  (((*_mdot_soln)(node_in_i) + (*_mdot_soln)(node_in_j)) / (Si_in + Sj_in) +
707  ((*_mdot_soln)(node_out_i) + (*_mdot_soln)(node_out_j)) / (Si_out + Sj_out)) *
708  Sij;
709  }
710  }
711  }
712  else
713  {
714  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
715  {
716  auto dz = _z_grid[iz] - _z_grid[iz - 1];
717  auto iz_ind = iz - first_node;
718  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
719  {
720  auto chans = _subchannel_mesh.getGapChannels(i_gap);
721  unsigned int i_ch = chans.first;
722  unsigned int j_ch = chans.second;
723  auto * node_in_i = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
724  auto * node_out_i = _subchannel_mesh.getChannelNode(i_ch, iz);
725  auto * node_in_j = _subchannel_mesh.getChannelNode(j_ch, iz - 1);
726  auto * node_out_j = _subchannel_mesh.getChannelNode(j_ch, iz);
727  auto Si_in = (*_S_flow_soln)(node_in_i);
728  auto Sj_in = (*_S_flow_soln)(node_in_j);
729  auto Si_out = (*_S_flow_soln)(node_out_i);
730  auto Sj_out = (*_S_flow_soln)(node_out_j);
731  // crossflow area between channels i,j (dz*gap_width)
732  auto Sij = dz * _subchannel_mesh.getGapWidth(i_gap);
733 
734  // Base value - I don't want to write it every time
735  PetscScalar base_value = _beta * 0.5 * Sij;
736 
737  // Bottom values
738  if (iz == first_node)
739  {
740  PetscScalar value_tl = -1.0 * base_value / (Si_in + Sj_in) *
741  ((*_mdot_soln)(node_in_i) + (*_mdot_soln)(node_in_j));
742  PetscInt row = i_gap + _n_gaps * iz_ind;
743  LibmeshPetscCall(
744  VecSetValues(_amc_turbulent_cross_flows_rhs, 1, &row, &value_tl, INSERT_VALUES));
745  }
746  else
747  {
748  PetscScalar value_tl = base_value / (Si_in + Sj_in);
749  PetscInt row = i_gap + _n_gaps * iz_ind;
750 
751  PetscInt col_ich = i_ch + _n_channels * (iz_ind - 1);
752  LibmeshPetscCall(MatSetValues(
753  _amc_turbulent_cross_flows_mat, 1, &row, 1, &col_ich, &value_tl, INSERT_VALUES));
754 
755  PetscInt col_jch = j_ch + _n_channels * (iz_ind - 1);
756  LibmeshPetscCall(MatSetValues(
757  _amc_turbulent_cross_flows_mat, 1, &row, 1, &col_jch, &value_tl, INSERT_VALUES));
758  }
759 
760  // Top values
761  PetscScalar value_bl = base_value / (Si_out + Sj_out);
762  PetscInt row = i_gap + _n_gaps * iz_ind;
763 
764  PetscInt col_ich = i_ch + _n_channels * iz_ind;
765  LibmeshPetscCall(MatSetValues(
766  _amc_turbulent_cross_flows_mat, 1, &row, 1, &col_ich, &value_bl, INSERT_VALUES));
767 
768  PetscInt col_jch = j_ch + _n_channels * iz_ind;
769  LibmeshPetscCall(MatSetValues(
770  _amc_turbulent_cross_flows_mat, 1, &row, 1, &col_jch, &value_bl, INSERT_VALUES));
771  }
772  }
773  LibmeshPetscCall(MatAssemblyBegin(_amc_turbulent_cross_flows_mat, MAT_FINAL_ASSEMBLY));
774  LibmeshPetscCall(MatAssemblyEnd(_amc_turbulent_cross_flows_mat, MAT_FINAL_ASSEMBLY));
775 
777  Vec loc_prod;
778  Vec loc_Wij;
779  LibmeshPetscCall(VecDuplicate(_amc_sys_mdot_rhs, &loc_prod));
780  LibmeshPetscCall(VecDuplicate(_Wij_vec, &loc_Wij));
781  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
782  loc_prod, *_mdot_soln, first_node, last_node, _n_channels));
783  LibmeshPetscCall(MatMult(_amc_turbulent_cross_flows_mat, loc_prod, loc_Wij));
784  LibmeshPetscCall(VecAXPY(loc_Wij, -1.0, _amc_turbulent_cross_flows_rhs));
786  loc_Wij, _WijPrime, first_node, last_node, _n_gaps));
787  LibmeshPetscCall(VecDestroy(&loc_prod));
788  LibmeshPetscCall(VecDestroy(&loc_Wij));
789  }
790 }
Mat _amc_turbulent_cross_flows_mat
Mass conservation - density time derivative No implicit matrix.
virtual const std::pair< unsigned int, unsigned int > & getGapChannels(unsigned int i_gap) const =0
Return a pair of inter-wrapper indices for a given gap index.
const bool _implicit_bool
Flag to define the usage of a implicit or explicit solution.
std::unique_ptr< SolutionHandle > _mdot_soln
PetscErrorCode populateDenseFromVector(const Vec &x, T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
const Real & _beta
Thermal diffusion coefficient used in turbulent crossflow.
std::vector< Real > _z_grid
axial location of nodes
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
virtual Real getGapWidth(unsigned int gap_index) const =0
Return gap width for a given gap index.
libMesh::DenseMatrix< Real > _WijPrime

◆ createPetscMatrix()

PetscErrorCode InterWrapper1PhaseProblem::createPetscMatrix ( Mat &  M,
PetscInt  n,
PetscInt  m 
)
protectedvirtualinherited

Definition at line 364 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

365 {
367  LibmeshPetscCall(MatCreate(PETSC_COMM_WORLD, &M));
368  LibmeshPetscCall(MatSetSizes(M, PETSC_DECIDE, PETSC_DECIDE, n, m));
369  LibmeshPetscCall(MatSetFromOptions(M));
370  LibmeshPetscCall(MatSetUp(M));
371  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
372 }
PetscFunctionBegin
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ createPetscVector()

PetscErrorCode InterWrapper1PhaseProblem::createPetscVector ( Vec &  v,
PetscInt  n 
)
protectedvirtualinherited

Petsc Functions.

Definition at line 352 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::computeWij(), InterWrapper1PhaseProblem::implicitPetscSolve(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

353 {
355  LibmeshPetscCall(VecCreate(PETSC_COMM_WORLD, &v));
356  LibmeshPetscCall(PetscObjectSetName((PetscObject)v, "Solution"));
357  LibmeshPetscCall(VecSetSizes(v, PETSC_DECIDE, n));
358  LibmeshPetscCall(VecSetFromOptions(v));
359  LibmeshPetscCall(VecZeroEntries(v));
360  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
361 }
PetscFunctionBegin
static const std::string v
Definition: NS.h:84
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ externalSolve()

void InterWrapper1PhaseProblem::externalSolve ( )
overridevirtualinherited

Implements ExternalProblem.

Reimplemented in TriInterWrapper1PhaseProblem.

Definition at line 2927 of file InterWrapper1PhaseProblem.C.

2928 {
2929  _console << "Executing subchannel solver\n";
2931  _console << "Solution initialized" << std::endl;
2932  Real P_error = 1.0;
2933  unsigned int P_it = 0;
2934  unsigned int P_it_max;
2935 
2936  if (_segregated_bool)
2937  P_it_max = 20 * _n_blocks;
2938  else
2939  P_it_max = 100;
2940 
2941  if ((_n_blocks == 1) && (_segregated_bool))
2942  P_it_max = 5;
2943 
2944  while ((P_error > _P_tol && P_it < P_it_max))
2945  {
2946  P_it += 1;
2947  if (P_it == P_it_max && _n_blocks != 1)
2948  {
2949  _console << "Reached maximum number of axial pressure iterations" << std::endl;
2950  _converged = false;
2951  }
2952  _console << "Solving Outer Iteration : " << P_it << std::endl;
2953  auto P_L2norm_old_axial = _P_soln->L2norm();
2954  for (unsigned int iblock = 0; iblock < _n_blocks; iblock++)
2955  {
2956  int last_level = (iblock + 1) * _block_size;
2957  int first_level = iblock * _block_size + 1;
2958  Real T_block_error = 1.0;
2959  auto T_it = 0;
2960  _console << "Solving Block: " << iblock << " From first level: " << first_level
2961  << " to last level: " << last_level << std::endl;
2962  while (T_block_error > _T_tol && T_it < _T_maxit)
2963  {
2964  T_it += 1;
2965  if (T_it == _T_maxit)
2966  {
2967  _console << "Reached maximum number of temperature iterations for block: " << iblock
2968  << std::endl;
2969  _converged = false;
2970  }
2971  auto T_L2norm_old_block = _T_soln->L2norm();
2972 
2973  if (_segregated_bool)
2974  {
2975  computeWijFromSolve(iblock);
2976 
2977  if (_compute_power)
2978  {
2979  computeh(iblock);
2980  _console << "Done with h solve" << std::endl;
2981  computeT(iblock);
2982  _console << "Done with T solve" << std::endl;
2983  }
2984  }
2985  else
2986  {
2987  LibmeshPetscCall(implicitPetscSolve(iblock));
2988  computeWijPrime(iblock);
2989 
2990  _console << "Done with main solve." << std::endl;
2992  {
2993  // Enthalpy is already solved from the monolithic solve
2994  computeT(iblock);
2995  }
2996  else
2997  {
2998  if (_compute_power)
2999  {
3000  computeh(iblock);
3001  computeT(iblock);
3002  }
3003  _console << "Done with thermal solve." << std::endl;
3004  }
3005  }
3006 
3007  if (_compute_density)
3008  computeRho(iblock);
3009 
3010  if (_compute_viscosity)
3011  computeMu(iblock);
3012 
3013  _console << "Done updating thermophysical properties." << std::endl;
3014 
3015  // We must do a global assembly to make sure data is parallel consistent before we do things
3016  // like compute L2 norms
3017  _aux->solution().close();
3018 
3019  auto T_L2norm_new = _T_soln->L2norm();
3020  T_block_error =
3021  std::abs((T_L2norm_new - T_L2norm_old_block) / (T_L2norm_old_block + 1E-14));
3022  _console << "T_block_error: " << T_block_error << std::endl;
3023  }
3024  }
3025  auto P_L2norm_new_axial = _P_soln->L2norm();
3026  P_error =
3027  std::abs((P_L2norm_new_axial - P_L2norm_old_axial) / (P_L2norm_old_axial + _P_out + 1E-14));
3028  _console << "P_error :" << P_error << std::endl;
3029  }
3030  // update old crossflow matrix
3031  _Wij_old = _Wij;
3032 
3033  auto power_in = 0.0;
3034  auto power_out = 0.0;
3035  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
3036  {
3037  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, 0);
3038  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, _n_cells);
3039  power_in += (*_mdot_soln)(node_in) * (*_h_soln)(node_in);
3040  power_out += (*_mdot_soln)(node_out) * (*_h_soln)(node_out);
3041  }
3042  _console << "Power added to coolant is: " << power_out - power_in << " Watt" << std::endl;
3043  if (_pin_mesh_exist)
3044  {
3045  _console << "Commencing calculation of Pin surface temperature \n";
3046  for (unsigned int i_pin = 0; i_pin < _n_assemblies; i_pin++)
3047  {
3048  for (unsigned int iz = 0; iz < _n_cells + 1; ++iz)
3049  {
3050  auto * pin_node = _subchannel_mesh.getPinNode(i_pin, iz);
3051  Real sumTemp = 0.0;
3052  // Calculate sum of pin surface temperatures that the channels around the pin see
3053  for (auto i_ch : _subchannel_mesh.getPinChannels(i_pin))
3054  {
3055  auto * node = _subchannel_mesh.getChannelNode(i_ch, iz);
3056 
3057  auto mu = (*_mu_soln)(node);
3058  auto S = (*_S_flow_soln)(node);
3059  auto w_perim = (*_w_perim_soln)(node);
3060  auto Dh_i = 4.0 * S / w_perim;
3061  auto Re = (((*_mdot_soln)(node) / S) * Dh_i / mu);
3062 
3063  auto k = _fp->k_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node));
3064  auto cp = _fp->cp_from_p_T((*_P_soln)(node) + _P_out, (*_T_soln)(node));
3065  auto Pr = (*_mu_soln)(node)*cp / k;
3066 
3067  auto Nu = 0.023 * std::pow(Re, 0.8) * std::pow(Pr, 0.4);
3068  auto hw = Nu * k / Dh_i;
3069 
3070  auto perimeter = 2.0 * (_subchannel_mesh.getSideX() + _subchannel_mesh.getSideY());
3071  sumTemp += (*_q_prime_soln)(pin_node) / (perimeter * hw) + (*_T_soln)(node);
3072  }
3073  _Tpin_soln->set(pin_node, 0.25 * sumTemp);
3074  }
3075  }
3076  }
3077  _aux->solution().close();
3078  _aux->update();
3079 }
const bool _compute_density
Flag that activates or deactivates the calculation of density.
const bool _segregated_bool
Segregated solve.
std::unique_ptr< SolutionHandle > _P_soln
const int & _T_maxit
Maximum iterations for the inner temperature loop.
const Real & _P_tol
Convergence tolerance for the pressure loop in external solve.
virtual void computeh(int iblock)
Computes Enthalpy per channel for block iblock.
std::unique_ptr< SolutionHandle > _Tpin_soln
virtual const Real & getSideX() const
Return side lengths of the assembly.
virtual void computeRho(int iblock)
Computes Density per channel for block iblock.
virtual void computeWijPrime(int iblock)
Computes turbulent crossflow per gap for block iblock.
virtual void computeT(int iblock)
Computes Temperature per channel for block iblock.
const Real & _P_out
Outlet Pressure.
virtual PetscErrorCode implicitPetscSolve(int iblock)
Computes implicit solve using PetSc.
const SinglePhaseFluidProperties * _fp
Solutions handles and link to TH tables properties.
virtual const std::vector< unsigned int > & getPinChannels(unsigned int i_pin) const =0
Return a vector of channel indices for a given Pin index.
static const std::string cp
Definition: NS.h:121
virtual const Real & getSideY() const
std::unique_ptr< SolutionHandle > _T_soln
static const std::string S
Definition: NS.h:163
virtual Node * getPinNode(unsigned int i_pin, unsigned iz) const =0
Get the pin mesh node for a given pin index and elevation index.
static const std::string mu
Definition: NS.h:123
std::shared_ptr< AuxiliarySystem > _aux
const bool _monolithic_thermal_bool
Thermal monolithic bool.
unsigned int _n_blocks
number of axial blocks
const InterWrapperMesh & _subchannel_mesh
const bool _compute_viscosity
Flag that activates or deactivates the calculation of viscosity.
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
virtual void computeWijFromSolve(int iblock)
Computes diversion crossflow per gap for block with index iblock Block is a partition of the whole do...
const Real & _T_tol
Convergence tolerance for the temperature loop in external solve.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void initializeSolution()
Function to initialize the solution fields.
const bool _pin_mesh_exist
Flag that informs if there is a pin mesh or not.
virtual void computeMu(int iblock)
Computes Viscosity per channel for block iblock.
std::unique_ptr< SolutionHandle > _h_soln
libMesh::DenseMatrix< Real > _Wij_old
const bool _compute_power
Flag that informs if we need to solve the Enthalpy/Temperature equations or not.
const ConsoleStream _console
libMesh::DenseMatrix< Real > & _Wij
bool _converged
Variable that informs whether we exited external solve with a converged solution or not...
MooseUnits pow(const MooseUnits &, int)
static const std::string k
Definition: NS.h:130

◆ implicitPetscSolve()

PetscErrorCode InterWrapper1PhaseProblem::implicitPetscSolve ( int  iblock)
protectedvirtualinherited

Computes implicit solve using PetSc.

Initializing flags

Assembling matrices

Setting up linear solver

Solving

Destroying solver elements

Recovering the solutions

Assigning the solutions to arrays

Populating Mass flow

Populating Pressure

Populating Crossflow

Populating Enthalpy

Populating sum_Wij

Destroying arrays

Definition at line 2393 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve().

2394 {
2395  Vec b_nest, x_nest; /* approx solution, RHS, exact solution */
2396  Mat A_nest; /* linear system matrix */
2397  KSP ksp; /* linear solver context */
2398  PC pc; /* preconditioner context */
2399 
2401  PetscInt Q = _monolithic_thermal_bool ? 4 : 3;
2402  std::vector<Mat> mat_array(Q * Q);
2403  std::vector<Vec> vec_array(Q);
2404 
2406  bool _axial_mass_flow_tight_coupling = true;
2407  bool _pressure_axial_momentum_tight_coupling = true;
2408  bool _pressure_cross_momentum_tight_coupling = true;
2409  unsigned int first_node = iblock * _block_size + 1;
2410  unsigned int last_node = (iblock + 1) * _block_size;
2411 
2413  // Computing sum of crossflows with previous iteration
2414  computeSumWij(iblock);
2415  // Assembling axial flux matrix
2416  computeMdot(iblock);
2417  // Computing turbulent crossflow with previous step axial mass flows
2418  computeWijPrime(iblock);
2419  // Assembling for Pressure Drop matrix
2420  computeDP(iblock);
2421  // Assembling pressure matrix
2422  computeP(iblock);
2423  // Assembling cross fluxes matrix
2424  computeWij(iblock);
2425  // If monolithic solve - Assembling enthalpy matrix
2427  computeh(iblock);
2428 
2429  _console << "Starting nested system." << std::endl;
2430  _console << "Number of simultaneous variables: " << Q << std::endl;
2431  // Mass conservation
2432  PetscInt field_num = 0;
2433  LibmeshPetscCall(
2434  MatDuplicate(_mc_axial_convection_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 0]));
2435  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 0], MAT_FINAL_ASSEMBLY));
2436  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 0], MAT_FINAL_ASSEMBLY));
2437  mat_array[Q * field_num + 1] = NULL;
2438  if (_axial_mass_flow_tight_coupling)
2439  {
2440  LibmeshPetscCall(MatDuplicate(_mc_sumWij_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 2]));
2441  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 2], MAT_FINAL_ASSEMBLY));
2442  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 2], MAT_FINAL_ASSEMBLY));
2443  }
2444  else
2445  {
2446  mat_array[Q * field_num + 2] = NULL;
2447  }
2449  {
2450  mat_array[Q * field_num + 3] = NULL;
2451  }
2452  LibmeshPetscCall(VecDuplicate(_mc_axial_convection_rhs, &vec_array[field_num]));
2453  LibmeshPetscCall(VecCopy(_mc_axial_convection_rhs, vec_array[field_num]));
2454  if (!_axial_mass_flow_tight_coupling)
2455  {
2456  Vec sumWij_loc;
2457  LibmeshPetscCall(VecDuplicate(_mc_axial_convection_rhs, &sumWij_loc));
2458  LibmeshPetscCall(VecSet(sumWij_loc, 0.0));
2459  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
2460  {
2461  auto iz_ind = iz - first_node;
2462  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
2463  {
2464  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
2465 
2466  PetscScalar value_vec_2 = -1.0 * (*_SumWij_soln)(node_out);
2467  PetscInt row_vec_2 = i_ch + _n_channels * iz_ind;
2468  LibmeshPetscCall(VecSetValues(sumWij_loc, 1, &row_vec_2, &value_vec_2, ADD_VALUES));
2469  }
2470  }
2471  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, sumWij_loc));
2472  LibmeshPetscCall(VecDestroy(&sumWij_loc));
2473  }
2474 
2475  _console << "Mass ok." << std::endl;
2476  // Axial momentum conservation
2477  field_num = 1;
2478  if (_pressure_axial_momentum_tight_coupling)
2479  {
2480  LibmeshPetscCall(
2481  MatDuplicate(_amc_sys_mdot_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 0]));
2482  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 0], MAT_FINAL_ASSEMBLY));
2483  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 0], MAT_FINAL_ASSEMBLY));
2484  }
2485  else
2486  {
2487  mat_array[Q * field_num + 0] = NULL;
2488  }
2489  LibmeshPetscCall(
2490  MatDuplicate(_amc_pressure_force_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 1]));
2491  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 1], MAT_FINAL_ASSEMBLY));
2492  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 1], MAT_FINAL_ASSEMBLY));
2493  mat_array[Q * field_num + 2] = NULL;
2495  {
2496  mat_array[Q * field_num + 3] = NULL;
2497  }
2498  LibmeshPetscCall(VecDuplicate(_amc_pressure_force_rhs, &vec_array[field_num]));
2499  LibmeshPetscCall(VecCopy(_amc_pressure_force_rhs, vec_array[field_num]));
2500  if (_pressure_axial_momentum_tight_coupling)
2501  {
2502  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, _amc_sys_mdot_rhs));
2503  }
2504  else
2505  {
2506  unsigned int last_node = (iblock + 1) * _block_size;
2507  unsigned int first_node = iblock * _block_size + 1;
2508  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2509  _prod, *_mdot_soln, first_node, last_node, _n_channels));
2510  Vec ls;
2511  LibmeshPetscCall(VecDuplicate(_amc_sys_mdot_rhs, &ls));
2512  LibmeshPetscCall(MatMult(_amc_sys_mdot_mat, _prod, ls));
2513  LibmeshPetscCall(VecAXPY(ls, -1.0, _amc_sys_mdot_rhs));
2514  LibmeshPetscCall(VecAXPY(vec_array[field_num], -1.0, ls));
2515  LibmeshPetscCall(VecDestroy(&ls));
2516  }
2517 
2518  _console << "Lin mom OK." << std::endl;
2519 
2520  // Cross momentum conservation
2521  field_num = 2;
2522  mat_array[Q * field_num + 0] = NULL;
2523  if (_pressure_cross_momentum_tight_coupling)
2524  {
2525  LibmeshPetscCall(
2526  MatDuplicate(_cmc_pressure_force_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 1]));
2527  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 1], MAT_FINAL_ASSEMBLY));
2528  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 1], MAT_FINAL_ASSEMBLY));
2529  }
2530  else
2531  {
2532  mat_array[Q * field_num + 1] = NULL;
2533  }
2534 
2535  LibmeshPetscCall(MatDuplicate(_cmc_sys_Wij_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 2]));
2536  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 2], MAT_FINAL_ASSEMBLY));
2537  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 2], MAT_FINAL_ASSEMBLY));
2539  {
2540  mat_array[Q * field_num + 3] = NULL;
2541  }
2542 
2543  LibmeshPetscCall(VecDuplicate(_cmc_sys_Wij_rhs, &vec_array[field_num]));
2544  LibmeshPetscCall(VecCopy(_cmc_sys_Wij_rhs, vec_array[field_num]));
2545  if (_pressure_cross_momentum_tight_coupling)
2546  {
2547  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, _cmc_pressure_force_rhs));
2548  }
2549  else
2550  {
2551  Vec sol_holder_P;
2552  LibmeshPetscCall(createPetscVector(sol_holder_P, _block_size * _n_gaps));
2553  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2554  _prodp, *_P_soln, iblock * _block_size, (iblock + 1) * _block_size - 1, _n_channels));
2555 
2556  LibmeshPetscCall(MatMult(_cmc_pressure_force_mat, _prodp, sol_holder_P));
2557  LibmeshPetscCall(VecAXPY(sol_holder_P, -1.0, _cmc_pressure_force_rhs));
2558  LibmeshPetscCall(VecScale(sol_holder_P, 1.0));
2559  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, sol_holder_P));
2560  LibmeshPetscCall(VecDestroy(&sol_holder_P));
2561  }
2562 
2563  _console << "Cross mom ok." << std::endl;
2564 
2565  // Energy conservation
2567  {
2568  field_num = 3;
2569  mat_array[Q * field_num + 0] = NULL;
2570  mat_array[Q * field_num + 1] = NULL;
2571  mat_array[Q * field_num + 2] = NULL;
2572  LibmeshPetscCall(MatDuplicate(_hc_sys_h_mat, MAT_COPY_VALUES, &mat_array[Q * field_num + 3]));
2573  LibmeshPetscCall(MatAssemblyBegin(mat_array[Q * field_num + 3], MAT_FINAL_ASSEMBLY));
2574  LibmeshPetscCall(MatAssemblyEnd(mat_array[Q * field_num + 3], MAT_FINAL_ASSEMBLY));
2575  LibmeshPetscCall(VecDuplicate(_hc_sys_h_rhs, &vec_array[field_num]));
2576  LibmeshPetscCall(VecCopy(_hc_sys_h_rhs, vec_array[field_num]));
2577  }
2578  _console << "Energy ok." << std::endl;
2579 
2580  // Relaxing linear system
2581  // Weaker relaxation
2582  if (true)
2583  {
2584  // Estimating cross-flow resistances to achieve realizable solves
2585  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2586  _prod, *_mdot_soln, first_node, last_node, _n_channels));
2587  Vec mdot_estimate;
2588  LibmeshPetscCall(createPetscVector(mdot_estimate, _block_size * _n_channels));
2589  Vec pmat_diag;
2590  LibmeshPetscCall(createPetscVector(pmat_diag, _block_size * _n_channels));
2591  Vec p_estimate;
2592  LibmeshPetscCall(createPetscVector(p_estimate, _block_size * _n_channels));
2593  Vec unity_vec;
2594  LibmeshPetscCall(createPetscVector(unity_vec, _block_size * _n_channels));
2595  LibmeshPetscCall(VecSet(unity_vec, 1.0));
2596  Vec sol_holder_P;
2597  LibmeshPetscCall(createPetscVector(sol_holder_P, _block_size * _n_gaps));
2598  Vec diag_Wij_loc;
2599  LibmeshPetscCall(createPetscVector(diag_Wij_loc, _block_size * _n_gaps));
2600  Vec Wij_estimate;
2601  LibmeshPetscCall(createPetscVector(Wij_estimate, _block_size * _n_gaps));
2602  Vec unity_vec_Wij;
2603  LibmeshPetscCall(createPetscVector(unity_vec_Wij, _block_size * _n_gaps));
2604  LibmeshPetscCall(VecSet(unity_vec_Wij, 1.0));
2605  Vec _Wij_loc_vec;
2606  LibmeshPetscCall(createPetscVector(_Wij_loc_vec, _block_size * _n_gaps));
2607  Vec _Wij_old_loc_vec;
2608  LibmeshPetscCall(createPetscVector(_Wij_old_loc_vec, _block_size * _n_gaps));
2609  LibmeshPetscCall(MatMult(mat_array[Q], _prod, mdot_estimate));
2610  LibmeshPetscCall(MatGetDiagonal(mat_array[Q + 1], pmat_diag));
2611  LibmeshPetscCall(VecAXPY(pmat_diag, 1e-10, unity_vec));
2612  LibmeshPetscCall(VecPointwiseDivide(p_estimate, mdot_estimate, pmat_diag));
2613  LibmeshPetscCall(MatMult(mat_array[2 * Q + 1], p_estimate, sol_holder_P));
2614  LibmeshPetscCall(VecAXPY(sol_holder_P, -1.0, _cmc_pressure_force_rhs));
2615  LibmeshPetscCall(MatGetDiagonal(mat_array[2 * Q + 2], diag_Wij_loc));
2616  LibmeshPetscCall(VecAXPY(diag_Wij_loc, 1e-10, unity_vec_Wij));
2617  LibmeshPetscCall(VecPointwiseDivide(Wij_estimate, sol_holder_P, diag_Wij_loc));
2618  Vec sumWij_loc;
2619  LibmeshPetscCall(createPetscVector(sumWij_loc, _block_size * _n_channels));
2620  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
2621  {
2622  auto iz_ind = iz - first_node;
2623  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
2624  {
2625  PetscScalar sumWij = 0.0;
2626  unsigned int counter = 0;
2627  for (auto i_gap : _subchannel_mesh.getChannelGaps(i_ch))
2628  {
2629  auto chans = _subchannel_mesh.getGapChannels(i_gap);
2630  unsigned int i_ch_loc = chans.first;
2631  PetscInt row_vec = i_ch_loc + _n_channels * iz_ind;
2632  PetscScalar loc_Wij_value;
2633  LibmeshPetscCall(VecGetValues(sol_holder_P, 1, &row_vec, &loc_Wij_value));
2634  sumWij += _subchannel_mesh.getCrossflowSign(i_ch, counter) * loc_Wij_value;
2635  counter++;
2636  }
2637  PetscInt row_vec = i_ch + _n_channels * iz_ind;
2638  LibmeshPetscCall(VecSetValues(sumWij_loc, 1, &row_vec, &sumWij, INSERT_VALUES));
2639  }
2640  }
2641 
2642  PetscScalar min_mdot;
2643  LibmeshPetscCall(VecAbs(_prod));
2644  LibmeshPetscCall(VecMin(_prod, NULL, &min_mdot));
2645  _console << "Minimum estimated mdot: " << min_mdot << std::endl;
2646 
2647  LibmeshPetscCall(VecAbs(sumWij_loc));
2648  LibmeshPetscCall(VecMax(sumWij_loc, NULL, &_max_sumWij));
2649  _max_sumWij = std::max(1e-10, _max_sumWij);
2650  _console << "Maximum estimated Wij: " << _max_sumWij << std::endl;
2651 
2653  _Wij_loc_vec, _Wij, first_node, last_node, _n_gaps));
2654  LibmeshPetscCall(VecAbs(_Wij_loc_vec));
2656  _Wij_old_loc_vec, _Wij_old, first_node, last_node, _n_gaps));
2657  LibmeshPetscCall(VecAbs(_Wij_old_loc_vec));
2658  LibmeshPetscCall(VecAXPY(_Wij_loc_vec, -1.0, _Wij_old_loc_vec));
2659  PetscScalar relax_factor;
2660  LibmeshPetscCall(VecAbs(_Wij_loc_vec));
2661 #if !PETSC_VERSION_LESS_THAN(3, 16, 0)
2662  LibmeshPetscCall(VecMean(_Wij_loc_vec, &relax_factor));
2663 #else
2664  LibmeshPetscCall(VecSum(_Wij_loc_vec, &relax_factor));
2665  relax_factor /= _block_size * _n_gaps;
2666 #endif
2667  relax_factor = relax_factor / _max_sumWij + 0.5;
2668  _console << "Relax base value: " << relax_factor << std::endl;
2669 
2670  PetscScalar resistance_relaxation = 0.9;
2671  _added_K = _max_sumWij / min_mdot;
2672  _console << "New cross resistance: " << _added_K << std::endl;
2673  _added_K = (_added_K * resistance_relaxation + (1.0 - resistance_relaxation) * _added_K_old) *
2674  relax_factor;
2675  _console << "Relaxed cross resistance: " << _added_K << std::endl;
2676  if (_added_K < 10 && _added_K >= 1.0)
2677  _added_K = 1.0; //(1.0 - resistance_relaxation);
2678  if (_added_K < 1.0 && _added_K >= 0.1)
2679  _added_K = 0.5;
2680  if (_added_K < 0.1 && _added_K >= 0.01)
2681  _added_K = 1. / 3.;
2682  if (_added_K < 1e-2 && _added_K >= 1e-3)
2683  _added_K = 0.1;
2684  if (_added_K < 1e-3)
2685  _added_K = 1.0 * _added_K;
2686  _console << "Actual added cross resistance: " << _added_K << std::endl;
2687  LibmeshPetscCall(VecScale(unity_vec_Wij, _added_K));
2689 
2690  // Adding cross resistances
2691  LibmeshPetscCall(MatDiagonalSet(mat_array[2 * Q + 2], unity_vec_Wij, ADD_VALUES));
2692  LibmeshPetscCall(VecDestroy(&mdot_estimate));
2693  LibmeshPetscCall(VecDestroy(&pmat_diag));
2694  LibmeshPetscCall(VecDestroy(&unity_vec));
2695  LibmeshPetscCall(VecDestroy(&p_estimate));
2696  LibmeshPetscCall(VecDestroy(&sol_holder_P));
2697  LibmeshPetscCall(VecDestroy(&diag_Wij_loc));
2698  LibmeshPetscCall(VecDestroy(&unity_vec_Wij));
2699  LibmeshPetscCall(VecDestroy(&Wij_estimate));
2700  LibmeshPetscCall(VecDestroy(&sumWij_loc));
2701  LibmeshPetscCall(VecDestroy(&_Wij_loc_vec));
2702  LibmeshPetscCall(VecDestroy(&_Wij_old_loc_vec));
2703 
2704  // Auto-computing relaxation factors
2705  PetscScalar relaxation_factor_mdot, relaxation_factor_P, relaxation_factor_Wij;
2706  relaxation_factor_mdot = 1.0;
2707  relaxation_factor_P = 1.0; // std::exp(-5.0);
2708  relaxation_factor_Wij = 0.1;
2709 
2710  _console << "Relax mdot: " << relaxation_factor_mdot << std::endl;
2711  _console << "Relax P: " << relaxation_factor_P << std::endl;
2712  _console << "Relax Wij: " << relaxation_factor_Wij << std::endl;
2713 
2714  PetscInt field_num = 0;
2715  Vec diag_mdot;
2716  LibmeshPetscCall(VecDuplicate(vec_array[field_num], &diag_mdot));
2717  LibmeshPetscCall(MatGetDiagonal(mat_array[Q * field_num + field_num], diag_mdot));
2718  LibmeshPetscCall(VecScale(diag_mdot, 1.0 / relaxation_factor_mdot));
2719  LibmeshPetscCall(
2720  MatDiagonalSet(mat_array[Q * field_num + field_num], diag_mdot, INSERT_VALUES));
2721  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2722  _prod, *_mdot_soln, first_node, last_node, _n_channels));
2723  LibmeshPetscCall(VecScale(diag_mdot, (1.0 - relaxation_factor_mdot)));
2724  LibmeshPetscCall(VecPointwiseMult(_prod, _prod, diag_mdot));
2725  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, _prod));
2726  LibmeshPetscCall(VecDestroy(&diag_mdot));
2727 
2728  _console << "mdot relaxed" << std::endl;
2729 
2730  field_num = 1;
2731  Vec diag_P;
2732  LibmeshPetscCall(VecDuplicate(vec_array[field_num], &diag_P));
2733  LibmeshPetscCall(MatGetDiagonal(mat_array[Q * field_num + field_num], diag_P));
2734  LibmeshPetscCall(VecScale(diag_P, 1.0 / relaxation_factor_P));
2735  LibmeshPetscCall(MatDiagonalSet(mat_array[Q * field_num + field_num], diag_P, INSERT_VALUES));
2736  _console << "Mat assembled" << std::endl;
2737  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2738  _prod, *_P_soln, first_node, last_node, _n_channels));
2739  LibmeshPetscCall(VecScale(diag_P, (1.0 - relaxation_factor_P)));
2740  LibmeshPetscCall(VecPointwiseMult(_prod, _prod, diag_P));
2741  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, _prod));
2742  LibmeshPetscCall(VecDestroy(&diag_P));
2743 
2744  _console << "P relaxed" << std::endl;
2745 
2746  field_num = 2;
2747  Vec diag_Wij;
2748  LibmeshPetscCall(VecDuplicate(vec_array[field_num], &diag_Wij));
2749  LibmeshPetscCall(MatGetDiagonal(mat_array[Q * field_num + field_num], diag_Wij));
2750  LibmeshPetscCall(VecScale(diag_Wij, 1.0 / relaxation_factor_Wij));
2751  LibmeshPetscCall(MatDiagonalSet(mat_array[Q * field_num + field_num], diag_Wij, INSERT_VALUES));
2753  _Wij_vec, _Wij, first_node, last_node, _n_gaps));
2754  LibmeshPetscCall(VecScale(diag_Wij, (1.0 - relaxation_factor_Wij)));
2755  LibmeshPetscCall(VecPointwiseMult(_Wij_vec, _Wij_vec, diag_Wij));
2756  LibmeshPetscCall(VecAXPY(vec_array[field_num], 1.0, _Wij_vec));
2757  LibmeshPetscCall(VecDestroy(&diag_Wij));
2758 
2759  _console << "Wij relaxed" << std::endl;
2760  }
2761  _console << "Linear solver relaxed." << std::endl;
2762 
2763  // Creating nested matrices
2764  LibmeshPetscCall(MatCreateNest(PETSC_COMM_WORLD, Q, NULL, Q, NULL, mat_array.data(), &A_nest));
2765  LibmeshPetscCall(VecCreateNest(PETSC_COMM_WORLD, Q, NULL, vec_array.data(), &b_nest));
2766  _console << "Nested system created." << std::endl;
2767 
2769  // Creating linear solver
2770  LibmeshPetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
2771  LibmeshPetscCall(KSPSetType(ksp, KSPFGMRES));
2772  // Setting KSP operators
2773  LibmeshPetscCall(KSPSetOperators(ksp, A_nest, A_nest));
2774  // Set KSP and PC options
2775  LibmeshPetscCall(KSPGetPC(ksp, &pc));
2776  LibmeshPetscCall(PCSetType(pc, PCFIELDSPLIT));
2777  LibmeshPetscCall(KSPSetTolerances(ksp, _rtol, _atol, _dtol, _maxit));
2778  // Splitting fields
2779  std::vector<IS> rows(Q);
2780  // IS rows[Q];
2781  PetscInt M = 0;
2782  LibmeshPetscCall(MatNestGetISs(A_nest, rows.data(), NULL));
2783  for (PetscInt j = 0; j < Q; ++j)
2784  {
2785  IS expand1;
2786  LibmeshPetscCall(ISDuplicate(rows[M], &expand1));
2787  M += 1;
2788  LibmeshPetscCall(PCFieldSplitSetIS(pc, NULL, expand1));
2789  LibmeshPetscCall(ISDestroy(&expand1));
2790  }
2791  _console << "Linear solver assembled." << std::endl;
2792 
2794  LibmeshPetscCall(VecDuplicate(b_nest, &x_nest));
2795  LibmeshPetscCall(VecSet(x_nest, 0.0));
2796  LibmeshPetscCall(KSPSolve(ksp, b_nest, x_nest));
2797 
2799  LibmeshPetscCall(VecDestroy(&b_nest));
2800  LibmeshPetscCall(MatDestroy(&A_nest));
2801  LibmeshPetscCall(KSPDestroy(&ksp));
2802  for (PetscInt i = 0; i < Q * Q; i++)
2803  {
2804  LibmeshPetscCall(MatDestroy(&mat_array[i]));
2805  }
2806  for (PetscInt i = 0; i < Q; i++)
2807  {
2808  LibmeshPetscCall(VecDestroy(&vec_array[i]));
2809  }
2810 
2812  Vec sol_mdot, sol_p, sol_Wij;
2813  _console << "Vectors created." << std::endl;
2814  PetscInt num_vecs;
2815  Vec * loc_vecs;
2816  LibmeshPetscCall(VecNestGetSubVecs(x_nest, &num_vecs, &loc_vecs));
2817  _console << "Starting extraction." << std::endl;
2818  LibmeshPetscCall(VecDuplicate(_mc_axial_convection_rhs, &sol_mdot));
2819  LibmeshPetscCall(VecCopy(loc_vecs[0], sol_mdot));
2820  LibmeshPetscCall(VecDuplicate(_amc_sys_mdot_rhs, &sol_p));
2821  LibmeshPetscCall(VecCopy(loc_vecs[1], sol_p));
2822  LibmeshPetscCall(VecDuplicate(_cmc_sys_Wij_rhs, &sol_Wij));
2823  LibmeshPetscCall(VecCopy(loc_vecs[2], sol_Wij));
2824  _console << "Getting individual field solutions from coupled solver." << std::endl;
2825 
2827  PetscScalar * sol_p_array;
2828  LibmeshPetscCall(VecGetArray(sol_p, &sol_p_array));
2829  PetscScalar * sol_Wij_array;
2830  LibmeshPetscCall(VecGetArray(sol_Wij, &sol_Wij_array));
2831 
2833  LibmeshPetscCall(populateSolutionChan<SolutionHandle>(
2834  sol_mdot, *_mdot_soln, first_node, last_node, _n_channels));
2835 
2837  for (unsigned int iz = last_node; iz > first_node - 1; iz--)
2838  {
2839  auto iz_ind = iz - first_node;
2840  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
2841  {
2842  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
2843  PetscScalar value = sol_p_array[iz_ind * _n_channels + i_ch];
2844  _P_soln->set(node_in, value);
2845  }
2846  }
2847 
2850  sol_Wij, _Wij, first_node, last_node - 1, _n_gaps));
2851 
2854  {
2855  Vec sol_h;
2856  LibmeshPetscCall(VecDuplicate(_hc_sys_h_rhs, &sol_h));
2857  LibmeshPetscCall(VecCopy(loc_vecs[3], sol_h));
2858  PetscScalar * sol_h_array;
2859  LibmeshPetscCall(VecGetArray(sol_h, &sol_h_array));
2860  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
2861  {
2862  auto iz_ind = iz - first_node;
2863  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
2864  {
2865  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
2866  auto h_out = sol_h_array[iz_ind * _n_channels + i_ch];
2867  if (h_out < 0)
2868  {
2869  mooseError(name(),
2870  " : Calculation of negative Enthalpy h_out = : ",
2871  h_out,
2872  " Axial Level= : ",
2873  iz);
2874  }
2875  _h_soln->set(node_out, h_out);
2876  }
2877  }
2878  LibmeshPetscCall(VecDestroy(&sol_h));
2879  }
2880 
2882  LibmeshPetscCall(MatMult(_mc_sumWij_mat, sol_Wij, _prod));
2883  LibmeshPetscCall(populateSolutionChan<SolutionHandle>(
2884  _prod, *_SumWij_soln, first_node, last_node, _n_channels));
2885 
2886  Vec sumWij_loc;
2887  LibmeshPetscCall(createPetscVector(sumWij_loc, _block_size * _n_channels));
2888  LibmeshPetscCall(populateVectorFromHandle<SolutionHandle>(
2889  _prod, *_SumWij_soln, first_node, last_node, _n_channels));
2890 
2891  LibmeshPetscCall(VecAbs(_prod));
2892  LibmeshPetscCall(VecMax(_prod, NULL, &_max_sumWij_new));
2893  _console << "Maximum estimated Wij new: " << _max_sumWij_new << std::endl;
2895  _console << "Correction factor: " << _correction_factor << std::endl;
2896  _console << "Solutions assigned to MOOSE variables." << std::endl;
2897 
2899  LibmeshPetscCall(VecDestroy(&x_nest));
2900  LibmeshPetscCall(VecDestroy(&sol_mdot));
2901  LibmeshPetscCall(VecDestroy(&sol_p));
2902  LibmeshPetscCall(VecDestroy(&sol_Wij));
2903  LibmeshPetscCall(VecDestroy(&sumWij_loc));
2904  _console << "Solutions destroyed." << std::endl;
2905 
2906  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
2907 }
std::unique_ptr< SolutionHandle > _P_soln
virtual void computeSumWij(int iblock)
Computes net diversion crossflow per channel for block iblock.
virtual void computeh(int iblock)
Computes Enthalpy per channel for block iblock.
virtual void computeDP(int iblock)
Computes Pressure Drop per channel for block iblock.
PetscScalar _added_K
Added resistances for monolithic convergence.
virtual const std::pair< unsigned int, unsigned int > & getGapChannels(unsigned int i_gap) const =0
Return a pair of inter-wrapper indices for a given gap index.
virtual void computeWijPrime(int iblock)
Computes turbulent crossflow per gap for block iblock.
std::unique_ptr< SolutionHandle > _SumWij_soln
Mat _cmc_pressure_force_mat
Cross momentum conservation - pressure force.
Mat _mc_sumWij_mat
Mass conservation Mass conservation - sum of cross fluxes.
const PetscReal & _dtol
The divergence tolerance for the ksp linear solver.
PetscFunctionBegin
virtual PetscErrorCode createPetscVector(Vec &v, PetscInt n)
Petsc Functions.
const PetscReal & _rtol
The relative convergence tolerance, (relative decrease) for the ksp linear solver.
Mat _mc_axial_convection_mat
Mass conservation - axial convection.
std::unique_ptr< SolutionHandle > _mdot_soln
virtual void computeMdot(int iblock)
Computes mass flow per channel for block iblock.
PetscErrorCode populateDenseFromVector(const Vec &x, T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
virtual const std::string & name() const
const PetscInt & _maxit
The maximum number of iterations to use for the ksp linear solver.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
virtual void computeP(int iblock)
Computes Pressure per channel for block iblock.
const bool _monolithic_thermal_bool
Thermal monolithic bool.
virtual void computeWij(int iblock)
Computes cross fluxes for block iblock.
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
Mat _cmc_sys_Wij_mat
Lateral momentum system matrix.
Mat _amc_pressure_force_mat
Axial momentum conservation - pressure force.
void mooseError(Args &&... args) const
PetscErrorCode populateVectorFromDense(Vec &x, const T &solution, const unsigned int first_axial_level, const unsigned int last_axial_level, const unsigned int cross_dimension)
std::unique_ptr< SolutionHandle > _h_soln
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
libMesh::DenseMatrix< Real > _Wij_old
const ConsoleStream _console
libMesh::DenseMatrix< Real > & _Wij
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)
virtual const std::vector< unsigned int > & getChannelGaps(unsigned int i_chan) const =0
Return a vector of gap indices for a given channel index.
Mat _amc_sys_mdot_mat
Axial momentum system matrix.
const PetscReal & _atol
The absolute convergence tolerance for the ksp linear solver.
virtual const Real & getCrossflowSign(unsigned int i_chan, unsigned int i_local) const =0
Return a signs for the cross flow given a inter-wrapper index and local neighbor index.

◆ initializeSolution()

void InterWrapper1PhaseProblem::initializeSolution ( )
protectedvirtualinherited

Function to initialize the solution fields.

Definition at line 2910 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

2911 {
2912  unsigned int last_node = _n_cells;
2913  unsigned int first_node = 1;
2914  for (unsigned int iz = first_node; iz < last_node + 1; iz++)
2915  {
2916  for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
2917  {
2918  auto * node_out = _subchannel_mesh.getChannelNode(i_ch, iz);
2919  auto * node_in = _subchannel_mesh.getChannelNode(i_ch, iz - 1);
2920  _mdot_soln->set(node_out, (*_mdot_soln)(node_in));
2921  }
2922  }
2923  _mdot_soln->close();
2924 }
std::unique_ptr< SolutionHandle > _mdot_soln
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.

◆ initialSetup()

void InterWrapper1PhaseProblem::initialSetup ( )
overridevirtualinherited

Reimplemented from ExternalProblem.

Definition at line 226 of file InterWrapper1PhaseProblem.C.

227 {
229  _fp = &getUserObject<SinglePhaseFluidProperties>(getParam<UserObjectName>("fp"));
230  _mdot_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::MASS_FLOW_RATE));
231  _SumWij_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::SUM_CROSSFLOW));
232  _P_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::PRESSURE));
233  _DP_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::PRESSURE_DROP));
234  _h_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::ENTHALPY));
235  _T_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::TEMPERATURE));
236  if (_pin_mesh_exist)
237  _Tpin_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::PIN_TEMPERATURE));
238  _rho_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::DENSITY));
239  _mu_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::VISCOSITY));
240  _S_flow_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::SURFACE_AREA));
241  _w_perim_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::WETTED_PERIMETER));
242  _q_prime_soln = std::make_unique<SolutionHandle>(getVariable(0, SubChannelApp::LINEAR_HEAT_RATE));
243 }
static const std::string PRESSURE_DROP
pressure drop
Definition: SubChannelApp.h:37
static const std::string MASS_FLOW_RATE
mass flow rate
Definition: SubChannelApp.h:29
std::unique_ptr< SolutionHandle > _w_perim_soln
static const std::string DENSITY
density
Definition: SubChannelApp.h:47
std::unique_ptr< SolutionHandle > _S_flow_soln
std::unique_ptr< SolutionHandle > _P_soln
std::unique_ptr< SolutionHandle > _Tpin_soln
std::unique_ptr< SolutionHandle > _SumWij_soln
std::unique_ptr< SolutionHandle > _rho_soln
std::unique_ptr< SolutionHandle > _mdot_soln
const SinglePhaseFluidProperties * _fp
Solutions handles and link to TH tables properties.
virtual const MooseVariableFieldBase & getVariable(const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const override
static const std::string WETTED_PERIMETER
wetted perimeter
Definition: SubChannelApp.h:51
static const std::string VISCOSITY
viscosity
Definition: SubChannelApp.h:49
static const std::string PIN_TEMPERATURE
pin temperature
Definition: SubChannelApp.h:43
static const std::string LINEAR_HEAT_RATE
linear heat rate
Definition: SubChannelApp.h:53
std::unique_ptr< SolutionHandle > _T_soln
static const std::string ENTHALPY
enthalpy
Definition: SubChannelApp.h:39
void initialSetup() override
std::unique_ptr< SolutionHandle > _mu_soln
std::unique_ptr< SolutionHandle > _DP_soln
static const std::string SUM_CROSSFLOW
sum of diversion crossflow
Definition: SubChannelApp.h:33
static const std::string PRESSURE
pressure
Definition: SubChannelApp.h:35
std::unique_ptr< SolutionHandle > _q_prime_soln
static const std::string SURFACE_AREA
surface area
Definition: SubChannelApp.h:31
const bool _pin_mesh_exist
Flag that informs if there is a pin mesh or not.
std::unique_ptr< SolutionHandle > _h_soln
static const std::string TEMPERATURE
temperature
Definition: SubChannelApp.h:41

◆ petscSnesSolver()

PetscErrorCode InterWrapper1PhaseProblem::petscSnesSolver ( int  iblock,
const libMesh::DenseVector< Real > &  solution,
libMesh::DenseVector< Real > &  root 
)
protectedvirtualinherited

Computes solution of nonlinear equation using snes.

Definition at line 2337 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::computeWijFromSolve().

2340 {
2341  SNES snes;
2342  KSP ksp;
2343  PC pc;
2344  Vec x, r;
2345  PetscMPIInt size;
2346  PetscScalar * xx;
2347 
2349  PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
2350  if (size > 1)
2351  SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Example is only for sequential runs");
2352  LibmeshPetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
2353  LibmeshPetscCall(VecCreate(PETSC_COMM_WORLD, &x));
2354  LibmeshPetscCall(VecSetSizes(x, PETSC_DECIDE, _block_size * _n_gaps));
2355  LibmeshPetscCall(VecSetFromOptions(x));
2356  LibmeshPetscCall(VecDuplicate(x, &r));
2357 
2358 #if PETSC_VERSION_LESS_THAN(3, 13, 0)
2359  LibmeshPetscCall(PetscOptionsSetValue(PETSC_NULL, "-snes_mf", PETSC_NULL));
2360 #else
2361  LibmeshPetscCall(SNESSetUseMatrixFree(snes, PETSC_FALSE, PETSC_TRUE));
2362 #endif
2363  problemInfo info;
2364  info.iblock = iblock;
2365  info.schp = this;
2366  LibmeshPetscCall(SNESSetFunction(snes, r, formFunctionIW, &info));
2367 
2368  LibmeshPetscCall(SNESGetKSP(snes, &ksp));
2369  LibmeshPetscCall(KSPGetPC(ksp, &pc));
2370  LibmeshPetscCall(PCSetType(pc, PCNONE));
2371  LibmeshPetscCall(KSPSetTolerances(ksp, _rtol, _atol, _dtol, _maxit));
2372  LibmeshPetscCall(SNESSetFromOptions(snes));
2373  LibmeshPetscCall(VecGetArray(x, &xx));
2374  for (unsigned int i = 0; i < _block_size * _n_gaps; i++)
2375  {
2376  xx[i] = solution(i);
2377  }
2378  LibmeshPetscCall(VecRestoreArray(x, &xx));
2379 
2380  LibmeshPetscCall(SNESSolve(snes, NULL, x));
2381  LibmeshPetscCall(VecGetArray(x, &xx));
2382  for (unsigned int i = 0; i < _block_size * _n_gaps; i++)
2383  root(i) = xx[i];
2384 
2385  LibmeshPetscCall(VecRestoreArray(x, &xx));
2386  LibmeshPetscCall(VecDestroy(&x));
2387  LibmeshPetscCall(VecDestroy(&r));
2388  LibmeshPetscCall(SNESDestroy(&snes));
2389  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
2390 }
MPI_Info info
const PetscReal & _dtol
The divergence tolerance for the ksp linear solver.
PetscFunctionBegin
const PetscReal & _rtol
The relative convergence tolerance, (relative decrease) for the ksp linear solver.
const PetscInt & _maxit
The maximum number of iterations to use for the ksp linear solver.
const std::vector< double > x
Real root(std::function< Real(Real)> const &f, Real x1, Real x2, Real tol=1.0e-12)
Finds the root of a function using Brent&#39;s method.
Definition: BrentsMethod.C:66
problem information to be used in the PETSc problem
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)
friend PetscErrorCode formFunctionIW(SNES snes, Vec x, Vec f, void *info)
creates the residual function to be used in the PETCs snes context
const PetscReal & _atol
The absolute convergence tolerance for the ksp linear solver.

◆ populateDenseFromVector()

template<class T >
PetscErrorCode InterWrapper1PhaseProblem::populateDenseFromVector ( const Vec &  x,
T &  solution,
const unsigned int  first_axial_level,
const unsigned int  last_axial_level,
const unsigned int  cross_dimension 
)
protectedinherited

Definition at line 294 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::computeWijPrime(), and InterWrapper1PhaseProblem::implicitPetscSolve().

299 {
300  PetscErrorCode ierr;
301  PetscScalar * xx;
302 
304  ierr = VecGetArray(x, &xx);
305  CHKERRQ(ierr);
306  for (unsigned int iz = first_axial_level; iz < last_axial_level + 1; iz++)
307  {
308  unsigned int iz_ind = iz - first_axial_level;
309  for (unsigned int i_l = 0; i_l < cross_dimension; i_l++)
310  {
311  loc_solution(i_l, iz) = xx[iz_ind * cross_dimension + i_l];
312  }
313  }
314  ierr = VecRestoreArray(x, &xx);
315  CHKERRQ(ierr);
316  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
317 }
PetscFunctionBegin
const std::vector< double > x
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ populateSolutionChan()

template<class T >
PetscErrorCode InterWrapper1PhaseProblem::populateSolutionChan ( const Vec &  x,
T &  solution,
const unsigned int  first_axial_level,
const unsigned int  last_axial_level,
const unsigned int  cross_dimension 
)
protectedinherited

Definition at line 423 of file InterWrapper1PhaseProblem.C.

428 {
429  PetscScalar * xx;
431  LibmeshPetscCall(VecGetArray(x, &xx));
432  Node * loc_node;
433  for (unsigned int iz = first_axial_level; iz < last_axial_level + 1; iz++)
434  {
435  unsigned int iz_ind = iz - first_axial_level;
436  for (unsigned int i_l = 0; i_l < cross_dimension; i_l++)
437  {
438  loc_node = _subchannel_mesh.getChannelNode(i_l, iz);
439  loc_solution.set(loc_node, xx[iz_ind * cross_dimension + i_l]);
440  }
441  }
442  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
443 }
PetscFunctionBegin
const std::vector< double > x
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ populateSolutionGap()

template<class T >
PetscErrorCode InterWrapper1PhaseProblem::populateSolutionGap ( const Vec &  x,
T &  solution,
const unsigned int  first_axial_level,
const unsigned int  last_axial_level,
const unsigned int  cross_dimension 
)
protectedinherited

Definition at line 447 of file InterWrapper1PhaseProblem.C.

452 {
453  PetscScalar * xx;
455  LibmeshPetscCall(VecGetArray(x, &xx));
456  for (unsigned int iz = first_axial_level; iz < last_axial_level + 1; iz++)
457  {
458  unsigned int iz_ind = iz - first_axial_level;
459  for (unsigned int i_l = 0; i_l < cross_dimension; i_l++)
460  {
461  loc_solution(iz * cross_dimension + i_l) = xx[iz_ind * cross_dimension + i_l];
462  }
463  }
464  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
465 }
PetscFunctionBegin
const std::vector< double > x
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ populateVectorFromDense()

template<class T >
PetscErrorCode InterWrapper1PhaseProblem::populateVectorFromDense ( Vec &  x,
const T &  solution,
const unsigned int  first_axial_level,
const unsigned int  last_axial_level,
const unsigned int  cross_dimension 
)
protectedinherited

Definition at line 376 of file InterWrapper1PhaseProblem.C.

Referenced by InterWrapper1PhaseProblem::computeSumWij(), InterWrapper1PhaseProblem::computeWij(), and InterWrapper1PhaseProblem::implicitPetscSolve().

381 {
382  PetscScalar * xx;
384  LibmeshPetscCall(VecGetArray(x, &xx));
385  for (unsigned int iz = first_axial_level; iz < last_axial_level; iz++)
386  {
387  unsigned int iz_ind = iz - first_axial_level;
388  for (unsigned int i_l = 0; i_l < cross_dimension; i_l++)
389  {
390  xx[iz_ind * cross_dimension + i_l] = loc_solution(i_l, iz);
391  }
392  }
393  LibmeshPetscCall(VecRestoreArray(x, &xx));
394  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
395 }
PetscFunctionBegin
const std::vector< double > x
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ populateVectorFromHandle()

template<class T >
PetscErrorCode InterWrapper1PhaseProblem::populateVectorFromHandle ( Vec &  x,
const T &  solution,
const unsigned int  first_axial_level,
const unsigned int  last_axial_level,
const unsigned int  cross_dimension 
)
protectedinherited

Definition at line 399 of file InterWrapper1PhaseProblem.C.

404 {
405  PetscScalar * xx;
407  LibmeshPetscCall(VecGetArray(x, &xx));
408  for (unsigned int iz = first_axial_level; iz < last_axial_level + 1; iz++)
409  {
410  unsigned int iz_ind = iz - first_axial_level;
411  for (unsigned int i_l = 0; i_l < cross_dimension; i_l++)
412  {
413  auto * loc_node = _subchannel_mesh.getChannelNode(i_l, iz);
414  xx[iz_ind * cross_dimension + i_l] = loc_solution(loc_node);
415  }
416  }
417  LibmeshPetscCall(VecRestoreArray(x, &xx));
418  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
419 }
PetscFunctionBegin
const std::vector< double > x
const InterWrapperMesh & _subchannel_mesh
virtual Node * getChannelNode(unsigned int i_chan, unsigned iz) const =0
Get the inter-wrapper mesh node for a given channel index and elevation index.
PetscFunctionReturn(LIBMESH_PETSC_SUCCESS)

◆ residualFunction()

libMesh::DenseVector< Real > InterWrapper1PhaseProblem::residualFunction ( int  iblock,
libMesh::DenseVector< Real solution 
)
protectedvirtualinherited

Computes Residual per gap for block iblock.

Definition at line 2295 of file InterWrapper1PhaseProblem.C.

Referenced by formFunctionIW().

2296 {
2297  unsigned int last_node = (iblock + 1) * _block_size;
2298  unsigned int first_node = iblock * _block_size;
2299  libMesh::DenseVector<Real> Wij_residual_vector(_n_gaps * _block_size, 0.0);
2300  // Assign the solution to the cross-flow matrix
2301  int i = 0;
2302  for (unsigned int iz = first_node + 1; iz < last_node + 1; iz++)
2303  {
2304  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
2305  {
2306  _Wij(i_gap, iz) = solution(i);
2307  i++;
2308  }
2309  }
2310 
2311  // Calculating sum of crossflows
2312  computeSumWij(iblock);
2313  // Solving axial flux
2314  computeMdot(iblock);
2315  // Calculation of turbulent Crossflow
2316  computeWijPrime(iblock);
2317  // Solving for Pressure Drop
2318  computeDP(iblock);
2319  // Solving for pressure
2320  computeP(iblock);
2321  // Solving cross fluxes
2322  computeWij(iblock);
2323 
2324  // Turn the residual matrix into a residual vector
2325  for (unsigned int iz = 0; iz < _block_size; iz++)
2326  {
2327  for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
2328  {
2329  int i = _n_gaps * iz + i_gap; // column wise transfer
2330  Wij_residual_vector(i) = _Wij_residual_matrix(i_gap, iz);
2331  }
2332  }
2333  return Wij_residual_vector;
2334 }
virtual void computeSumWij(int iblock)
Computes net diversion crossflow per channel for block iblock.
virtual void computeDP(int iblock)
Computes Pressure Drop per channel for block iblock.
virtual void computeWijPrime(int iblock)
Computes turbulent crossflow per gap for block iblock.
virtual void computeMdot(int iblock)
Computes mass flow per channel for block iblock.
virtual void computeP(int iblock)
Computes Pressure per channel for block iblock.
virtual void computeWij(int iblock)
Computes cross fluxes for block iblock.
libMesh::DenseMatrix< Real > & _Wij
libMesh::DenseMatrix< Real > _Wij_residual_matrix

◆ solverSystemConverged()

bool InterWrapper1PhaseProblem::solverSystemConverged ( const unsigned int  )
overridevirtualinherited

Reimplemented from ExternalProblem.

Definition at line 310 of file InterWrapper1PhaseProblem.C.

311 {
312  return _converged;
313 }
bool _converged
Variable that informs whether we exited external solve with a converged solution or not...

◆ syncSolutions()

void InterWrapper1PhaseProblem::syncSolutions ( Direction  direction)
overridevirtualinherited

Implements ExternalProblem.

Definition at line 3082 of file InterWrapper1PhaseProblem.C.

3083 {
3084 }

◆ validParams()

InputParameters QuadInterWrapper1PhaseProblem::validParams ( )
static

Definition at line 15 of file QuadInterWrapper1PhaseProblem.C.

16 {
18  params.addClassDescription(
19  "Solver class for interwrapper of assemblies in a square-lattice arrangement");
20  return params;
21 }
static InputParameters validParams()
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _added_K

PetscScalar InterWrapper1PhaseProblem::_added_K = 0.0
protectedinherited

Added resistances for monolithic convergence.

Definition at line 282 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::implicitPetscSolve().

◆ _added_K_old

PetscScalar InterWrapper1PhaseProblem::_added_K_old = 1000.0
protectedinherited

◆ _amc_advective_derivative_mat

Mat InterWrapper1PhaseProblem::_amc_advective_derivative_mat
protectedinherited

Axial momentum conservation - advective (Eulerian) derivative.

Definition at line 227 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeDP(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _amc_advective_derivative_rhs

Vec InterWrapper1PhaseProblem::_amc_advective_derivative_rhs
protectedinherited

◆ _amc_cross_derivative_mat

Mat InterWrapper1PhaseProblem::_amc_cross_derivative_mat
protectedinherited

Axial momentum conservation - cross flux derivative.

Definition at line 230 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeDP(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _amc_cross_derivative_rhs

Vec InterWrapper1PhaseProblem::_amc_cross_derivative_rhs
protectedinherited

◆ _amc_friction_force_mat

Mat InterWrapper1PhaseProblem::_amc_friction_force_mat
protectedinherited

◆ _amc_friction_force_rhs

Vec InterWrapper1PhaseProblem::_amc_friction_force_rhs
protectedinherited

◆ _amc_gravity_rhs

Vec InterWrapper1PhaseProblem::_amc_gravity_rhs
protectedinherited

Axial momentum conservation - buoyancy force No implicit matrix.

Definition at line 237 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeDP(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _amc_pressure_force_mat

Mat InterWrapper1PhaseProblem::_amc_pressure_force_mat
protectedinherited

◆ _amc_pressure_force_rhs

Vec InterWrapper1PhaseProblem::_amc_pressure_force_rhs
protectedinherited

◆ _amc_sys_mdot_mat

Mat InterWrapper1PhaseProblem::_amc_sys_mdot_mat
protectedinherited

◆ _amc_sys_mdot_rhs

Vec InterWrapper1PhaseProblem::_amc_sys_mdot_rhs
protectedinherited

◆ _amc_time_derivative_mat

Mat InterWrapper1PhaseProblem::_amc_time_derivative_mat
protectedinherited

◆ _amc_time_derivative_rhs

Vec InterWrapper1PhaseProblem::_amc_time_derivative_rhs
protectedinherited

◆ _amc_turbulent_cross_flows_mat

Mat InterWrapper1PhaseProblem::_amc_turbulent_cross_flows_mat
protectedinherited

Mass conservation - density time derivative No implicit matrix.

Axial momentum Axial momentum conservation - compute turbulent cross fluxes

Definition at line 221 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeWijPrime(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _amc_turbulent_cross_flows_rhs

Vec InterWrapper1PhaseProblem::_amc_turbulent_cross_flows_rhs
protectedinherited

◆ _atol

const PetscReal& InterWrapper1PhaseProblem::_atol
protectedinherited

◆ _beta

const Real& InterWrapper1PhaseProblem::_beta
protectedinherited

Thermal diffusion coefficient used in turbulent crossflow.

Definition at line 127 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::computeWijPrime().

◆ _block_size

const unsigned int InterWrapper1PhaseProblem::_block_size
protectedinherited

◆ _cmc_advective_derivative_mat

Mat InterWrapper1PhaseProblem::_cmc_advective_derivative_mat
protectedinherited

Cross momentum conservation - advective (Eulerian) derivative.

Definition at line 250 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeWij(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _cmc_advective_derivative_rhs

Vec InterWrapper1PhaseProblem::_cmc_advective_derivative_rhs
protectedinherited

◆ _cmc_friction_force_mat

Mat InterWrapper1PhaseProblem::_cmc_friction_force_mat
protectedinherited

◆ _cmc_friction_force_rhs

Vec InterWrapper1PhaseProblem::_cmc_friction_force_rhs
protectedinherited

◆ _cmc_pressure_force_mat

Mat InterWrapper1PhaseProblem::_cmc_pressure_force_mat
protectedinherited

◆ _cmc_pressure_force_rhs

Vec InterWrapper1PhaseProblem::_cmc_pressure_force_rhs
protectedinherited

◆ _cmc_sys_Wij_mat

Mat InterWrapper1PhaseProblem::_cmc_sys_Wij_mat
protectedinherited

◆ _cmc_sys_Wij_rhs

Vec InterWrapper1PhaseProblem::_cmc_sys_Wij_rhs
protectedinherited

◆ _cmc_time_derivative_mat

Mat InterWrapper1PhaseProblem::_cmc_time_derivative_mat
protectedinherited

Cross momentum Cross momentum conservation - time derivative.

Definition at line 247 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeWij(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _cmc_time_derivative_rhs

Vec InterWrapper1PhaseProblem::_cmc_time_derivative_rhs
protectedinherited

◆ _cmc_Wij_channel_dummy

Vec InterWrapper1PhaseProblem::_cmc_Wij_channel_dummy
protectedinherited

◆ _compute_density

const bool InterWrapper1PhaseProblem::_compute_density
protectedinherited

Flag that activates or deactivates the calculation of density.

Definition at line 113 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

◆ _compute_power

const bool InterWrapper1PhaseProblem::_compute_power
protectedinherited

Flag that informs if we need to solve the Enthalpy/Temperature equations or not.

Definition at line 117 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

◆ _compute_viscosity

const bool InterWrapper1PhaseProblem::_compute_viscosity
protectedinherited

Flag that activates or deactivates the calculation of viscosity.

Definition at line 115 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

◆ _converged

bool InterWrapper1PhaseProblem::_converged
protectedinherited

◆ _correction_factor

PetscScalar InterWrapper1PhaseProblem::_correction_factor = 1.0
protectedinherited

◆ _CT

const Real& InterWrapper1PhaseProblem::_CT
protectedinherited

Turbulent modeling parameter used in axial momentum equation.

Definition at line 129 of file InterWrapper1PhaseProblem.h.

Referenced by TriInterWrapper1PhaseProblem::computeDP(), and InterWrapper1PhaseProblem::computeDP().

◆ _DP_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_DP_soln
protectedinherited

◆ _dt

const Real& InterWrapper1PhaseProblem::_dt
protectedinherited

◆ _dtol

const PetscReal& InterWrapper1PhaseProblem::_dtol
protectedinherited

◆ _fp

const SinglePhaseFluidProperties* InterWrapper1PhaseProblem::_fp
protectedinherited

◆ _g_grav

const Real InterWrapper1PhaseProblem::_g_grav
protectedinherited

◆ _global_counter

PetscInt InterWrapper1PhaseProblem::_global_counter = 0
protectedinherited

No implicit matrix.

Definition at line 279 of file InterWrapper1PhaseProblem.h.

◆ _h_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_h_soln
protectedinherited

◆ _hc_added_heat_rhs

Vec InterWrapper1PhaseProblem::_hc_added_heat_rhs
protectedinherited

◆ _hc_advective_derivative_mat

Mat InterWrapper1PhaseProblem::_hc_advective_derivative_mat
protectedinherited

Enthalpy conservation - advective (Eulerian) derivative;.

Definition at line 268 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeh(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _hc_advective_derivative_rhs

Vec InterWrapper1PhaseProblem::_hc_advective_derivative_rhs
protectedinherited

◆ _hc_cross_derivative_mat

Mat InterWrapper1PhaseProblem::_hc_cross_derivative_mat
protectedinherited

◆ _hc_cross_derivative_rhs

Vec InterWrapper1PhaseProblem::_hc_cross_derivative_rhs
protectedinherited

◆ _hc_sys_h_mat

Mat InterWrapper1PhaseProblem::_hc_sys_h_mat
protectedinherited

◆ _hc_sys_h_rhs

Vec InterWrapper1PhaseProblem::_hc_sys_h_rhs
protectedinherited

◆ _hc_time_derivative_mat

Mat InterWrapper1PhaseProblem::_hc_time_derivative_mat
protectedinherited

Enthalpy Enthalpy conservation - time derivative.

Definition at line 265 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::cleanUp(), InterWrapper1PhaseProblem::computeh(), and InterWrapper1PhaseProblem::InterWrapper1PhaseProblem().

◆ _hc_time_derivative_rhs

Vec InterWrapper1PhaseProblem::_hc_time_derivative_rhs
protectedinherited

◆ _implicit_bool

const bool InterWrapper1PhaseProblem::_implicit_bool
protectedinherited

◆ _interpolation_scheme

const MooseEnum InterWrapper1PhaseProblem::_interpolation_scheme
protectedinherited

The interpolation method used in constructing the systems.

Definition at line 145 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::computeInterpolationCoefficients().

◆ _kij

const Real& InterWrapper1PhaseProblem::_kij
protectedinherited

Definition at line 101 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::computeWij().

◆ _max_sumWij

PetscScalar InterWrapper1PhaseProblem::_max_sumWij
protectedinherited

◆ _max_sumWij_new

PetscScalar InterWrapper1PhaseProblem::_max_sumWij_new
protectedinherited

◆ _maxit

const PetscInt& InterWrapper1PhaseProblem::_maxit
protectedinherited

◆ _mc_axial_convection_mat

Mat InterWrapper1PhaseProblem::_mc_axial_convection_mat
protectedinherited

◆ _mc_axial_convection_rhs

Vec InterWrapper1PhaseProblem::_mc_axial_convection_rhs
protectedinherited

◆ _mc_sumWij_mat

Mat InterWrapper1PhaseProblem::_mc_sumWij_mat
protectedinherited

◆ _mdot_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_mdot_soln
protectedinherited

◆ _monolithic_thermal_bool

const bool InterWrapper1PhaseProblem::_monolithic_thermal_bool
protectedinherited

◆ _mu_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_mu_soln
protectedinherited

◆ _n_assemblies

const unsigned int InterWrapper1PhaseProblem::_n_assemblies
protectedinherited

◆ _n_blocks

unsigned int InterWrapper1PhaseProblem::_n_blocks
protectedinherited

◆ _n_cells

const unsigned int InterWrapper1PhaseProblem::_n_cells
protectedinherited

◆ _n_channels

const unsigned int InterWrapper1PhaseProblem::_n_channels
protectedinherited

◆ _n_gaps

const unsigned int InterWrapper1PhaseProblem::_n_gaps
protectedinherited

◆ _one

Real InterWrapper1PhaseProblem::_one
protectedinherited

Definition at line 109 of file InterWrapper1PhaseProblem.h.

◆ _P_out

const Real& InterWrapper1PhaseProblem::_P_out
protectedinherited

◆ _P_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_P_soln
protectedinherited

◆ _P_tol

const Real& InterWrapper1PhaseProblem::_P_tol
protectedinherited

Convergence tolerance for the pressure loop in external solve.

Definition at line 131 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

◆ _pin_mesh_exist

const bool InterWrapper1PhaseProblem::_pin_mesh_exist
protectedinherited

◆ _prod

Vec InterWrapper1PhaseProblem::_prod
protectedinherited

◆ _prodp

Vec InterWrapper1PhaseProblem::_prodp
protectedinherited

◆ _q_prime_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_q_prime_soln
protectedinherited

◆ _rho_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_rho_soln
protectedinherited

◆ _rtol

const PetscReal& InterWrapper1PhaseProblem::_rtol
protectedinherited

◆ _S_flow_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_S_flow_soln
protectedinherited

◆ _segregated_bool

const bool InterWrapper1PhaseProblem::_segregated_bool
protectedinherited

◆ _staggered_pressure_bool

const bool InterWrapper1PhaseProblem::_staggered_pressure_bool
protectedinherited

Flag to define the usage of staggered or collocated pressure.

Definition at line 149 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::computeP(), and InterWrapper1PhaseProblem::computeWij().

◆ _subchannel_mesh

const InterWrapperMesh& InterWrapper1PhaseProblem::_subchannel_mesh
protectedinherited

◆ _SumWij_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_SumWij_soln
protectedinherited

◆ _T_maxit

const int& InterWrapper1PhaseProblem::_T_maxit
protectedinherited

Maximum iterations for the inner temperature loop.

Definition at line 135 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

◆ _T_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_T_soln
protectedinherited

◆ _T_tol

const Real& InterWrapper1PhaseProblem::_T_tol
protectedinherited

Convergence tolerance for the temperature loop in external solve.

Definition at line 133 of file InterWrapper1PhaseProblem.h.

Referenced by InterWrapper1PhaseProblem::externalSolve(), and TriInterWrapper1PhaseProblem::externalSolve().

◆ _Tpin_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_Tpin_soln
protectedinherited

◆ _TR

Real InterWrapper1PhaseProblem::_TR
protectedinherited

◆ _w_perim_soln

std::unique_ptr<SolutionHandle> InterWrapper1PhaseProblem::_w_perim_soln
protectedinherited

◆ _Wij

libMesh::DenseMatrix<Real>& InterWrapper1PhaseProblem::_Wij
protectedinherited

◆ _Wij_old

libMesh::DenseMatrix<Real> InterWrapper1PhaseProblem::_Wij_old
protectedinherited

◆ _Wij_residual_matrix

libMesh::DenseMatrix<Real> InterWrapper1PhaseProblem::_Wij_residual_matrix
protectedinherited

◆ _Wij_vec

Vec InterWrapper1PhaseProblem::_Wij_vec
protectedinherited

◆ _WijPrime

libMesh::DenseMatrix<Real> InterWrapper1PhaseProblem::_WijPrime
protectedinherited

◆ _z_grid

std::vector<Real> InterWrapper1PhaseProblem::_z_grid
protectedinherited

The documentation for this class was generated from the following files: