www.mooseframework.org
MooseVariable.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "MooseVariable.h"
16 #include "SubProblem.h"
17 #include "SystemBase.h"
18 #include "Assembly.h"
19 #include "MooseMesh.h"
20 
21 // libMesh
22 #include "libmesh/numeric_vector.h"
23 #include "libmesh/dof_map.h"
24 #include "libmesh/elem.h"
25 #include "libmesh/quadrature.h"
26 #include "libmesh/dense_vector.h"
27 
28 MooseVariable::MooseVariable(unsigned int var_num,
29  const FEType & fe_type,
30  SystemBase & sys,
31  Assembly & assembly,
32  Moose::VarKindType var_kind)
33  : MooseVariableBase(var_num, fe_type, sys, assembly, var_kind),
34 
35  _qrule(_assembly.qRule()),
36  _qrule_face(_assembly.qRuleFace()),
37  _qrule_neighbor(_assembly.qRuleNeighbor()),
38  _elem(_assembly.elem()),
39  _current_side(_assembly.side()),
40  _neighbor(_assembly.neighbor()),
41 
42  _need_u_old(false),
43  _need_u_older(false),
44  _need_u_previous_nl(false),
45 
46  _need_grad_old(false),
47  _need_grad_older(false),
48  _need_grad_previous_nl(false),
49 
50  _need_second(false),
51  _need_second_old(false),
52  _need_second_older(false),
53  _need_second_previous_nl(false),
54 
55  _need_u_old_neighbor(false),
56  _need_u_older_neighbor(false),
57  _need_u_previous_nl_neighbor(false),
58 
59  _need_grad_old_neighbor(false),
60  _need_grad_older_neighbor(false),
61  _need_grad_previous_nl_neighbor(false),
62 
63  _need_second_neighbor(false),
64  _need_second_old_neighbor(false),
65  _need_second_older_neighbor(false),
66  _need_second_previous_nl_neighbor(false),
67 
68  _need_nodal_u(false),
69  _need_nodal_u_old(false),
70  _need_nodal_u_older(false),
71  _need_nodal_u_previous_nl(false),
72  _need_nodal_u_dot(false),
73 
74  _need_nodal_u_neighbor(false),
75  _need_nodal_u_old_neighbor(false),
76  _need_nodal_u_older_neighbor(false),
77  _need_nodal_u_previous_nl_neighbor(false),
78  _need_nodal_u_dot_neighbor(false),
79 
80  _need_solution_dofs(false),
81  _need_solution_dofs_old(false),
82  _need_solution_dofs_older(false),
83  _need_solution_dofs_neighbor(false),
84  _need_solution_dofs_old_neighbor(false),
85  _need_solution_dofs_older_neighbor(false),
86 
87  _phi(_assembly.fePhi(_fe_type)),
88  _grad_phi(_assembly.feGradPhi(_fe_type)),
89 
90  _phi_face(_assembly.fePhiFace(_fe_type)),
91  _grad_phi_face(_assembly.feGradPhiFace(_fe_type)),
92 
93  _phi_neighbor(_assembly.fePhiNeighbor(_fe_type)),
94  _grad_phi_neighbor(_assembly.feGradPhiNeighbor(_fe_type)),
95 
96  _phi_face_neighbor(_assembly.fePhiFaceNeighbor(_fe_type)),
97  _grad_phi_face_neighbor(_assembly.feGradPhiFaceNeighbor(_fe_type)),
98 
99  _normals(_assembly.normals()),
100 
101  _is_defined(false),
102  _has_nodal_value(false),
103  _has_nodal_value_neighbor(false),
104 
105  _node(_assembly.node()),
106  _is_defined_neighbor(false),
107  _node_neighbor(_assembly.nodeNeighbor())
108 {
110 
111  // FIXME: continuity of FE type seems equivalent with the definition of nodal variables.
112  // Continuity does not depend on the FE dimension, so we just pass in a valid dimension.
113  _is_nodal = _assembly.getFE(feType(), _sys.mesh().dimension())->get_continuity() != DISCONTINUOUS;
114 }
115 
117 {
118  _u.release();
119  _u_bak.release();
120  _u_old.release();
122  _u_older.release();
125 
126  _grad_u.release();
133 
134  _second_u.release();
141 
142  _u_dot.release();
146 
151 
152  _nodal_u.release();
158 
165 
167 
172 
177 
182 }
183 
184 const std::set<SubdomainID> &
186 {
187  return _sys.system().variable(_var_num).active_subdomains();
188 }
189 
190 bool
192 {
193  return _sys.system().variable(_var_num).active_on_subdomain(subdomain);
194 }
195 
196 bool
198 {
199  return _is_nodal;
200 }
201 
202 void
204 {
205  _dof_indices.clear();
206 }
207 
208 void
210 {
211  _dof_map.dof_indices(_elem, _dof_indices, _var_num);
212  _has_nodal_value = false;
214 
215  // FIXME: remove this when the Richard's module is migrated to use the new NodalCoupleable
216  // interface.
217  if (_dof_indices.size() > 0)
218  _is_defined = true;
219  else
220  _is_defined = false;
221 }
222 
223 void
225 {
227  _has_nodal_value = false;
229 }
230 
231 void
233 {
234  _has_nodal_value = false;
236 }
237 
238 void
240 {
241  _dof_map.dof_indices(_elem, _dof_indices, _var_num);
242  _nodal_u.resize(_dof_indices.size());
243 
244  unsigned int nqp = _qrule->n_points();
245  _u.resize(nqp);
246 }
247 
248 void
250 {
251  if (_node->n_dofs(_sys.number(), _var_num) > 0)
252  {
253  _nodal_dof_index = _node->dof_number(_sys.number(), _var_num, 0);
254  _dof_indices.resize(1);
256  _is_defined = true;
257  }
258  else
259  {
260  _dof_indices.clear(); // Clear these so Assembly doesn't think there's dofs here
261  _is_defined = false;
262  }
263 }
264 
265 void
267 {
268  if (_node_neighbor->n_dofs(_sys.number(), _var_num) > 0)
269  {
271  _dof_indices_neighbor.resize(1);
273  _is_defined_neighbor = true;
274  }
275  else
276  {
277  _dof_indices_neighbor.clear(); // Clear these so Assembly doesn't think there's dofs here
278  _is_defined_neighbor = false;
279  }
280 }
281 
282 void
284 {
285  /* FIXME: this method is only for elemental auxiliary variables, so
286  * we may want to rename it */
287  _dof_map.dof_indices(_elem, _dof_indices, _var_num);
288  if (_elem->n_dofs(_sys.number(), _var_num) > 0)
289  {
290  // FIXME: check if the following is equivalent with '_nodal_dof_index = _dof_indices[0];'?
291  _nodal_dof_index = _elem->dof_number(_sys.number(), _var_num, 0);
292  libmesh_assert(_dof_indices.size());
293  _nodal_u.resize(_dof_indices.size());
295 
296  _is_defined = true;
297  }
298  else
299  _is_defined = false;
300 }
301 
302 void
304 {
305  if (_neighbor)
306  {
308  if (_neighbor->n_dofs(_sys.number(), _var_num) > 0)
309  {
311 
312  libmesh_assert(_dof_indices_neighbor.size());
315 
316  _is_defined_neighbor = true;
317  }
318  else
319  _is_defined_neighbor = false;
320  }
321  else
322  _is_defined_neighbor = false;
323 }
324 
325 void
326 MooseVariable::reinitNodes(const std::vector<dof_id_type> & nodes)
327 {
328  _dof_indices.clear();
329  for (const auto & node_id : nodes)
330  {
331  Node * nd = _subproblem.mesh().getMesh().query_node_ptr(node_id);
332  if (nd && (_subproblem.mesh().isSemiLocal(nd)))
333  {
334  if (nd->n_dofs(_sys.number(), _var_num) > 0)
335  {
336  dof_id_type dof = nd->dof_number(_sys.number(), _var_num, 0);
337  _dof_indices.push_back(dof);
338  }
339  }
340  }
341 
342  if (_dof_indices.size() > 0)
343  _is_defined = true;
344  else
345  _is_defined = false;
346 }
347 
348 void
349 MooseVariable::reinitNodesNeighbor(const std::vector<dof_id_type> & nodes)
350 {
351  _dof_indices_neighbor.clear();
352  for (const auto & node_id : nodes)
353  {
354  Node * nd = _subproblem.mesh().getMesh().query_node_ptr(node_id);
355  if (nd && (_subproblem.mesh().isSemiLocal(nd)))
356  {
357  if (nd->n_dofs(_sys.number(), _var_num) > 0)
358  {
359  dof_id_type dof = nd->dof_number(_sys.number(), _var_num, 0);
360  _dof_indices_neighbor.push_back(dof);
361  }
362  }
363  }
364 
365  if (_dof_indices_neighbor.size() > 0)
366  _is_defined_neighbor = true;
367  else
368  _is_defined_neighbor = false;
369 }
370 
371 void
372 MooseVariable::getDofIndices(const Elem * elem, std::vector<dof_id_type> & dof_indices)
373 {
374  _dof_map.dof_indices(elem, dof_indices, _var_num);
375 }
376 
377 void
378 MooseVariable::insert(NumericVector<Number> & residual)
379 {
380  if (_has_nodal_value)
381  residual.insert(&_nodal_u[0], _dof_indices);
382 
384  residual.insert(&_nodal_u_neighbor[0], _dof_indices_neighbor);
385 }
386 
387 void
388 MooseVariable::add(NumericVector<Number> & residual)
389 {
390  if (_has_nodal_value)
391  residual.add_vector(&_nodal_u[0], _dof_indices);
392 
394  residual.add_vector(&_nodal_u_neighbor[0], _dof_indices_neighbor);
395 }
396 
397 const VariablePhiValue &
399 {
400  return _phi;
401 }
402 
403 const VariablePhiGradient &
405 {
406  return _grad_phi;
407 }
408 
409 const VariablePhiSecond &
411 {
413  return *_second_phi;
414 }
415 
416 const VariablePhiValue &
418 {
419  return _phi_face;
420 }
421 
422 const VariablePhiGradient &
424 {
425  return _grad_phi_face;
426 }
427 
428 const VariablePhiSecond &
430 {
432  return *_second_phi_face;
433 }
434 
435 const VariablePhiValue &
437 {
438  return _phi_neighbor;
439 }
440 
441 const VariablePhiGradient &
443 {
444  return _grad_phi_neighbor;
445 }
446 
447 const VariablePhiSecond &
449 {
451  return *_second_phi_neighbor;
452 }
453 
454 const VariablePhiValue &
456 {
457  return _phi_face_neighbor;
458 }
459 
460 const VariablePhiGradient &
462 {
464 }
465 
466 const VariablePhiSecond &
468 {
471 }
472 
473 const VariableValue &
475 {
476  if (isNodal())
477  {
478  _need_nodal_u = true;
479  return _nodal_u;
480  }
481  else
482  mooseError("Nodal values can be requested only on nodal variables, variable '",
483  name(),
484  "' is not nodal.");
485 }
486 
487 const VariableValue &
489 {
490  if (isNodal())
491  {
492  _need_nodal_u_old = true;
493  return _nodal_u_old;
494  }
495  else
496  mooseError("Nodal values can be requested only on nodal variables, variable '",
497  name(),
498  "' is not nodal.");
499 }
500 
501 const VariableValue &
503 {
504  if (isNodal())
505  {
506  _need_nodal_u_older = true;
507  return _nodal_u_older;
508  }
509  else
510  mooseError("Nodal values can be requested only on nodal variables, variable '",
511  name(),
512  "' is not nodal.");
513 }
514 
515 const VariableValue &
517 {
518  if (isNodal())
519  {
521  return _nodal_u_previous_nl;
522  }
523  else
524  mooseError("Nodal values can be requested only on nodal variables, variable '",
525  name(),
526  "' is not nodal.");
527 }
528 
529 const VariableValue &
531 {
532  if (isNodal())
533  {
534  _need_nodal_u_dot = true;
535  return _nodal_u_dot;
536  }
537  else
538  mooseError("Nodal values can be requested only on nodal variables, variable '",
539  name(),
540  "' is not nodal.");
541 }
542 
543 const VariableValue &
545 {
546  if (isNodal())
547  {
548  _need_nodal_u_neighbor = true;
549  return _nodal_u_neighbor;
550  }
551  else
552  mooseError("Nodal values can be requested only on nodal variables, variable '",
553  name(),
554  "' is not nodal.");
555 }
556 
557 const VariableValue &
559 {
560  if (isNodal())
561  {
563  return _nodal_u_old_neighbor;
564  }
565  else
566  mooseError("Nodal values can be requested only on nodal variables, variable '",
567  name(),
568  "' is not nodal.");
569 }
570 
571 const VariableValue &
573 {
574  if (isNodal())
575  {
578  }
579  else
580  mooseError("Nodal values can be requested only on nodal variables, variable '",
581  name(),
582  "' is not nodal.");
583 }
584 
585 const VariableValue &
587 {
588  if (isNodal())
589  {
592  }
593  else
594  mooseError("Nodal values can be requested only on nodal variables, variable '",
595  name(),
596  "' is not nodal.");
597 }
598 
599 const VariableValue &
601 {
602  if (isNodal())
603  {
605  return _nodal_u_dot_neighbor;
606  }
607  else
608  mooseError("Nodal values can be requested only on nodal variables, variable '",
609  name(),
610  "' is not nodal.");
611 }
612 
613 // FIXME: this and computeElemeValues() could be refactored to reuse most of
614 // the common code, instead of duplicating it.
615 void
616 MooseVariable::computePerturbedElemValues(unsigned int perturbation_idx,
617  Real perturbation_scale,
618  Real & perturbation)
619 {
620 
621  bool is_transient = _subproblem.isTransient();
622  unsigned int nqp = _qrule->n_points();
623 
624  _u_bak = _u;
626  _u.resize(nqp);
627  _grad_u.resize(nqp);
628 
629  if (_need_second)
631  _second_u.resize(nqp);
632 
633  if (is_transient)
634  {
635  _u_dot_bak = _u_dot;
636  _u_dot.resize(nqp);
638  _du_dot_du.resize(nqp);
639 
640  if (_need_u_old)
641  _u_old_bak = _u_old;
642  _u_old.resize(nqp);
643 
644  if (_need_u_older)
646  _u_older.resize(nqp);
647 
648  if (_need_grad_old)
650  _grad_u_old.resize(nqp);
651 
652  if (_need_grad_older)
654  _grad_u_older.resize(nqp);
655 
656  if (_need_second_old)
658  _second_u_old.resize(nqp);
659 
660  if (_need_second_older)
662  _second_u_older.resize(nqp);
663  }
664 
665  for (unsigned int i = 0; i < nqp; ++i)
666  {
667  _u[i] = 0;
668  _grad_u[i] = 0;
669 
670  if (_need_second)
671  _second_u[i] = 0;
672 
673  if (is_transient)
674  {
675  _u_dot[i] = 0;
676  _du_dot_du[i] = 0;
677 
678  if (_need_u_old)
679  _u_old[i] = 0;
680 
681  if (_need_u_older)
682  _u_older[i] = 0;
683 
684  if (_need_grad_old)
685  _grad_u_old[i] = 0;
686 
687  if (_need_grad_older)
688  _grad_u_older[i] = 0;
689 
690  if (_need_second_old)
691  _second_u_old[i] = 0;
692 
693  if (_need_second_older)
694  _second_u_older[i] = 0;
695  }
696  }
697 
698  unsigned int num_dofs = _dof_indices.size();
699 
700  const NumericVector<Real> & current_solution = *_sys.currentSolution();
701  const NumericVector<Real> & solution_old = _sys.solutionOld();
702  const NumericVector<Real> & solution_older = _sys.solutionOlder();
703  const NumericVector<Real> & u_dot = _sys.solutionUDot();
704  const Real & du_dot_du = _sys.duDotDu();
705 
706  dof_id_type idx = 0;
707  Real soln_local = 0;
708  Real soln_old_local = 0;
709  Real soln_older_local = 0;
710  Real u_dot_local = 0;
711 
712  Real phi_local = 0;
713  const RealGradient * dphi_qp = NULL;
714  const RealTensor * d2phi_local = NULL;
715 
716  RealGradient * grad_u_qp = NULL;
717 
718  RealGradient * grad_u_old_qp = NULL;
719  RealGradient * grad_u_older_qp = NULL;
720 
721  RealTensor * second_u_qp = NULL;
722 
723  RealTensor * second_u_old_qp = NULL;
724  RealTensor * second_u_older_qp = NULL;
725 
726  for (unsigned int i = 0; i < num_dofs; i++)
727  {
728  idx = _dof_indices[i];
729  soln_local = current_solution(idx);
730  if (i == perturbation_idx)
731  {
732  // Compute the size of the perturbation.
733  // For the PETSc DS differencing method we use the magnitude of the variable at the "node"
734  // to determine the differencing parameters. The WP method could use the element L2 norm of
735  // the variable instead.
736  perturbation = soln_local;
737  // HACK: the use of fabs() and < assume Real is double or similar. Otherwise need to use
738  // PetscAbsScalar, PetscRealPart, etc.
739  if (fabs(perturbation) < 1.0e-16)
740  perturbation = (perturbation < 0. ? -1.0 : 1.0) * 0.1;
741  perturbation *= perturbation_scale;
742  soln_local += perturbation;
743  }
744  if (is_transient)
745  {
747  soln_old_local = solution_old(idx);
748 
750  soln_older_local = solution_older(idx);
751 
752  u_dot_local = u_dot(idx);
753  }
754 
755  for (unsigned int qp = 0; qp < nqp; qp++)
756  {
757  phi_local = _phi[i][qp];
758  dphi_qp = &_grad_phi[i][qp];
759 
760  grad_u_qp = &_grad_u[qp];
761 
762  if (is_transient)
763  {
764  if (_need_grad_old)
765  grad_u_old_qp = &_grad_u_old[qp];
766 
767  if (_need_grad_older)
768  grad_u_older_qp = &_grad_u_older[qp];
769  }
770 
772  {
773  d2phi_local = &(*_second_phi)[i][qp];
774 
775  if (_need_second)
776  {
777  second_u_qp = &_second_u[qp];
778  second_u_qp->add_scaled(*d2phi_local, soln_local);
779  }
780 
781  if (is_transient)
782  {
783  if (_need_second_old)
784  second_u_old_qp = &_second_u_old[qp];
785 
786  if (_need_second_older)
787  second_u_older_qp = &_second_u_older[qp];
788  }
789  }
790 
791  _u[qp] += phi_local * soln_local;
792 
793  grad_u_qp->add_scaled(*dphi_qp, soln_local);
794 
795  if (is_transient)
796  {
797  _u_dot[qp] += phi_local * u_dot_local;
798  _du_dot_du[qp] = du_dot_du;
799 
800  if (_need_u_old)
801  _u_old[qp] += phi_local * soln_old_local;
802 
803  if (_need_u_older)
804  _u_older[qp] += phi_local * soln_older_local;
805 
806  if (_need_grad_old)
807  grad_u_old_qp->add_scaled(*dphi_qp, soln_old_local);
808 
809  if (_need_grad_older)
810  grad_u_older_qp->add_scaled(*dphi_qp, soln_older_local);
811 
812  if (_need_second_old)
813  second_u_old_qp->add_scaled(*d2phi_local, soln_old_local);
814 
815  if (_need_second_older)
816  second_u_older_qp->add_scaled(*d2phi_local, soln_older_local);
817  }
818  }
819  }
820 }
821 
822 void
824 {
825  _u = _u_bak;
827  if (_need_second)
829 
830  if (_subproblem.isTransient())
831  {
832  _u_dot = _u_dot_bak;
834 
835  if (_need_u_old)
836  _u_old = _u_old_bak;
837 
838  if (_need_u_older)
840 
841  if (_need_grad_old)
843 
844  if (_need_grad_older)
846 
847  if (_need_second_old)
849 
850  if (_need_second_older)
852  }
853 }
854 
855 void
857 {
858 
859  bool is_transient = _subproblem.isTransient();
860  unsigned int nqp = _qrule->n_points();
861 
862  _u.resize(nqp);
863  _grad_u.resize(nqp);
864 
865  if (_need_second)
866  _second_u.resize(nqp);
867 
869  _u_previous_nl.resize(nqp);
870 
873 
876 
877  if (is_transient)
878  {
879  _u_dot.resize(nqp);
880  _du_dot_du.resize(nqp);
881 
882  if (_need_u_old)
883  _u_old.resize(nqp);
884 
885  if (_need_u_older)
886  _u_older.resize(nqp);
887 
888  if (_need_grad_old)
889  _grad_u_old.resize(nqp);
890 
891  if (_need_grad_older)
892  _grad_u_older.resize(nqp);
893 
894  if (_need_second_old)
895  _second_u_old.resize(nqp);
896 
897  if (_need_second_older)
898  _second_u_older.resize(nqp);
899  }
900 
901  for (unsigned int i = 0; i < nqp; ++i)
902  {
903  _u[i] = 0;
904  _grad_u[i] = 0;
905 
906  if (_need_second)
907  _second_u[i] = 0;
908 
910  _u_previous_nl[i] = 0;
911 
913  _grad_u_previous_nl[i] = 0;
914 
916  _second_u_previous_nl[i] = 0;
917 
918  if (is_transient)
919  {
920  _u_dot[i] = 0;
921  _du_dot_du[i] = 0;
922 
923  if (_need_u_old)
924  _u_old[i] = 0;
925 
926  if (_need_u_older)
927  _u_older[i] = 0;
928 
929  if (_need_grad_old)
930  _grad_u_old[i] = 0;
931 
932  if (_need_grad_older)
933  _grad_u_older[i] = 0;
934 
935  if (_need_second_old)
936  _second_u_old[i] = 0;
937 
938  if (_need_second_older)
939  _second_u_older[i] = 0;
940  }
941  }
942 
943  unsigned int num_dofs = _dof_indices.size();
944 
945  if (_need_nodal_u)
946  _nodal_u.resize(num_dofs);
947 
949  _nodal_u_previous_nl.resize(num_dofs);
950 
951  if (is_transient)
952  {
953  if (_need_nodal_u_old)
954  _nodal_u_old.resize(num_dofs);
956  _nodal_u_older.resize(num_dofs);
957  if (_need_nodal_u_dot)
958  _nodal_u_dot.resize(num_dofs);
959  }
960 
962  _solution_dofs.resize(num_dofs);
963 
965  _solution_dofs_old.resize(num_dofs);
966 
968  _solution_dofs_older.resize(num_dofs);
969 
970  const NumericVector<Real> & current_solution = *_sys.currentSolution();
971  const NumericVector<Real> & solution_old = _sys.solutionOld();
972  const NumericVector<Real> & solution_older = _sys.solutionOlder();
973  const NumericVector<Real> * solution_prev_nl = _sys.solutionPreviousNewton();
974  const NumericVector<Real> & u_dot = _sys.solutionUDot();
975  const Real & du_dot_du = _sys.duDotDu();
976 
977  dof_id_type idx = 0;
978  Real soln_local = 0;
979  Real soln_old_local = 0;
980  Real soln_older_local = 0;
981  Real soln_previous_nl_local = 0;
982  Real u_dot_local = 0;
983 
984  Real phi_local = 0;
985  const RealGradient * dphi_qp = NULL;
986  const RealTensor * d2phi_local = NULL;
987 
988  RealGradient * grad_u_qp = NULL;
989 
990  RealGradient * grad_u_old_qp = NULL;
991  RealGradient * grad_u_older_qp = NULL;
992  RealGradient * grad_u_previous_nl_qp = NULL;
993 
994  RealTensor * second_u_qp = NULL;
995 
996  RealTensor * second_u_old_qp = NULL;
997  RealTensor * second_u_older_qp = NULL;
998  RealTensor * second_u_previous_nl_qp = NULL;
999 
1000  for (unsigned int i = 0; i < num_dofs; i++)
1001  {
1002  idx = _dof_indices[i];
1003  soln_local = current_solution(idx);
1004 
1005  if (_need_nodal_u)
1006  _nodal_u[i] = soln_local;
1007 
1010  soln_previous_nl_local = (*solution_prev_nl)(idx);
1011 
1013  _nodal_u_previous_nl[i] = soln_previous_nl_local;
1014 
1015  if (_need_solution_dofs)
1016  _solution_dofs(i) = soln_local;
1017 
1018  if (is_transient)
1019  {
1021  soln_old_local = solution_old(idx);
1022 
1024  soln_older_local = solution_older(idx);
1025 
1026  if (_need_nodal_u_old)
1027  _nodal_u_old[i] = soln_old_local;
1028  if (_need_nodal_u_older)
1029  _nodal_u_older[i] = soln_older_local;
1030 
1031  u_dot_local = u_dot(idx);
1032  if (_need_nodal_u_dot)
1033  _nodal_u_dot[i] = u_dot_local;
1034 
1036  _solution_dofs_old(i) = solution_old(idx);
1037 
1039  _solution_dofs_older(i) = solution_older(idx);
1040  }
1041 
1042  for (unsigned int qp = 0; qp < nqp; qp++)
1043  {
1044  phi_local = _phi[i][qp];
1045  dphi_qp = &_grad_phi[i][qp];
1046 
1047  grad_u_qp = &_grad_u[qp];
1048 
1050  grad_u_previous_nl_qp = &_grad_u_previous_nl[qp];
1051 
1052  if (is_transient)
1053  {
1054  if (_need_grad_old)
1055  grad_u_old_qp = &_grad_u_old[qp];
1056 
1057  if (_need_grad_older)
1058  grad_u_older_qp = &_grad_u_older[qp];
1059  }
1060 
1062  {
1063  d2phi_local = &(*_second_phi)[i][qp];
1064 
1065  if (_need_second)
1066  {
1067  second_u_qp = &_second_u[qp];
1068  second_u_qp->add_scaled(*d2phi_local, soln_local);
1069  }
1070 
1072  {
1073  second_u_previous_nl_qp = &_second_u_previous_nl[qp];
1074  second_u_previous_nl_qp->add_scaled(*d2phi_local, soln_previous_nl_local);
1075  }
1076 
1077  if (is_transient)
1078  {
1079  if (_need_second_old)
1080  second_u_old_qp = &_second_u_old[qp];
1081 
1082  if (_need_second_older)
1083  second_u_older_qp = &_second_u_older[qp];
1084  }
1085  }
1086 
1087  _u[qp] += phi_local * soln_local;
1088 
1089  grad_u_qp->add_scaled(*dphi_qp, soln_local);
1090 
1091  if (_need_u_previous_nl)
1092  _u_previous_nl[qp] += phi_local * soln_previous_nl_local;
1094  grad_u_previous_nl_qp->add_scaled(*dphi_qp, soln_previous_nl_local);
1095 
1096  if (is_transient)
1097  {
1098  _u_dot[qp] += phi_local * u_dot_local;
1099  _du_dot_du[qp] = du_dot_du;
1100 
1101  if (_need_u_old)
1102  _u_old[qp] += phi_local * soln_old_local;
1103 
1104  if (_need_u_older)
1105  _u_older[qp] += phi_local * soln_older_local;
1106 
1107  if (_need_grad_old)
1108  grad_u_old_qp->add_scaled(*dphi_qp, soln_old_local);
1109 
1110  if (_need_grad_older)
1111  grad_u_older_qp->add_scaled(*dphi_qp, soln_older_local);
1112 
1113  if (_need_second_old)
1114  second_u_old_qp->add_scaled(*d2phi_local, soln_old_local);
1115 
1116  if (_need_second_older)
1117  second_u_older_qp->add_scaled(*d2phi_local, soln_older_local);
1118  }
1119  }
1120  }
1121 }
1122 
1123 void
1125 {
1126  bool is_transient = _subproblem.isTransient();
1127  unsigned int nqp = _qrule_face->n_points();
1128  _u.resize(nqp);
1129  _grad_u.resize(nqp);
1130 
1131  if (_need_second)
1132  _second_u.resize(nqp);
1133 
1134  if (_need_u_previous_nl)
1135  _u_previous_nl.resize(nqp);
1136 
1139 
1142 
1143  if (is_transient)
1144  {
1145  _u_dot.resize(nqp);
1146  _du_dot_du.resize(nqp);
1147 
1148  if (_need_u_old)
1149  _u_old.resize(nqp);
1150 
1151  if (_need_u_older)
1152  _u_older.resize(nqp);
1153 
1154  if (_need_grad_old)
1155  _grad_u_old.resize(nqp);
1156 
1157  if (_need_grad_older)
1158  _grad_u_older.resize(nqp);
1159 
1160  if (_need_second_old)
1161  _second_u_old.resize(nqp);
1162 
1163  if (_need_second_older)
1164  _second_u_older.resize(nqp);
1165  }
1166 
1167  for (unsigned int i = 0; i < nqp; ++i)
1168  {
1169  _u[i] = 0;
1170  _grad_u[i] = 0;
1171 
1172  if (_need_second)
1173  _second_u[i] = 0;
1174 
1175  if (_need_u_previous_nl)
1176  _u_previous_nl[i] = 0;
1177 
1179  _grad_u_previous_nl[i] = 0;
1180 
1182  _second_u_previous_nl[i] = 0;
1183 
1184  if (_subproblem.isTransient())
1185  {
1186  _u_dot[i] = 0;
1187  _du_dot_du[i] = 0;
1188 
1189  if (_need_u_old)
1190  _u_old[i] = 0;
1191 
1192  if (_need_u_older)
1193  _u_older[i] = 0;
1194 
1195  if (_need_grad_old)
1196  _grad_u_old[i] = 0;
1197 
1198  if (_need_grad_older)
1199  _grad_u_older[i] = 0;
1200 
1201  if (_need_second_old)
1202  _second_u_old[i] = 0;
1203 
1204  if (_need_second_older)
1205  _second_u_older[i] = 0;
1206  }
1207  }
1208 
1209  unsigned int num_dofs = _dof_indices.size();
1210 
1211  if (_need_nodal_u)
1212  _nodal_u.resize(num_dofs);
1214  _nodal_u_previous_nl.resize(num_dofs);
1215  if (is_transient)
1216  {
1217  if (_need_nodal_u_old)
1218  _nodal_u_old.resize(num_dofs);
1219  if (_need_nodal_u_older)
1220  _nodal_u_older.resize(num_dofs);
1221  if (_need_nodal_u_dot)
1222  _nodal_u_dot.resize(num_dofs);
1223  }
1224 
1225  if (_need_solution_dofs)
1226  _solution_dofs.resize(num_dofs);
1227 
1229  _solution_dofs_old.resize(num_dofs);
1230 
1232  _solution_dofs_older.resize(num_dofs);
1233 
1234  const NumericVector<Real> & current_solution = *_sys.currentSolution();
1235  const NumericVector<Real> & solution_old = _sys.solutionOld();
1236  const NumericVector<Real> & solution_older = _sys.solutionOlder();
1237  const NumericVector<Real> * solution_prev_nl = _sys.solutionPreviousNewton();
1238  const NumericVector<Real> & u_dot = _sys.solutionUDot();
1239  const Real & du_dot_du = _sys.duDotDu();
1240 
1241  dof_id_type idx;
1242  Real soln_local;
1243  Real soln_old_local = 0;
1244  Real soln_older_local = 0;
1245  Real soln_previous_nl_local = 0;
1246  Real u_dot_local = 0;
1247 
1248  Real phi_local;
1249  RealGradient dphi_local;
1250  RealTensor d2phi_local;
1251 
1252  for (unsigned int i = 0; i < num_dofs; ++i)
1253  {
1254  idx = _dof_indices[i];
1255  soln_local = current_solution(idx);
1256 
1257  if (_need_nodal_u)
1258  _nodal_u[i] = soln_local;
1259 
1262  soln_previous_nl_local = (*solution_prev_nl)(idx);
1264  _nodal_u_previous_nl[i] = soln_previous_nl_local;
1265 
1266  if (_need_solution_dofs)
1267  _solution_dofs(i) = soln_local;
1268 
1269  if (is_transient)
1270  {
1272  soln_old_local = solution_old(idx);
1273 
1275  soln_older_local = solution_older(idx);
1276 
1277  if (_need_nodal_u_old)
1278  _nodal_u_old[i] = soln_old_local;
1279  if (_need_nodal_u_older)
1280  _nodal_u_older[i] = soln_older_local;
1281 
1282  u_dot_local = u_dot(idx);
1283  if (_need_nodal_u_dot)
1284  _nodal_u_dot[i] = u_dot_local;
1285 
1287  _solution_dofs_old(i) = solution_old(idx);
1288 
1290  _solution_dofs_older(i) = solution_older(idx);
1291  }
1292 
1293  for (unsigned int qp = 0; qp < nqp; ++qp)
1294  {
1295  phi_local = _phi_face[i][qp];
1296  dphi_local = _grad_phi_face[i][qp];
1297 
1299  d2phi_local = (*_second_phi_face)[i][qp];
1300 
1301  _u[qp] += phi_local * soln_local;
1302  _grad_u[qp] += dphi_local * soln_local;
1303 
1304  if (_need_second)
1305  _second_u[qp] += d2phi_local * soln_local;
1306 
1307  if (_need_u_previous_nl)
1308  _u_previous_nl[qp] += phi_local * soln_previous_nl_local;
1309 
1311  _grad_u_previous_nl[qp] += dphi_local * soln_previous_nl_local;
1312 
1314  _second_u_previous_nl[qp] += d2phi_local * soln_previous_nl_local;
1315 
1316  if (is_transient)
1317  {
1318  _u_dot[qp] += phi_local * u_dot_local;
1319  _du_dot_du[qp] = du_dot_du;
1320 
1321  if (_need_u_old)
1322  _u_old[qp] += phi_local * soln_old_local;
1323 
1324  if (_need_u_older)
1325  _u_older[qp] += phi_local * soln_older_local;
1326 
1327  if (_need_grad_old)
1328  _grad_u_old[qp] += dphi_local * soln_old_local;
1329 
1330  if (_need_grad_older)
1331  _grad_u_older[qp] += dphi_local * soln_older_local;
1332 
1333  if (_need_second_old)
1334  _second_u_old[qp] += d2phi_local * soln_old_local;
1335 
1336  if (_need_second_older)
1337  _second_u_older[qp] += d2phi_local * soln_older_local;
1338  }
1339  }
1340  }
1341 }
1342 
1343 void
1345 {
1346  bool is_transient = _subproblem.isTransient();
1347  unsigned int nqp = _qrule_neighbor->n_points();
1348 
1349  _u_neighbor.resize(nqp);
1350  _grad_u_neighbor.resize(nqp);
1351 
1354 
1355  if (is_transient)
1356  {
1357  _u_dot_neighbor.resize(nqp);
1359 
1361  _u_old_neighbor.resize(nqp);
1362 
1365 
1368 
1371 
1374 
1377  }
1378 
1379  for (unsigned int i = 0; i < nqp; ++i)
1380  {
1381  _u_neighbor[i] = 0;
1382  _grad_u_neighbor[i] = 0;
1383 
1385  _second_u_neighbor[i] = 0;
1386 
1387  if (_subproblem.isTransient())
1388  {
1389  _u_dot_neighbor[i] = 0;
1390  _du_dot_du_neighbor[i] = 0;
1391 
1393  _u_old_neighbor[i] = 0;
1394 
1396  _u_older_neighbor[i] = 0;
1397 
1399  _grad_u_old_neighbor[i] = 0;
1400 
1402  _grad_u_older_neighbor[i] = 0;
1403 
1405  _second_u_old_neighbor[i] = 0;
1406 
1408  _second_u_older_neighbor[i] = 0;
1409  }
1410  }
1411 
1412  unsigned int num_dofs = _dof_indices_neighbor.size();
1413 
1415  _nodal_u_neighbor.resize(num_dofs);
1416  if (is_transient)
1417  {
1419  _nodal_u_old_neighbor.resize(num_dofs);
1421  _nodal_u_older_neighbor.resize(num_dofs);
1423  _nodal_u_dot_neighbor.resize(num_dofs);
1424  }
1425 
1427  _solution_dofs_neighbor.resize(num_dofs);
1428 
1430  _solution_dofs_old_neighbor.resize(num_dofs);
1431 
1433  _solution_dofs_older_neighbor.resize(num_dofs);
1434 
1435  const NumericVector<Real> & current_solution = *_sys.currentSolution();
1436  const NumericVector<Real> & solution_old = _sys.solutionOld();
1437  const NumericVector<Real> & solution_older = _sys.solutionOlder();
1438  const NumericVector<Real> & u_dot = _sys.solutionUDot();
1439  const Real & du_dot_du = _sys.duDotDu();
1440 
1441  dof_id_type idx;
1442  Real soln_local;
1443  Real soln_old_local = 0;
1444  Real soln_older_local = 0;
1445  Real u_dot_local = 0;
1446 
1447  Real phi_local;
1448  RealGradient dphi_local;
1449  RealTensor d2phi_local;
1450 
1451  for (unsigned int i = 0; i < num_dofs; ++i)
1452  {
1453  idx = _dof_indices_neighbor[i];
1454  soln_local = current_solution(idx);
1455 
1457  _nodal_u_neighbor[i] = soln_local;
1458 
1460  _solution_dofs_neighbor(i) = soln_local;
1461 
1462  if (is_transient)
1463  {
1465  soln_old_local = solution_old(idx);
1466 
1468  soln_older_local = solution_older(idx);
1469 
1471  _nodal_u_old_neighbor[i] = soln_old_local;
1473  _nodal_u_older_neighbor[i] = soln_older_local;
1474 
1475  u_dot_local = u_dot(idx);
1477  _nodal_u_dot_neighbor[i] = u_dot_local;
1478 
1480  _solution_dofs_old_neighbor(i) = solution_old(idx);
1481 
1483  _solution_dofs_older_neighbor(i) = solution_older(idx);
1484  }
1485 
1486  for (unsigned int qp = 0; qp < nqp; ++qp)
1487  {
1488  phi_local = _phi_face_neighbor[i][qp];
1489  dphi_local = _grad_phi_face_neighbor[i][qp];
1490 
1492  d2phi_local = (*_second_phi_face_neighbor)[i][qp];
1493 
1494  _u_neighbor[qp] += phi_local * soln_local;
1495  _grad_u_neighbor[qp] += dphi_local * soln_local;
1496 
1498  _second_u_neighbor[qp] += d2phi_local * soln_local;
1499 
1500  if (is_transient)
1501  {
1502  _u_dot_neighbor[qp] += phi_local * u_dot_local;
1503  _du_dot_du_neighbor[qp] = du_dot_du;
1504 
1506  _u_old_neighbor[qp] += phi_local * soln_old_local;
1507 
1509  _u_older_neighbor[qp] += phi_local * soln_older_local;
1510 
1512  _grad_u_old_neighbor[qp] += dphi_local * soln_old_local;
1513 
1515  _grad_u_older_neighbor[qp] += dphi_local * soln_older_local;
1516 
1518  _second_u_old_neighbor[qp] += d2phi_local * soln_old_local;
1519 
1521  _second_u_older_neighbor[qp] += d2phi_local * soln_older_local;
1522  }
1523  }
1524  }
1525 }
1526 
1527 void
1529 {
1530  bool is_transient = _subproblem.isTransient();
1531  unsigned int nqp = _qrule_neighbor->n_points();
1532 
1533  _u_neighbor.resize(nqp);
1534  _grad_u_neighbor.resize(nqp);
1535 
1538 
1539  if (is_transient)
1540  {
1542  _u_old_neighbor.resize(nqp);
1543 
1546 
1549 
1552 
1555 
1558  }
1559 
1560  for (unsigned int i = 0; i < nqp; ++i)
1561  {
1562  _u_neighbor[i] = 0;
1563  _grad_u_neighbor[i] = 0;
1564 
1566  _second_u_neighbor[i] = 0;
1567 
1568  if (_subproblem.isTransient())
1569  {
1571  _u_old_neighbor[i] = 0;
1572 
1574  _u_older_neighbor[i] = 0;
1575 
1577  _grad_u_old_neighbor[i] = 0;
1578 
1580  _grad_u_older_neighbor[i] = 0;
1581 
1583  _second_u_old_neighbor[i] = 0;
1584 
1586  _second_u_older_neighbor[i] = 0;
1587  }
1588  }
1589 
1590  unsigned int num_dofs = _dof_indices_neighbor.size();
1591 
1593  _nodal_u_neighbor.resize(num_dofs);
1594  if (is_transient)
1595  {
1597  _nodal_u_old_neighbor.resize(num_dofs);
1599  _nodal_u_older_neighbor.resize(num_dofs);
1601  _nodal_u_dot_neighbor.resize(num_dofs);
1602  }
1603 
1605  _solution_dofs_neighbor.resize(num_dofs);
1606 
1608  _solution_dofs_old_neighbor.resize(num_dofs);
1609 
1611  _solution_dofs_older_neighbor.resize(num_dofs);
1612 
1613  const NumericVector<Real> & current_solution = *_sys.currentSolution();
1614  const NumericVector<Real> & solution_old = _sys.solutionOld();
1615  const NumericVector<Real> & solution_older = _sys.solutionOlder();
1616  const NumericVector<Real> & u_dot = _sys.solutionUDot();
1617 
1618  dof_id_type idx;
1619  Real soln_local;
1620  Real soln_old_local = 0;
1621  Real soln_older_local = 0;
1622  Real u_dot_local = 0;
1623 
1624  Real phi_local;
1625  RealGradient dphi_local;
1626  RealTensor d2phi_local;
1627 
1628  for (unsigned int i = 0; i < num_dofs; ++i)
1629  {
1630  idx = _dof_indices_neighbor[i];
1631  soln_local = current_solution(idx);
1632 
1634  _nodal_u_neighbor[i] = soln_local;
1635 
1637  _solution_dofs_neighbor(i) = soln_local;
1638 
1639  if (is_transient)
1640  {
1642  soln_old_local = solution_old(idx);
1643 
1645  soln_older_local = solution_older(idx);
1646 
1648  _nodal_u_old_neighbor[i] = soln_old_local;
1650  _nodal_u_older_neighbor[i] = soln_older_local;
1651 
1652  u_dot_local = u_dot(idx);
1654  _nodal_u_dot_neighbor[i] = u_dot_local;
1655 
1657  _solution_dofs_old_neighbor(i) = solution_old(idx);
1658 
1660  _solution_dofs_older_neighbor(i) = solution_older(idx);
1661  }
1662 
1663  for (unsigned int qp = 0; qp < nqp; ++qp)
1664  {
1665  phi_local = _phi_neighbor[i][qp];
1666  dphi_local = _grad_phi_neighbor[i][qp];
1667 
1669  d2phi_local = (*_second_phi_neighbor)[i][qp];
1670 
1671  _u_neighbor[qp] += phi_local * soln_local;
1672  _grad_u_neighbor[qp] += dphi_local * soln_local;
1673 
1675  _second_u_neighbor[qp] += d2phi_local * soln_local;
1676 
1677  if (is_transient)
1678  {
1680  _u_old_neighbor[qp] += phi_local * soln_old_local;
1681 
1683  _u_older_neighbor[qp] += phi_local * soln_older_local;
1684 
1686  _grad_u_old_neighbor[qp] += dphi_local * soln_old_local;
1687 
1689  _grad_u_older_neighbor[qp] += dphi_local * soln_older_local;
1690 
1692  _second_u_old_neighbor[qp] += d2phi_local * soln_old_local;
1693 
1695  _second_u_older_neighbor[qp] += d2phi_local * soln_older_local;
1696  }
1697  }
1698  }
1699 }
1700 
1701 void
1703 {
1704  if (_is_defined)
1705  {
1706  const unsigned int n = _dof_indices.size();
1707  mooseAssert(n, "Defined but empty?");
1708  _nodal_u.resize(n);
1710 
1712  {
1715  }
1716 
1717  if (_subproblem.isTransient())
1718  {
1719  _nodal_u_old.resize(n);
1723 
1724  _nodal_u_dot.resize(n);
1726  for (unsigned int i = 0; i < n; i++)
1727  {
1729  _nodal_du_dot_du[i] = _sys.duDotDu();
1730  }
1731  }
1732  }
1733  else
1734  {
1735  _nodal_u.resize(0);
1738  if (_subproblem.isTransient())
1739  {
1740  _nodal_u_old.resize(0);
1742  _nodal_u_dot.resize(0);
1744  }
1745  }
1746 }
1747 
1748 void
1750 {
1752  {
1753  const unsigned int n = _dof_indices_neighbor.size();
1754  mooseAssert(n, "Defined but empty?");
1757 
1758  if (_subproblem.isTransient())
1759  {
1764 
1767  for (unsigned int i = 0; i < n; i++)
1768  {
1771  }
1772  }
1773  }
1774  else
1775  {
1777  if (_subproblem.isTransient())
1778  {
1783  }
1784  }
1785 }
1786 
1787 void
1788 MooseVariable::setNodalValue(const DenseVector<Number> & values)
1789 {
1790  for (unsigned int i = 0; i < values.size(); i++)
1791  _nodal_u[i] = values(i);
1792 
1793  _has_nodal_value = true;
1794 
1795  for (unsigned int qp = 0; qp < _u.size(); qp++)
1796  {
1797  _u[qp] = 0;
1798  for (unsigned int i = 0; i < _nodal_u.size(); i++)
1799  _u[qp] += _phi[i][qp] * _nodal_u[i];
1800  }
1801 }
1802 
1803 void
1804 MooseVariable::setNodalValueNeighbor(const DenseVector<Number> & values)
1805 {
1806  for (unsigned int i = 0; i < values.size(); i++)
1807  _nodal_u_neighbor[i] = values(i);
1809 
1810  if (isNodal())
1811  mooseError("Variable " + name() + " has to be nodal!");
1812  else
1813  {
1814  for (unsigned int qp = 0; qp < _u_neighbor.size(); qp++)
1815  {
1816  _u_neighbor[qp] = 0;
1817  for (unsigned int i = 0; i < _nodal_u_neighbor.size(); i++)
1819  }
1820  }
1821 }
1822 
1823 void
1824 MooseVariable::setNodalValue(Number value, unsigned int idx /* = 0*/)
1825 {
1826  _nodal_u[idx] = value; // update variable nodal value
1827  _has_nodal_value = true;
1828 
1829  // Update the qp values as well
1830  for (unsigned int qp = 0; qp < _u.size(); qp++)
1831  _u[qp] = value;
1832 }
1833 
1834 void
1836 {
1837  _nodal_u_neighbor[0] = value; // update variable nodal value
1839 
1840  if (!isNodal()) // If this is an elemental variable, then update the qp values as well
1841  {
1842  for (unsigned int qp = 0; qp < _u_neighbor.size(); qp++)
1843  _u_neighbor[qp] = value;
1844  }
1845 }
1846 
1847 void
1848 MooseVariable::computeIncrementAtQps(const NumericVector<Number> & increment_vec)
1849 {
1850  unsigned int nqp = _qrule->n_points();
1851 
1852  _increment.resize(nqp);
1853  // Compute the increment at each quadrature point
1854  unsigned int num_dofs = _dof_indices.size();
1855  for (unsigned int qp = 0; qp < nqp; qp++)
1856  {
1857  _increment[qp] = 0;
1858  for (unsigned int i = 0; i < num_dofs; i++)
1859  _increment[qp] += _phi[i][qp] * increment_vec(_dof_indices[i]);
1860  }
1861 }
1862 
1863 void
1864 MooseVariable::computeIncrementAtNode(const NumericVector<Number> & increment_vec)
1865 {
1866  if (!isNodal())
1867  mooseError("computeIncrementAtNode can only be called for nodal variables");
1868 
1869  _increment.resize(1);
1870 
1871  // Compute the increment for the current DOF
1872  _increment[0] = increment_vec(_dof_indices[0]);
1873 }
1874 
1875 Number
1877 {
1878  mooseAssert(_subproblem.mesh().isSemiLocal(const_cast<Node *>(&node)), "Node is not Semilocal");
1879 
1880  // Make sure that the node has DOFs
1881  /* Note, this is a reproduction of an assert within libMesh::Node::dof_number, this is done to
1882  * produce a
1883  * better error (see misc/check_error.node_value_off_block) */
1884  mooseAssert(node.n_dofs(_sys.number(), _var_num) > 0,
1885  "Node " << node.id() << " does not contain any dofs for the "
1886  << _sys.system().variable_name(_var_num)
1887  << " variable");
1888 
1889  dof_id_type dof = node.dof_number(_sys.number(), _var_num, 0);
1890 
1891  return (*_sys.currentSolution())(dof);
1892 }
1893 
1894 Number
1896 {
1897  mooseAssert(_subproblem.mesh().isSemiLocal(const_cast<Node *>(&node)), "Node is not Semilocal");
1898 
1899  // Make sure that the node has DOFs
1900  /* Note, this is a reproduction of an assert within libMesh::Node::dof_number, this is done to
1901  * produce a
1902  * better error (see misc/check_error.node_value_off_block) */
1903  mooseAssert(node.n_dofs(_sys.number(), _var_num) > 0,
1904  "Node " << node.id() << " does not contain any dofs for the "
1905  << _sys.system().variable_name(_var_num)
1906  << " variable");
1907 
1908  dof_id_type dof = node.dof_number(_sys.number(), _var_num, 0);
1909  return _sys.solutionOld()(dof);
1910 }
1911 
1912 Number
1914 {
1915  mooseAssert(_subproblem.mesh().isSemiLocal(const_cast<Node *>(&node)), "Node is not Semilocal");
1916 
1917  // Make sure that the node has DOFs
1918  /* Note, this is a reproduction of an assert within libMesh::Node::dof_number, this is done to
1919  * produce a
1920  * better error (see misc/check_error.node_value_off_block) */
1921  mooseAssert(node.n_dofs(_sys.number(), _var_num) > 0,
1922  "Node " << node.id() << " does not contain any dofs for the "
1923  << _sys.system().variable_name(_var_num)
1924  << " variable");
1925 
1926  dof_id_type dof = node.dof_number(_sys.number(), _var_num, 0);
1927  return _sys.solutionOlder()(dof);
1928 }
1929 
1930 Real
1931 MooseVariable::getValue(const Elem * elem, const std::vector<std::vector<Real>> & phi) const
1932 {
1933  std::vector<dof_id_type> dof_indices;
1934  _dof_map.dof_indices(elem, dof_indices, _var_num);
1935 
1936  Real value = 0;
1937  if (isNodal())
1938  {
1939  for (unsigned int i = 0; i < dof_indices.size(); ++i)
1940  {
1941  // The zero index is because we only have one point that the phis are evaluated at
1942  value += phi[i][0] * (*_sys.currentSolution())(dof_indices[i]);
1943  }
1944  }
1945  else
1946  {
1947  mooseAssert(dof_indices.size() == 1, "Wrong size for dof indices");
1948  value = (*_sys.currentSolution())(dof_indices[0]);
1949  }
1950 
1951  return value;
1952 }
1953 
1955 MooseVariable::getGradient(const Elem * elem,
1956  const std::vector<std::vector<RealGradient>> & grad_phi) const
1957 {
1958  std::vector<dof_id_type> dof_indices;
1959  _dof_map.dof_indices(elem, dof_indices, _var_num);
1960 
1961  RealGradient value;
1962  if (isNodal())
1963  {
1964  for (unsigned int i = 0; i < dof_indices.size(); ++i)
1965  {
1966  // The zero index is because we only have one point that the phis are evaluated at
1967  value += grad_phi[i][0] * (*_sys.currentSolution())(dof_indices[i]);
1968  }
1969  }
1970  else
1971  {
1972  mooseAssert(dof_indices.size() == 1, "Wrong size for dof indices");
1973  value = 0.0;
1974  }
1975 
1976  return value;
1977 }
1978 
1979 Real
1980 MooseVariable::getElementalValue(const Elem * elem, unsigned int idx) const
1981 {
1982  std::vector<dof_id_type> dof_indices;
1983  _dof_map.dof_indices(elem, dof_indices, _var_num);
1984 
1985  return (*_sys.currentSolution())(dof_indices[idx]);
1986 }
bool _need_grad_older_neighbor
const VariablePhiValue & phiFace()
virtual MooseMesh & mesh()=0
VariableValue _increment
const VariablePhiGradient & _grad_phi
RealVectorValue RealGradient
Definition: Assembly.h:43
const VariablePhiSecond & feSecondPhiFace(FEType type)
Definition: Assembly.C:301
void restoreUnperturbedElemValues()
Restore the values the variable had before a call to computePerturbedElemValues().
const VariablePhiSecond & secondPhiFace()
bool _need_second_old_neighbor
bool _need_u_old_neighbor
VariableSecond _second_u_previous_nl_neighbor
bool _need_solution_dofs_older_neighbor
FEType _fe_type
The FEType associated with this variable.
const VariableValue & nodalValue()
const VariablePhiSecond * _second_phi_neighbor
bool _need_second_older
bool _need_nodal_u_older
VariableValue _nodal_u_old
dof_id_type _nodal_dof_index_neighbor
Number getNodalValueOld(const Node &node)
Get the old value of this variable at given node.
DenseVector< Number > _solution_dofs_older
VariableValue _u_older_bak
const FEType & feType() const
Get the type of finite element object.
Keeps track of stuff related to assembling.
Definition: Assembly.h:63
bool _need_nodal_u_previous_nl
virtual unsigned int dimension() const
Returns MeshBase::mesh_dimsension(), (not MeshBase::spatial_dimension()!) of the underlying libMesh m...
Definition: MooseMesh.C:1945
VariableGradient _grad_u_old_bak
void prepareNeighbor()
const VariablePhiSecond * _second_phi
const VariableValue & nodalValueNeighbor()
const DofMap & _dof_map
DOF map.
VariableSecond _second_u_older
subdomain_id_type SubdomainID
Definition: MooseTypes.h:77
void resize(const unsigned int size)
Change the number of elements the array can store.
Definition: MooseArray.h:205
const VariablePhiValue & phiNeighbor()
VariableSecond _second_u_old
virtual NumericVector< Number > & solutionOld()=0
bool _need_grad_previous_nl
VariableGradient _grad_u_neighbor
VariableGradient _grad_u_old_neighbor
void mooseError(Args &&...args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:182
bool _is_defined_neighbor
If the variable is defined at the neighbor node (used in compute nodal values)
virtual void computeElemValuesFace()
Compute values at facial quadrature points.
DenseVector< Number > _solution_dofs_neighbor
bool _is_defined
If the variable is defined at the node (used in compute nodal values)
const VariablePhiSecond & feSecondPhiFaceNeighbor(FEType type)
Definition: Assembly.C:345
const VariablePhiSecond & feSecondPhi(FEType type)
Definition: Assembly.C:279
bool _need_grad_old_neighbor
QBase *& _qrule_face
Quadrature rule for the face.
virtual bool isTransient() const =0
const VariablePhiValue & _phi_neighbor
bool _has_nodal_value
If true, the nodal value gets inserted on calling insert()
const Node *& _node
VariableValue _nodal_u_older_neighbor
bool _need_solution_dofs_old
RealTensorValue RealTensor
Definition: Assembly.h:46
const VariablePhiGradient & _grad_phi_face
VariableValue _u_old_neighbor
void computeNodalNeighborValues()
Compute nodal values of this variable in the neighbor.
VariableValue _du_dot_du
derivative of u_dot wrt u
VariableValue _u_previous_nl_neighbor
VariableGradient _grad_u_older
bool _need_second_previous_nl
void computePerturbedElemValues(unsigned i, Real scale, Real &h)
Compute values at interior quadrature points when this variable&#39;s elem dof i is perturbed by h in the...
void computeIncrementAtQps(const NumericVector< Number > &increment_vec)
Compute and store incremental change in solution at QPs based on increment_vec.
Base class for a system (of equations)
Definition: SystemBase.h:91
VariableValue _u_older_neighbor
const VariableValue & nodalValuePreviousNL()
VariableGradient _grad_u_old
FEBase *& getFE(FEType type, unsigned int dim)
Get a reference to a pointer that will contain the current volume FE.
Definition: Assembly.C:244
const std::string & name() const
Get the variable name.
virtual NumericVector< Number > * solutionPreviousNewton()=0
const VariablePhiSecond * _second_phi_face_neighbor
bool _need_nodal_u_older_neighbor
std::vector< dof_id_type > _dof_indices
DOF indices.
const VariablePhiSecond & secondPhiFaceNeighbor()
const VariableValue & nodalValueDotNeighbor()
void reinitAuxNeighbor()
const std::set< SubdomainID > & activeSubdomains() const
void reinitNodesNeighbor(const std::vector< dof_id_type > &nodes)
void reinitNodes(const std::vector< dof_id_type > &nodes)
VariableValue _nodal_u
VariableSecond _second_u
VariableValue _nodal_u_neighbor
virtual bool isNodal() const override
Is this variable nodal.
virtual Number & duDotDu()
Definition: SystemBase.h:156
VariableValue _nodal_du_dot_du
nodal values of derivative of u_dot wrt u
VariableSecond _second_u_old_neighbor
void computeNodalValues()
Compute nodal values of this variable.
VariableValue _u
QBase *& _qrule
Quadrature rule for interior.
VariableGradient _grad_u
const VariableValue & nodalValueOldNeighbor()
const VariablePhiSecond & secondPhi()
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2355
SubProblem & _subproblem
Problem this variable is part of.
virtual NumericVector< Number > & solutionUDot()
Definition: SystemBase.h:157
VariableValue _u_dot_neighbor
const VariablePhiSecond * _second_phi_face
const VariablePhiGradient & _grad_phi_neighbor
const Node *& _node_neighbor
const VariablePhiValue & phi()
bool _need_nodal_u_previous_nl_neighbor
SystemBase & _sys
System this variable is part of.
dof_id_type _nodal_dof_index
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:156
const VariablePhiGradient & gradPhi()
const Node *& node()
void setNodalValue(Number value, unsigned int idx=0)
Set the nodal value for this variable to keep everything up to date.
bool _need_nodal_u_old_neighbor
VariableValue _u_dot_bak_neighbor
bool _need_nodal_u_old
VariableValue _nodal_u_previous_nl
bool activeOnSubdomain(SubdomainID subdomain) const
Is the variable active on the subdomain?
Number getElementalValue(const Elem *elem, unsigned int idx=0) const
Retrieve the Elemental DOF.
VariableGradient _grad_u_previous_nl
VariableValue _du_dot_du_bak_neighbor
const VariablePhiValue & _phi_face_neighbor
bool _has_nodal_value_neighbor
void buildFE(FEType type)
Build FEs with a type.
Definition: Assembly.C:168
VariableValue _u_previous_nl
RealGradient getGradient(const Elem *elem, const std::vector< std::vector< RealGradient >> &phi) const
VariableValue _u_old
VariableValue _nodal_u_previous_nl_neighbor
QBase *& _qrule_neighbor
Quadrature rule for the neighbor.
VariableValue _u_neighbor
virtual System & system()=0
Get the reference to the libMesh system.
VariableSecond _second_u_old_bak
const VariablePhiValue & _phi_face
bool _need_solution_dofs_neighbor
VariableValue _u_bak
DenseVector< Number > _solution_dofs_old_neighbor
virtual void computeNeighborValues()
Compute values at quadrature points for the neighbor.
const VariablePhiGradient & gradPhiFaceNeighbor()
Real getValue(const Elem *elem, const std::vector< std::vector< Real >> &phi) const
Compute the variable value at a point on an element.
DenseVector< Number > _solution_dofs_old
bool _need_nodal_u_dot_neighbor
const VariablePhiValue & _phi
bool _need_second_neighbor
const VariableValue & nodalValueOld()
virtual unsigned int number()
Gets the number of this system.
Definition: SystemBase.C:604
bool _need_nodal_u_dot
virtual void computeElemValues()
Compute values at interior quadrature points.
bool _is_nodal
if variable is nodal
bool _need_solution_dofs_older
VariableValue _nodal_du_dot_du_neighbor
VariableValue _nodal_u_dot
nodal values of u_dot
VariableValue _u_dot_bak
void clearDofIndices()
Clear out the dof indices.
const VariablePhiSecond & secondPhiNeighbor()
VariableSecond _second_u_older_bak
Assembly & _assembly
Assembly data.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
DenseVector< Number > _solution_dofs_older_neighbor
const Elem *& _elem
current element
VariableGradient _grad_u_bak
Number getNodalValueOlder(const Node &node)
Get the t-2 value of this variable at given node.
PetscInt n
void release()
Manually deallocates the data pointer.
Definition: MooseArray.h:56
VariableSecond _second_u_previous_nl
VariableValue _du_dot_du_bak
const VariableValue & nodalValueOlderNeighbor()
const VariablePhiSecond & feSecondPhiNeighbor(FEType type)
Definition: Assembly.C:323
VariableGradient _grad_u_previous_nl_neighbor
bool isSemiLocal(Node *node)
Returns true if the node is semi-local.
Definition: MooseMesh.C:556
unsigned int _var_num
variable number (from libMesh)
VariableValue _nodal_u_old_neighbor
virtual MooseMesh & mesh()
Definition: SystemBase.h:102
virtual void computeNeighborValuesFace()
Compute values at facial quadrature points for the neighbor.
VariableSecond _second_u_older_neighbor
virtual NumericVector< Number > & solutionOlder()=0
const VariablePhiGradient & gradPhiFace()
void computeIncrementAtNode(const NumericVector< Number > &increment_vec)
Compute and store incremental change at the current node based on increment_vec.
bool _need_solution_dofs
VariableValue _nodal_u_older
const Elem *& _neighbor
neighboring element
const VariablePhiGradient & _grad_phi_face_neighbor
const VariableValue & nodalValuePreviousNLNeighbor()
DenseVector< Number > _solution_dofs
local elemental DoFs
VariableSecond _second_u_bak
void getDofIndices(const Elem *elem, std::vector< dof_id_type > &dof_indices)
Get dof indices for the variable.
void add(NumericVector< Number > &residual)
const VariablePhiValue & phiFaceNeighbor()
MooseVariable(unsigned int var_num, const FEType &fe_type, SystemBase &sys, Assembly &assembly, Moose::VarKindType var_kind)
Definition: MooseVariable.C:28
VariableValue _u_dot
u_dot (time derivative)
std::vector< dof_id_type > _dof_indices_neighbor
DOF indices (neighbor)
void setNodalValueNeighbor(Number value)
Set the neighbor nodal value for this variable.
void insert(NumericVector< Number > &residual)
bool _need_solution_dofs_old_neighbor
bool _need_second_older_neighbor
virtual ~MooseVariable()
bool _need_u_older_neighbor
VariableSecond _second_u_neighbor
VariableGradient _grad_u_older_bak
const VariableValue & nodalValueDot()
VariableValue _du_dot_du_neighbor
Number getNodalValue(const Node &node)
Get the value of this variable at given node.
const VariableValue & nodalValueOlder()
bool _need_u_previous_nl
VariableGradient _grad_u_older_neighbor
VariableValue _u_old_bak
void reinitNodeNeighbor()
VariableValue _u_older
const VariablePhiGradient & gradPhiNeighbor()
VariableValue _nodal_u_dot_neighbor
virtual const NumericVector< Number > *& currentSolution()=0
The solution vector that is currently being operated on.
bool _need_nodal_u_neighbor