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

#include <TaggingInterface.h>

Inheritance diagram for TaggingInterface:
[legend]

Classes

class  MatrixTagsKey
 Class that is used as a parameter to some of our matrix tag APIs that allows only blessed framework classes to call them. More...
 
class  VectorTagsKey
 Class that is used as a parameter to some of our vector tag APIs that allows only blessed framework classes to call them. More...
 

Public Types

enum  ResidualTagType { ResidualTagType::NonReference, ResidualTagType::Reference }
 Enumerate whether a (residual) vector tag is to be of a non-reference or reference tag type. More...
 

Public Member Functions

 TaggingInterface (const MooseObject *moose_object)
 
virtual ~TaggingInterface ()
 
void useVectorTag (const TagName &tag_name, VectorTagsKey)
 
void useMatrixTag (const TagName &tag_name, MatrixTagsKey)
 
void useVectorTag (TagID tag_id, VectorTagsKey)
 
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
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Member Functions

void prepareVectorTag (Assembly &assembly, unsigned int ivar)
 Prepare data for computing element residual according to active tags. More...
 
void prepareVectorTag (Assembly &assembly, unsigned int ivar, ResidualTagType tag_type)
 Prepare vector tags in a reference residual problem context. More...
 
void prepareVectorTagNeighbor (Assembly &assembly, unsigned int ivar)
 Prepare data for computing element residual the according to active tags for DG and interface kernels. More...
 
void prepareVectorTagLower (Assembly &assembly, unsigned int ivar)
 Prepare data for computing the residual according to active tags for mortar constraints. More...
 
void prepareMatrixTag (Assembly &assembly, unsigned int ivar, unsigned int jvar)
 Prepare data for computing element jacobian according to the active tags. More...
 
void prepareMatrixTag (Assembly &assembly, unsigned int ivar, unsigned int jvar, DenseMatrix< Number > &k) const
 
void prepareMatrixTagNonlocal (Assembly &assembly, unsigned int ivar, unsigned int jvar)
 Prepare data for computing nonlocal element jacobian according to the active tags. More...
 
void prepareMatrixTagNeighbor (Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type)
 Prepare data for computing element jacobian according to the active tags for DG and interface kernels. More...
 
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)
 Prepare data for computing the jacobian according to the active tags for mortar. More...
 
void accumulateTaggedLocalResidual ()
 Local residual blocks will be appended by adding the current local kernel residual. More...
 
void assignTaggedLocalResidual ()
 Local residual blocks will assigned as the current local kernel residual. More...
 
void accumulateTaggedLocalMatrix ()
 Local Jacobian blocks will be appended by adding the current local kernel Jacobian. More...
 
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 ()
 Nonlocal Jacobian blocks will be appended by adding the current nonlocal kernel Jacobian. More...
 
void assignTaggedLocalMatrix ()
 Local Jacobian blocks will assigned as the current local kernel Jacobian. More...
 
