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

Kernel that implements the stress tensor and advection terms for the momentum equation. More...

#include <LinearWCNSFVMomentumFlux.h>

Inheritance diagram for LinearWCNSFVMomentumFlux:
[legend]

Public Types

enum  ResidualTagType { ResidualTagType::NonReference, ResidualTagType::Reference }
 
typedef DataFileName DataFileParameterType
 

Public Member Functions

 LinearWCNSFVMomentumFlux (const InputParameters &params)
 Class constructor. More...
 
virtual Real computeElemMatrixContribution () override
 
virtual Real computeNeighborMatrixContribution () override
 
virtual Real computeElemRightHandSideContribution () override
 
virtual Real computeNeighborRightHandSideContribution () override
 
virtual Real computeBoundaryMatrixContribution (const LinearFVBoundaryCondition &bc) override
 
virtual Real computeBoundaryRHSContribution (const LinearFVBoundaryCondition &bc) override
 
virtual void setupFaceData (const FaceInfo *face_info) override
 Set the current FaceInfo object. More...
 
virtual void addMatrixContribution () override
 
virtual void addRightHandSideContribution () override
 
virtual bool hasFaceSide (const FaceInfo &fi, bool fi_elem_side) const override
 
void setCurrentFaceArea (const Real area)
 
virtual const MooseLinearVariableFV< Real > & variable () const override
 
void linkTaggedVectorsAndMatrices (const std::set< TagID > &vector_tags, const std::set< TagID > &matrix_tags)
 
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
 
virtual void initialSetup ()
 
virtual void timestepSetup ()
 
virtual void jacobianSetup ()
 
virtual void residualSetup ()
 
virtual void subdomainSetup ()
 
virtual void customSetup (const ExecFlagType &)
 
const ExecFlagEnumgetExecuteOnEnum () const
 
const FunctiongetFunction (const std::string &name) const
 
const FunctiongetFunctionByName (const FunctionName &name) const
 
bool hasFunction (const std::string &param_name) const
 
bool hasFunctionByName (const FunctionName &name) const
 
UserObjectName getUserObjectName (const std::string &param_name) const
 
const T & getUserObject (const std::string &param_name, bool is_dependency=true) const
 
const T & getUserObjectByName (const UserObjectName &object_name, bool is_dependency=true) const
 
const UserObjectgetUserObjectBase (const std::string &param_name, bool is_dependency=true) const
 
const UserObjectgetUserObjectBaseByName (const UserObjectName &object_name, bool is_dependency=true) const
 
bool isImplicit ()
 
Moose::StateArg determineState () const
 
bool isDefaultPostprocessorValue (const std::string &param_name, const unsigned int index=0) const
 
bool hasPostprocessor (const std::string &param_name, const unsigned int index=0) const
 
bool hasPostprocessorByName (const PostprocessorName &name) const
 
std::size_t coupledPostprocessors (const std::string &param_name) const
 
const PostprocessorName & getPostprocessorName (const std::string &param_name, const unsigned int index=0) const
 
