https://mooseframework.inl.gov
TaggingInterface.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "TaggingInterface.h"
11 #include "Conversion.h"
12 #include "FEProblem.h"
13 #include "Assembly.h"
16 
17 #include "libmesh/dense_vector.h"
18 
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 }
55 
57  : _subproblem(*moose_object->parameters().getCheckedPointerParam<SubProblem *>("_subproblem")),
58  _moose_object(*moose_object),
59  _tag_params(_moose_object.parameters())
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 }
158 
159 #ifdef MOOSE_KOKKOS_ENABLED
162  : _subproblem(object._subproblem),
163  _moose_object(object._moose_object),
164  _tag_params(object._tag_params)
165 {
166 }
167 #endif
168 
169 void
171 {
172  if (!_subproblem.vectorTagExists(tag_name))
173  mooseError("Vector tag ", tag_name, " does not exist in system");
174 
175  _vector_tags.insert(_subproblem.getVectorTagID(tag_name));
176 }
177 
178 void
180 {
181  if (!_subproblem.matrixTagExists(tag_name))
182  mooseError("Matrix tag ", tag_name, " does not exist in system");
183 
184  _matrix_tags.insert(_subproblem.getMatrixTagID(tag_name));
185 }
186 
187 void
189 {
190  if (!_subproblem.vectorTagExists(tag_id))
191  mooseError("Vector tag ", tag_id, " does not exist in system");
192 
193  _vector_tags.insert(tag_id);
194 }
195 
196 void
198 {
199  if (!_subproblem.matrixTagExists(tag_id))
200  mooseError("Matrix tag ", tag_id, " does not exist in system");
201 
202  _matrix_tags.insert(tag_id);
203 }
204 
205 void
206 TaggingInterface::prepareVectorTag(Assembly & assembly, const unsigned int ivar)
207 {
209 }
210 
211 void
213  const unsigned int ivar,
214  const ResidualTagType tag_type)
215 {
216  if (tag_type == ResidualTagType::NonReference)
218  else
220 }
221 
222 void
224  const unsigned int ivar,
225  const std::set<TagID> & vector_tags,
226  const std::set<TagID> & absolute_value_vector_tags)
227 {
228  auto prepare = [this, ivar, &assembly](auto & re_blocks, const auto & tags)
229  {
230  re_blocks.clear();
231  re_blocks.reserve(tags.size());
232  for (const auto tag_id : tags)
233  {
234  const auto & tag = _subproblem.getVectorTag(tag_id);
235  re_blocks.push_back(&assembly.residualBlock(ivar, Assembly::LocalDataKey{}, tag._type_id));
236  }
237  };
238 
239  prepare(_re_blocks, vector_tags);
240  prepare(_absre_blocks, absolute_value_vector_tags);
241 
242  _local_re.resize(_re_blocks.empty()
243  ? (_absre_blocks.empty() ? std::size_t(0) : _absre_blocks[0]->size())
244  : _re_blocks[0]->size());
245 }
246 
247 void
249 {
250  _re_blocks.resize(_vector_tags.size());
251  mooseAssert(_vector_tags.size() >= 1, "we need at least one active tag");
252  auto vector_tag = _vector_tags.begin();
253  for (MooseIndex(_vector_tags) i = 0; i < _vector_tags.size(); i++, ++vector_tag)
254  {
255  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
257  }
258  _local_re.resize(_re_blocks[0]->size());
259 
260  _absre_blocks.resize(_abs_vector_tags.size());
261  vector_tag = _abs_vector_tags.begin();
262  for (MooseIndex(_abs_vector_tags) i = 0; i < _abs_vector_tags.size(); i++, ++vector_tag)
263  {
264  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
265  _absre_blocks[i] =
266  &assembly.residualBlockNeighbor(ivar, Assembly::LocalDataKey{}, tag._type_id);
267  }
268 }
269 
270 void
271 TaggingInterface::prepareVectorTagLower(Assembly & assembly, unsigned int ivar)
272 {
273  _re_blocks.resize(_vector_tags.size());
274  mooseAssert(_vector_tags.size() >= 1, "we need at least one active tag");
275  auto vector_tag = _vector_tags.begin();
276  for (MooseIndex(_vector_tags) i = 0; i < _vector_tags.size(); i++, ++vector_tag)
277  {
278  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
279  _re_blocks[i] = &assembly.residualBlockLower(ivar, Assembly::LocalDataKey{}, tag._type_id);
280  }
281  _local_re.resize(_re_blocks[0]->size());
282 
283  _absre_blocks.resize(_abs_vector_tags.size());
284  vector_tag = _abs_vector_tags.begin();
285  for (MooseIndex(_abs_vector_tags) i = 0; i < _abs_vector_tags.size(); i++, ++vector_tag)
286  {
287  const VectorTag & tag = _subproblem.getVectorTag(*vector_tag);
289  }
290 }
291 
292 void
293 TaggingInterface::prepareMatrixTag(Assembly & assembly, unsigned int ivar, unsigned int jvar)
294 {
295  _ke_blocks.resize(_matrix_tags.size());
296  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
297  auto mat_vector = _matrix_tags.begin();
298  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
299  _ke_blocks[i] = &assembly.jacobianBlock(ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
300 
301  _local_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
302 }
303 
304 void
306  unsigned int ivar,
307  unsigned int jvar,
308  DenseMatrix<Number> & k) const
309 {
310  mooseAssert(!_matrix_tags.empty(), "No matrix tags exist");
311  const auto & ij_mat =
312  assembly.jacobianBlock(ivar, jvar, Assembly::LocalDataKey{}, *_matrix_tags.begin());
313  k.resize(ij_mat.m(), ij_mat.n());
314 }
315 
316 void
318  unsigned int ivar,
319  unsigned int jvar)
320 {
321  _ke_blocks.resize(_matrix_tags.size());
322  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
323  auto mat_vector = _matrix_tags.begin();
324  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
325  _ke_blocks[i] =
326  &assembly.jacobianBlockNonlocal(ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
327 
328  _nonlocal_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
329 }
330 
331 void
333  unsigned int ivar,
334  unsigned int jvar,
336 {
337  _ke_blocks.resize(_matrix_tags.size());
338  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
339  auto mat_vector = _matrix_tags.begin();
340  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
341  _ke_blocks[i] =
342  &assembly.jacobianBlockNeighbor(type, ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
343 
344  _local_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
345 }
346 
347 void
349  unsigned int ivar,
350  unsigned int jvar,
352  DenseMatrix<Number> & k) const
353 {
354  mooseAssert(!_matrix_tags.empty(), "No matrix tags exist");
355  const auto & ij_mat = assembly.jacobianBlockNeighbor(
356  type, ivar, jvar, Assembly::LocalDataKey{}, *_matrix_tags.begin());
357  k.resize(ij_mat.m(), ij_mat.n());
358 }
359 
360 void
362  unsigned int ivar,
363  unsigned int jvar,
365 {
366  _ke_blocks.resize(_matrix_tags.size());
367  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
368  auto mat_vector = _matrix_tags.begin();
369  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
370  _ke_blocks[i] =
371  &assembly.jacobianBlockMortar(type, ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
372 
373  _local_ke.resize(_ke_blocks[0]->m(), _ke_blocks[0]->n());
374 }
375 
376 void
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 }
385 
386 void
388 {
389  for (auto & re : _re_blocks)
390  *re = _local_re;
391  for (auto & absre : _absre_blocks)
392  for (const auto i : index_range(_local_re))
393  (*absre)(i) = std::abs(_local_re(i));
394 }
395 
396 void
398 {
399  for (auto & ke : _ke_blocks)
400  *ke += _local_ke;
401 }
402 
403 void
405  const unsigned int ivar,
406  const unsigned int jvar,
407  const DenseMatrix<Number> & k)
408 {
409  _ke_blocks.resize(_matrix_tags.size());
410  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
411  auto mat_vector = _matrix_tags.begin();
412  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
413  _ke_blocks[i] = &assembly.jacobianBlock(ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
414  mooseAssert(_ke_blocks[0]->m() == k.m() && _ke_blocks[0]->n() == k.n(),
415  "Passed-in k must match the blocks we are about to sum into");
416  for (auto & ke : _ke_blocks)
417  *ke += k;
418 }
419 
420 void
422  const unsigned int ivar,
423  const unsigned int jvar,
424  const Moose::DGJacobianType type,
425  const DenseMatrix<Number> & k)
426 {
427  _ke_blocks.resize(_matrix_tags.size());
428  mooseAssert(_matrix_tags.size() >= 1, "we need at least one active tag");
429  auto mat_vector = _matrix_tags.begin();
430  for (MooseIndex(_matrix_tags) i = 0; i < _matrix_tags.size(); i++, ++mat_vector)
431  _ke_blocks[i] =
432  &assembly.jacobianBlockNeighbor(type, ivar, jvar, Assembly::LocalDataKey{}, *mat_vector);
433  mooseAssert(_ke_blocks[0]->m() == k.m() && _ke_blocks[0]->n() == k.n(),
434  "Passed-in k must match the blocks we are about to sum into");
435  for (auto & ke : _ke_blocks)
436  *ke += k;
437 }
438 
439 void
441 {
442  for (auto & ke : _ke_blocks)
443  *ke += _nonlocal_ke;
444 }
445 
446 void
448 {
449  for (auto & ke : _ke_blocks)
450  *ke = _local_ke;
451 }
452 
virtual TagID getVectorTagID(const TagName &tag_name) const
Get a TagID from a TagName.
Definition: SubProblem.C:203
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
void accumulateTaggedNonlocalMatrix()
Nonlocal Jacobian blocks will be appended by adding the current nonlocal kernel Jacobian.
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.
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
Keeps track of stuff related to assembling.
Definition: Assembly.h:109
unsigned int TagID
Definition: MooseTypes.h:210
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:1153
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
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.
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:3112
void assignTaggedLocalMatrix()
Local Jacobian blocks will assigned as the current local kernel Jacobian.
Class that is used as a parameter to some of our matrix tag APIs that allows only blessed framework c...
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseBase.h:131
void resize(const unsigned int n)
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...
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
void assignTaggedLocalResidual()
Local residual blocks will assigned as the current local kernel residual.
std::vector< DenseVector< Number > * > _absre_blocks
Residual blocks for absolute value residual tags.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int m() const
ResidualTagType
Enumerate whether a (residual) vector tag is to be of a non-reference or reference tag type...
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.
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:1142
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
std::set< TagID > _abs_vector_tags
The absolute value residual tag ids.
InputParameters emptyInputParameters()
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:3153
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...
TaggingInterface(const MooseObject *moose_object)
std::set< TagID > _non_ref_abs_vector_tags
A set to hold absolute value vector tags excluding the reference residual tag.
void prepareMatrixTagNonlocal(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing nonlocal element jacobian according to the active tags.
const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:103
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:27
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
DenseVector< Number > & residualBlockLower(unsigned int var_num, LocalDataKey, TagID tag_id)
Get residual block for lower.
Definition: Assembly.h:1133
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.
virtual bool vectorTagExists(const TagID tag_id) const
Check to see if a particular Tag exists.
Definition: SubProblem.h:201
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
void useMatrixTag(const TagName &tag_name, MatrixTagsKey)
Uses a reference residual to define relative convergence criteria.
DGJacobianType
Definition: MooseTypes.h:750
void prepareVectorTagNeighbor(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags for DG and interface kernels...
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:78
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:1124
DenseVector< Number > & residualBlock(unsigned int var_num, LocalDataKey, TagID tag_id)
Get local residual block for a variable and a tag.
Definition: Assembly.h:1115
ConstraintJacobianType
Definition: MooseTypes.h:797
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.
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
const MooseObject & _moose_object
Moose objct this tag works on.
void resize(const unsigned int new_m, const unsigned int new_n)
virtual ~TaggingInterface()
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...
std::vector< DenseMatrix< Number > * > _ke_blocks
Kernel blocks Vectors For each Tag.
void prepareVectorTagLower(Assembly &assembly, unsigned int ivar)
Prepare data for computing the residual according to active tags for mortar constraints.
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...
unsigned int n() const
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
Storage for all of the information pretaining to a vector tag.
Definition: VectorTag.h:17
DenseMatrix< Number > _nonlocal_ke
Holds nonlocal Jacobian entries as they are accumulated by this Kernel.
static InputParameters validParams()
virtual const VectorTag & getVectorTag(const TagID tag_id) const
Get a VectorTag from a TagID.
Definition: SubProblem.C:161
auto index_range(const T &sizable)
Key structure for APIs adding/caching local element residuals/Jacobians.
Definition: Assembly.h:862
void useVectorTag(const TagName &tag_name, VectorTagsKey)
Class that is used as a parameter to some of our vector tag APIs that allows only blessed framework c...
virtual bool matrixTagExists(const TagName &tag_name) const
Check to see if a particular Tag exists.
Definition: SubProblem.C:328
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...