www.mooseframework.org
MooseVariableData.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 "MooseVariableData.h"
11 #include "MooseVariableField.h"
12 #include "Assembly.h"
13 #include "MooseError.h"
14 #include "DisplacedSystem.h"
15 #include "TimeIntegrator.h"
16 #include "MooseVariableFE.h"
17 #include "MooseTypes.h"
18 #include "DualRealOps.h"
19 
20 #include "libmesh/quadrature.h"
21 #include "libmesh/fe_base.h"
22 #include "libmesh/system.h"
23 #include "libmesh/type_n_tensor.h"
24 
25 template <typename OutputType>
27  SystemBase & sys,
28  THREAD_ID tid,
29  Moose::ElementType element_type,
30  const QBase * const & qrule_in,
31  const QBase * const & qrule_face_in,
32  const Node * const & node,
33  const Elem * const & elem)
34 
35  : MooseVariableDataBase<OutputType>(var, sys, tid),
36  _fe_type(var.feType()),
37  _var_num(var.number()),
38  _assembly(_subproblem.assembly(_tid, var.kind() == Moose::VAR_NONLINEAR ? sys.number() : 0)),
39  _element_type(element_type),
40  _ad_zero(0),
41  _need_ad_u_dot(false),
42  _need_ad_u_dotdot(false),
43  _need_second(false),
44  _need_second_old(false),
45  _need_second_older(false),
46  _need_second_previous_nl(false),
47  _need_curl(false),
48  _need_curl_old(false),
49  _need_curl_older(false),
50  _need_div(false),
51  _need_div_old(false),
52  _need_div_older(false),
53  _need_ad(false),
54  _need_ad_u(false),
55  _need_ad_grad_u(false),
56  _need_ad_grad_u_dot(false),
57  _need_ad_second_u(false),
58  _has_dof_indices(false),
59  _qrule(qrule_in),
60  _qrule_face(qrule_face_in),
61  _use_dual(var.useDual()),
62  _second_phi_assembly_method(nullptr),
63  _second_phi_face_assembly_method(nullptr),
64  _curl_phi_assembly_method(nullptr),
65  _curl_phi_face_assembly_method(nullptr),
66  _div_phi_assembly_method(nullptr),
67  _div_phi_face_assembly_method(nullptr),
68  _ad_grad_phi_assembly_method(nullptr),
69  _ad_grad_phi_face_assembly_method(nullptr),
70  _time_integrator(nullptr),
71  _node(node),
72  _elem(elem),
73  _displaced(dynamic_cast<const DisplacedSystem *>(&_sys) ? true : false),
74  _current_side(_assembly.side())
75 {
76  // FIXME: continuity of FE type seems equivalent with the definition of nodal variables.
77  // Continuity does not depend on the FE dimension, so we just pass in a valid dimension.
78  if (_fe_type.family == NEDELEC_ONE || _fe_type.family == LAGRANGE_VEC ||
79  _fe_type.family == MONOMIAL_VEC || _fe_type.family == RAVIART_THOMAS)
80  _continuity = _assembly.getVectorFE(_fe_type, _sys.mesh().dimension())->get_continuity();
81  else
82  _continuity = _assembly.getFE(_fe_type, _sys.mesh().dimension())->get_continuity();
83 
85 
87 
88  switch (_element_type)
89  {
91  {
92  _phi_assembly_method = &Assembly::fePhi<OutputShape>;
93  _phi_face_assembly_method = &Assembly::fePhiFace<OutputShape>;
94  _grad_phi_assembly_method = &Assembly::feGradPhi<OutputShape>;
95  _grad_phi_face_assembly_method = &Assembly::feGradPhiFace<OutputShape>;
96  _second_phi_assembly_method = &Assembly::feSecondPhi<OutputShape>;
97  _second_phi_face_assembly_method = &Assembly::feSecondPhiFace<OutputShape>;
98  _curl_phi_assembly_method = &Assembly::feCurlPhi<OutputShape>;
99  _curl_phi_face_assembly_method = &Assembly::feCurlPhiFace<OutputShape>;
100  _div_phi_assembly_method = &Assembly::feDivPhi<OutputShape>;
101  _div_phi_face_assembly_method = &Assembly::feDivPhiFace<OutputShape>;
102  _ad_grad_phi_assembly_method = &Assembly::feADGradPhi<OutputShape>;
103  _ad_grad_phi_face_assembly_method = &Assembly::feADGradPhiFace<OutputShape>;
104 
107  break;
108  }
110  {
111  _phi_assembly_method = &Assembly::fePhiNeighbor<OutputShape>;
112  _phi_face_assembly_method = &Assembly::fePhiFaceNeighbor<OutputShape>;
113  _grad_phi_assembly_method = &Assembly::feGradPhiNeighbor<OutputShape>;
114  _grad_phi_face_assembly_method = &Assembly::feGradPhiFaceNeighbor<OutputShape>;
115  _second_phi_assembly_method = &Assembly::feSecondPhiNeighbor<OutputShape>;
116  _second_phi_face_assembly_method = &Assembly::feSecondPhiFaceNeighbor<OutputShape>;
117  _curl_phi_assembly_method = &Assembly::feCurlPhiNeighbor<OutputShape>;
118  _curl_phi_face_assembly_method = &Assembly::feCurlPhiFaceNeighbor<OutputShape>;
119  _div_phi_assembly_method = &Assembly::feDivPhiNeighbor<OutputShape>;
120  _div_phi_face_assembly_method = &Assembly::feDivPhiFaceNeighbor<OutputShape>;
121 
122  _ad_grad_phi = nullptr;
123  _ad_grad_phi_face = nullptr;
124  break;
125  }
127  {
128  if (_use_dual)
129  {
130  _phi_assembly_method = &Assembly::feDualPhiLower<OutputType>;
131  _phi_face_assembly_method = &Assembly::feDualPhiLower<OutputType>; // Place holder
132  _grad_phi_assembly_method = &Assembly::feGradDualPhiLower<OutputType>;
133  _grad_phi_face_assembly_method = &Assembly::feGradDualPhiLower<OutputType>; // Place holder
134  }
135  else
136  {
137  _phi_assembly_method = &Assembly::fePhiLower<OutputType>;
138  _phi_face_assembly_method = &Assembly::fePhiLower<OutputType>; // Place holder
139  _grad_phi_assembly_method = &Assembly::feGradPhiLower<OutputType>;
140  _grad_phi_face_assembly_method = &Assembly::feGradPhiLower<OutputType>; // Place holder
141  }
142 
143  _ad_grad_phi = nullptr;
144  _ad_grad_phi_face = nullptr;
145  break;
146  }
147  }
152 }
153 
154 template <typename OutputType>
155 void
157 {
158  switch (gm_type)
159  {
160  case Moose::Volume:
161  {
162  _current_qrule = _qrule;
163  _current_phi = _phi;
164  _current_grad_phi = _grad_phi;
165  _current_second_phi = _second_phi;
166  _current_curl_phi = _curl_phi;
167  _current_div_phi = _div_phi;
168  _current_ad_grad_phi = _ad_grad_phi;
169  break;
170  }
171  case Moose::Face:
172  {
173  _current_qrule = _qrule_face;
174  _current_phi = _phi_face;
175  _current_grad_phi = _grad_phi_face;
176  _current_second_phi = _second_phi_face;
177  _current_curl_phi = _curl_phi_face;
178  _current_div_phi = _div_phi_face;
179  _current_ad_grad_phi = _ad_grad_phi_face;
180  break;
181  }
182  }
183 }
184 
185 template <typename OutputType>
188 {
189  if (_sys.solutionUDot())
190  {
191  _need_u_dot = true;
192  return _u_dot;
193  }
194  else
195  mooseError("MooseVariableFE: Time derivative of solution (`u_dot`) is not stored. Please set "
196  "uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
197 }
198 
199 template <typename OutputType>
202 {
203  if (_sys.solutionUDotDot())
204  {
205  _need_u_dotdot = true;
206  return _u_dotdot;
207  }
208  else
209  mooseError("MooseVariableFE: Second time derivative of solution (`u_dotdot`) is not stored. "
210  "Please set uDotDotRequested() to true in FEProblemBase before requesting "
211  "`u_dotdot`.");
212 }
213 
214 template <typename OutputType>
217 {
218  if (_sys.solutionUDotOld())
219  {
220  _need_u_dot_old = true;
221  return _u_dot_old;
222  }
223  else
224  mooseError("MooseVariableFE: Old time derivative of solution (`u_dot_old`) is not stored. "
225  "Please set uDotOldRequested() to true in FEProblemBase before requesting "
226  "`u_dot_old`.");
227 }
228 
229 template <typename OutputType>
232 {
233  if (_sys.solutionUDotDotOld())
234  {
235  _need_u_dotdot_old = true;
236  return _u_dotdot_old;
237  }
238  else
239  mooseError("MooseVariableFE: Old second time derivative of solution (`u_dotdot_old`) is not "
240  "stored. Please set uDotDotOldRequested() to true in FEProblemBase before "
241  "requesting `u_dotdot_old`");
242 }
243 
244 template <typename OutputType>
247 {
248  if (_sys.solutionUDot())
249  {
250  _need_grad_dot = true;
251  return _grad_u_dot;
252  }
253  else
254  mooseError("MooseVariableFE: Time derivative of solution (`u_dot`) is not stored. Please set "
255  "uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
256 }
257 
258 template <typename OutputType>
261 {
262  if (_sys.solutionUDotDot())
263  {
264  _need_grad_dotdot = true;
265  return _grad_u_dotdot;
266  }
267  else
268  mooseError("MooseVariableFE: Second time derivative of solution (`u_dotdot`) is not stored. "
269  "Please set uDotDotRequested() to true in FEProblemBase before requesting "
270  "`u_dotdot`.");
271 }
272 
273 template <typename OutputType>
276 {
277  secondPhi();
278  secondPhiFace();
279  switch (state)
280  {
281  case Moose::Current:
282  {
283  _need_second = true;
284  return _second_u;
285  }
286 
287  case Moose::Old:
288  {
289  _need_second_old = true;
290  return _second_u_old;
291  }
292 
293  case Moose::Older:
294  {
295  _need_second_older = true;
296  return _second_u_older;
297  }
298 
299  case Moose::PreviousNL:
300  {
301  _need_second_previous_nl = true;
302  return _second_u_previous_nl;
303  }
304 
305  default:
306  // We should never get here but gcc requires that we have a default. See
307  // htpps://stackoverflow.com/questions/18680378/after-defining-case-for-all-enum-values-compiler-still-says-control-reaches-e
308  mooseError("Unknown SolutionState!");
309  }
310 }
311 
312 template <typename OutputType>
315 {
316  curlPhi();
317  curlPhiFace();
318  switch (state)
319  {
320  case Moose::Current:
321  {
322  _need_curl = true;
323  return _curl_u;
324  }
325 
326  case Moose::Old:
327  {
328  _need_curl_old = true;
329  return _curl_u_old;
330  }
331 
332  case Moose::Older:
333  {
334  _need_curl_older = true;
335  return _curl_u_older;
336  }
337 
338  default:
339  mooseError("We don't currently support curl from the previous non-linear iteration");
340  }
341 }
342 
343 template <typename OutputType>
346 {
347  divPhi();
348  divPhiFace();
349  switch (state)
350  {
351  case Moose::Current:
352  {
353  _need_div = true;
354  return _div_u;
355  }
356 
357  case Moose::Old:
358  {
359  _need_div_old = true;
360  return _div_u_old;
361  }
362 
363  case Moose::Older:
364  {
365  _need_div_older = true;
366  return _div_u_older;
367  }
368 
369  default:
370  mooseError("We don't currently support divergence from the previous non-linear iteration");
371  }
372 }
373 
374 template <typename OutputType>
377 {
378  _second_phi = &_second_phi_assembly_method(_assembly, _fe_type);
379  return *_second_phi;
380 }
381 
382 template <typename OutputType>
385 {
386  _second_phi_face = &_second_phi_face_assembly_method(_assembly, _fe_type);
387  return *_second_phi_face;
388 }
389 
390 template <typename OutputType>
393 {
394  _curl_phi = &_curl_phi_assembly_method(_assembly, _fe_type);
395  return *_curl_phi;
396 }
397 
398 template <typename OutputType>
401 {
402  _curl_phi_face = &_curl_phi_face_assembly_method(_assembly, _fe_type);
403  return *_curl_phi_face;
404 }
405 
406 template <typename OutputType>
409 {
410  _div_phi = &_div_phi_assembly_method(_assembly, _fe_type);
411  return *_div_phi;
412 }
413 
414 template <typename OutputType>
417 {
418  _div_phi_face = &_div_phi_face_assembly_method(_assembly, _fe_type);
419  return *_div_phi_face;
420 }
421 
422 template <typename OutputType>
423 void
425 {
426  unsigned int num_dofs = _dof_indices.size();
427 
428  if (num_dofs > 0)
429  fetchDoFValues();
430 
431  bool is_transient = _subproblem.isTransient();
432  unsigned int nqp = _current_qrule->n_points();
433  auto && active_coupleable_matrix_tags = _subproblem.getActiveFEVariableCoupleableMatrixTags(_tid);
434 
435  for (auto tag : _required_vector_tags)
436  {
437  if (_need_vector_tag_u[tag])
438  _vector_tag_u[tag].resize(nqp);
439  if (_need_vector_tag_grad[tag])
440  _vector_tag_grad[tag].resize(nqp);
441  }
442 
443  for (auto tag : active_coupleable_matrix_tags)
444  if (_need_matrix_tag_u[tag])
445  _matrix_tag_u[tag].resize(nqp);
446 
447  if (_need_second)
448  _second_u.resize(nqp);
449 
450  if (_need_curl)
451  _curl_u.resize(nqp);
452 
453  if (_need_div)
454  _div_u.resize(nqp);
455 
456  if (_need_second_previous_nl)
457  _second_u_previous_nl.resize(nqp);
458 
459  if (is_transient)
460  {
461  if (_need_u_dot)
462  _u_dot.resize(nqp);
463 
464  if (_need_u_dotdot)
465  _u_dotdot.resize(nqp);
466 
467  if (_need_u_dot_old)
468  _u_dot_old.resize(nqp);
469 
470  if (_need_u_dotdot_old)
471  _u_dotdot_old.resize(nqp);
472 
473  if (_need_du_dot_du)
474  _du_dot_du.resize(nqp);
475 
476  if (_need_du_dotdot_du)
477  _du_dotdot_du.resize(nqp);
478 
479  if (_need_grad_dot)
480  _grad_u_dot.resize(nqp);
481 
482  if (_need_grad_dotdot)
483  _grad_u_dotdot.resize(nqp);
484 
485  if (_need_second_old)
486  _second_u_old.resize(nqp);
487 
488  if (_need_curl_old)
489  _curl_u_old.resize(nqp);
490 
491  if (_need_div_old)
492  _div_u_old.resize(nqp);
493 
494  if (_need_second_older)
495  _second_u_older.resize(nqp);
496  }
497 
498  for (unsigned int i = 0; i < nqp; ++i)
499  {
500  for (auto tag : _required_vector_tags)
501  {
502  if (_need_vector_tag_u[tag])
503  _vector_tag_u[tag][i] = 0;
504  if (_need_vector_tag_grad[tag])
505  _vector_tag_grad[tag][i] = 0;
506  }
507 
508  for (auto tag : active_coupleable_matrix_tags)
509  if (_need_matrix_tag_u[tag])
510  _matrix_tag_u[tag][i] = 0;
511 
512  if (_need_second)
513  _second_u[i] = 0;
514 
515  if (_need_curl)
516  _curl_u[i] = 0;
517 
518  if (_need_div)
519  _div_u[i] = 0;
520 
521  if (_need_second_previous_nl)
522  _second_u_previous_nl[i] = 0;
523 
524  if (is_transient)
525  {
526  if (_need_u_dot)
527  _u_dot[i] = 0;
528 
529  if (_need_u_dotdot)
530  _u_dotdot[i] = 0;
531 
532  if (_need_u_dot_old)
533  _u_dot_old[i] = 0;
534 
535  if (_need_u_dotdot_old)
536  _u_dotdot_old[i] = 0;
537 
538  if (_need_du_dot_du)
539  _du_dot_du[i] = 0;
540 
541  if (_need_du_dotdot_du)
542  _du_dotdot_du[i] = 0;
543 
544  if (_need_grad_dot)
545  _grad_u_dot[i] = 0;
546 
547  if (_need_grad_dotdot)
548  _grad_u_dotdot[i] = 0;
549 
550  if (_need_second_old)
551  _second_u_old[i] = 0;
552 
553  if (_need_second_older)
554  _second_u_older[i] = 0;
555 
556  if (_need_curl_old)
557  _curl_u_old[i] = 0;
558 
559  if (_need_div_old)
560  _div_u_old[i] = 0;
561  }
562  }
563 
564  bool second_required =
565  _need_second || _need_second_old || _need_second_older || _need_second_previous_nl;
566  bool curl_required = _need_curl || _need_curl_old;
567  bool div_required = _need_div || _need_div_old;
568 
569  for (unsigned int i = 0; i < num_dofs; i++)
570  {
571  for (unsigned int qp = 0; qp < nqp; qp++)
572  {
573  const OutputType phi_local = (*_current_phi)[i][qp];
574  const typename OutputTools<OutputType>::OutputGradient dphi_qp = (*_current_grad_phi)[i][qp];
575 
576  if (is_transient)
577  {
578  if (_need_u_dot)
579  _u_dot[qp] += phi_local * _dof_values_dot[i];
580 
581  if (_need_u_dotdot)
582  _u_dotdot[qp] += phi_local * _dof_values_dotdot[i];
583 
584  if (_need_u_dot_old)
585  _u_dot_old[qp] += phi_local * _dof_values_dot_old[i];
586 
587  if (_need_u_dotdot_old)
588  _u_dotdot_old[qp] += phi_local * _dof_values_dotdot_old[i];
589 
590  if (_need_grad_dot)
591  _grad_u_dot[qp].add_scaled(dphi_qp, _dof_values_dot[i]);
592 
593  if (_need_grad_dotdot)
594  _grad_u_dotdot[qp].add_scaled(dphi_qp, _dof_values_dotdot[i]);
595 
596  if (_need_du_dot_du)
597  _du_dot_du[qp] = _dof_du_dot_du[i];
598 
599  if (_need_du_dotdot_du)
600  _du_dotdot_du[qp] = _dof_du_dotdot_du[i];
601  }
602 
603  if (second_required)
604  {
605  mooseAssert(
606  _current_second_phi,
607  "We're requiring a second calculation but have not set a second shape function!");
608  const typename OutputTools<OutputType>::OutputSecond d2phi_local =
609  (*_current_second_phi)[i][qp];
610 
611  if (_need_second)
612  _second_u[qp].add_scaled(d2phi_local, _vector_tags_dof_u[_solution_tag][i]);
613 
614  if (_need_second_previous_nl)
615  _second_u_previous_nl[qp].add_scaled(d2phi_local,
616  _vector_tags_dof_u[_previous_nl_solution_tag][i]);
617 
618  if (is_transient)
619  {
620  if (_need_second_old)
621  _second_u_old[qp].add_scaled(d2phi_local, _vector_tags_dof_u[_old_solution_tag][i]);
622 
623  if (_need_second_older)
624  _second_u_older[qp].add_scaled(d2phi_local, _vector_tags_dof_u[_older_solution_tag][i]);
625  }
626  }
627 
628  if (curl_required)
629  {
630  mooseAssert(_current_curl_phi,
631  "We're requiring a curl calculation but have not set a curl shape function!");
632  const OutputType curl_phi_local = (*_current_curl_phi)[i][qp];
633 
634  if (_need_curl)
635  _curl_u[qp] += curl_phi_local * _vector_tags_dof_u[_solution_tag][i];
636 
637  if (is_transient && _need_curl_old)
638  _curl_u_old[qp] += curl_phi_local * _vector_tags_dof_u[_old_solution_tag][i];
639  }
640 
641  if (div_required)
642  {
643  mooseAssert(
644  _current_div_phi,
645  "We're requiring a divergence calculation but have not set a div shape function!");
646  const OutputShapeDivergence div_phi_local = (*_current_div_phi)[i][qp];
647 
648  if (_need_div)
649  _div_u[qp] += div_phi_local * _vector_tags_dof_u[_solution_tag][i];
650 
651  if (is_transient && _need_div_old)
652  _div_u_old[qp] += div_phi_local * _vector_tags_dof_u[_old_solution_tag][i];
653  }
654 
655  for (auto tag : _required_vector_tags)
656  {
657  if (_sys.hasVector(tag) && _sys.getVector(tag).closed())
658  {
659  if (_need_vector_tag_u[tag])
660  _vector_tag_u[tag][qp] += phi_local * _vector_tags_dof_u[tag][i];
661  if (_need_vector_tag_grad[tag])
662  _vector_tag_grad[tag][qp].add_scaled(dphi_qp, _vector_tags_dof_u[tag][i]);
663  }
664  }
665 
666  for (auto tag : active_coupleable_matrix_tags)
667  if (_need_matrix_tag_u[tag])
668  _matrix_tag_u[tag][qp] += phi_local * _matrix_tags_dof_u[tag][i];
669  }
670  }
671 
672  // Automatic differentiation
673  if (_need_ad)
674  computeAD(num_dofs, nqp);
675 }
676 
677 template <>
678 void
680 {
681  unsigned int num_dofs = _dof_indices.size();
682 
683  if (num_dofs > 0)
684  fetchDoFValues();
685 
686  bool is_transient = _subproblem.isTransient();
687  unsigned int nqp = _current_qrule->n_points();
688  auto && active_coupleable_matrix_tags = _subproblem.getActiveFEVariableCoupleableMatrixTags(_tid);
689 
690  // Map grad_phi using Eigen so that we can perform array operations easier
691  if (_qrule == _current_qrule)
692  {
693  _mapped_grad_phi.resize(num_dofs);
694  for (unsigned int i = 0; i < num_dofs; i++)
695  {
696  _mapped_grad_phi[i].resize(nqp, Eigen::Map<RealDIMValue>(nullptr));
697  for (unsigned int qp = 0; qp < nqp; qp++)
698  // Note: this does NOT do any allocation. It is "reconstructing" the object in place
699  new (&_mapped_grad_phi[i][qp])
700  Eigen::Map<RealDIMValue>(const_cast<Real *>(&(*_current_grad_phi)[i][qp](0)));
701  }
702  }
703  else
704  {
705  _mapped_grad_phi_face.resize(num_dofs);
706  for (unsigned int i = 0; i < num_dofs; i++)
707  {
708  _mapped_grad_phi_face[i].resize(nqp, Eigen::Map<RealDIMValue>(nullptr));
709  for (unsigned int qp = 0; qp < nqp; qp++)
710  // Note: this does NOT do any allocation. It is "reconstructing" the object in place
711  new (&_mapped_grad_phi_face[i][qp])
712  Eigen::Map<RealDIMValue>(const_cast<Real *>(&(*_current_grad_phi)[i][qp](0)));
713  }
714  }
715 
716  for (auto tag : _required_vector_tags)
717  {
718  if (_need_vector_tag_u[tag])
719  _vector_tag_u[tag].resize(nqp);
720  if (_need_vector_tag_grad[tag])
721  _vector_tag_grad[tag].resize(nqp);
722  }
723 
724  for (auto tag : active_coupleable_matrix_tags)
725  if (_need_matrix_tag_u[tag])
726  _matrix_tag_u[tag].resize(nqp);
727 
728  if (_need_second)
729  _second_u.resize(nqp);
730 
731  if (_need_curl)
732  _curl_u.resize(nqp);
733 
734  if (_need_div)
735  _div_u.resize(nqp);
736 
737  if (_need_second_previous_nl)
738  _second_u_previous_nl.resize(nqp);
739 
740  if (is_transient)
741  {
742  if (_need_u_dot)
743  _u_dot.resize(nqp);
744 
745  if (_need_u_dotdot)
746  _u_dotdot.resize(nqp);
747 
748  if (_need_u_dot_old)
749  _u_dot_old.resize(nqp);
750 
751  if (_need_u_dotdot_old)
752  _u_dotdot_old.resize(nqp);
753 
754  if (_need_du_dot_du)
755  _du_dot_du.resize(nqp);
756 
757  if (_need_du_dotdot_du)
758  _du_dotdot_du.resize(nqp);
759 
760  if (_need_grad_dot)
761  _grad_u_dot.resize(nqp);
762 
763  if (_need_grad_dotdot)
764  _grad_u_dotdot.resize(nqp);
765 
766  if (_need_second_old)
767  _second_u_old.resize(nqp);
768 
769  if (_need_curl_old)
770  _curl_u_old.resize(nqp);
771 
772  if (_need_div_old)
773  _div_u_old.resize(nqp);
774 
775  if (_need_second_older)
776  _second_u_older.resize(nqp);
777  }
778 
779  for (unsigned int i = 0; i < nqp; ++i)
780  {
781  for (auto tag : _required_vector_tags)
782  {
783  if (_need_vector_tag_u[tag])
784  _vector_tag_u[tag][i].setZero(_count);
785  if (_need_vector_tag_grad[tag])
786  _vector_tag_grad[tag][i].setZero(_count, LIBMESH_DIM);
787  }
788 
789  for (auto tag : active_coupleable_matrix_tags)
790  if (_need_matrix_tag_u[tag])
791  _matrix_tag_u[tag][i].setZero(_count);
792 
793  if (_need_second)
794  _second_u[i].setZero(_count, LIBMESH_DIM * LIBMESH_DIM);
795 
796  if (_need_curl)
797  _curl_u[i].setZero(_count);
798 
799  if (_need_div)
800  _div_u[i].setZero(_count);
801 
802  if (_need_second_previous_nl)
803  _second_u_previous_nl[i].setZero(_count, LIBMESH_DIM * LIBMESH_DIM);
804 
805  if (is_transient)
806  {
807  if (_need_u_dot)
808  _u_dot[i].setZero(_count);
809 
810  if (_need_u_dotdot)
811  _u_dotdot[i].setZero(_count);
812 
813  if (_need_u_dot_old)
814  _u_dot_old[i].setZero(_count);
815 
816  if (_need_u_dotdot_old)
817  _u_dotdot_old[i].setZero(_count);
818 
819  if (_need_du_dot_du)
820  _du_dot_du[i] = 0;
821 
822  if (_need_du_dotdot_du)
823  _du_dotdot_du[i] = 0;
824 
825  if (_need_grad_dot)
826  _grad_u_dot[i].setZero(_count, LIBMESH_DIM);
827 
828  if (_need_grad_dotdot)
829  _grad_u_dotdot[i].setZero(_count, LIBMESH_DIM);
830 
831  if (_need_second_old)
832  _second_u_old[i].setZero(_count, LIBMESH_DIM * LIBMESH_DIM);
833 
834  if (_need_second_older)
835  _second_u_older[i].setZero(_count, LIBMESH_DIM * LIBMESH_DIM);
836 
837  if (_need_curl_old)
838  _curl_u_old[i].setZero(_count);
839 
840  if (_need_div_old)
841  _div_u_old[i].setZero(_count);
842  }
843  }
844 
845  bool second_required =
846  _need_second || _need_second_old || _need_second_older || _need_second_previous_nl;
847  bool curl_required = _need_curl || _need_curl_old;
848  bool div_required = _need_div || _need_div_old;
849 
850  for (unsigned int i = 0; i < num_dofs; i++)
851  {
852  for (unsigned int qp = 0; qp < nqp; qp++)
853  {
854  const OutputShape phi_local = (*_current_phi)[i][qp];
855  const OutputShapeGradient dphi_qp = (*_current_grad_phi)[i][qp];
856 
857  if (is_transient)
858  {
859  if (_need_u_dot)
860  _u_dot[qp] += phi_local * _dof_values_dot[i];
861 
862  if (_need_u_dotdot)
863  _u_dotdot[qp] += phi_local * _dof_values_dotdot[i];
864 
865  if (_need_u_dot_old)
866  _u_dot_old[qp] += phi_local * _dof_values_dot_old[i];
867 
868  if (_need_u_dotdot_old)
869  _u_dotdot_old[qp] += phi_local * _dof_values_dotdot_old[i];
870 
871  if (_need_grad_dot)
872  for (const auto d : make_range(Moose::dim))
873  _grad_u_dot[qp].col(d) += dphi_qp(d) * _dof_values_dot[i];
874 
875  if (_need_grad_dotdot)
876  for (const auto d : make_range(Moose::dim))
877  _grad_u_dotdot[qp].col(d) += dphi_qp(d) * _dof_values_dotdot[i];
878 
879  if (_need_du_dot_du)
880  _du_dot_du[qp] = _dof_du_dot_du[i];
881 
882  if (_need_du_dotdot_du)
883  _du_dotdot_du[qp] = _dof_du_dotdot_du[i];
884  }
885 
886  if (second_required)
887  {
888  mooseAssert(
889  _current_second_phi,
890  "We're requiring a second calculation but have not set a second shape function!");
891  const RealTensorValue d2phi_local = (*_current_second_phi)[i][qp];
892 
893  if (_need_second)
894  for (unsigned int d = 0, d1 = 0; d1 < LIBMESH_DIM; ++d1)
895  for (unsigned int d2 = 0; d2 < LIBMESH_DIM; ++d2)
896  _second_u[qp].col(d++) += d2phi_local(d1, d2) * _vector_tags_dof_u[_solution_tag][i];
897 
898  if (_need_second_previous_nl)
899  for (unsigned int d = 0, d1 = 0; d1 < LIBMESH_DIM; ++d1)
900  for (unsigned int d2 = 0; d2 < LIBMESH_DIM; ++d2)
901  _second_u_previous_nl[qp].col(d++) +=
902  d2phi_local(d1, d2) * _vector_tags_dof_u[_previous_nl_solution_tag][i];
903 
904  if (is_transient)
905  {
906  if (_need_second_old)
907  for (unsigned int d = 0, d1 = 0; d1 < LIBMESH_DIM; ++d1)
908  for (unsigned int d2 = 0; d2 < LIBMESH_DIM; ++d2)
909  _second_u_old[qp].col(d++) +=
910  d2phi_local(d1, d2) * _vector_tags_dof_u[_old_solution_tag][i];
911 
912  if (_need_second_older)
913  for (unsigned int d = 0, d1 = 0; d1 < LIBMESH_DIM; ++d1)
914  for (unsigned int d2 = 0; d2 < LIBMESH_DIM; ++d2)
915  _second_u_older[qp].col(d++) +=
916  d2phi_local(d1, d2) * _vector_tags_dof_u[_older_solution_tag][i];
917  }
918  }
919 
920  if (curl_required)
921  {
922  mooseAssert(_current_curl_phi,
923  "We're requiring a curl calculation but have not set a curl shape function!");
924  const auto curl_phi_local = (*_current_curl_phi)[i][qp];
925 
926  if (_need_curl)
927  _curl_u[qp] += curl_phi_local * _vector_tags_dof_u[_solution_tag][i];
928 
929  if (is_transient && _need_curl_old)
930  _curl_u_old[qp] += curl_phi_local * _vector_tags_dof_u[_old_solution_tag][i];
931  }
932 
933  if (div_required)
934  {
935  mooseAssert(_current_div_phi,
936  "We're requiring a divergence calculation but have not set a divergence shape "
937  "function!");
938  const auto div_phi_local = (*_current_div_phi)[i][qp];
939 
940  if (_need_div)
941  _div_u[qp] += div_phi_local * _vector_tags_dof_u[_solution_tag][i];
942 
943  if (is_transient && _need_div_old)
944  _div_u_old[qp] += div_phi_local * _vector_tags_dof_u[_old_solution_tag][i];
945  }
946 
947  for (auto tag : _required_vector_tags)
948  {
949  if (_need_vector_tag_u[tag])
950  _vector_tag_u[tag][qp] += phi_local * _vector_tags_dof_u[tag][i];
951  if (_need_vector_tag_grad[tag])
952  for (const auto d : make_range(Moose::dim))
953  _vector_tag_grad[tag][qp].col(d) += dphi_qp(d) * _vector_tags_dof_u[tag][i];
954  }
955 
956  for (auto tag : active_coupleable_matrix_tags)
957  if (_need_matrix_tag_u[tag])
958  _matrix_tag_u[tag][qp] += phi_local * _matrix_tags_dof_u[tag][i];
959  }
960  }
961  // No AD support for array variable yet.
962 }
963 
964 template <typename OutputType>
965 void
967 {
968  if (_dof_indices.size() == 0)
969  return;
970 
971  if (_elem->p_level())
972  {
973  // The optimizations in this routine are not appropriate after p-refinement
974  computeValues();
975  return;
976  }
977 
978  bool is_transient = _subproblem.isTransient();
979  unsigned int nqp = _current_qrule->n_points();
980 
981  if (_need_second)
982  _second_u.resize(nqp);
983 
984  if (_need_second_previous_nl)
985  _second_u_previous_nl.resize(nqp);
986 
987  if (is_transient)
988  {
989  if (_need_u_dot)
990  _u_dot.resize(nqp);
991 
992  if (_need_u_dotdot)
993  _u_dotdot.resize(nqp);
994 
995  if (_need_u_dot_old)
996  _u_dot_old.resize(nqp);
997 
998  if (_need_u_dotdot_old)
999  _u_dotdot_old.resize(nqp);
1000 
1001  if (_need_du_dot_du)
1002  _du_dot_du.resize(nqp);
1003 
1004  if (_need_du_dotdot_du)
1005  _du_dotdot_du.resize(nqp);
1006 
1007  if (_need_second_old)
1008  _second_u_old.resize(nqp);
1009 
1010  if (_need_second_older)
1011  _second_u_older.resize(nqp);
1012  }
1013 
1014  if (is_transient)
1015  {
1016  if (_need_dof_values_dot)
1017  _dof_values_dot.resize(1);
1018  if (_need_dof_values_dotdot)
1019  _dof_values_dotdot.resize(1);
1020  if (_need_dof_values_dot_old)
1021  _dof_values_dot_old.resize(1);
1022  if (_need_dof_values_dotdot_old)
1023  _dof_values_dotdot_old.resize(1);
1024  }
1025 
1026  const dof_id_type & idx = _dof_indices[0];
1027  Real u_dot = 0;
1028  Real u_dotdot = 0;
1029  Real u_dot_old = 0;
1030  Real u_dotdot_old = 0;
1031  const Real & du_dot_du = _sys.duDotDu();
1032  const Real & du_dotdot_du = _sys.duDotDotDu();
1033 
1034  if (is_transient)
1035  {
1036  if (_sys.solutionUDot())
1037  u_dot = (*_sys.solutionUDot())(idx);
1038  if (_sys.solutionUDotDot())
1039  u_dotdot = (*_sys.solutionUDotDot())(idx);
1040  if (_sys.solutionUDotOld())
1041  u_dot_old = (*_sys.solutionUDotOld())(idx);
1042  if (_sys.solutionUDotDotOld())
1043  u_dotdot_old = (*_sys.solutionUDotDotOld())(idx);
1044 
1045  if (_need_dof_values_dot)
1046  _dof_values_dot[0] = u_dot;
1047 
1048  if (_need_dof_values_dotdot)
1049  _dof_values_dotdot[0] = u_dotdot;
1050  }
1051 
1052  auto phi = (*_current_phi)[0][0];
1053 
1054  if (is_transient)
1055  {
1056  if (_need_u_dot)
1057  _u_dot[0] = phi * u_dot;
1058 
1059  if (_need_u_dotdot)
1060  _u_dotdot[0] = phi * u_dotdot;
1061 
1062  if (_need_u_dot_old)
1063  _u_dot_old[0] = phi * u_dot_old;
1064 
1065  if (_need_u_dotdot_old)
1066  _u_dotdot_old[0] = phi * u_dotdot_old;
1067 
1068  if (_need_du_dot_du)
1069  _du_dot_du[0] = du_dot_du;
1070 
1071  if (_need_du_dotdot_du)
1072  _du_dotdot_du[0] = du_dotdot_du;
1073  }
1074 
1075  for (unsigned qp = 1; qp < nqp; ++qp)
1076  {
1077  if (is_transient)
1078  {
1079  if (_need_u_dot)
1080  _u_dot[qp] = _u_dot[0];
1081 
1082  if (_need_u_dotdot)
1083  _u_dotdot[qp] = _u_dotdot[0];
1084 
1085  if (_need_u_dot_old)
1086  _u_dot_old[qp] = _u_dot_old[0];
1087 
1088  if (_need_u_dotdot_old)
1089  _u_dotdot_old[qp] = _u_dotdot_old[0];
1090 
1091  if (_need_du_dot_du)
1092  _du_dot_du[qp] = _du_dot_du[0];
1093 
1094  if (_need_du_dotdot_du)
1095  _du_dotdot_du[qp] = _du_dotdot_du[0];
1096  }
1097  }
1098 
1099  auto && active_coupleable_matrix_tags = _subproblem.getActiveFEVariableCoupleableMatrixTags(_tid);
1100 
1101  for (auto tag : _required_vector_tags)
1102  {
1103  if (_need_vector_tag_u[tag] || _need_vector_tag_grad[tag] || _need_vector_tag_dof_u[tag])
1104  if ((_subproblem.vectorTagType(tag) == Moose::VECTOR_TAG_RESIDUAL &&
1105  _subproblem.safeAccessTaggedVectors()) ||
1106  _subproblem.vectorTagType(tag) == Moose::VECTOR_TAG_SOLUTION)
1107  // tag is defined on problem but may not be used by a system
1108  if (_sys.hasVector(tag) && _sys.getVector(tag).closed())
1109  {
1110  auto & vec = _sys.getVector(tag);
1111  _vector_tags_dof_u[tag].resize(1);
1112  _vector_tags_dof_u[tag][0] = vec(_dof_indices[0]);
1113  }
1114 
1115  if (_need_vector_tag_u[tag])
1116  {
1117  _vector_tag_u[tag].resize(nqp);
1118  auto v = phi * _vector_tags_dof_u[tag][0];
1119  for (unsigned int qp = 0; qp < nqp; ++qp)
1120  _vector_tag_u[tag][qp] = v;
1121  }
1122  if (_need_vector_tag_grad[tag])
1123  _vector_tag_grad[tag].resize(nqp);
1124  }
1125 
1126  if (_subproblem.safeAccessTaggedMatrices())
1127  {
1128  auto & active_coupleable_matrix_tags =
1129  _subproblem.getActiveFEVariableCoupleableMatrixTags(_tid);
1130  for (auto tag : active_coupleable_matrix_tags)
1131  {
1132  _matrix_tags_dof_u[tag].resize(1);
1133  if (_need_matrix_tag_dof_u[tag] || _need_matrix_tag_u[tag])
1134  if (_sys.hasMatrix(tag) && _sys.matrixTagActive(tag) && _sys.getMatrix(tag).closed())
1135  {
1136  auto & mat = _sys.getMatrix(tag);
1137  {
1138  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1139  _matrix_tags_dof_u[tag][0] = mat(_dof_indices[0], _dof_indices[0]);
1140  }
1141  }
1142  }
1143  }
1144  for (auto tag : active_coupleable_matrix_tags)
1145  if (_need_matrix_tag_u[tag])
1146  {
1147  _matrix_tag_u[tag].resize(nqp);
1148  auto v = phi * _matrix_tags_dof_u[tag][0];
1149  for (unsigned int qp = 0; qp < nqp; ++qp)
1150  _matrix_tag_u[tag][qp] = v;
1151  }
1152 }
1153 
1154 template <>
1155 void
1157 {
1158  // Fixeme: will think of optimization later
1159  computeValues();
1160 }
1161 
1162 template <typename OutputType>
1163 void
1164 MooseVariableData<OutputType>::computeAD(const unsigned int num_dofs, const unsigned int nqp)
1165 {
1166  // Have to do this because upon construction this won't initialize any of the derivatives
1167  // (because DualNumber::do_derivatives is false at that time).
1168  _ad_zero = 0;
1169 
1170  _ad_dof_values.resize(num_dofs);
1171  if (_need_ad_u)
1172  _ad_u.resize(nqp);
1173 
1174  if (_need_ad_grad_u)
1175  _ad_grad_u.resize(nqp);
1176 
1177  if (_need_ad_second_u)
1178  _ad_second_u.resize(nqp);
1179 
1180  if (_need_ad_u_dot)
1181  {
1182  _ad_dofs_dot.resize(num_dofs);
1183  _ad_u_dot.resize(nqp);
1184  }
1185  if (_need_ad_grad_u_dot)
1186  _ad_grad_u_dot.resize(nqp);
1187 
1188  if (_need_ad_u_dotdot)
1189  {
1190  _ad_dofs_dotdot.resize(num_dofs);
1191  _ad_u_dotdot.resize(nqp);
1192  }
1193 
1194  const bool do_derivatives = Moose::doDerivatives(_subproblem, _sys);
1195 
1196  for (unsigned int qp = 0; qp < nqp; qp++)
1197  {
1198  if (_need_ad_u)
1199  _ad_u[qp] = _ad_zero;
1200 
1201  if (_need_ad_grad_u)
1202  _ad_grad_u[qp] = _ad_zero;
1203 
1204  if (_need_ad_second_u)
1205  _ad_second_u[qp] = _ad_zero;
1206 
1207  if (_need_ad_u_dot)
1208  _ad_u_dot[qp] = _ad_zero;
1209 
1210  if (_need_ad_u_dotdot)
1211  _ad_u_dotdot[qp] = _ad_zero;
1212 
1213  if (_need_ad_grad_u_dot)
1214  _ad_grad_u_dot[qp] = _ad_zero;
1215  }
1216 
1217  for (unsigned int i = 0; i < num_dofs; i++)
1218  {
1219  _ad_dof_values[i] = (*_sys.currentSolution())(_dof_indices[i]);
1220 
1221  // NOTE! You have to do this AFTER setting the value!
1222  if (do_derivatives)
1223  Moose::derivInsert(_ad_dof_values[i].derivatives(), _dof_indices[i], 1.);
1224 
1225  if (_need_ad_u_dot && _time_integrator && _time_integrator->dt())
1226  {
1227  _ad_dofs_dot[i] = _ad_dof_values[i];
1228  _time_integrator->computeADTimeDerivatives(_ad_dofs_dot[i],
1229  _dof_indices[i],
1230  _need_ad_u_dotdot ? _ad_dofs_dotdot[i]
1231  : _ad_real_dummy);
1232  }
1233  }
1234 
1235  // Now build up the solution at each quadrature point:
1236  for (unsigned int i = 0; i < num_dofs; i++)
1237  {
1238  for (unsigned int qp = 0; qp < nqp; qp++)
1239  {
1240  if (_need_ad_u)
1241  _ad_u[qp] += _ad_dof_values[i] * (*_current_phi)[i][qp];
1242 
1243  if (_need_ad_grad_u)
1244  {
1245  // The latter check here is for handling the fact that we have not yet implemented
1246  // calculation of ad_grad_phi for neighbor and neighbor-face, so if we are in that
1247  // situation we need to default to using the non-ad grad_phi
1248  if (_displaced && _current_ad_grad_phi)
1249  _ad_grad_u[qp] += _ad_dof_values[i] * (*_current_ad_grad_phi)[i][qp];
1250  else
1251  _ad_grad_u[qp] += _ad_dof_values[i] * (*_current_grad_phi)[i][qp];
1252  }
1253 
1254  if (_need_ad_second_u)
1255  // Note that this will not carry any derivatives with respect to displacements because
1256  // those calculations have not yet been implemented in Assembly
1257  _ad_second_u[qp] += _ad_dof_values[i] * (*_current_second_phi)[i][qp];
1258 
1259  if (_need_ad_u_dot && _time_integrator && _time_integrator->dt())
1260  {
1261  _ad_u_dot[qp] += (*_current_phi)[i][qp] * _ad_dofs_dot[i];
1262  if (_need_ad_u_dotdot)
1263  _ad_u_dotdot[qp] += (*_current_phi)[i][qp] * _ad_dofs_dotdot[i];
1264  }
1265 
1266  if (_need_ad_grad_u_dot && _time_integrator && _time_integrator->dt())
1267  {
1268  // The latter check here is for handling the fact that we have not yet implemented
1269  // calculation of ad_grad_phi for neighbor and neighbor-face, so if we are in that
1270  // situation we need to default to using the non-ad grad_phi
1271  if (_displaced && _current_ad_grad_phi)
1272  _ad_grad_u_dot[qp] += _ad_dofs_dot[i] * (*_current_ad_grad_phi)[i][qp];
1273  else
1274  _ad_grad_u_dot[qp] += _ad_dofs_dot[i] * (*_current_grad_phi)[i][qp];
1275  }
1276  }
1277  }
1278 
1279  if (_need_ad_u_dot && !_time_integrator)
1280  for (MooseIndex(nqp) qp = 0; qp < nqp; ++qp)
1281  {
1282  _ad_u_dot[qp] = _u_dot[qp];
1283  if (_need_ad_u_dotdot)
1284  _ad_u_dotdot[qp] = _u_dotdot[qp];
1285  }
1286 
1287  if (_need_ad_grad_u_dot && !_time_integrator)
1288  for (MooseIndex(nqp) qp = 0; qp < nqp; ++qp)
1289  _ad_grad_u_dot[qp] = _grad_u_dot[qp];
1290 }
1291 
1292 template <>
1293 void
1294 MooseVariableData<RealEigenVector>::computeAD(const unsigned int /*num_dofs*/,
1295  const unsigned int /*nqp*/)
1296 {
1297  mooseError("AD for array variable has not been implemented");
1298 }
1299 
1300 template <typename OutputType>
1301 void
1302 MooseVariableData<OutputType>::setDofValue(const OutputData & value, unsigned int index)
1303 {
1304  auto & dof_values = _vector_tags_dof_u[_solution_tag];
1305  dof_values[index] = value;
1306  _has_dof_values = true;
1307 
1308  auto & u = _vector_tag_u[_solution_tag];
1309  for (unsigned int qp = 0; qp < u.size(); qp++)
1310  {
1311  u[qp] = (*_phi)[0][qp] * dof_values[0];
1312 
1313  for (unsigned int i = 1; i < dof_values.size(); i++)
1314  u[qp] += (*_phi)[i][qp] * dof_values[i];
1315  }
1316 }
1317 
1318 template <typename OutputType>
1319 void
1320 MooseVariableData<OutputType>::setDofValues(const DenseVector<OutputData> & values)
1321 {
1322  auto & dof_values = _vector_tags_dof_u[_solution_tag];
1323  for (unsigned int i = 0; i < values.size(); i++)
1324  dof_values[i] = values(i);
1325 
1326  _has_dof_values = true;
1327 
1328  auto & u = _vector_tag_u[_solution_tag];
1329  for (unsigned int qp = 0; qp < u.size(); qp++)
1330  {
1331  u[qp] = (*_phi)[0][qp] * dof_values[0];
1332  for (unsigned int i = 1; i < dof_values.size(); i++)
1333  u[qp] += (*_phi)[i][qp] * dof_values[i];
1334  }
1335 }
1336 
1337 template <typename OutputType>
1338 void
1340  const OutputData & v)
1341 {
1342  residual.set(_nodal_dof_index, v);
1343 }
1344 
1345 template <>
1346 void
1348  const RealEigenVector & v)
1349 {
1350  for (unsigned int j = 0; j < _count; ++j)
1351  residual.set(_nodal_dof_index + j, v(j));
1352 }
1353 
1354 template <typename OutputType>
1357 {
1358  mooseAssert(_subproblem.mesh().isSemiLocal(const_cast<Node *>(&node)), "Node is not Semilocal");
1359 
1360  // Make sure that the node has DOFs
1361  /* Note, this is a reproduction of an assert within libMesh::Node::dof_number, this is done to
1362  * produce a better error (see misc/check_error.node_value_off_block) */
1363  mooseAssert(node.n_dofs(_sys.number(), _var_num) > 0,
1364  "Node " << node.id() << " does not contain any dofs for the "
1365  << _sys.system().variable_name(_var_num) << " variable");
1366 
1367  dof_id_type dof = node.dof_number(_sys.number(), _var_num, 0);
1368 
1369  switch (state)
1370  {
1371  case Moose::Current:
1372  return (*_sys.currentSolution())(dof);
1373 
1374  case Moose::Old:
1375  return _sys.solutionOld()(dof);
1376 
1377  case Moose::Older:
1378  return _sys.solutionOlder()(dof);
1379 
1380  default:
1381  mooseError("PreviousNL not currently supported for getNodalValue");
1382  }
1383 }
1384 
1385 template <>
1388  Moose::SolutionState state) const
1389 {
1390  mooseAssert(_subproblem.mesh().isSemiLocal(const_cast<Node *>(&node)), "Node is not Semilocal");
1391 
1392  // Make sure that the node has DOFs
1393  /* Note, this is a reproduction of an assert within libMesh::Node::dof_number, this is done to
1394  * produce a better error (see misc/check_error.node_value_off_block) */
1395  mooseAssert(node.n_dofs(_sys.number(), _var_num) > 0,
1396  "Node " << node.id() << " does not contain any dofs for the "
1397  << _sys.system().variable_name(_var_num) << " variable");
1398 
1399  dof_id_type dof = node.dof_number(_sys.number(), _var_num, 0);
1400 
1401  RealEigenVector v(_count);
1402  switch (state)
1403  {
1404  case Moose::Current:
1405  for (unsigned int i = 0; i < _count; ++i)
1406  v(i) = (*_sys.currentSolution())(dof++);
1407  break;
1408 
1409  case Moose::Old:
1410  for (unsigned int i = 0; i < _count; ++i)
1411  v(i) = _sys.solutionOld()(dof++);
1412  break;
1413 
1414  case Moose::Older:
1415  for (unsigned int i = 0; i < _count; ++i)
1416  v(i) = _sys.solutionOlder()(dof++);
1417  break;
1418 
1419  default:
1420  mooseError("PreviousNL not currently supported for getNodalValue");
1421  }
1422  return v;
1423 }
1424 
1425 template <typename OutputType>
1428  Moose::SolutionState state,
1429  unsigned int idx) const
1430 {
1431  std::vector<dof_id_type> dof_indices;
1432  _dof_map.dof_indices(elem, dof_indices, _var_num);
1433 
1434  switch (state)
1435  {
1436  case Moose::Current:
1437  return (*_sys.currentSolution())(dof_indices[idx]);
1438 
1439  case Moose::Old:
1440  return _sys.solutionOld()(dof_indices[idx]);
1441 
1442  case Moose::Older:
1443  return _sys.solutionOlder()(dof_indices[idx]);
1444 
1445  default:
1446  mooseError("PreviousNL not currently supported for getElementalValue");
1447  }
1448 }
1449 
1450 template <>
1453  Moose::SolutionState state,
1454  unsigned int idx) const
1455 {
1456  std::vector<dof_id_type> dof_indices;
1457  _dof_map.dof_indices(elem, dof_indices, _var_num);
1458 
1459  dof_id_type dof = dof_indices[idx];
1460 
1461  RealEigenVector v(_count);
1462 
1463  switch (state)
1464  {
1465  case Moose::Current:
1466  for (unsigned int i = 0; i < _count; ++i)
1467  v(i) = (*_sys.currentSolution())(dof++);
1468  break;
1469 
1470  case Moose::Old:
1471  for (unsigned int i = 0; i < _count; ++i)
1472  v(i) = _sys.solutionOld()(dof++);
1473  break;
1474 
1475  case Moose::Older:
1476  for (unsigned int i = 0; i < _count; ++i)
1477  v(i) = _sys.solutionOlder()(dof++);
1478  break;
1479 
1480  default:
1481  mooseError("PreviousNL not currently supported for getElementalValue");
1482  }
1483  return v;
1484 }
1485 
1486 template <typename OutputType>
1487 void
1489  std::vector<dof_id_type> & dof_indices) const
1490 {
1491  _dof_map.dof_indices(elem, dof_indices, _var_num);
1492 }
1493 
1494 template <typename OutputType>
1495 void
1497  const DenseVector<Number> & v) const
1498 {
1499  sol.add_vector(v, _dof_indices);
1500 }
1501 
1502 template <>
1503 void
1505  const DenseVector<Number> & v) const
1506 {
1507  unsigned int p = 0;
1508  for (unsigned int j = 0; j < _count; ++j)
1509  {
1510  unsigned int inc = (isNodal() ? j : j * _dof_indices.size());
1511  for (unsigned int i = 0; i < _dof_indices.size(); ++i)
1512  sol.add(_dof_indices[i] + inc, v(p++));
1513  }
1514 }
1515 
1516 template <typename OutputType>
1519 {
1520  if (_sys.solutionUDot())
1521  {
1522  _need_dof_values_dot = true;
1523  return _dof_values_dot;
1524  }
1525  else
1526  mooseError("MooseVariableData: Time derivative of solution (`u_dot`) is not stored. Please set "
1527  "uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
1528 }
1529 
1530 template <typename OutputType>
1533 {
1534  if (_sys.solutionUDotDot())
1535  {
1536  _need_dof_values_dotdot = true;
1537  return _dof_values_dotdot;
1538  }
1539  else
1540  mooseError("MooseVariableData: Second time derivative of solution (`u_dotdot`) is not stored. "
1541  "Please set uDotDotRequested() to true in FEProblemBase before requesting "
1542  "`u_dotdot`.");
1543 }
1544 
1545 template <typename OutputType>
1548 {
1549  if (_sys.solutionUDotOld())
1550  {
1551  _need_dof_values_dot_old = true;
1552  return _dof_values_dot_old;
1553  }
1554  else
1555  mooseError("MooseVariableData: Old time derivative of solution (`u_dot_old`) is not stored. "
1556  "Please set uDotOldRequested() to true in FEProblemBase before requesting "
1557  "`u_dot_old`.");
1558 }
1559 
1560 template <typename OutputType>
1563 {
1564  if (_sys.solutionUDotDotOld())
1565  {
1566  _need_dof_values_dotdot_old = true;
1567  return _dof_values_dotdot_old;
1568  }
1569  else
1570  mooseError("MooseVariableData: Old second time derivative of solution (`u_dotdot_old`) is not "
1571  "stored. Please set uDotDotOldRequested() to true in FEProblemBase before "
1572  "requesting `u_dotdot_old`.");
1573 }
1574 
1575 template <typename OutputType>
1576 const MooseArray<Number> &
1578 {
1579  _need_dof_du_dot_du = true;
1580  return _dof_du_dot_du;
1581 }
1582 
1583 template <typename OutputType>
1584 const MooseArray<Number> &
1586 {
1587  _need_dof_du_dotdot_du = true;
1588  return _dof_du_dotdot_du;
1589 }
1590 
1591 template <typename OutputType>
1592 void
1594 {
1595  unsigned int nqp = _qrule->n_points();
1596 
1597  _increment.resize(nqp);
1598  // Compute the increment at each quadrature point
1599  unsigned int num_dofs = _dof_indices.size();
1600  for (unsigned int qp = 0; qp < nqp; qp++)
1601  {
1602  _increment[qp] = 0;
1603  for (unsigned int i = 0; i < num_dofs; i++)
1604  _increment[qp] += (*_phi)[i][qp] * increment_vec(_dof_indices[i]);
1605  }
1606 }
1607 
1608 template <>
1609 void
1611  const NumericVector<Number> & increment_vec)
1612 {
1613  unsigned int nqp = _qrule->n_points();
1614 
1615  _increment.resize(nqp);
1616  // Compute the increment at each quadrature point
1617  unsigned int num_dofs = _dof_indices.size();
1618  if (isNodal())
1619  {
1620  for (unsigned int qp = 0; qp < nqp; qp++)
1621  {
1622  for (unsigned int i = 0; i < num_dofs; i++)
1623  for (unsigned int j = 0; j < _count; j++)
1624  _increment[qp](j) += (*_phi)[i][qp] * increment_vec(_dof_indices[i] + j);
1625  }
1626  }
1627  else
1628  {
1629  for (unsigned int qp = 0; qp < nqp; qp++)
1630  {
1631  unsigned int n = 0;
1632  for (unsigned int j = 0; j < _count; j++)
1633  for (unsigned int i = 0; i < num_dofs; i++)
1634  {
1635  _increment[qp](j) += (*_phi)[i][qp] * increment_vec(_dof_indices[i] + n);
1636  n += num_dofs;
1637  }
1638  }
1639  }
1640 }
1641 
1642 template <typename OutputType>
1643 void
1645 {
1646  if (!isNodal())
1647  mooseError("computeIncrementAtNode can only be called for nodal variables");
1648 
1649  _increment.resize(1);
1650 
1651  // Compute the increment for the current DOF
1652  _increment[0] = increment_vec(_dof_indices[0]);
1653 }
1654 
1655 template <>
1656 void
1658  const NumericVector<Number> & increment_vec)
1659 {
1660  if (!isNodal())
1661  mooseError("computeIncrementAtNode can only be called for nodal variables");
1662 
1663  _increment.resize(1);
1664 
1665  // Compute the increment for the current DOF
1666  if (isNodal())
1667  for (unsigned int j = 0; j < _count; j++)
1668  _increment[0](j) = increment_vec(_dof_indices[0] + j);
1669  else
1670  {
1671  unsigned int n = 0;
1672  for (unsigned int j = 0; j < _count; j++)
1673  {
1674  _increment[0](j) = increment_vec(_dof_indices[0] + n);
1675  n += _dof_indices.size();
1676  }
1677  }
1678 }
1679 
1680 template <typename OutputType>
1681 const OutputType &
1683 {
1684  if (isNodal())
1685  {
1686  if (_sys.solutionUDot())
1687  {
1688  _need_dof_values_dot = true;
1689  return _nodal_value_dot;
1690  }
1691  else
1692  mooseError(
1693  "MooseVariableData: Time derivative of solution (`u_dot`) is not stored. Please set "
1694  "uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
1695  }
1696  else
1697  mooseError("Nodal values can be requested only on nodal variables, variable '",
1698  var().name(),
1699  "' is not nodal.");
1700 }
1701 
1702 template <typename OutputType>
1703 const OutputType &
1705 {
1706  if (isNodal())
1707  {
1708  if (_sys.solutionUDotDot())
1709  {
1710  _need_dof_values_dotdot = true;
1711  return _nodal_value_dotdot;
1712  }
1713  else
1714  mooseError(
1715  "MooseVariableData: Second time derivative of solution (`u_dotdot`) is not stored. "
1716  "Please set uDotDotRequested() to true in FEProblemBase before requesting "
1717  "`u_dotdot`.");
1718  }
1719  else
1720  mooseError("Nodal values can be requested only on nodal variables, variable '",
1721  var().name(),
1722  "' is not nodal.");
1723 }
1724 
1725 template <typename OutputType>
1726 const OutputType &
1728 {
1729  if (isNodal())
1730  {
1731  if (_sys.solutionUDotOld())
1732  {
1733  _need_dof_values_dot_old = true;
1734  return _nodal_value_dot_old;
1735  }
1736  else
1737  mooseError("MooseVariableData: Old time derivative of solution (`u_dot_old`) is not stored. "
1738  "Please set uDotOldRequested() to true in FEProblemBase before requesting "
1739  "`u_dot_old`.");
1740  }
1741  else
1742  mooseError("Nodal values can be requested only on nodal variables, variable '",
1743  var().name(),
1744  "' is not nodal.");
1745 }
1746 
1747 template <typename OutputType>
1748 const OutputType &
1750 {
1751  if (isNodal())
1752  {
1753  if (_sys.solutionUDotDotOld())
1754  {
1755  _need_dof_values_dotdot_old = true;
1756  return _nodal_value_dotdot_old;
1757  }
1758  else
1759  mooseError(
1760  "MooseVariableData: Old second time derivative of solution (`u_dotdot_old`) is not "
1761  "stored. Please set uDotDotOldRequested() to true in FEProblemBase before "
1762  "requesting `u_dotdot_old`.");
1763  }
1764  else
1765  mooseError("Nodal values can be requested only on nodal variables, variable '",
1766  var().name(),
1767  "' is not nodal.");
1768 }
1769 
1770 template <typename OutputType>
1771 void
1773 {
1774  if (_has_dof_indices)
1775  {
1776  fetchDoFValues();
1777  assignNodalValue();
1778 
1779  if (_need_ad)
1780  fetchADDoFValues();
1781  }
1782  else
1783  zeroSizeDofValues();
1784 }
1785 
1786 template <typename OutputType>
1787 void
1789 {
1790  auto n = _dof_indices.size();
1791  libmesh_assert(n);
1792  _ad_dof_values.resize(n);
1793 
1794  const bool do_derivatives =
1795  ADReal::do_derivatives && _sys.number() == _subproblem.currentNlSysNum();
1796 
1797  for (decltype(n) i = 0; i < n; ++i)
1798  {
1799  _ad_dof_values[i] = _vector_tags_dof_u[_solution_tag][i];
1800  if (do_derivatives)
1801  Moose::derivInsert(_ad_dof_values[i].derivatives(), _dof_indices[i], 1.);
1802  assignADNodalValue(_ad_dof_values[i], i);
1803  }
1804 }
1805 
1806 template <>
1807 void
1809 {
1810  mooseError("I do not know how to support AD with array variables");
1811 }
1812 
1813 template <>
1814 void
1815 MooseVariableData<Real>::assignADNodalValue(const DualReal & value, const unsigned int &)
1816 {
1817  _ad_nodal_value = value;
1818 }
1819 
1820 template <>
1821 void
1823  const unsigned int & component)
1824 {
1825  _ad_nodal_value(component) = value;
1826 }
1827 
1828 template <typename OutputType>
1829 void
1831 {
1832  _dof_map.dof_indices(_elem, _dof_indices, _var_num);
1833  _vector_tags_dof_u[_solution_tag].resize(_dof_indices.size());
1834 
1835  unsigned int nqp = _qrule->n_points();
1836  _vector_tag_u[_solution_tag].resize(nqp);
1837 }
1838 
1839 template <typename OutputType>
1840 void
1842 {
1843  _dof_map.dof_indices(_elem, _dof_indices, _var_num);
1844  _has_dof_values = false;
1845 
1846  // FIXME: remove this when the Richard's module is migrated to use the new NodalCoupleable
1847  // interface.
1848  if (_dof_indices.size() > 0)
1849  _has_dof_indices = true;
1850  else
1851  _has_dof_indices = false;
1852 }
1853 
1854 template <typename OutputType>
1855 void
1857 {
1858  if (std::size_t n_dofs = _node->n_dofs(_sys.number(), _var_num))
1859  {
1860  _dof_indices.resize(n_dofs);
1861  for (std::size_t i = 0; i < n_dofs; ++i)
1862  _dof_indices[i] = _node->dof_number(_sys.number(), _var_num, i);
1863  // For standard variables. _nodal_dof_index is retrieved by nodalDofIndex() which is used in
1864  // NodalBC for example
1865  _nodal_dof_index = _dof_indices[0];
1866  _has_dof_indices = true;
1867  }
1868  else
1869  {
1870  _dof_indices.clear(); // Clear these so Assembly doesn't think there's dofs here
1871  _has_dof_indices = false;
1872  }
1873 }
1874 
1875 template <typename OutputType>
1876 void
1878 {
1879  /* FIXME: this method is only for elemental auxiliary variables, so
1880  * we may want to rename it */
1881  if (_elem)
1882  {
1883  _dof_map.dof_indices(_elem, _dof_indices, _var_num);
1884  if (_elem->n_dofs(_sys.number(), _var_num) > 0)
1885  {
1886  // FIXME: check if the following is equivalent with '_nodal_dof_index = _dof_indices[0];'?
1887  _nodal_dof_index = _elem->dof_number(_sys.number(), _var_num, 0);
1888 
1889  fetchDoFValues();
1890 
1891  for (auto & dof_u : _vector_tags_dof_u)
1892  dof_u.resize(_dof_indices.size());
1893 
1894  for (auto & dof_u : _matrix_tags_dof_u)
1895  dof_u.resize(_dof_indices.size());
1896 
1897  _has_dof_indices = true;
1898  }
1899  else
1900  _has_dof_indices = false;
1901  }
1902  else
1903  _has_dof_indices = false;
1904 }
1905 
1906 template <typename OutputType>
1907 void
1908 MooseVariableData<OutputType>::reinitNodes(const std::vector<dof_id_type> & nodes)
1909 {
1910  _dof_indices.clear();
1911  for (const auto & node_id : nodes)
1912  {
1913  auto && nd = _subproblem.mesh().getMesh().query_node_ptr(node_id);
1914  if (nd && (_subproblem.mesh().isSemiLocal(const_cast<Node *>(nd))))
1915  {
1916  if (nd->n_dofs(_sys.number(), _var_num) > 0)
1917  {
1918  dof_id_type dof = nd->dof_number(_sys.number(), _var_num, 0);
1919  _dof_indices.push_back(dof);
1920  }
1921  }
1922  }
1923 
1924  if (_dof_indices.size() > 0)
1925  _has_dof_indices = true;
1926  else
1927  _has_dof_indices = false;
1928 }
1929 
1930 template class MooseVariableData<Real>;
1931 template class MooseVariableData<RealVectorValue>;
1932 template class MooseVariableData<RealEigenVector>;
const FieldVariableDivergence & divSln(Moose::SolutionState state) const
Local solution divergence getter.
std::string name(const ElemQuality q)
const Assembly & _assembly
const FEVectorBase *const & getVectorFE(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current volume FEVector.
Definition: Assembly.h:164
std::function< const typename OutputTools< OutputShape >::VariablePhiSecond &(const Assembly &, FEType)> _second_phi_face_assembly_method
OutputData getNodalValue(const Node &node, Moose::SolutionState state) const
void assignADNodalValue(const DualReal &value, const unsigned int &component)
Helper methods for assigning nodal values from their corresponding solution values (dof values as the...
const FieldVariableValue & uDotDotOld() const
void computeAD(const unsigned int num_dofs, const unsigned int nqp)
compute AD things
void insertNodalValue(NumericVector< Number > &residual, const OutputData &v)
Write a nodal value to the passed-in solution vector.
std::function< const typename OutputTools< OutputShape >::VariablePhiDivergence &(const Assembly &, FEType)> _div_phi_face_assembly_method
const FieldVariablePhiCurl & curlPhi() const
curl_phi getter
void reinitNodes(const std::vector< dof_id_type > &nodes)
Set _dof_indices to the degrees of freedom existing on the passed-in nodes.
LAGRANGE_VEC
void prepare()
Get the dof indices corresponding to the current element.
TensorTools::IncrementRank< OutputGradient >::type OutputSecond
Definition: MooseTypes.h:269
Class for stuff related to variables.
Definition: Coupleable.h:37
void setDofValue(const OutputData &value, unsigned int index)
dof value setters
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:284
std::function< const ADTemplateVariablePhiGradient< OutputShape > &(const Assembly &, FEType)> _ad_grad_phi_assembly_method
void computeValues()
compute the variable values
void computeIncrementAtQps(const NumericVector< Number > &increment_vec)
Compute and store incremental change in solution at QPs based on increment_vec.
virtual void add_vector(const Number *v, const std::vector< numeric_index_type > &dof_indices)
const FieldVariableValue & uDotDot() const
void getDofIndices(const Elem *elem, std::vector< dof_id_type > &dof_indices) const
DualNumber< Real, DNDerivativeType, true > DualReal
Definition: DualReal.h:49
Moose::DOFType< OutputType >::type OutputData
TimeIntegrator * getTimeIntegrator()
Definition: SystemBase.h:869
const OutputType & nodalValueDotOld() const
TensorTools::IncrementRank< OutputShape >::type OutputShapeGradient
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
Definition: Moose.h:148
const FieldVariableValue & uDot() const
std::function< const typename OutputTools< OutputShape >::VariablePhiGradient &(const Assembly &, FEType)> _grad_phi_assembly_method
FEContinuity _continuity
Continuity type of the variable.
void computeMonomialValues()
compute the values for const monomial variables
MooseVariableData(const MooseVariableField< OutputType > &var, SystemBase &sys, THREAD_ID tid, Moose::ElementType element_type, const QBase *const &qrule_in, const QBase *const &qrule_face_in, const Node *const &node, const Elem *const &elem)
void computeNodalValues()
compute nodal things
const MooseArray< Number > & dofValuesDuDotDotDu() const
ElementType
Definition: MooseTypes.h:676
Base class for a system (of equations)
Definition: SystemBase.h:84
void prepareIC()
prepare the initial condition
MONOMIAL_VEC
const FieldVariableGradient & gradSlnDot() const
Local time derivative of solution gradient getter.
bool doDerivatives(const SubProblem &subproblem, const SystemBase &sys)
Definition: ADUtils.C:79
void setDofValues(const DenseVector< OutputData > &values)
Set local DOF values and evaluate the values on quadrature points.
OutputData getElementalValue(const Elem *elem, Moose::SolutionState state, unsigned int idx=0) const
const DoFValue & dofValuesDotDotOld() const
std::function< const typename OutputTools< OutputShape >::VariablePhiDivergence &(const Assembly &, FEType)> _div_phi_assembly_method
const MooseArray< Number > & dofValuesDuDotDu() const
const DoFValue & dofValuesDot() const
const FieldVariablePhiDivergence & divPhi() const
divergence_phi getter
std::function< const ADTemplateVariablePhiGradient< OutputShape > &(const Assembly &, FEType)> _ad_grad_phi_face_assembly_method
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
void addSolution(NumericVector< Number > &sol, const DenseVector< Number > &v) const
Add passed in local DOF values to a solution vector.
const OutputType & nodalValueDot() const
std::function< const typename OutputTools< OutputShape >::VariablePhiCurl &(const Assembly &, FEType)> _curl_phi_face_assembly_method
virtual unsigned int dimension() const
Returns MeshBase::mesh_dimension(), (not MeshBase::spatial_dimension()!) of the underlying libMesh me...
Definition: MooseMesh.C:2678
const FieldVariablePhiSecond & secondPhiFace() const
second_phi_face getter
C_ZERO
TensorTools::DecrementRank< OutputShape >::type OutputShapeDivergence
const FieldVariableGradient & gradSlnDotDot() const
Local second time derivative of solution gradient getter.
TensorTools::IncrementRank< OutputType >::type OutputGradient
Definition: MooseTypes.h:268
libmesh_assert(ctx)
void reinitNode()
Prepare degrees of freedom for the current node.
const FEType & _fe_type
std::function< const typename OutputTools< OutputShape >::VariablePhiSecond &(const Assembly &, FEType)> _second_phi_assembly_method
const TimeIntegrator * _time_integrator
Pointer to time integrator.
std::function< const typename OutputTools< OutputShape >::VariablePhiValue &(const Assembly &, FEType)> _phi_face_assembly_method
const DoFValue & dofValuesDotDot() const
Moose::ElementType _element_type
The element type this object is storing data for. This is either Element, Neighbor, or Lower.
const ADTemplateVariablePhiGradient< OutputShape > * _ad_grad_phi
const OutputType & nodalValueDotDotOld() const
const FieldVariablePhiValue * _phi_face
RAVIART_THOMAS
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
GeometryType
Definition: MooseTypes.h:236
std::function< const typename OutputTools< OutputShape >::VariablePhiCurl &(const Assembly &, FEType)> _curl_phi_assembly_method
const FieldVariableCurl & curlSln(Moose::SolutionState state) const
Local solution curl getter.
C_ONE
const FieldVariableValue & uDotOld() const
SolutionState
Definition: MooseTypes.h:221
Moose::ShapeType< OutputType >::type OutputShape
const ADTemplateVariablePhiGradient< OutputShape > * _ad_grad_phi_face
const FieldVariablePhiSecond & secondPhi() const
second_phi getter
const OutputType & nodalValueDotDot() const
virtual MooseMesh & mesh()
Definition: SystemBase.h:96
SystemBase & _sys
The MOOSE system which ultimately holds the vectors and matrices relevant to this variable data...
IntRange< T > make_range(T beg, T end)
const FieldVariablePhiGradient * _grad_phi_face
void derivInsert(NumberArray< N, Real > &derivs, dof_id_type index, Real value)
const FieldVariablePhiCurl & curlPhiFace() const
curl_phi_face getter
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:138
std::function< const typename OutputTools< OutputType >::VariablePhiValue &(const Assembly &, FEType)> _phi_assembly_method
virtual void set(const numeric_index_type i, const Number value)=0
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
const DoFValue & dofValuesDotOld() const
NEDELEC_ONE
const FEBase *const & getFE(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current volume FE.
Definition: Assembly.h:116
std::function< const typename OutputTools< OutputShape >::VariablePhiGradient &(const Assembly &, FEType)> _grad_phi_face_assembly_method
void setGeometry(Moose::GeometryType gm_type)
Set the geometry type before calculating variables values.
void reinitAux()
Prepare dof indices and solution values for elemental auxiliary variables.
virtual void add(const numeric_index_type i, const Number value)=0
const FieldVariablePhiDivergence & divPhiFace() const
divergence_phi_face getter
const FieldVariableSecond & secondSln(Moose::SolutionState state) const
Local solution second spatial derivative getter.
void computeIncrementAtNode(const NumericVector< Number > &increment_vec)
Compute and store incremental change at the current node based on increment_vec.
bool _is_nodal
if variable is nodal
const FieldVariablePhiGradient * _grad_phi
const FieldVariablePhiValue * _phi
unsigned int THREAD_ID
Definition: MooseTypes.h:198
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)
uint8_t dof_id_type