const VectorPostprocessorValuegetVectorPostprocessorValue (const std::string &param_name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValue (const std::string &param_name, const std::string &vector_name, bool needs_broadcast) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueByName (const VectorPostprocessorName &name, const std::string &vector_name, bool needs_broadcast) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOld (const std::string &param_name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOld (const std::string &param_name, const std::string &vector_name, bool needs_broadcast) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOldByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOldByName (const VectorPostprocessorName &name, const std::string &vector_name, bool needs_broadcast) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValue (const std::string &param_name, const std::string &vector_name) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValueByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValueOld (const std::string &param_name, const std::string &vector_name) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValueOldByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
bool hasVectorPostprocessor (const std::string &param_name, const std::string &vector_name) const
 
bool hasVectorPostprocessor (const std::string &param_name) const
 
bool hasVectorPostprocessorByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
bool hasVectorPostprocessorByName (const VectorPostprocessorName &name) const
 
const VectorPostprocessorName & getVectorPostprocessorName (const std::string &param_name) const
 
void setRandomResetFrequency (ExecFlagType exec_flag)
 
unsigned long getRandomLong () const
 
Real getRandomReal () const
 
unsigned int getSeed (std::size_t id)
 
unsigned int getMasterSeed () const
 
bool isNodal () const
 
ExecFlagType getResetOnTime () const
 
void setRandomDataPointer (RandomData *random_data)
 
virtual void meshChanged ()
 
void useVectorTag (const TagName &tag_name, VectorTagsKey)
 
void useVectorTag (TagID tag_id, VectorTagsKey)
 
void useMatrixTag (const TagName &tag_name, MatrixTagsKey)
 
void useMatrixTag (TagID tag_id, MatrixTagsKey)
 
bool isVectorTagged ()
 
bool isMatrixTagged ()
 
bool hasVectorTags () const
 
const std::set< TagID > & getVectorTags (VectorTagsKey) const
 
const std::set< TagID > & getMatrixTags (MatrixTagsKey) const
 
const std::vector< SubdomainName > & blocks () const
 
unsigned int numBlocks () const
 
virtual const std::set< SubdomainID > & blockIDs () const
 
unsigned int blocksMaxDimension () const
 
bool hasBlocks (const SubdomainName &name) const
 
bool hasBlocks (const std::vector< SubdomainName > &names) const
 
bool hasBlocks (const std::set< SubdomainName > &names) const
 
bool hasBlocks (SubdomainID id) const
 
bool hasBlocks (const std::vector< SubdomainID > &ids) const
 
bool hasBlocks (const std::set< SubdomainID > &ids) const
 
bool isBlockSubset (const std::set< SubdomainID > &ids) const
 
bool isBlockSubset (const std::vector< SubdomainID > &ids) const
 
bool hasBlockMaterialProperty (const std::string &prop_name)
 
const std::set< SubdomainID > & meshBlockIDs () const
 
virtual bool blockRestricted () const
 
virtual void checkVariable (const MooseVariableFieldBase &variable) const
 
MooseVariableBasemooseVariableBase () const
 
MooseVariableField< Real > & mooseVariableField ()
 
MooseVariableFE< Real > * mooseVariable () const
 
MooseVariableFV< Real > * mooseVariableFV () const
 
MooseLinearVariableFV< Real > * mooseLinearVariableFV () const
 
const std::set< MooseVariableFieldBase *> & getMooseVariableDependencies () const
 
std::set< MooseVariableFieldBase *> checkAllVariables (const DofObjectType &dof_object, const std::set< MooseVariableFieldBase * > &vars_to_omit={})
 
std::set< MooseVariableFieldBase *> checkVariables (const DofObjectType &dof_object, const std::set< MooseVariableFieldBase * > &vars_to_check)
 
void addMooseVariableDependency (MooseVariableFieldBase *var)
 
void addMooseVariableDependency (const std::vector< MooseVariableFieldBase * > &vars)
 
Moose::FaceArg makeFace (const FaceInfo &fi, const Moose::FV::LimiterType limiter_type, const bool elem_is_upwind, const bool correct_skewness=false, const Moose::StateArg *state_limiter=nullptr) const
 
Moose::FaceArg makeCDFace (const FaceInfo &fi, const bool correct_skewness=false) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
const PostprocessorValuegetPostprocessorValue (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValue (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOld (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOld (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOlder (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOlder (const std::string &param_name, const unsigned int index=0) const
 
virtual const PostprocessorValuegetPostprocessorValueByName (const PostprocessorName &name) const
 
virtual const PostprocessorValuegetPostprocessorValueByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOldByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOldByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOlderByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOlderByName (const PostprocessorName &name) const
 
bool isVectorPostprocessorDistributed (const std::string &param_name) const
 
bool isVectorPostprocessorDistributed (const std::string &param_name) const
 
bool isVectorPostprocessorDistributedByName (const VectorPostprocessorName &name) const
 
bool isVectorPostprocessorDistributedByName (const VectorPostprocessorName &name) const
 
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 std::string deduceFunctorName (const std::string &name, const InputParameters &params)
 
static void setRMParamsAdvection (const InputParameters &obj_params, InputParameters &rm_params, const unsigned short conditional_extended_layers)
 
static void setRMParamsDiffusion (const InputParameters &obj_params, InputParameters &rm_params, const unsigned short conditional_extended_layers)
 

Public Attributes

const ConsoleStream _console
 

Protected Member Functions

Real computeInternalAdvectionElemMatrixContribution ()
 Computes the matrix contribution of the advective flux on the element side of current face when the face is an internal face (doesn't have associated boundary conditions). More...
 
Real computeInternalAdvectionNeighborMatrixContribution ()
 Computes the matrix contribution of the advective flux on the neighbor side of current face when the face is an internal face (doesn't have associated boundary conditions). More...
 
Real computeInternalStressMatrixContribution ()
 Computes the matrix contribution of the stress term on the current face when the face is an internal face (doesn't have associated boundary conditions). More...
 
Real computeInternalStressRHSContribution ()
 Computes the right hand side contribution of the stress term on the current face when the face is an internal face (doesn't have associated boundary conditions). More...
 
Real computeStressBoundaryMatrixContribution (const LinearFVAdvectionDiffusionBC *bc)
 Computes the matrix contributions of the boundary conditions resulting from the stress tensor. More...
 
Real computeStressBoundaryRHSContribution (const LinearFVAdvectionDiffusionBC *bc)
 Computes the right hand side contributions of the boundary conditions resulting from the stress tensor. More...
 
Real computeAdvectionBoundaryMatrixContribution (const LinearFVAdvectionDiffusionBC *bc)
 Computes the matrix contributions of the boundary conditions resulting from the advection term. More...
 
Real computeAdvectionBoundaryRHSContribution (const LinearFVAdvectionDiffusionBC *bc)
 Computes the right hand side contributions of the boundary conditions resulting from the advection term. More...
 
std::string deduceFunctorName (const std::string &name) const
 
Moose::FaceArg singleSidedFaceArg (const FaceInfo *fi, Moose::FV::LimiterType limiter_type=Moose::FV::LimiterType::CentralDifference, bool correct_skewness=false) const
 
void requestVariableCellGradient (const std::string &variable_name)
 
virtual void addUserObjectDependencyHelper (const UserObject &) const
 
virtual void addPostprocessorDependencyHelper (const PostprocessorName &) const
 
virtual void addVectorPostprocessorDependencyHelper (const VectorPostprocessorName &) 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
 
void prepareVectorTag (Assembly &assembly, unsigned int ivar)
 
void prepareVectorTag (Assembly &assembly, unsigned int ivar, ResidualTagType tag_type)
 
void prepareVectorTag (Assembly &assembly, unsigned int ivar, ResidualTagType tag_type)
 
void prepareVectorTag (Assembly &assembly, unsigned int ivar, ResidualTagType tag_type)
 
void prepareVectorTag (Assembly &assembly, unsigned int ivar, ResidualTagType tag_type)
 
void prepareVectorTagNeighbor (Assembly &assembly, unsigned int ivar)
 
void prepareVectorTagLower (Assembly &assembly, unsigned int ivar)
 
void prepareMatrixTag (Assembly &assembly, unsigned int ivar, unsigned int jvar)
 
void prepareMatrixTag (Assembly &assembly, unsigned int ivar, unsigned int jvar, DenseMatrix< Number > &k) const
 
void prepareMatrixTagNonlocal (Assembly &assembly, unsigned int ivar, unsigned int jvar)
 
void prepareMatrixTagNeighbor (Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type)
 
void prepareMatrixTagNeighbor (Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type, DenseMatrix< Number > &k) const
 
void prepareMatrixTagLower (Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::ConstraintJacobianType type)
 
void accumulateTaggedLocalResidual ()
 
void assignTaggedLocalResidual ()
 
void accumulateTaggedLocalMatrix ()
 
void accumulateTaggedLocalMatrix (Assembly &assembly, unsigned int ivar, unsigned int jvar, const DenseMatrix< Number > &k)
 
void accumulateTaggedLocalMatrix (Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type, const DenseMatrix< Number > &k)
 
void accumulateTaggedNonlocalMatrix ()
 
void assignTaggedLocalMatrix ()
 
void addResiduals (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addResiduals (Assembly &assembly, const DenseVector< T > &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addResiduals (Assembly &assembly, const ADResidualsPacket &packet)
 
void addResidualsAndJacobian (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addResidualsAndJacobian (Assembly &assembly, const ADResidualsPacket &packet)
 
void addJacobian (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addJacobian (Assembly &assembly, const ADResidualsPacket &packet)
 
void addJacobian (Assembly &assembly, DenseMatrix< Real > &local_k, const std::vector< dof_id_type > &row_indices, const std::vector< dof_id_type > &column_indices, Real scaling_factor)
 
void addResidualsWithoutConstraints (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addResidualsAndJacobianWithoutConstraints (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addJacobianWithoutConstraints (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 
void addJacobianElement (Assembly &assembly, Real value, dof_id_type row_index, dof_id_type column_index, Real scaling_factor)
 
void setResidual (SystemBase &sys, const T &residual, MooseVariableFE< T > &var)
 
void setResidual (SystemBase &sys, Real residual, dof_id_type dof_index)
 
void setResidual (SystemBase &sys, SetResidualFunctor set_residual_functor)
 
virtual bool hasBlockMaterialPropertyHelper (const std::string &prop_name)
 
void initializeBlockRestrictable (const MooseObject *moose_object)
 
Moose::CoordinateSystemType getBlockCoordSystem ()
 
const Moose::Functor< T > & getFunctor (const std::string &name)
 
const Moose::Functor< T > & getFunctor (const std::string &name, THREAD_ID tid)
 
const Moose::Functor< T > & getFunctor (const std::string &name, SubProblem &subproblem)
 
const Moose::Functor< T > & getFunctor (const std::string &name, SubProblem &subproblem, THREAD_ID tid)
 
bool isFunctor (const std::string &name) const
 
bool isFunctor (const std::string &name, const SubProblem &subproblem) const
 
Moose::ElemArg makeElemArg (const Elem *elem, bool correct_skewnewss=false) const
 
void checkFunctorSupportsSideIntegration (const std::string &name, bool qp_integration)
 
virtual const OutputTools< Real >::VariableValuevalue ()
 
virtual const OutputTools< Real >::VariableValuevalueOld ()
 
virtual const OutputTools< Real >::VariableValuevalueOlder ()
 
virtual const OutputTools< Real >::VariableValuedot ()
 
virtual const OutputTools< Real >::VariableValuedotDot ()
 
virtual const OutputTools< Real >::VariableValuedotOld ()
 
virtual const OutputTools< Real >::VariableValuedotDotOld ()
 
virtual const VariableValuedotDu ()
 
virtual const VariableValuedotDotDu ()
 
virtual const OutputTools< Real >::VariableGradientgradient ()
 
virtual const OutputTools< Real >::VariableGradientgradientOld ()
 
virtual const OutputTools< Real >::VariableGradientgradientOlder ()
 
virtual const OutputTools< Real >::VariableSecondsecond ()
 
virtual const OutputTools< Real >::VariableSecondsecondOld ()
 
virtual const OutputTools< Real >::VariableSecondsecondOlder ()
 
virtual const OutputTools< Real >::VariableTestSecondsecondTest ()
 
virtual const OutputTools< Real >::VariableTestSecondsecondTestFace ()
 
virtual const OutputTools< Real >::VariablePhiSecondsecondPhi ()
 
virtual const OutputTools< Real >::VariablePhiSecondsecondPhiFace ()
 

Protected Attributes

const unsigned int _dim
 The dimension of the mesh. More...
 
const RhieChowMassFlux_mass_flux_provider
 The Rhie-Chow user object that provides us with the face velocity. More...
 
const Moose::Functor< Real > & _mu
 The functor for the dynamic viscosity. More...
 
const bool _use_nonorthogonal_correction
 Switch to enable/disable nonorthogonal correction in the stress term. More...
 
const bool _use_deviatoric_terms
 Switch to enable/disable deviatoric parts in the stress term. More...
 
std::pair< Real, Real_advected_interp_coeffs
 Container for the current advected interpolation coefficients on the face to make sure we don't compute it multiple times for different terms. More...
 
Real _face_mass_flux
 Container for the mass flux on the face which will be reused in the advection term's matrix and right hand side contribution. More...
 
Real _boundary_normal_factor
 Multiplier that ensures the normal of the boundary always points outwards, even in cases when the boundary is within the mesh. More...
 
Real _stress_matrix_contribution
 The cached matrix contribution. More...
 
Real _stress_rhs_contribution
 The cached right hand side contribution. More...
 
Moose::FV::InterpMethod _advected_interp_method
 The interpolation method to use for the advected quantity. More...
 
const unsigned int _index
 Index x|y|z, this is mainly to handle the deviatoric parts correctly in in the stress term. More...
 
const MooseLinearVariableFVReal *const _u_var
 Velocity in direction x. More...
 
const MooseLinearVariableFVReal *const _v_var
 Velocity in direction y. More...
 
const MooseLinearVariableFVReal *const _w_var
 Velocity in direction z. More...
 
const FaceInfo_current_face_info
 
Real _current_face_area
 
FaceInfo::VarFaceNeighbors _current_face_type
 
bool _cached_matrix_contribution
 
bool _cached_rhs_contribution
 
const bool _force_boundary_execution
 
DenseVector< dof_id_type_dof_indices
 
DenseMatrix< Real_matrix_contribution
 
DenseVector< Real_rhs_contribution
 
MooseLinearVariableFV< Real > & _var
 
const unsigned int _var_num
 
const unsigned int _sys_num
 
FEProblemBase_fe_problem
 
SystemBase_sys
 
libMesh::LinearImplicitSystem_linear_system
 
const THREAD_ID _tid
 
MooseMesh_mesh
 
std::vector< NumericVector< Number > *> _vectors
 
std::vector< SparseMatrix< Number > *> _matrices
 
const bool & _enabled
 
MooseApp_app
 
const std::string _type
 
const std::string _name
 
const InputParameters_pars
 
Factory_factory
 
ActionFactory_action_factory
 
const ExecFlagEnum_execute_enum
 
const ExecFlagType_current_execute_flag
 
const InputParameters_ti_params
 
FEProblemBase_ti_feproblem
 
bool _is_implicit
 
Real_t
 
const Real_t_old
 
int_t_step
 
Real_dt
 
Real_dt_old
 
bool _is_transient
 
MooseApp_restartable_app
 
const std::string _restartable_system_name
 
const THREAD_ID _restartable_tid
 
const bool _restartable_read_only
 
FEProblemBase_mci_feproblem
 
SubProblem_subproblem
 
DenseVector< Number_local_re
 
DenseMatrix< Number_local_ke
 
DenseMatrix< Number_nonlocal_ke
 
const MaterialData_blk_material_data
 
bool _nodal
 
MooseVariableFE< Real > * _variable
 
MooseVariableFV< Real > * _fv_variable
 
MooseLinearVariableFV< Real > * _linear_fv_variable
 
MooseVariableField< Real > * _field_variable
 
Assembly_mvi_assembly
 
const Parallel::Communicator & _communicator
 

Detailed Description

Kernel that implements the stress tensor and advection terms for the momentum equation.

Definition at line 22 of file LinearWCNSFVMomentumFlux.h.

Constructor & Destructor Documentation

◆ LinearWCNSFVMomentumFlux()

LinearWCNSFVMomentumFlux::LinearWCNSFVMomentumFlux ( const InputParameters params)

Class constructor.

Parameters
paramsThe InputParameters for the kernel.

Definition at line 49 of file LinearWCNSFVMomentumFlux.C.

50  : LinearFVFluxKernel(params),
52  _mass_flux_provider(getUserObject<RhieChowMassFlux>("rhie_chow_user_object")),
53  _mu(getFunctor<Real>(getParam<MooseFunctorName>(NS::mu))),
54  _use_nonorthogonal_correction(getParam<bool>("use_nonorthogonal_correction")),
55  _use_deviatoric_terms(getParam<bool>("use_deviatoric_terms")),
56  _advected_interp_coeffs(std::make_pair<Real, Real>(0, 0)),
57  _face_mass_flux(0.0),
61  _index(getParam<MooseEnum>("momentum_component")),
62  _u_var(dynamic_cast<const MooseLinearVariableFVReal *>(
63  &_fe_problem.getVariable(_tid, getParam<SolverVariableName>("u")))),
64  _v_var(params.isParamValid("v")
65  ? dynamic_cast<const MooseLinearVariableFVReal *>(
66  &_fe_problem.getVariable(_tid, getParam<SolverVariableName>("v")))
67  : nullptr),
68  _w_var(params.isParamValid("w")
69  ? dynamic_cast<const MooseLinearVariableFVReal *>(
70  &_fe_problem.getVariable(_tid, getParam<SolverVariableName>("w")))
71  : nullptr)
72 {
73  // We only need gradients if the nonorthogonal correction is enabled or when we request the
74  // computation of the deviatoric parts of the stress tensor.
77 
78  Moose::FV::setInterpolationMethod(*this, _advected_interp_method, "advected_interp_method");
79 
80  if (!_u_var)
81  paramError("u", "the u velocity must be a MooseLinearVariableFVReal.");
82 
83  if (_dim >= 2 && !_v_var)
84  paramError("v",
85  "In two or more dimensions, the v velocity must be supplied and it must be a "
86  "MooseLinearVariableFVReal.");
87 
88  if (_dim >= 3 && !_w_var)
89  paramError("w",
90  "In three-dimensions, the w velocity must be supplied and it must be a "
91  "MooseLinearVariableFVReal.");
92 }
virtual MooseMesh & mesh()=0
SubProblem & _subproblem
const unsigned int _index
Index x|y|z, this is mainly to handle the deviatoric parts correctly in in the stress term...
const MooseLinearVariableFVReal *const _w_var
Velocity in direction z.
MooseLinearVariableFV< Real > & _var
const bool _use_nonorthogonal_correction
Switch to enable/disable nonorthogonal correction in the stress term.
const RhieChowMassFlux & _mass_flux_provider
The Rhie-Chow user object that provides us with the face velocity.
Real _boundary_normal_factor
Multiplier that ensures the normal of the boundary always points outwards, even in cases when the bou...
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
const Moose::Functor< Real > & _mu
The functor for the dynamic viscosity.
const THREAD_ID _tid
virtual unsigned int dimension() const
static const std::string mu
Definition: NS.h:123
std::pair< Real, Real > _advected_interp_coeffs
Container for the current advected interpolation coefficients on the face to make sure we don&#39;t compu...
Real _stress_rhs_contribution
The cached right hand side contribution.
void paramError(const std::string &param, Args... args) const
const unsigned int _dim
The dimension of the mesh.
Real _stress_matrix_contribution
The cached matrix contribution.
LinearFVFluxKernel(const InputParameters &params)
Real _face_mass_flux
Container for the mass flux on the face which will be reused in the advection term&#39;s matrix and right...
const MooseLinearVariableFVReal *const _v_var
Velocity in direction y.
Moose::FV::InterpMethod _advected_interp_method
The interpolation method to use for the advected quantity.
FEProblemBase & _fe_problem
bool setInterpolationMethod(const MooseObject &obj, Moose::FV::InterpMethod &interp_method, const std::string &param_name)
const bool _use_deviatoric_terms
Switch to enable/disable deviatoric parts in the stress term.
const MooseLinearVariableFVReal *const _u_var
Velocity in direction x.
bool isParamValid(const std::string &name) const

Member Function Documentation

◆ computeAdvectionBoundaryMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeAdvectionBoundaryMatrixContribution ( const LinearFVAdvectionDiffusionBC bc)
protected

Computes the matrix contributions of the boundary conditions resulting from the advection term.

Parameters
bcThe boundary condition whose contributions should be used

Definition at line 360 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeBoundaryMatrixContribution().

362 {
363  const auto boundary_value_matrix_contrib = bc->computeBoundaryValueMatrixContribution();
364  return boundary_value_matrix_contrib * _face_mass_flux;
365 }
virtual Real computeBoundaryValueMatrixContribution() const=0
Real _face_mass_flux
Container for the mass flux on the face which will be reused in the advection term&#39;s matrix and right...

◆ computeAdvectionBoundaryRHSContribution()

Real LinearWCNSFVMomentumFlux::computeAdvectionBoundaryRHSContribution ( const LinearFVAdvectionDiffusionBC bc)
protected

Computes the right hand side contributions of the boundary conditions resulting from the advection term.

Parameters
bcThe boundary condition whose contributions should be used

Definition at line 368 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeBoundaryRHSContribution().

370 {
371  const auto boundary_value_rhs_contrib = bc->computeBoundaryValueRHSContribution();
372  return -boundary_value_rhs_contrib * _face_mass_flux;
373 }
Real _face_mass_flux
Container for the mass flux on the face which will be reused in the advection term&#39;s matrix and right...
virtual Real computeBoundaryValueRHSContribution() const=0

◆ computeBoundaryMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeBoundaryMatrixContribution ( const LinearFVBoundaryCondition bc)
overridevirtual

Implements LinearFVFluxKernel.

Definition at line 123 of file LinearWCNSFVMomentumFlux.C.

124 {
125  const auto * const adv_diff_bc = static_cast<const LinearFVAdvectionDiffusionBC *>(&bc);
126  mooseAssert(adv_diff_bc, "This should be a valid BC!");
127  return (computeStressBoundaryMatrixContribution(adv_diff_bc) +
130 }
Real computeAdvectionBoundaryMatrixContribution(const LinearFVAdvectionDiffusionBC *bc)
Computes the matrix contributions of the boundary conditions resulting from the advection term...
Real computeStressBoundaryMatrixContribution(const LinearFVAdvectionDiffusionBC *bc)
Computes the matrix contributions of the boundary conditions resulting from the stress tensor...

◆ computeBoundaryRHSContribution()

Real LinearWCNSFVMomentumFlux::computeBoundaryRHSContribution ( const LinearFVBoundaryCondition bc)
overridevirtual

Implements LinearFVFluxKernel.

Definition at line 133 of file LinearWCNSFVMomentumFlux.C.

134 {
135  const auto * const adv_diff_bc = static_cast<const LinearFVAdvectionDiffusionBC *>(&bc);
136 
137  mooseAssert(adv_diff_bc, "This should be a valid BC!");
138  return (computeStressBoundaryRHSContribution(adv_diff_bc) +
141 }
Real computeStressBoundaryRHSContribution(const LinearFVAdvectionDiffusionBC *bc)
Computes the right hand side contributions of the boundary conditions resulting from the stress tenso...
Real computeAdvectionBoundaryRHSContribution(const LinearFVAdvectionDiffusionBC *bc)
Computes the right hand side contributions of the boundary conditions resulting from the advection te...

◆ computeElemMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeElemMatrixContribution ( )
overridevirtual

Implements LinearFVFluxKernel.

Definition at line 95 of file LinearWCNSFVMomentumFlux.C.

96 {
100 }
Real computeInternalAdvectionElemMatrixContribution()
Computes the matrix contribution of the advective flux on the element side of current face when the f...
Real computeInternalStressMatrixContribution()
Computes the matrix contribution of the stress term on the current face when the face is an internal ...

◆ computeElemRightHandSideContribution()

Real LinearWCNSFVMomentumFlux::computeElemRightHandSideContribution ( )
overridevirtual

Implements LinearFVFluxKernel.

Definition at line 111 of file LinearWCNSFVMomentumFlux.C.

112 {
114 }
Real computeInternalStressRHSContribution()
Computes the right hand side contribution of the stress term on the current face when the face is an ...

◆ computeInternalAdvectionElemMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeInternalAdvectionElemMatrixContribution ( )
protected

Computes the matrix contribution of the advective flux on the element side of current face when the face is an internal face (doesn't have associated boundary conditions).

Definition at line 144 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeElemMatrixContribution().

145 {
147 }
std::pair< Real, Real > _advected_interp_coeffs
Container for the current advected interpolation coefficients on the face to make sure we don&#39;t compu...
Real _face_mass_flux
Container for the mass flux on the face which will be reused in the advection term&#39;s matrix and right...

◆ computeInternalAdvectionNeighborMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeInternalAdvectionNeighborMatrixContribution ( )
protected

Computes the matrix contribution of the advective flux on the neighbor side of current face when the face is an internal face (doesn't have associated boundary conditions).

Definition at line 150 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeNeighborMatrixContribution().

151 {
153 }
std::pair< Real, Real > _advected_interp_coeffs
Container for the current advected interpolation coefficients on the face to make sure we don&#39;t compu...
Real _face_mass_flux
Container for the mass flux on the face which will be reused in the advection term&#39;s matrix and right...

◆ computeInternalStressMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeInternalStressMatrixContribution ( )
protected

Computes the matrix contribution of the stress term on the current face when the face is an internal face (doesn't have associated boundary conditions).

Definition at line 156 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeElemMatrixContribution(), and computeNeighborMatrixContribution().

157 {
158  // If we don't have the value yet, we compute it
160  {
161  const auto face_arg = makeCDFace(*_current_face_info);
162 
163  // If we requested nonorthogonal correction, we use the normal component of the
164  // cell to face vector.
165  const auto d = _use_nonorthogonal_correction
166  ? std::abs(_current_face_info->dCN() * _current_face_info->normal())
168 
169  // Cache the matrix contribution
172  }
173 
175 }
Moose::StateArg determineState() const
const bool _use_nonorthogonal_correction
Switch to enable/disable nonorthogonal correction in the stress term.
const Moose::Functor< Real > & _mu
The functor for the dynamic viscosity.
const FaceInfo * _current_face_info
const Point & normal() const
Real dCNMag() const
Real _stress_matrix_contribution
The cached matrix contribution.
bool _cached_matrix_contribution
Moose::FaceArg makeCDFace(const FaceInfo &fi, const bool correct_skewness=false) const
const Point & dCN() const

◆ computeInternalStressRHSContribution()

Real LinearWCNSFVMomentumFlux::computeInternalStressRHSContribution ( )
protected

Computes the right hand side contribution of the stress term on the current face when the face is an internal face (doesn't have associated boundary conditions).

Definition at line 178 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeElemRightHandSideContribution(), and computeNeighborRightHandSideContribution().

179 {
180  // We can have contributions to the right hand side in two occasions:
181  // (1) when we use nonorthogonal correction for the normal gradients
182  // (2) when we request the deviatoric parts of the stress tensor. (needed for space-dependent
183  // viscosities for example)
185  {
186  // scenario (1), we need to add the nonorthogonal correction. In 1D, we don't have
187  // any correction so we just skip this part
189  {
190  const auto face_arg = makeCDFace(*_current_face_info);
191  const auto state_arg = determineState();
192 
193  // Get the gradients from the adjacent cells
194  const auto grad_elem = _var.gradSln(*_current_face_info->elemInfo());
195  const auto & grad_neighbor = _var.gradSln(*_current_face_info->neighborInfo());
196 
197  // Interpolate the two gradients to the face
198  const auto interp_coeffs =
200 
201  const auto correction_vector =
205 
206  // Cache the matrix contribution
208  _mu(face_arg, state_arg) *
209  (interp_coeffs.first * grad_elem + interp_coeffs.second * grad_neighbor) *
210  correction_vector;
211  }
212  // scenario (2), we will have to account for the deviatoric parts of the stress tensor.
214  {
215  // Interpolate the two gradients to the face
216  const auto interp_coeffs =
218 
219  const auto u_grad_elem = _u_var->gradSln(*_current_face_info->elemInfo());
220  const auto u_grad_neighbor = _u_var->gradSln(*_current_face_info->neighborInfo());
221 
222  Real trace_elem = 0;
223  Real trace_neighbor = 0;
224  RealVectorValue deviatoric_vector_elem;
225  RealVectorValue deviatoric_vector_neighbor;
226 
227  deviatoric_vector_elem(0) = u_grad_elem(_index);
228  deviatoric_vector_neighbor(0) = u_grad_neighbor(_index);
229  trace_elem += u_grad_elem(0);
230  trace_neighbor += u_grad_neighbor(0);
231 
232  const auto face_arg = makeCDFace(*_current_face_info);
233  const auto state_arg = determineState();
234 
235  if (_dim > 1)
236  {
237  const auto v_grad_elem = _v_var->gradSln(*_current_face_info->elemInfo());
238  const auto v_grad_neighbor = _v_var->gradSln(*_current_face_info->neighborInfo());
239 
240  deviatoric_vector_elem(1) = v_grad_elem(_index);
241  deviatoric_vector_neighbor(1) = v_grad_neighbor(_index);
242  trace_elem += v_grad_elem(1);
243  trace_neighbor += v_grad_neighbor(1);
244  if (_dim > 2)
245  {
246  const auto w_grad_elem = _w_var->gradSln(*_current_face_info->elemInfo());
247  const auto w_grad_neighbor = _w_var->gradSln(*_current_face_info->neighborInfo());
248 
249  deviatoric_vector_elem(2) = w_grad_elem(_index);
250  deviatoric_vector_neighbor(2) = w_grad_neighbor(_index);
251  trace_elem += w_grad_elem(2);
252  trace_neighbor += w_grad_neighbor(2);
253  }
254  }
255  deviatoric_vector_elem(_index) -= 2 / 3 * trace_elem;
256  deviatoric_vector_neighbor(_index) -= 2 / 3 * trace_neighbor;
257 
258  _stress_rhs_contribution += _mu(face_arg, state_arg) *
259  (interp_coeffs.first * deviatoric_vector_elem +
260  interp_coeffs.second * deviatoric_vector_neighbor) *
262  }
264  }
265 
267 }
const unsigned int _index
Index x|y|z, this is mainly to handle the deviatoric parts correctly in in the stress term...
const MooseLinearVariableFVReal *const _w_var
Velocity in direction z.
std::pair< Real, Real > interpCoeffs(const InterpMethod m, const FaceInfo &fi, const bool one_is_elem, const T &face_flux=0.0)
Moose::StateArg determineState() const
const ElemInfo * neighborInfo() const
MooseLinearVariableFV< Real > & _var
const bool _use_nonorthogonal_correction
Switch to enable/disable nonorthogonal correction in the stress term.
const ElemInfo * elemInfo() const
const Moose::Functor< Real > & _mu
The functor for the dynamic viscosity.
const FaceInfo * _current_face_info
Real _stress_rhs_contribution
The cached right hand side contribution.
const Point & normal() const
const unsigned int _dim
The dimension of the mesh.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Point & eCN() const
const MooseLinearVariableFVReal *const _v_var
Velocity in direction y.
const VectorValue< Real > gradSln(const ElemInfo &elem_info) const
Moose::FaceArg makeCDFace(const FaceInfo &fi, const bool correct_skewness=false) const
const bool _use_deviatoric_terms
Switch to enable/disable deviatoric parts in the stress term.
const MooseLinearVariableFVReal *const _u_var
Velocity in direction x.

◆ computeNeighborMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeNeighborMatrixContribution ( )
overridevirtual

Implements LinearFVFluxKernel.

Definition at line 103 of file LinearWCNSFVMomentumFlux.C.

104 {
108 }
Real computeInternalAdvectionNeighborMatrixContribution()
Computes the matrix contribution of the advective flux on the neighbor side of current face when the ...
Real computeInternalStressMatrixContribution()
Computes the matrix contribution of the stress term on the current face when the face is an internal ...

◆ computeNeighborRightHandSideContribution()

Real LinearWCNSFVMomentumFlux::computeNeighborRightHandSideContribution ( )
overridevirtual

Implements LinearFVFluxKernel.

Definition at line 117 of file LinearWCNSFVMomentumFlux.C.

118 {
120 }
Real computeInternalStressRHSContribution()
Computes the right hand side contribution of the stress term on the current face when the face is an ...

◆ computeStressBoundaryMatrixContribution()

Real LinearWCNSFVMomentumFlux::computeStressBoundaryMatrixContribution ( const LinearFVAdvectionDiffusionBC bc)
protected

Computes the matrix contributions of the boundary conditions resulting from the stress tensor.

Parameters
bcThe boundary condition whose contributions should be used

Definition at line 270 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeBoundaryMatrixContribution().

272 {
273  auto grad_contrib = bc->computeBoundaryGradientMatrixContribution();
274  // If the boundary condition does not include the diffusivity contribution then
275  // add it here.
277  {
278  const auto face_arg = singleSidedFaceArg(_current_face_info);
279  grad_contrib *= _mu(face_arg, determineState());
280  }
281 
282  return grad_contrib;
283 }
Moose::FaceArg singleSidedFaceArg(const FaceInfo *fi, Moose::FV::LimiterType limiter_type=Moose::FV::LimiterType::CentralDifference, bool correct_skewness=false) const
Moose::StateArg determineState() const
virtual Real computeBoundaryGradientMatrixContribution() const=0
const Moose::Functor< Real > & _mu
The functor for the dynamic viscosity.
const FaceInfo * _current_face_info
virtual bool includesMaterialPropertyMultiplier() const

◆ computeStressBoundaryRHSContribution()

Real LinearWCNSFVMomentumFlux::computeStressBoundaryRHSContribution ( const LinearFVAdvectionDiffusionBC bc)
protected

Computes the right hand side contributions of the boundary conditions resulting from the stress tensor.

Parameters
bcThe boundary condition whose contributions should be used

Definition at line 286 of file LinearWCNSFVMomentumFlux.C.

Referenced by computeBoundaryRHSContribution().

288 {
289  const auto face_arg = singleSidedFaceArg(_current_face_info);
290  auto grad_contrib = bc->computeBoundaryGradientRHSContribution();
291 
292  // If the boundary condition does not include the diffusivity contribution then
293  // add it here.
295  grad_contrib *= _mu(face_arg, determineState());
296 
297  // We add the nonorthogonal corrector for the face here. Potential idea: we could do
298  // this in the boundary condition too. For now, however, we keep it like this.
300  {
301  // We support internal boundaries as well. In that case we have to decide on which side
302  // of the boundary we are on.
303  const auto elem_info = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM)
306 
307  // Unit vector to the boundary. Unfortunately, we have to recompute it because the value
308  // stored in the face info is only correct for external boundaries
309  const auto e_Cf = _current_face_info->faceCentroid() - elem_info->centroid();
310  const auto correction_vector =
311  _current_face_info->normal() - 1 / (_current_face_info->normal() * e_Cf) * e_Cf;
312 
313  grad_contrib += _mu(face_arg, determineState()) * _var.gradSln(*elem_info) *
314  _boundary_normal_factor * correction_vector;
315  }
316 
318  {
319  // We might be on a face which is an internal boundary so we want to make sure we
320  // get the gradient from the right side.
321  const auto elem_info = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM)
324 
325  const auto u_grad_elem = _u_var->gradSln(*elem_info);
326 
327  Real trace_elem = 0;
328  RealVectorValue frace_grad_approx;
329 
330  trace_elem += u_grad_elem(0);
331  frace_grad_approx(0) = u_grad_elem(_index);
332 
333  const auto state_arg = determineState();
334 
335  if (_dim > 1)
336  {
337  const auto v_grad_elem = _v_var->gradSln(*elem_info);
338  trace_elem += v_grad_elem(1);
339  frace_grad_approx(1) = v_grad_elem(_index);
340 
341  if (_dim > 2)
342  {
343  const auto w_grad_elem = _w_var->gradSln(*elem_info);
344  trace_elem += w_grad_elem(2);
345  frace_grad_approx(2) = w_grad_elem(_index);
346  }
347  }
348 
349  frace_grad_approx(_index) -= 2 / 3 * trace_elem;
350 
351  // We support internal boundaries too so we have to make sure the normal points always outward
352  grad_contrib += _mu(face_arg, state_arg) * frace_grad_approx * _boundary_normal_factor *
354  }
355 
356  return grad_contrib;
357 }
const unsigned int _index
Index x|y|z, this is mainly to handle the deviatoric parts correctly in in the stress term...
const MooseLinearVariableFVReal *const _w_var
Velocity in direction z.
virtual Real computeBoundaryGradientRHSContribution() const=0
Moose::FaceArg singleSidedFaceArg(const FaceInfo *fi, Moose::FV::LimiterType limiter_type=Moose::FV::LimiterType::CentralDifference, bool correct_skewness=false) const
Moose::StateArg determineState() const
const ElemInfo * neighborInfo() const
const Point & faceCentroid() const
MooseLinearVariableFV< Real > & _var
const bool _use_nonorthogonal_correction
Switch to enable/disable nonorthogonal correction in the stress term.
const ElemInfo * elemInfo() const
FaceInfo::VarFaceNeighbors _current_face_type
Real _boundary_normal_factor
Multiplier that ensures the normal of the boundary always points outwards, even in cases when the bou...
const Moose::Functor< Real > & _mu
The functor for the dynamic viscosity.
const FaceInfo * _current_face_info
const Point & normal() const
virtual bool includesMaterialPropertyMultiplier() const
const unsigned int _dim
The dimension of the mesh.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const MooseLinearVariableFVReal *const _v_var
Velocity in direction y.
const VectorValue< Real > gradSln(const ElemInfo &elem_info) const
const bool _use_deviatoric_terms
Switch to enable/disable deviatoric parts in the stress term.
const MooseLinearVariableFVReal *const _u_var
Velocity in direction x.

◆ setupFaceData()

void LinearWCNSFVMomentumFlux::setupFaceData ( const FaceInfo face_info)
overridevirtual

Set the current FaceInfo object.

We override this here to make sure the face velocity evaluation happens only once and that it can be reused for the matrix and right hand side contributions.

Parameters
face_infoThe face info which will be used as current face info

Reimplemented from LinearFVFluxKernel.

Definition at line 376 of file LinearWCNSFVMomentumFlux.C.

377 {
379 
380  // Multiplier that ensures the normal of the boundary always points outwards, even in cases
381  // when the boundary is within the mesh.
383 
384  // Caching the mass flux on the face which will be reused in the advection term's matrix and
385  // right hand side contributions
387 
388  // Caching the interpolation coefficients so they will be reused for the matrix and right hand
389  // side terms
392 
393  // We'll have to set this to zero to make sure that we don't accumulate values over multiple
394  // faces. The matrix contribution should be fine.
396 }
std::pair< Real, Real > interpCoeffs(const InterpMethod m, const FaceInfo &fi, const bool one_is_elem, const T &face_flux=0.0)
Real getMassFlux(const FaceInfo &fi) const
Get the face velocity times density (used in advection terms)
virtual void setupFaceData(const FaceInfo *face_info)
const RhieChowMassFlux & _mass_flux_provider
The Rhie-Chow user object that provides us with the face velocity.
FaceInfo::VarFaceNeighbors _current_face_type
Real _boundary_normal_factor
Multiplier that ensures the normal of the boundary always points outwards, even in cases when the bou...
const FaceInfo * _current_face_info
std::pair< Real, Real > _advected_interp_coeffs
Container for the current advected interpolation coefficients on the face to make sure we don&#39;t compu...
Real _stress_rhs_contribution
The cached right hand side contribution.
Real _face_mass_flux
Container for the mass flux on the face which will be reused in the advection term&#39;s matrix and right...
Moose::FV::InterpMethod _advected_interp_method
The interpolation method to use for the advected quantity.

◆ validParams()

InputParameters LinearWCNSFVMomentumFlux::validParams ( )
static

Definition at line 21 of file LinearWCNSFVMomentumFlux.C.

22 {
24  params.addClassDescription("Represents the matrix and right hand side contributions of the "
25  "stress and advection terms of the momentum equation.");
26  params.addRequiredParam<SolverVariableName>("u", "The velocity in the x direction.");
27  params.addParam<SolverVariableName>("v", "The velocity in the y direction.");
28  params.addParam<SolverVariableName>("w", "The velocity in the z direction.");
29  params.addRequiredParam<UserObjectName>(
30  "rhie_chow_user_object",
31  "The rhie-chow user-object which is used to determine the face velocity.");
32  params.addRequiredParam<MooseFunctorName>(NS::mu, "The diffusion coefficient.");
33  MooseEnum momentum_component("x=0 y=1 z=2");
35  "momentum_component",
36  momentum_component,
37  "The component of the momentum equation that this kernel applies to.");
38  params.addParam<bool>(
39  "use_nonorthogonal_correction",
40  true,
41  "If the nonorthogonal correction should be used when computing the normal gradient.");
42  params.addParam<bool>(
43  "use_deviatoric_terms", false, "If deviatoric terms in the stress terms need to be used.");
44 
46  return params;
47 }
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
void addRequiredParam(const std::string &name, const std::string &doc_string)
static InputParameters validParams()
InputParameters advectedInterpolationParameter()
static const std::string mu
Definition: NS.h:123
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _advected_interp_coeffs

std::pair<Real, Real> LinearWCNSFVMomentumFlux::_advected_interp_coeffs
protected

Container for the current advected interpolation coefficients on the face to make sure we don't compute it multiple times for different terms.

Definition at line 106 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeInternalAdvectionElemMatrixContribution(), computeInternalAdvectionNeighborMatrixContribution(), and setupFaceData().

◆ _advected_interp_method

Moose::FV::InterpMethod LinearWCNSFVMomentumFlux::_advected_interp_method
protected

The interpolation method to use for the advected quantity.

Definition at line 123 of file LinearWCNSFVMomentumFlux.h.

Referenced by LinearWCNSFVMomentumFlux(), and setupFaceData().

◆ _boundary_normal_factor

Real LinearWCNSFVMomentumFlux::_boundary_normal_factor
protected

Multiplier that ensures the normal of the boundary always points outwards, even in cases when the boundary is within the mesh.

Definition at line 114 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeStressBoundaryRHSContribution(), and setupFaceData().

◆ _dim

const unsigned int LinearWCNSFVMomentumFlux::_dim
protected

◆ _face_mass_flux

Real LinearWCNSFVMomentumFlux::_face_mass_flux
protected

Container for the mass flux on the face which will be reused in the advection term's matrix and right hand side contribution.

Definition at line 110 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeAdvectionBoundaryMatrixContribution(), computeAdvectionBoundaryRHSContribution(), computeInternalAdvectionElemMatrixContribution(), computeInternalAdvectionNeighborMatrixContribution(), and setupFaceData().

◆ _index

const unsigned int LinearWCNSFVMomentumFlux::_index
protected

Index x|y|z, this is mainly to handle the deviatoric parts correctly in in the stress term.

Definition at line 127 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeInternalStressRHSContribution(), and computeStressBoundaryRHSContribution().

◆ _mass_flux_provider

const RhieChowMassFlux& LinearWCNSFVMomentumFlux::_mass_flux_provider
protected

The Rhie-Chow user object that provides us with the face velocity.

Definition at line 93 of file LinearWCNSFVMomentumFlux.h.

Referenced by setupFaceData().

◆ _mu

const Moose::Functor<Real>& LinearWCNSFVMomentumFlux::_mu
protected

◆ _stress_matrix_contribution

Real LinearWCNSFVMomentumFlux::_stress_matrix_contribution
protected

The cached matrix contribution.

Definition at line 117 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeInternalStressMatrixContribution().

◆ _stress_rhs_contribution

Real LinearWCNSFVMomentumFlux::_stress_rhs_contribution
protected

The cached right hand side contribution.

Definition at line 120 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeInternalStressRHSContribution(), and setupFaceData().

◆ _u_var

const MooseLinearVariableFVReal* const LinearWCNSFVMomentumFlux::_u_var
protected

◆ _use_deviatoric_terms

const bool LinearWCNSFVMomentumFlux::_use_deviatoric_terms
protected

Switch to enable/disable deviatoric parts in the stress term.

Definition at line 102 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeInternalStressRHSContribution(), computeStressBoundaryRHSContribution(), and LinearWCNSFVMomentumFlux().

◆ _use_nonorthogonal_correction

const bool LinearWCNSFVMomentumFlux::_use_nonorthogonal_correction
protected

Switch to enable/disable nonorthogonal correction in the stress term.

Definition at line 99 of file LinearWCNSFVMomentumFlux.h.

Referenced by computeInternalStressMatrixContribution(), computeInternalStressRHSContribution(), computeStressBoundaryRHSContribution(), and LinearWCNSFVMomentumFlux().

◆ _v_var

const MooseLinearVariableFVReal* const LinearWCNSFVMomentumFlux::_v_var
protected

◆ _w_var

const MooseLinearVariableFVReal* const LinearWCNSFVMomentumFlux::_w_var
protected

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