template<typename Residuals , typename Indices >
void addResiduals (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided incoming residuals corresponding to the provided dof indices. More...
 
template<typename T , typename Indices >
void addResiduals (Assembly &assembly, const DenseVector< T > &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided incoming residuals corresponding to the provided dof indices. More...
 
template<typename Residuals , typename Indices >
void addResidualsAndJacobian (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided dof indices. More...
 
template<typename Residuals , typename Indices >
void addJacobian (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided residual derivatives into the Jacobian for the provided dof indices. More...
 
void addResiduals (Assembly &assembly, const ADResidualsPacket &packet)
 Add the provided incoming residuals corresponding to the provided dof indices. More...
 
void addResidualsAndJacobian (Assembly &assembly, const ADResidualsPacket &packet)
 Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided dof indices. More...
 
void addJacobian (Assembly &assembly, const ADResidualsPacket &packet)
 Add the provided residual derivatives into the Jacobian for the provided dof indices. More...
 
template<typename Residuals , typename Indices >
void addResidualsWithoutConstraints (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided incoming residuals corresponding to the provided dof indices. More...
 
template<typename Residuals , typename Indices >
void addResidualsAndJacobianWithoutConstraints (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided dof indices. More...
 
template<typename Residuals , typename Indices >
void addJacobianWithoutConstraints (Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
 Add the provided residual derivatives into the Jacobian for the provided dof indices. More...
 
void addJacobianElement (Assembly &assembly, Real value, dof_id_type row_index, dof_id_type column_index, Real scaling_factor)
 Add into a single Jacobian element. More...
 
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)
 Add a local Jacobian matrix. More...
 
template<typename T >
void setResidual (SystemBase &sys, const T &residual, MooseVariableFE< T > &var)
 Set residual using the variables' insertion API. More...
 
void setResidual (SystemBase &sys, Real residual, dof_id_type dof_index)
 Set residual at a specified degree of freedom index. More...
 
template<typename SetResidualFunctor >
void setResidual (SystemBase &sys, SetResidualFunctor set_residual_functor)
 Set residuals using the provided functor. More...
 

Protected Attributes

SubProblem_subproblem
 SubProblem that contains tag info. More...
 
DenseVector< Number_local_re
 Holds local residual entries as they are accumulated by this Kernel. More...
 
DenseMatrix< Number_local_ke
 Holds local Jacobian entries as they are accumulated by this Kernel. More...
 
DenseMatrix< Number_nonlocal_ke
 Holds nonlocal Jacobian entries as they are accumulated by this Kernel. More...
 

Private Member Functions

void prepareVectorTagInternal (Assembly &assembly, unsigned int ivar, const std::set< TagID > &vector_tags, const std::set< TagID > &absolute_value_vector_tags)
 Prepare data for computing element residual according to the specified tags Residual blocks for different tags will be extracted from Assembly. More...
 

Private Attributes

std::set< TagID_vector_tags
 The vector tag ids this Kernel will contribute to. More...
 
std::set< TagID_abs_vector_tags
 The absolute value residual tag ids. More...
 
std::set< TagID_matrix_tags
 The matrices this Kernel will contribute to. More...
 
std::set< TagID_non_ref_vector_tags
 A set to hold vector tags excluding the reference residual tag. More...
 
std::set< TagID_non_ref_abs_vector_tags
 A set to hold absolute value vector tags excluding the reference residual tag. More...
 
std::set< TagID_ref_vector_tags
 A set of either size 1 or 0. More...
 
std::set< TagID_ref_abs_vector_tags
 A set of either size 1 or 0. More...
 
const MooseObject_moose_object
 Moose objct this tag works on. More...
 
const InputParameters_tag_params
 Parameters from moose object. More...
 
std::vector< DenseVector< Number > * > _re_blocks
 Residual blocks Vectors For each Tag. More...
 
std::vector< DenseVector< Number > * > _absre_blocks
 Residual blocks for absolute value residual tags. More...
 
std::vector< DenseMatrix< Number > * > _ke_blocks
 Kernel blocks Vectors For each Tag. More...
 
std::vector< Real_absolute_residuals
 A container to hold absolute values of residuals passed into addResiduals. More...
 

Friends

class NonlinearSystemBase
 

Detailed Description

Definition at line 48 of file TaggingInterface.h.

Member Enumeration Documentation

◆ ResidualTagType

Enumerate whether a (residual) vector tag is to be of a non-reference or reference tag type.

Enumerator
NonReference 
Reference 

Definition at line 91 of file TaggingInterface.h.

92  {
93  NonReference,
94  Reference
95  };

Constructor & Destructor Documentation

◆ TaggingInterface()

TaggingInterface::TaggingInterface ( const MooseObject moose_object)

Definition at line 56 of file TaggingInterface.C.

57  : _subproblem(*moose_object->parameters().getCheckedPointerParam<SubProblem *>("_subproblem")),
58  _moose_object(*moose_object),
60 {
61  auto & vector_tag_names = _tag_params.get<MultiMooseEnum>("vector_tags");
62 
63  if (!vector_tag_names.isValid())
64  {
65  if (!_tag_params.get<bool>("matrix_only"))
66  mooseError("MUST provide at least one vector_tag for Kernel: ", _moose_object.name());
67  }
68  else
69  {
70  for (auto & vector_tag_name : vector_tag_names)
71  {
72  const TagID vector_tag_id = _subproblem.getVectorTagID(vector_tag_name.name());
74  mooseError("Vector tag '",
75  vector_tag_name.name(),
76  "' for Kernel '",
78  "' is not a residual vector tag");
79  _vector_tags.insert(vector_tag_id);
80  }
81  }
82 
83  // Add extra vector tags. These tags should be created in the System already, otherwise
84  // we can not add the extra tags
85  auto & extra_vector_tags = _tag_params.get<std::vector<TagName>>("extra_vector_tags");
86 
87  for (auto & vector_tag_name : extra_vector_tags)
88  {
89  const TagID vector_tag_id = _subproblem.getVectorTagID(vector_tag_name);
91  mooseError("Extra vector tag '",
92  vector_tag_name,
93  "' for Kernel '",
95  "' is not a residual vector tag");
96  _vector_tags.insert(vector_tag_id);
97  }
98 
99  // Add absolue value vector tags. These tags should be created in the System already, otherwise
100  // we can not add the extra tags
101  auto & abs_vector_tags = _tag_params.get<std::vector<TagName>>("absolute_value_vector_tags");
102 
103  for (auto & vector_tag_name : abs_vector_tags)
104  {
105  const TagID vector_tag_id = _subproblem.getVectorTagID(vector_tag_name);
107  mooseError("Absolute value vector tag '",
108  vector_tag_name,
109  "' for Kernel '",
111  "' is not a residual vector tag");
112  _abs_vector_tags.insert(vector_tag_id);
113  }
114 
115  auto & matrix_tag_names = _tag_params.get<MultiMooseEnum>("matrix_tags");
116 
117  if (matrix_tag_names.isValid())
118  for (auto & matrix_tag_name : matrix_tag_names)
119  _matrix_tags.insert(_subproblem.getMatrixTagID(matrix_tag_name.name()));
120 
121  auto & extra_matrix_tags = _tag_params.get<std::vector<TagName>>("extra_matrix_tags");
122 
123  for (auto & matrix_tag_name : extra_matrix_tags)
124  _matrix_tags.insert(_subproblem.getMatrixTagID(matrix_tag_name));
125 
126  _re_blocks.resize(_vector_tags.size());
127  _absre_blocks.resize(_abs_vector_tags.size());
128  _ke_blocks.resize(_matrix_tags.size());
129 
130  const auto * const fe_problem =
131  moose_object->parameters().getCheckedPointerParam<FEProblemBase *>("_fe_problem_base");
132 
133  for (const auto & conv : fe_problem->getConvergenceObjects())
134  {
135  const auto * const ref_conv = dynamic_cast<const ReferenceResidualConvergence *>(conv.get());
136  if (ref_conv)
137  {
138  const auto reference_tag = ref_conv->referenceVectorTagID({});
139  auto create_tags_split =
140  [reference_tag](const auto & tags, auto & non_ref_tags, auto & ref_tags)
141  {
142  for (const auto tag : tags)
143  if (tag == reference_tag)
144  ref_tags.insert(tag);
145  else
146  non_ref_tags.insert(tag);
147  };
150  }
151  else
152  {
155  }
156  }
157 }
virtual TagID getVectorTagID(const TagName &tag_name) const
Get a TagID from a TagName.
Definition: SubProblem.C:203
std::set< TagID > _ref_abs_vector_tags
A set of either size 1 or 0.
SubProblem & _subproblem
SubProblem that contains tag info.
const InputParameters & _tag_params
Parameters from moose object.
unsigned int TagID
Definition: MooseTypes.h:210
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:302
std::vector< std::pair< R1, R2 > > get(const std::string &param1, const std::string &param2) const
Combine two vector parameters into a single vector of pairs.
T getCheckedPointerParam(const std::string &name, const std::string &error_string="") const
Verifies that the requested parameter exists and is not NULL and returns it to the caller...
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
std::vector< DenseVector< Number > * > _re_blocks
Residual blocks Vectors For each Tag.
std::set< TagID > _ref_vector_tags
A set of either size 1 or 0.
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:57
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
std::set< TagID > _non_ref_abs_vector_tags
A set to hold absolute value vector tags excluding the reference residual tag.
virtual TagID getMatrixTagID(const TagName &tag_name) const
Get a TagID from a TagName.
Definition: SubProblem.C:342
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
virtual Moose::VectorTagType vectorTagType(const TagID tag_id) const
Definition: SubProblem.C:231
TagID referenceVectorTagID(ReferenceVectorTagIDKey) const
Returns the tag ID associated with the reference vector tag ID key.
std::set< TagID > _non_ref_vector_tags
A set to hold vector tags excluding the reference residual tag.
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
Uses a reference residual to define relative convergence criteria.
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:78
const MooseObject & _moose_object
Moose objct this tag works on.
const InputParameters & parameters() const
Get the parameters of the object.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type...
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.

◆ ~TaggingInterface()

TaggingInterface::~TaggingInterface ( )
virtual

Definition at line 443 of file TaggingInterface.C.

443 {}

Member Function Documentation

◆ accumulateTaggedLocalMatrix() [1/3]

void TaggingInterface::accumulateTaggedLocalMatrix ( )
protected

Local Jacobian blocks will be appended by adding the current local kernel Jacobian.

It should be called after the local element matrix has been computed.

Definition at line 387 of file TaggingInterface.C.

Referenced by DGKernel::computeElemNeighJacobian(), ElemElemConstraint::computeElemNeighJacobian(), ArrayDGKernel::computeElemNeighJacobian(), MassLumpedTimeDerivative::computeJacobian(), TimeDerivative::computeJacobian(), VectorTimeDerivative::computeJacobian(), ScalarKernel::computeJacobian(), Kernel::computeJacobian(), ODEKernel::computeJacobian(), VectorKernel::computeJacobian(), ArrayKernel::computeJacobian(), IntegratedBC::computeJacobian(), VectorIntegratedBC::computeJacobian(), EigenKernel::computeJacobian(), ArrayIntegratedBC::computeJacobian(), NodeElemConstraint::computeJacobian(), NonlocalIntegratedBC::computeJacobian(), KernelGrad::computeJacobian(), KernelValue::computeJacobian(), NonlocalKernel::computeJacobian(), MortarConstraint::computeJacobian(), NodeFaceConstraint::computeJacobian(), LowerDIntegratedBC::computeLowerDJacobian(), ArrayLowerDIntegratedBC::computeLowerDJacobian(), DGLowerDKernel::computeLowerDJacobian(), ArrayDGLowerDKernel::computeLowerDJacobian(), LowerDIntegratedBC::computeLowerDOffDiagJacobian(), ArrayLowerDIntegratedBC::computeLowerDOffDiagJacobian(), DGKernel::computeOffDiagElemNeighJacobian(), ArrayDGKernel::computeOffDiagElemNeighJacobian(), Kernel::computeOffDiagJacobian(), VectorKernel::computeOffDiagJacobian(), EigenKernel::computeOffDiagJacobian(), ArrayKernel::computeOffDiagJacobian(), VectorIntegratedBC::computeOffDiagJacobian(), IntegratedBC::computeOffDiagJacobian(), ArrayIntegratedBC::computeOffDiagJacobian(), NodeElemConstraint::computeOffDiagJacobian(), NonlocalIntegratedBC::computeOffDiagJacobian(), NonlocalKernel::computeOffDiagJacobian(), KernelValue::computeOffDiagJacobian(), KernelGrad::computeOffDiagJacobian(), NodeFaceConstraint::computeOffDiagJacobian(), ODEKernel::computeOffDiagJacobianScalar(), VectorKernel::computeOffDiagJacobianScalar(), ArrayKernel::computeOffDiagJacobianScalar(), IntegratedBC::computeOffDiagJacobianScalar(), VectorIntegratedBC::computeOffDiagJacobianScalar(), Kernel::computeOffDiagJacobianScalar(), ArrayIntegratedBC::computeOffDiagJacobianScalar(), ScalarLagrangeMultiplier::computeOffDiagJacobianScalar(), DGLowerDKernel::computeOffDiagLowerDJacobian(), and ArrayDGLowerDKernel::computeOffDiagLowerDJacobian().

388 {
389  for (auto & ke : _ke_blocks)
390  *ke += _local_ke;
391 }
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.

◆ accumulateTaggedLocalMatrix() [2/3]

void TaggingInterface::accumulateTaggedLocalMatrix ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar,
const DenseMatrix< Number > &  k 
)
protected

Definition at line 394 of file TaggingInterface.C.

398 {
399  _ke_blocks.resize(_matrix_tags.size());
400  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
401  auto mat_vector = _matrix_tags.begin();
402  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
403  _ke_blocks[i] = &assembly.jacobianBlock(ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
404  mooseAssert(_ke_blocks[0]->m() == k.m() && _ke_blocks[0]->n() == k.n(),
405  "Passed-in k must match the blocks we are about to sum into");
406  for (auto & ke : _ke_blocks)
407  *ke += k;
408 }
unsigned int m() const
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block for a pair of variables and a tag.
Definition: Assembly.h:1113
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
unsigned int n() const
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ accumulateTaggedLocalMatrix() [3/3]

void TaggingInterface::accumulateTaggedLocalMatrix ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar,
Moose::DGJacobianType  type,
const DenseMatrix< Number > &  k 
)
protected

Definition at line 411 of file TaggingInterface.C.

416 {
417  _ke_blocks.resize(_matrix_tags.size());
418  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
419  auto mat_vector = _matrix_tags.begin();
420  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
421  _ke_blocks[i] =
422  &assembly.jacobianBlockNeighbor(type, ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
423  mooseAssert(_ke_blocks[0]->m() == k.m() && _ke_blocks[0]->n() == k.n(),
424  "Passed-in k must match the blocks we are about to sum into");
425  for (auto & ke : _ke_blocks)
426  *ke += k;
427 }
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block of a DG Jacobian type for a pair of variables and a tag. ...
Definition: Assembly.C:3117
unsigned int m() const
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
unsigned int n() const
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ accumulateTaggedLocalResidual()

void TaggingInterface::accumulateTaggedLocalResidual ( )
protected

Local residual blocks will be appended by adding the current local kernel residual.

It should be called after the local element vector has been computed.

Definition at line 367 of file TaggingInterface.C.

Referenced by FVInterfaceKernel::addResidual(), ADDGKernel::computeElemNeighResidual(), DGKernel::computeElemNeighResidual(), ElemElemConstraint::computeElemNeighResidual(), ArrayDGKernel::computeElemNeighResidual(), DGLowerDKernel::computeLowerDResidual(), ArrayDGLowerDKernel::computeLowerDResidual(), ScalarKernel::computeResidual(), Kernel::computeResidual(), VectorKernel::computeResidual(), ArrayKernel::computeResidual(), LowerDIntegratedBC::computeResidual(), ADScalarKernel::computeResidual(), TimeKernel::computeResidual(), ODEKernel::computeResidual(), ODETimeKernel::computeResidual(), VectorTimeKernel::computeResidual(), VectorIntegratedBC::computeResidual(), IntegratedBC::computeResidual(), ArrayLowerDIntegratedBC::computeResidual(), ArrayIntegratedBC::computeResidual(), NodeElemConstraint::computeResidual(), EigenKernel::computeResidual(), ADMortarConstraint::computeResidual(), FVBoundaryScalarLagrangeMultiplierConstraint::computeResidual(), FVScalarLagrangeMultiplierConstraint::computeResidual(), FVFluxBC::computeResidual(), MortarConstraint::computeResidual(), KernelValue::computeResidual(), KernelGrad::computeResidual(), FVElementalKernel::computeResidual(), FVFluxKernel::computeResidual(), and NodeFaceConstraint::computeResidual().

368 {
369  for (auto & re : _re_blocks)
370  *re += _local_re;
371  for (auto & absre : _absre_blocks)
372  for (const auto i : index_range(_local_re))
373  (*absre)(i) += std::abs(_local_re(i));
374 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
std::vector< DenseVector< Number > * > _re_blocks
Residual blocks Vectors For each Tag.
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
auto index_range(const T &sizable)

◆ accumulateTaggedNonlocalMatrix()

void TaggingInterface::accumulateTaggedNonlocalMatrix ( )
protected

Nonlocal Jacobian blocks will be appended by adding the current nonlocal kernel Jacobian.

It should be called after the nonlocal element matrix has been computed.

Definition at line 430 of file TaggingInterface.C.

Referenced by NonlocalIntegratedBC::computeNonlocalJacobian(), NonlocalKernel::computeNonlocalJacobian(), NonlocalKernel::computeNonlocalOffDiagJacobian(), and NonlocalIntegratedBC::computeNonlocalOffDiagJacobian().

431 {
432  for (auto & ke : _ke_blocks)
433  *ke += _nonlocal_ke;
434 }
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
DenseMatrix< Number > _nonlocal_ke
Holds nonlocal Jacobian entries as they are accumulated by this Kernel.

◆ addJacobian() [1/3]

template<typename Residuals , typename Indices >
void TaggingInterface::addJacobian ( Assembly assembly,
const Residuals &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided residual derivatives into the Jacobian for the provided dof indices.

Definition at line 515 of file TaggingInterface.h.

Referenced by addJacobian(), addResidualsAndJacobian(), ADScalarKernel::computeADJacobian(), ADDGKernel::computeElemNeighJacobian(), DiffusionLHDGDirichletBC::computeJacobian(), DiffusionLHDGPrescribedGradientBC::computeJacobian(), IPHDGBC::computeJacobian(), DiffusionLHDGKernel::computeJacobian(), IPHDGKernel::computeJacobian(), ADNodeElemConstraint::computeJacobian(), FVElementalKernel::computeJacobian(), ADNodalKernel::computeJacobian(), NodalConstraint::computeJacobian(), DiffusionLHDGKernel::computeJacobianOnSide(), IPHDGKernel::computeJacobianOnSide(), ADDGKernel::computeOffDiagElemNeighJacobian(), FVScalarLagrangeMultiplierConstraint::computeOffDiagJacobian(), MortarScalarBase::computeOffDiagJacobianScalar(), KernelScalarBase::computeOffDiagJacobianScalarLocal(), MortarScalarBase::computeScalarJacobian(), KernelScalarBase::computeScalarJacobian(), MortarScalarBase::computeScalarOffDiagJacobian(), KernelScalarBase::computeScalarOffDiagJacobian(), KernelScalarBase::computeScalarOffDiagJacobianScalar(), and MortarScalarBase::computeScalarOffDiagJacobianScalar().

519 {
520  assembly.cacheJacobian(
521  residuals, dof_indices, scaling_factor, Assembly::LocalDataKey{}, _matrix_tags);
522 }
void cacheJacobian(GlobalDataKey)
Takes the values that are currently in _sub_Kee and appends them to the cached values.
Definition: Assembly.C:4065
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ addJacobian() [2/3]

void TaggingInterface::addJacobian ( Assembly assembly,
const ADResidualsPacket packet 
)
inlineprotected

Add the provided residual derivatives into the Jacobian for the provided dof indices.

Definition at line 608 of file TaggingInterface.h.

609 {
610  addJacobian(assembly, packet.residuals, packet.dof_indices, packet.scaling_factor);
611 }
const std::vector< dof_id_type > & dof_indices
const Real scaling_factor
const DenseVector< ADReal > & residuals
void addJacobian(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided residual derivatives into the Jacobian for the provided dof indices.

◆ addJacobian() [3/3]

void TaggingInterface::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 
)
inlineprotected

Add a local Jacobian matrix.

Definition at line 558 of file TaggingInterface.h.

563 {
564  for (const auto matrix_tag : _matrix_tags)
565  assembly.cacheJacobianBlock(
566  local_k, row_indices, column_indices, scaling_factor, Assembly::LocalDataKey{}, matrix_tag);
567 }
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void cacheJacobianBlock(DenseMatrix< Number > &jac_block, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices, Real scaling_factor, LocalDataKey, TagID tag)
Cache a local Jacobian block with the provided rows (idof_indices) and columns (jdof_indices) for eve...
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ addJacobianElement()

void TaggingInterface::addJacobianElement ( Assembly assembly,
Real  value,
dof_id_type  row_index,
dof_id_type  column_index,
Real  scaling_factor 
)
inlineprotected

Add into a single Jacobian element.

Definition at line 547 of file TaggingInterface.h.

Referenced by ArrayNodalBC::computeJacobian(), VectorNodalBC::computeJacobian(), NodalBC::computeJacobian(), NodalKernel::computeJacobian(), NodalConstraint::computeJacobian(), ArrayNodalBC::computeOffDiagJacobian(), VectorNodalBC::computeOffDiagJacobian(), NodalBC::computeOffDiagJacobian(), NodalKernel::computeOffDiagJacobian(), and MortarConstraintBase::zeroInactiveLMDofs().

552 {
553  assembly.cacheJacobian(
554  row_index, column_index, value * scaling_factor, Assembly::LocalDataKey{}, _matrix_tags);
555 }
void cacheJacobian(GlobalDataKey)
Takes the values that are currently in _sub_Kee and appends them to the cached values.
Definition: Assembly.C:4065
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ addJacobianWithoutConstraints()

template<typename Residuals , typename Indices >
void TaggingInterface::addJacobianWithoutConstraints ( Assembly assembly,
const Residuals &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided residual derivatives into the Jacobian for the provided dof indices.

This API should only be used if the caller knows that no libMesh-level constraints (hanging nodes or periodic boundary conditions) apply to the provided dof indices

Definition at line 537 of file TaggingInterface.h.

Referenced by addResidualsAndJacobianWithoutConstraints().

541 {
543  residuals, dof_indices, scaling_factor, Assembly::LocalDataKey{}, _matrix_tags);
544 }
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void cacheJacobianWithoutConstraints(const Residuals &residuals, const Indices &row_indices, Real scaling_factor, LocalDataKey, const std::set< TagID > &matrix_tags)
Process the derivatives() data of a vector of ADReals.
Definition: Assembly.h:3151
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ addResiduals() [1/3]

template<typename Residuals , typename Indices >
void TaggingInterface::addResiduals ( Assembly assembly,
const Residuals &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided incoming residuals corresponding to the provided dof indices.

Definition at line 448 of file TaggingInterface.h.

Referenced by addResiduals(), addResidualsAndJacobian(), FVScalarLagrangeMultiplierInterface::computeResidual(), DiffusionLHDGKernel::computeResidual(), IPHDGKernel::computeResidual(), DiffusionLHDGDirichletBC::computeResidual(), DiffusionLHDGPrescribedGradientBC::computeResidual(), IPHDGBC::computeResidual(), ADNodeElemConstraint::computeResidual(), TimeNodalKernel::computeResidual(), FVBoundaryScalarLagrangeMultiplierConstraint::computeResidual(), FVScalarLagrangeMultiplierConstraint::computeResidual(), NodalKernel::computeResidual(), ADNodalKernel::computeResidual(), ADKernelScalarBase::computeResidual(), NodalConstraint::computeResidual(), ADMortarScalarBase::computeResidual(), MortarScalarBase::computeResidual(), DiffusionLHDGKernel::computeResidualOnSide(), IPHDGKernel::computeResidualOnSide(), KernelScalarBase::computeScalarResidual(), and MortarConstraintBase::zeroInactiveLMDofs().

452 {
453  assembly.cacheResiduals(
454  residuals, dof_indices, scaling_factor, Assembly::LocalDataKey{}, _vector_tags);
455  if (!_abs_vector_tags.empty())
456  {
457  _absolute_residuals.resize(residuals.size());
458  for (const auto i : index_range(residuals))
460 
462  dof_indices,
463  scaling_factor,
466  }
467 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
auto raw_value(const Eigen::Map< T > &in)
Definition: EigenADReal.h:73
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
void cacheResiduals(const Residuals &residuals, const Indices &row_indices, Real scaling_factor, LocalDataKey, const std::set< TagID > &vector_tags)
Process the supplied residual values.
Definition: Assembly.h:3031
auto index_range(const T &sizable)
std::vector< Real > _absolute_residuals
A container to hold absolute values of residuals passed into addResiduals.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ addResiduals() [2/3]

template<typename T , typename Indices >
void TaggingInterface::addResiduals ( Assembly assembly,
const DenseVector< T > &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided incoming residuals corresponding to the provided dof indices.

Definition at line 471 of file TaggingInterface.h.

475 {
476  addResiduals(assembly, residuals.get_values(), dof_indices, scaling_factor);
477 }
void addResiduals(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.

◆ addResiduals() [3/3]

void TaggingInterface::addResiduals ( Assembly assembly,
const ADResidualsPacket packet 
)
inlineprotected

Add the provided incoming residuals corresponding to the provided dof indices.

Definition at line 596 of file TaggingInterface.h.

597 {
598  addResiduals(assembly, packet.residuals, packet.dof_indices, packet.scaling_factor);
599 }
const std::vector< dof_id_type > & dof_indices
void addResiduals(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.
const Real scaling_factor
const DenseVector< ADReal > & residuals

◆ addResidualsAndJacobian() [1/2]

template<typename Residuals , typename Indices >
void TaggingInterface::addResidualsAndJacobian ( Assembly assembly,
const Residuals &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided dof indices.

Definition at line 504 of file TaggingInterface.h.

Referenced by addResidualsAndJacobian(), FVScalarLagrangeMultiplierInterface::computeJacobian(), FVBoundaryScalarLagrangeMultiplierConstraint::computeJacobian(), FVFluxBC::computeJacobian(), ADKernelScalarBase::computeJacobian(), FVFluxKernel::computeJacobian(), FVInterfaceKernel::computeJacobian(), IPHDGBC::computeResidualAndJacobian(), IPHDGKernel::computeResidualAndJacobian(), FVScalarLagrangeMultiplierConstraint::computeResidualAndJacobian(), FVElementalKernel::computeResidualAndJacobian(), ADKernelScalarBase::computeResidualAndJacobian(), and IPHDGKernel::computeResidualAndJacobianOnSide().

508 {
509  addResiduals(assembly, residuals, dof_indices, scaling_factor);
510  addJacobian(assembly, residuals, dof_indices, scaling_factor);
511 }
void addResiduals(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.
void addJacobian(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided residual derivatives into the Jacobian for the provided dof indices.

◆ addResidualsAndJacobian() [2/2]

void TaggingInterface::addResidualsAndJacobian ( Assembly assembly,
const ADResidualsPacket packet 
)
inlineprotected

Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided dof indices.

Definition at line 602 of file TaggingInterface.h.

603 {
604  addResidualsAndJacobian(assembly, packet.residuals, packet.dof_indices, packet.scaling_factor);
605 }
void addResidualsAndJacobian(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided d...
const std::vector< dof_id_type > & dof_indices
const Real scaling_factor
const DenseVector< ADReal > & residuals

◆ addResidualsAndJacobianWithoutConstraints()

template<typename Residuals , typename Indices >
void TaggingInterface::addResidualsAndJacobianWithoutConstraints ( Assembly assembly,
const Residuals &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided incoming residuals and derivatives for the Jacobian, corresponding to the provided dof indices.

This API should only be used if the caller knows that no libMesh-level constraints (hanging nodes or periodic boundary conditions) apply to the provided dof indices

Definition at line 526 of file TaggingInterface.h.

Referenced by ADMortarConstraint::computeJacobian(), and ADMortarScalarBase::computeJacobian().

530 {
531  addResidualsWithoutConstraints(assembly, residuals, dof_indices, scaling_factor);
532  addJacobianWithoutConstraints(assembly, residuals, dof_indices, scaling_factor);
533 }
void addResidualsWithoutConstraints(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.
void addJacobianWithoutConstraints(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided residual derivatives into the Jacobian for the provided dof indices.

◆ addResidualsWithoutConstraints()

template<typename Residuals , typename Indices >
void TaggingInterface::addResidualsWithoutConstraints ( Assembly assembly,
const Residuals &  residuals,
const Indices &  dof_indices,
Real  scaling_factor 
)
protected

Add the provided incoming residuals corresponding to the provided dof indices.

This API should only be used if the caller knows that no libMesh-level constraints (hanging nodes or periodic boundary conditions) apply to the provided dof indices

Definition at line 481 of file TaggingInterface.h.

Referenced by addResidualsAndJacobianWithoutConstraints().

485 {
487  residuals, dof_indices, scaling_factor, Assembly::LocalDataKey{}, _vector_tags);
488  if (!_abs_vector_tags.empty())
489  {
490  _absolute_residuals.resize(residuals.size());
491  for (const auto i : index_range(residuals))
493 
495  dof_indices,
496  scaling_factor,
499  }
500 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
auto raw_value(const Eigen::Map< T > &in)
Definition: EigenADReal.h:73
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
auto index_range(const T &sizable)
std::vector< Real > _absolute_residuals
A container to hold absolute values of residuals passed into addResiduals.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833
void cacheResidualsWithoutConstraints(const Residuals &residuals, const Indices &row_indices, Real scaling_factor, LocalDataKey, const std::set< TagID > &vector_tags)
Process the supplied residual values.
Definition: Assembly.h:3071

◆ assignTaggedLocalMatrix()

void TaggingInterface::assignTaggedLocalMatrix ( )
protected

Local Jacobian blocks will assigned as the current local kernel Jacobian.

It should be called after the local element matrix has been computed.

Definition at line 437 of file TaggingInterface.C.

Referenced by NodalEqualValueConstraint::computeJacobian().

438 {
439  for (auto & ke : _ke_blocks)
440  *ke = _local_ke;
441 }
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.

◆ assignTaggedLocalResidual()

void TaggingInterface::assignTaggedLocalResidual ( )
protected

Local residual blocks will assigned as the current local kernel residual.

It should be called after the local element vector has been computed.

Definition at line 377 of file TaggingInterface.C.

Referenced by NodalEqualValueConstraint::computeResidual(), and NodeFaceConstraint::computeResidual().

378 {
379  for (auto & re : _re_blocks)
380  *re = _local_re;
381  for (auto & absre : _absre_blocks)
382  for (const auto i : index_range(_local_re))
383  (*absre)(i) = std::abs(_local_re(i));
384 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
std::vector< DenseVector< Number > * > _re_blocks
Residual blocks Vectors For each Tag.
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
auto index_range(const T &sizable)

◆ getMatrixTags()

const std::set<TagID>& TaggingInterface::getMatrixTags ( MatrixTagsKey  ) const
inline

Definition at line 113 of file TaggingInterface.h.

Referenced by LinearSystemContributionObject::linkTaggedVectorsAndMatrices().

113 { return _matrix_tags; }
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.

◆ getVectorTags()

const std::set<TagID>& TaggingInterface::getVectorTags ( VectorTagsKey  ) const
inline

Definition at line 111 of file TaggingInterface.h.

Referenced by LinearSystemContributionObject::linkTaggedVectorsAndMatrices().

111 { return _vector_tags; }
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.

◆ hasVectorTags()

bool TaggingInterface::hasVectorTags ( ) const
inline

Definition at line 109 of file TaggingInterface.h.

Referenced by Kernel::computeResidual().

109 { return !_vector_tags.empty(); }
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.

◆ isMatrixTagged()

bool TaggingInterface::isMatrixTagged ( )
inline

Definition at line 107 of file TaggingInterface.h.

107 { return _matrix_tags.size() > 0; }
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.

◆ isVectorTagged()

bool TaggingInterface::isVectorTagged ( )
inline

Definition at line 105 of file TaggingInterface.h.

105 { return _vector_tags.size() > 0; }
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.

◆ prepareMatrixTag() [1/2]

void TaggingInterface::prepareMatrixTag ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar 
)
protected

Prepare data for computing element jacobian according to the active tags.

Jacobian blocks for different tags will be extracted from Assembly. A local Jacobian will be zeroed. It should be called right before the local element matrix is computed.

Definition at line 283 of file TaggingInterface.C.

Referenced by DGKernel::computeElemNeighJacobian(), ArrayDGKernel::computeElemNeighJacobian(), ScalarKernel::computeJacobian(), TimeDerivative::computeJacobian(), VectorTimeDerivative::computeJacobian(), MassLumpedTimeDerivative::computeJacobian(), ODEKernel::computeJacobian(), VectorKernel::computeJacobian(), Kernel::computeJacobian(), ArrayKernel::computeJacobian(), VectorIntegratedBC::computeJacobian(), IntegratedBC::computeJacobian(), ArrayIntegratedBC::computeJacobian(), EigenKernel::computeJacobian(), NodalEqualValueConstraint::computeJacobian(), NonlocalIntegratedBC::computeJacobian(), KernelGrad::computeJacobian(), KernelValue::computeJacobian(), NonlocalKernel::computeJacobian(), DGKernel::computeOffDiagElemNeighJacobian(), ArrayDGKernel::computeOffDiagElemNeighJacobian(), VectorKernel::computeOffDiagJacobian(), Kernel::computeOffDiagJacobian(), EigenKernel::computeOffDiagJacobian(), ArrayKernel::computeOffDiagJacobian(), VectorIntegratedBC::computeOffDiagJacobian(), IntegratedBC::computeOffDiagJacobian(), ArrayIntegratedBC::computeOffDiagJacobian(), NonlocalIntegratedBC::computeOffDiagJacobian(), NonlocalKernel::computeOffDiagJacobian(), KernelValue::computeOffDiagJacobian(), KernelGrad::computeOffDiagJacobian(), ODEKernel::computeOffDiagJacobianScalar(), VectorKernel::computeOffDiagJacobianScalar(), ArrayKernel::computeOffDiagJacobianScalar(), VectorIntegratedBC::computeOffDiagJacobianScalar(), IntegratedBC::computeOffDiagJacobianScalar(), ArrayIntegratedBC::computeOffDiagJacobianScalar(), Kernel::computeOffDiagJacobianScalar(), and ScalarLagrangeMultiplier::computeOffDiagJacobianScalar().

284 {
285  _ke_blocks.resize(_matrix_tags.size());
286  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
287  auto mat_vector = _matrix_tags.begin();
288  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
289  _ke_blocks[i] = &assembly.jacobianBlock(ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
290 
291  _local_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
292 }
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block for a pair of variables and a tag.
Definition: Assembly.h:1113
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void resize(const unsigned int new_m, const unsigned int new_n)
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareMatrixTag() [2/2]

void TaggingInterface::prepareMatrixTag ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar,
DenseMatrix< Number > &  k 
) const
protected

Definition at line 295 of file TaggingInterface.C.

299 {
300  mooseAssert(!_matrix_tags.empty(), "No matrix tags exist");
301  const auto & ij_mat =
302  assembly.jacobianBlock(ivar, jvar, Assembly::LocalDataKey{}, *_matrix_tags.begin());
303  k.resize(ij_mat.m(), ij_mat.n());
304 }
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block for a pair of variables and a tag.
Definition: Assembly.h:1113
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void resize(const unsigned int new_m, const unsigned int new_n)
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareMatrixTagLower()

void TaggingInterface::prepareMatrixTagLower ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar,
Moose::ConstraintJacobianType  type 
)
protected

Prepare data for computing the jacobian according to the active tags for mortar.

Jacobian blocks for different tags will be extracted from Assembly. A local Jacobian will be zeroed. It should be called right before the local element matrix is computed.

Definition at line 351 of file TaggingInterface.C.

Referenced by MortarConstraint::computeJacobian(), LowerDIntegratedBC::computeLowerDJacobian(), ArrayLowerDIntegratedBC::computeLowerDJacobian(), DGLowerDKernel::computeLowerDJacobian(), ArrayDGLowerDKernel::computeLowerDJacobian(), LowerDIntegratedBC::computeLowerDOffDiagJacobian(), ArrayLowerDIntegratedBC::computeLowerDOffDiagJacobian(), DGLowerDKernel::computeOffDiagLowerDJacobian(), and ArrayDGLowerDKernel::computeOffDiagLowerDJacobian().

355 {
356  _ke_blocks.resize(_matrix_tags.size());
357  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
358  auto mat_vector = _matrix_tags.begin();
359  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
360  _ke_blocks[i] =
361  &assembly.jacobianBlockMortar(type, ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
362 
363  _local_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
364 }
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
DenseMatrix< Number > & jacobianBlockMortar(Moose::ConstraintJacobianType type, unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Returns the jacobian block for the given mortar Jacobian type.
Definition: Assembly.C:3158
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void resize(const unsigned int new_m, const unsigned int new_n)
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareMatrixTagNeighbor() [1/2]

void TaggingInterface::prepareMatrixTagNeighbor ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar,
Moose::DGJacobianType  type 
)
protected

Prepare data for computing element jacobian according to the active tags for DG and interface kernels.

Jacobian blocks for different tags will be extracted from Assembly. A local Jacobian will be zeroed. It should be called right before the local element matrix is computed.

Definition at line 322 of file TaggingInterface.C.

Referenced by DGKernel::computeElemNeighJacobian(), ElemElemConstraint::computeElemNeighJacobian(), ArrayDGKernel::computeElemNeighJacobian(), NodeElemConstraint::computeJacobian(), NodeFaceConstraint::computeJacobian(), DGKernel::computeOffDiagElemNeighJacobian(), ArrayDGKernel::computeOffDiagElemNeighJacobian(), NodeElemConstraint::computeOffDiagJacobian(), and NodeFaceConstraint::computeOffDiagJacobian().

326 {
327  _ke_blocks.resize(_matrix_tags.size());
328  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
329  auto mat_vector = _matrix_tags.begin();
330  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
331  _ke_blocks[i] =
332  &assembly.jacobianBlockNeighbor(type, ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
333 
334  _local_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
335 }
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block of a DG Jacobian type for a pair of variables and a tag. ...
Definition: Assembly.C:3117
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void resize(const unsigned int new_m, const unsigned int new_n)
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareMatrixTagNeighbor() [2/2]

void TaggingInterface::prepareMatrixTagNeighbor ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar,
Moose::DGJacobianType  type,
DenseMatrix< Number > &  k 
) const
protected

Definition at line 338 of file TaggingInterface.C.

343 {
344  mooseAssert(!_matrix_tags.empty(), "No matrix tags exist");
345  const auto & ij_mat = assembly.jacobianBlockNeighbor(
346  type, ivar, jvar, Assembly::LocalDataKey{}, *_matrix_tags.begin());
347  k.resize(ij_mat.m(), ij_mat.n());
348 }
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block of a DG Jacobian type for a pair of variables and a tag. ...
Definition: Assembly.C:3117
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void resize(const unsigned int new_m, const unsigned int new_n)
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareMatrixTagNonlocal()

void TaggingInterface::prepareMatrixTagNonlocal ( Assembly assembly,
unsigned int  ivar,
unsigned int  jvar 
)
protected

Prepare data for computing nonlocal element jacobian according to the active tags.

Jacobian blocks for different tags will be extracted from Assembly. A nonlocal Jacobian will be zeroed. It should be called right before the nonlocal element matrix is computed.

Definition at line 307 of file TaggingInterface.C.

Referenced by NonlocalIntegratedBC::computeNonlocalJacobian(), NonlocalKernel::computeNonlocalJacobian(), NonlocalKernel::computeNonlocalOffDiagJacobian(), and NonlocalIntegratedBC::computeNonlocalOffDiagJacobian().

310 {
311  _ke_blocks.resize(_matrix_tags.size());
312  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
313  auto mat_vector = _matrix_tags.begin();
314  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
315  _ke_blocks[i] =
316  &assembly.jacobianBlockNonlocal(ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
317 
318  _nonlocal_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
319 }
DenseMatrix< Number > & jacobianBlockNonlocal(unsigned int ivar, unsigned int jvar, LocalDataKey, TagID tag)
Get local Jacobian block from non-local contribution for a pair of variables and a tag...
Definition: Assembly.h:1124
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
void resize(const unsigned int new_m, const unsigned int new_n)
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
DenseMatrix< Number > _nonlocal_ke
Holds nonlocal Jacobian entries as they are accumulated by this Kernel.
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareVectorTag() [1/2]

void TaggingInterface::prepareVectorTag ( Assembly assembly,
unsigned int  ivar 
)
protected

Prepare data for computing element residual according to active tags.

Residual blocks for different tags will be extracted from Assembly. A local residual will be zeroed. It should be called right before the local element vector is computed.

Definition at line 196 of file TaggingInterface.C.

Referenced by FVInterfaceKernel::addResidual(), ADDGKernel::computeElemNeighResidual(), DGKernel::computeElemNeighResidual(), ElemElemConstraint::computeElemNeighResidual(), ArrayDGKernel::computeElemNeighResidual(), ScalarKernel::computeResidual(), Kernel::computeResidual(), VectorKernel::computeResidual(), ArrayKernel::computeResidual(), VectorTimeKernel::computeResidual(), ODEKernel::computeResidual(), TimeKernel::computeResidual(), ODETimeKernel::computeResidual(), ADScalarKernel::computeResidual(), IntegratedBC::computeResidual(), VectorIntegratedBC::computeResidual(), EigenKernel::computeResidual(), NodeElemConstraint::computeResidual(), ArrayIntegratedBC::computeResidual(), NodalEqualValueConstraint::computeResidual(), ADMortarConstraint::computeResidual(), FVScalarLagrangeMultiplierConstraint::computeResidual(), FVBoundaryScalarLagrangeMultiplierConstraint::computeResidual(), MortarConstraint::computeResidual(), FVFluxBC::computeResidual(), KernelValue::computeResidual(), KernelGrad::computeResidual(), FVElementalKernel::computeResidual(), FVFluxKernel::computeResidual(), and NodeFaceConstraint::computeResidual().

197 {
199 }
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
void prepareVectorTagInternal(Assembly &assembly, unsigned int ivar, const std::set< TagID > &vector_tags, const std::set< TagID > &absolute_value_vector_tags)
Prepare data for computing element residual according to the specified tags Residual blocks for diffe...

◆ prepareVectorTag() [2/2]

void TaggingInterface::prepareVectorTag ( Assembly assembly,
unsigned int  ivar,
ResidualTagType  tag_type 
)
protected

Prepare vector tags in a reference residual problem context.

Parameters
AssemblyThe assembly object that we obtain the local residual blocks from
ivarThe variable which we are retrieving the local residual blocks for
ref_problemA pointer to a reference residual problem. This can be a nullptr
tag_typeWhat type of tags to prepare

Definition at line 202 of file TaggingInterface.C.

205 {
206  if (tag_type == ResidualTagType::NonReference)
208  else
210 }
std::set< TagID > _ref_abs_vector_tags
A set of either size 1 or 0.
std::set< TagID > _ref_vector_tags
A set of either size 1 or 0.
std::set< TagID > _non_ref_abs_vector_tags
A set to hold absolute value vector tags excluding the reference residual tag.
std::set< TagID > _non_ref_vector_tags
A set to hold vector tags excluding the reference residual tag.
void prepareVectorTagInternal(Assembly &assembly, unsigned int ivar, const std::set< TagID > &vector_tags, const std::set< TagID > &absolute_value_vector_tags)
Prepare data for computing element residual according to the specified tags Residual blocks for diffe...

◆ prepareVectorTagInternal()

void TaggingInterface::prepareVectorTagInternal ( Assembly assembly,
unsigned int  ivar,
const std::set< TagID > &  vector_tags,
const std::set< TagID > &  absolute_value_vector_tags 
)
private

Prepare data for computing element residual according to the specified tags Residual blocks for different tags will be extracted from Assembly.

A local residual will be zeroed. It should be called right before the local element vector is computed.

Definition at line 213 of file TaggingInterface.C.

Referenced by prepareVectorTag().

217 {
218  auto prepare = [this, ivar, &assembly](auto & re_blocks, const auto & tags)
219  {
220  re_blocks.clear();
221  re_blocks.reserve(tags.size());
222  for (const auto tag_id : tags)
223  {
224  const auto & tag = _subproblem.getVectorTag(tag_id);
225  re_blocks.push_back(&assembly.residualBlock(ivar, Assembly::LocalDataKey{}, tag._type_id));
226  }
227  };
228 
229  prepare(_re_blocks, vector_tags);
230  prepare(_absre_blocks, absolute_value_vector_tags);
231 
232  _local_re.resize(_re_blocks.empty()
233  ? (_absre_blocks.empty() ? std::size_t(0) : _absre_blocks[0]->size())
234  : _re_blocks[0]->size());
235 }
SubProblem & _subproblem
SubProblem that contains tag info.
void resize(const unsigned int n)
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
std::vector< DenseVector< Number > * > _re_blocks
Residual blocks Vectors For each Tag.
DenseVector< Number > & residualBlock(unsigned int var_num, LocalDataKey, TagID tag_id)
Get local residual block for a variable and a tag.
Definition: Assembly.h:1086
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
virtual const VectorTag & getVectorTag(const TagID tag_id) const
Get a VectorTag from a TagID.
Definition: SubProblem.C:161
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareVectorTagLower()

void TaggingInterface::prepareVectorTagLower ( Assembly assembly,
unsigned int  ivar 
)
protected

Prepare data for computing the residual according to active tags for mortar constraints.

Residual blocks for different tags will be extracted from Assembly. A local residual will be zeroed. It should be called right before the local element vector is computed.

Definition at line 261 of file TaggingInterface.C.

Referenced by DGLowerDKernel::computeLowerDResidual(), ArrayDGLowerDKernel::computeLowerDResidual(), LowerDIntegratedBC::computeResidual(), ArrayLowerDIntegratedBC::computeResidual(), ADMortarConstraint::computeResidual(), and MortarConstraint::computeResidual().

262 {
263  _re_blocks.resize(_vector_tags.size());
264  mooseAssert(_vector_tags.size() >= 1, "we need at least one active tag");
265  auto vector_tag = _vector_tags.begin();
266  for (MooseIndex(_vector_tags) i = 0; i < _vector_tags.size(); i++, ++vector_tag)
267  {
268  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
269  _re_blocks[i] = &assembly.residualBlockLower(ivar, Assembly::LocalDataKey{}, tag._type_id);
270  }
271  _local_re.resize(_re_blocks[0]->size());
272 
273  _absre_blocks.resize(_abs_vector_tags.size());
274  vector_tag = _abs_vector_tags.begin();
275  for (MooseIndex(_abs_vector_tags) i = 0; i < _abs_vector_tags.size(); i++, ++vector_tag)
276  {
277  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
279  }
280 }
SubProblem & _subproblem
SubProblem that contains tag info.
void resize(const unsigned int n)
TagTypeID _type_id
The index for this tag into a vector that contains tags of only its type ordered by ID...
Definition: VectorTag.h:47
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
std::vector< DenseVector< Number > * > _re_blocks
Residual blocks Vectors For each Tag.
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
DenseVector< Number > & residualBlockLower(unsigned int var_num, LocalDataKey, TagID tag_id)
Get residual block for lower.
Definition: Assembly.h:1104
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
Storage for all of the information pretaining to a vector tag.
Definition: VectorTag.h:17
virtual const VectorTag & getVectorTag(const TagID tag_id) const
Get a VectorTag from a TagID.
Definition: SubProblem.C:161
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ prepareVectorTagNeighbor()

void TaggingInterface::prepareVectorTagNeighbor ( Assembly assembly,
unsigned int  ivar 
)
protected

Prepare data for computing element residual the according to active tags for DG and interface kernels.

Residual blocks for different tags will be extracted from Assembly. A local residual will be zeroed. It should be called right before the local element vector is computed.

Definition at line 238 of file TaggingInterface.C.

Referenced by FVInterfaceKernel::addResidual(), ADDGKernel::computeElemNeighResidual(), DGKernel::computeElemNeighResidual(), ElemElemConstraint::computeElemNeighResidual(), ArrayDGKernel::computeElemNeighResidual(), NodeElemConstraint::computeResidual(), ADMortarConstraint::computeResidual(), FVFluxBC::computeResidual(), MortarConstraint::computeResidual(), FVFluxKernel::computeResidual(), and NodeFaceConstraint::computeResidual().

239 {
240  _re_blocks.resize(_vector_tags.size());
241  mooseAssert(_vector_tags.size() >= 1, "we need at least one active tag");
242  auto vector_tag = _vector_tags.begin();
243  for (MooseIndex(_vector_tags) i = 0; i < _vector_tags.size(); i++, ++vector_tag)
244  {
245  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
247  }
248  _local_re.resize(_re_blocks[0]->size());
249 
250  _absre_blocks.resize(_abs_vector_tags.size());
251  vector_tag = _abs_vector_tags.begin();
252  for (MooseIndex(_abs_vector_tags) i = 0; i < _abs_vector_tags.size(); i++, ++vector_tag)
253  {
254  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
255  _absre_blocks[i] =
256  &assembly.residualBlockNeighbor(ivar, Assembly::LocalDataKey{}, tag._type_id);
257  }
258 }
SubProblem & _subproblem
SubProblem that contains tag info.
void resize(const unsigned int n)
TagTypeID _type_id
The index for this tag into a vector that contains tags of only its type ordered by ID...
Definition: VectorTag.h:47
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
std::vector< DenseVector< Number > * > _re_blocks
Residual blocks Vectors For each Tag.
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, LocalDataKey, TagID tag_id)
Get local neighbor residual block for a variable and a tag.
Definition: Assembly.h:1095
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
Storage for all of the information pretaining to a vector tag.
Definition: VectorTag.h:17
virtual const VectorTag & getVectorTag(const TagID tag_id) const
Get a VectorTag from a TagID.
Definition: SubProblem.C:161
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:833

◆ setResidual() [1/3]

template<typename T >
void TaggingInterface::setResidual ( SystemBase sys,
const T &  residual,
MooseVariableFE< T > &  var 
)
protected

Set residual using the variables' insertion API.

Definition at line 571 of file TaggingInterface.h.

Referenced by ArrayNodalBC::computeResidual(), VectorNodalBC::computeResidual(), and NodalBC::computeResidual().

572 {
573  for (const auto tag_id : _vector_tags)
574  if (sys.hasVector(tag_id))
575  var.insertNodalValue(sys.getVector(tag_id), residual);
576 }
bool hasVector(const std::string &tag_name) const
Check if the named vector exists in the system.
Definition: SystemBase.C:907
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
void insertNodalValue(libMesh::NumericVector< libMesh::Number > &residual, const OutputData &v)
Write a nodal value to the passed-in solution vector.
virtual NumericVector< Number > & getVector(const std::string &name)
Get a raw NumericVector by name.
Definition: SystemBase.C:916

◆ setResidual() [2/3]

void TaggingInterface::setResidual ( SystemBase sys,
Real  residual,
dof_id_type  dof_index 
)
inlineprotected

Set residual at a specified degree of freedom index.

Definition at line 579 of file TaggingInterface.h.

580 {
581  for (const auto tag_id : _vector_tags)
582  if (sys.hasVector(tag_id))
583  sys.getVector(tag_id).set(dof_index, residual);
584 }
bool hasVector(const std::string &tag_name) const
Check if the named vector exists in the system.
Definition: SystemBase.C:907
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
virtual void set(const numeric_index_type i, const Number value)=0
virtual NumericVector< Number > & getVector(const std::string &name)
Get a raw NumericVector by name.
Definition: SystemBase.C:916

◆ setResidual() [3/3]

template<typename SetResidualFunctor >
void TaggingInterface::setResidual ( SystemBase sys,
SetResidualFunctor  set_residual_functor 
)
protected

Set residuals using the provided functor.

Definition at line 588 of file TaggingInterface.h.

589 {
590  for (const auto tag_id : _vector_tags)
591  if (sys.hasVector(tag_id))
592  set_residual_functor(sys.getVector(tag_id));
593 }
bool hasVector(const std::string &tag_name) const
Check if the named vector exists in the system.
Definition: SystemBase.C:907
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
virtual NumericVector< Number > & getVector(const std::string &name)
Get a raw NumericVector by name.
Definition: SystemBase.C:916

◆ useMatrixTag() [1/2]

void TaggingInterface::useMatrixTag ( const TagName &  tag_name,
MatrixTagsKey   
)

Definition at line 169 of file TaggingInterface.C.

170 {
171  if (!_subproblem.matrixTagExists(tag_name))
172  mooseError("Matrix tag ", tag_name, " does not exist in system");
173 
174  _matrix_tags.insert(_subproblem.getMatrixTagID(tag_name));
175 }
SubProblem & _subproblem
SubProblem that contains tag info.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:302
virtual TagID getMatrixTagID(const TagName &tag_name) const
Get a TagID from a TagName.
Definition: SubProblem.C:342
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
virtual bool matrixTagExists(const TagName &tag_name) const
Check to see if a particular Tag exists.
Definition: SubProblem.C:328

◆ useMatrixTag() [2/2]

void TaggingInterface::useMatrixTag ( TagID  tag_id,
MatrixTagsKey   
)

Definition at line 187 of file TaggingInterface.C.

188 {
189  if (!_subproblem.matrixTagExists(tag_id))
190  mooseError("Matrix tag ", tag_id, " does not exist in system");
191 
192  _matrix_tags.insert(tag_id);
193 }
SubProblem & _subproblem
SubProblem that contains tag info.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:302
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
virtual bool matrixTagExists(const TagName &tag_name) const
Check to see if a particular Tag exists.
Definition: SubProblem.C:328

◆ useVectorTag() [1/2]

void TaggingInterface::useVectorTag ( const TagName &  tag_name,
VectorTagsKey   
)

Definition at line 160 of file TaggingInterface.C.

161 {
162  if (!_subproblem.vectorTagExists(tag_name))
163  mooseError("Vector tag ", tag_name, " does not exist in system");
164 
165  _vector_tags.insert(_subproblem.getVectorTagID(tag_name));
166 }
virtual TagID getVectorTagID(const TagName &tag_name) const
Get a TagID from a TagName.
Definition: SubProblem.C:203
SubProblem & _subproblem
SubProblem that contains tag info.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:302
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
virtual bool vectorTagExists(const TagID tag_id) const
Check to see if a particular Tag exists.
Definition: SubProblem.h:201

◆ useVectorTag() [2/2]

void TaggingInterface::useVectorTag ( TagID  tag_id,
VectorTagsKey   
)

Definition at line 178 of file TaggingInterface.C.

179 {
180  if (!_subproblem.vectorTagExists(tag_id))
181  mooseError("Vector tag ", tag_id, " does not exist in system");
182 
183  _vector_tags.insert(tag_id);
184 }
SubProblem & _subproblem
SubProblem that contains tag info.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:302
std::set< TagID > _vector_tags
The vector tag ids this Kernel will contribute to.
virtual bool vectorTagExists(const TagID tag_id) const
Check to see if a particular Tag exists.
Definition: SubProblem.h:201

◆ validParams()

InputParameters TaggingInterface::validParams ( )
static

Definition at line 20 of file TaggingInterface.C.

Referenced by ResidualObject::validParams(), LinearSystemContributionObject::validParams(), LinearFVBoundaryCondition::validParams(), FVInterfaceKernel::validParams(), and FVBoundaryCondition::validParams().

21 {
22 
24 
25  // These are the default names for tags, but users will be able to add their own
26  MultiMooseEnum vtags("nontime time", "nontime", true);
27  MultiMooseEnum mtags("nontime system", "system", true);
28 
29  params.addParam<bool>(
30  "matrix_only", false, "Whether this object is only doing assembly to matrices (no vectors)");
31 
32  params.addParam<MultiMooseEnum>(
33  "vector_tags", vtags, "The tag for the vectors this Kernel should fill");
34 
35  params.addParam<MultiMooseEnum>(
36  "matrix_tags", mtags, "The tag for the matrices this Kernel should fill");
37 
38  params.addParam<std::vector<TagName>>("extra_vector_tags",
39  "The extra tags for the vectors this Kernel should fill");
40 
41  params.addParam<std::vector<TagName>>(
42  "absolute_value_vector_tags",
43  "The tags for the vectors this residual object should fill with the "
44  "absolute value of the residual contribution");
45 
46  params.addParam<std::vector<TagName>>("extra_matrix_tags",
47  "The extra tags for the matrices this Kernel should fill");
48 
49  params.addParamNamesToGroup(
50  "vector_tags matrix_tags extra_vector_tags extra_matrix_tags absolute_value_vector_tags",
51  "Contribution to tagged field data");
52 
53  return params;
54 }
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
InputParameters emptyInputParameters()
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type...
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...

Friends And Related Function Documentation

◆ NonlinearSystemBase

friend class NonlinearSystemBase
friend

Definition at line 430 of file TaggingInterface.h.

Member Data Documentation

◆ _abs_vector_tags

std::set<TagID> TaggingInterface::_abs_vector_tags
private

◆ _absolute_residuals

std::vector<Real> TaggingInterface::_absolute_residuals
private

A container to hold absolute values of residuals passed into addResiduals.

We maintain this data member to avoid constant dynamic heap allocations

Definition at line 428 of file TaggingInterface.h.

Referenced by addResiduals(), and addResidualsWithoutConstraints().

◆ _absre_blocks

std::vector<DenseVector<Number> *> TaggingInterface::_absre_blocks
private

◆ _ke_blocks

std::vector<DenseMatrix<Number> *> TaggingInterface::_ke_blocks
private

◆ _local_ke

DenseMatrix<Number> TaggingInterface::_local_ke
protected

Holds local Jacobian entries as they are accumulated by this Kernel.

Definition at line 367 of file TaggingInterface.h.

Referenced by accumulateTaggedLocalMatrix(), assignTaggedLocalMatrix(), ADDGKernel::computeElemNeighJacobian(), DGKernel::computeElemNeighJacobian(), ElemElemConstraint::computeElemNeighJacobian(), ArrayDGKernel::computeElemNeighJacobian(), VectorTimeDerivative::computeJacobian(), ScalarKernel::computeJacobian(), MassLumpedTimeDerivative::computeJacobian(), TimeDerivative::computeJacobian(), Kernel::computeJacobian(), ODEKernel::computeJacobian(), VectorKernel::computeJacobian(), ArrayKernel::computeJacobian(), IntegratedBC::computeJacobian(), VectorIntegratedBC::computeJacobian(), EigenKernel::computeJacobian(), ArrayIntegratedBC::computeJacobian(), NodalEqualValueConstraint::computeJacobian(), NodeElemConstraint::computeJacobian(), NonlocalIntegratedBC::computeJacobian(), KernelGrad::computeJacobian(), KernelValue::computeJacobian(), NonlocalKernel::computeJacobian(), MortarConstraint::computeJacobian(), NodeFaceConstraint::computeJacobian(), LowerDIntegratedBC::computeLowerDJacobian(), ArrayLowerDIntegratedBC::computeLowerDJacobian(), DGLowerDKernel::computeLowerDJacobian(), ArrayDGLowerDKernel::computeLowerDJacobian(), LowerDIntegratedBC::computeLowerDOffDiagJacobian(), ArrayLowerDIntegratedBC::computeLowerDOffDiagJacobian(), DGKernel::computeOffDiagElemNeighJacobian(), ArrayDGKernel::computeOffDiagElemNeighJacobian(), VectorKernel::computeOffDiagJacobian(), Kernel::computeOffDiagJacobian(), ArrayKernel::computeOffDiagJacobian(), EigenKernel::computeOffDiagJacobian(), VectorIntegratedBC::computeOffDiagJacobian(), IntegratedBC::computeOffDiagJacobian(), ArrayIntegratedBC::computeOffDiagJacobian(), NodeElemConstraint::computeOffDiagJacobian(), NonlocalIntegratedBC::computeOffDiagJacobian(), NonlocalKernel::computeOffDiagJacobian(), KernelGrad::computeOffDiagJacobian(), KernelValue::computeOffDiagJacobian(), NodeFaceConstraint::computeOffDiagJacobian(), ODEKernel::computeOffDiagJacobianScalar(), VectorKernel::computeOffDiagJacobianScalar(), ArrayKernel::computeOffDiagJacobianScalar(), IntegratedBC::computeOffDiagJacobianScalar(), VectorIntegratedBC::computeOffDiagJacobianScalar(), Kernel::computeOffDiagJacobianScalar(), ArrayIntegratedBC::computeOffDiagJacobianScalar(), ScalarLagrangeMultiplier::computeOffDiagJacobianScalar(), MortarScalarBase::computeOffDiagJacobianScalar(), KernelScalarBase::computeOffDiagJacobianScalarLocal(), DGLowerDKernel::computeOffDiagLowerDJacobian(), ArrayDGLowerDKernel::computeOffDiagLowerDJacobian(), MortarScalarBase::computeScalarJacobian(), KernelScalarBase::computeScalarJacobian(), MortarScalarBase::computeScalarOffDiagJacobian(), KernelScalarBase::computeScalarOffDiagJacobian(), KernelScalarBase::computeScalarOffDiagJacobianScalar(), MortarScalarBase::computeScalarOffDiagJacobianScalar(), prepareMatrixTag(), prepareMatrixTagLower(), and prepareMatrixTagNeighbor().

◆ _local_re

DenseVector<Number> TaggingInterface::_local_re
protected

Holds local residual entries as they are accumulated by this Kernel.

Definition at line 364 of file TaggingInterface.h.

Referenced by accumulateTaggedLocalResidual(), FVInterfaceKernel::addResidual(), assignTaggedLocalResidual(), ADDGKernel::computeElemNeighResidual(), DGKernel::computeElemNeighResidual(), ElemElemConstraint::computeElemNeighResidual(), ArrayDGKernel::computeElemNeighResidual(), DGLowerDKernel::computeLowerDResidual(), ArrayDGLowerDKernel::computeLowerDResidual(), ScalarKernel::computeResidual(), VectorKernel::computeResidual(), Kernel::computeResidual(), LowerDIntegratedBC::computeResidual(), ArrayKernel::computeResidual(), ODEKernel::computeResidual(), ADScalarKernel::computeResidual(), VectorTimeKernel::computeResidual(), TimeKernel::computeResidual(), ODETimeKernel::computeResidual(), ArrayLowerDIntegratedBC::computeResidual(), VectorIntegratedBC::computeResidual(), IntegratedBC::computeResidual(), ArrayIntegratedBC::computeResidual(), EigenKernel::computeResidual(), NodeElemConstraint::computeResidual(), NodalEqualValueConstraint::computeResidual(), ADMortarConstraint::computeResidual(), FVScalarLagrangeMultiplierConstraint::computeResidual(), FVBoundaryScalarLagrangeMultiplierConstraint::computeResidual(), FVFluxBC::computeResidual(), MortarConstraint::computeResidual(), KernelValue::computeResidual(), KernelGrad::computeResidual(), FVElementalKernel::computeResidual(), FVFluxKernel::computeResidual(), NodeFaceConstraint::computeResidual(), prepareVectorTagInternal(), prepareVectorTagLower(), and prepareVectorTagNeighbor().

◆ _matrix_tags

std::set<TagID> TaggingInterface::_matrix_tags
private

◆ _moose_object

const MooseObject& TaggingInterface::_moose_object
private

Moose objct this tag works on.

Definition at line 412 of file TaggingInterface.h.

Referenced by TaggingInterface().

◆ _non_ref_abs_vector_tags

std::set<TagID> TaggingInterface::_non_ref_abs_vector_tags
private

A set to hold absolute value vector tags excluding the reference residual tag.

If there is no reference residual problem, this container is the same as _abs_vector_tags;

Definition at line 399 of file TaggingInterface.h.

Referenced by prepareVectorTag(), and TaggingInterface().

◆ _non_ref_vector_tags

std::set<TagID> TaggingInterface::_non_ref_vector_tags
private

A set to hold vector tags excluding the reference residual tag.

If there is no reference residual problem, this container is the same as _vector_tags;

Definition at line 395 of file TaggingInterface.h.

Referenced by prepareVectorTag(), and TaggingInterface().

◆ _nonlocal_ke

DenseMatrix<Number> TaggingInterface::_nonlocal_ke
protected

◆ _re_blocks

std::vector<DenseVector<Number> *> TaggingInterface::_re_blocks
private

◆ _ref_abs_vector_tags

std::set<TagID> TaggingInterface::_ref_abs_vector_tags
private

A set of either size 1 or 0.

If we have a reference residual problem and _abs_vector_tags holds the reference vector tag, then this set holds the reference vector tags, otherwise it holds nothing

Definition at line 409 of file TaggingInterface.h.

Referenced by prepareVectorTag(), and TaggingInterface().

◆ _ref_vector_tags

std::set<TagID> TaggingInterface::_ref_vector_tags
private

A set of either size 1 or 0.

If we have a reference residual problem and _vector_tags holds the reference vector tag, then this set holds the reference vector tags, otherwise it holds nothing

Definition at line 404 of file TaggingInterface.h.

Referenced by prepareVectorTag(), and TaggingInterface().

◆ _subproblem

SubProblem& TaggingInterface::_subproblem
protected

◆ _tag_params

const InputParameters& TaggingInterface::_tag_params
private

Parameters from moose object.

Definition at line 415 of file TaggingInterface.h.

Referenced by TaggingInterface().

◆ _vector_tags

std::set<TagID> TaggingInterface::_vector_tags
private

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