www.mooseframework.org
MooseVariable.h
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 #ifndef MOOSEVARIABLE_H
16 #define MOOSEVARIABLE_H
17 
18 #include "MooseTypes.h"
19 #include "MooseVariableBase.h"
20 
21 // Forward declarations
22 class Assembly;
23 class SubProblem;
24 class SystemBase;
25 
26 // libMesh forward declarations
27 namespace libMesh
28 {
29 class Elem;
30 class Node;
31 class QBase;
32 template <typename T>
33 class NumericVector;
34 }
35 
36 #include "libmesh/dense_vector.h"
37 
44 {
45 public:
46  MooseVariable(unsigned int var_num,
47  const FEType & fe_type,
48  SystemBase & sys,
49  Assembly & assembly,
50  Moose::VarKindType var_kind);
51  virtual ~MooseVariable();
52 
57  void clearDofIndices();
58 
59  void prepare();
60 
61  void prepareNeighbor();
62  void prepareAux();
63  void prepareIC();
64 
65  void reinitNode();
66  void reinitNodeNeighbor();
67  void reinitAux();
68  void reinitAuxNeighbor();
69 
70  void reinitNodes(const std::vector<dof_id_type> & nodes);
71  void reinitNodesNeighbor(const std::vector<dof_id_type> & nodes);
72 
73  const std::set<SubdomainID> & activeSubdomains() const;
74 
80  bool activeOnSubdomain(SubdomainID subdomain) const;
81 
86  virtual bool isNodal() const override;
87 
91  const Elem *& currentElem() { return _elem; }
92 
96  unsigned int & currentSide() { return _current_side; }
97 
101  const Elem *& neighbor() { return _neighbor; }
102 
107  {
108  return _need_second || _need_second_old || _need_second_older || _need_second_previous_nl;
109  }
110 
111  const VariablePhiValue & phi();
112  const VariablePhiGradient & gradPhi();
113  const VariablePhiSecond & secondPhi();
114 
115  const VariablePhiValue & phiFace();
116  const VariablePhiGradient & gradPhiFace();
117  const VariablePhiSecond & secondPhiFace();
118 
119  const VariablePhiValue & phiNeighbor();
120  const VariablePhiGradient & gradPhiNeighbor();
121  const VariablePhiSecond & secondPhiNeighbor();
122 
123  const VariablePhiValue & phiFaceNeighbor();
124  const VariablePhiGradient & gradPhiFaceNeighbor();
125  const VariablePhiSecond & secondPhiFaceNeighbor();
126 
127  const MooseArray<Point> & normals() { return _normals; }
128 
129  // damping
130  VariableValue & increment() { return _increment; }
131 
132  const VariableValue & sln() { return _u; }
134  {
135  _need_u_old = true;
136  return _u_old;
137  }
139  {
140  _need_u_older = true;
141  return _u_older;
142  }
144  {
145  _need_u_previous_nl = true;
146  return _u_previous_nl;
147  }
148  const VariableGradient & gradSln() { return _grad_u; }
150  {
151  _need_grad_old = true;
152  return _grad_u_old;
153  }
155  {
156  _need_grad_older = true;
157  return _grad_u_older;
158  }
160  {
161  _need_grad_previous_nl = true;
162  return _grad_u_previous_nl;
163  }
165  {
166  _need_second = true;
167  secondPhi();
168  secondPhiFace();
169  return _second_u;
170  }
172  {
173  _need_second_old = true;
174  secondPhi();
175  secondPhiFace();
176  return _second_u_old;
177  }
179  {
180  _need_second_older = true;
181  secondPhi();
182  secondPhiFace();
183  return _second_u_older;
184  }
186  {
187  _need_second_previous_nl = true;
188  secondPhi();
189  secondPhiFace();
190  return _second_u_previous_nl;
191  }
192 
193  const VariableValue & uDot() { return _u_dot; }
194  const VariableValue & duDotDu() { return _du_dot_du; }
195 
196  const Node *& node() { return _node; }
197  dof_id_type & nodalDofIndex() { return _nodal_dof_index; }
198  bool isNodalDefined() { return _is_defined; }
199  const VariableValue & nodalSln() { return _nodal_u; }
200  const VariableValue & nodalSlnOld() { return _nodal_u_old; }
201  const VariableValue & nodalSlnOlder() { return _nodal_u_older; }
203  {
204  _need_nodal_u_previous_nl = true;
205  return _nodal_u_previous_nl;
206  }
207  const VariableValue & nodalSlnDot() { return _nodal_u_dot; }
208  const VariableValue & nodalSlnDuDotDu() { return _nodal_du_dot_du; }
209 
210  const VariableValue & nodalValue();
211  const VariableValue & nodalValueOld();
212  const VariableValue & nodalValueOlder();
213  const VariableValue & nodalValuePreviousNL();
214  const VariableValue & nodalValueDot();
215 
216  const VariableValue & nodalValueNeighbor();
217  const VariableValue & nodalValueOldNeighbor();
218  const VariableValue & nodalValueOlderNeighbor();
219  const VariableValue & nodalValuePreviousNLNeighbor();
220  const VariableValue & nodalValueDotNeighbor();
221 
222  const VariableValue & slnNeighbor() { return _u_neighbor; }
224  {
225  _need_u_old_neighbor = true;
226  return _u_old_neighbor;
227  }
229  {
230  _need_u_older_neighbor = true;
231  return _u_older_neighbor;
232  }
234  {
235  _need_u_previous_nl_neighbor = true;
236  return _u_previous_nl_neighbor;
237  }
238  const VariableGradient & gradSlnNeighbor() { return _grad_u_neighbor; }
240  {
241  _need_grad_old_neighbor = true;
242  return _grad_u_old_neighbor;
243  }
245  {
246  _need_grad_older_neighbor = true;
247  return _grad_u_older_neighbor;
248  }
250  {
251  _need_grad_previous_nl_neighbor = true;
252  return _grad_u_previous_nl_neighbor;
253  }
255  {
256  _need_second_neighbor = true;
257  secondPhiFaceNeighbor();
258  return _second_u_neighbor;
259  }
261  {
262  _need_second_old_neighbor = true;
263  secondPhiFaceNeighbor();
264  return _second_u_old_neighbor;
265  }
267  {
268  _need_second_older_neighbor = true;
269  secondPhiFaceNeighbor();
270  return _second_u_older_neighbor;
271  }
273  {
274  _need_second_previous_nl_neighbor = true;
275  secondPhiFaceNeighbor();
276  return _second_u_previous_nl_neighbor;
277  }
278 
279  const VariableValue & uDotNeighbor() { return _u_dot_neighbor; }
280  const VariableValue & duDotDuNeighbor() { return _du_dot_du_neighbor; }
281 
282  const Node *& nodeNeighbor() { return _node_neighbor; }
283  dof_id_type & nodalDofIndexNeighbor() { return _nodal_dof_index_neighbor; }
284  bool isNodalNeighborDefined() { return _is_defined_neighbor; }
285  const VariableValue & nodalSlnNeighbor() { return _nodal_u_neighbor; }
286  const VariableValue & nodalSlnOldNeighbor() { return _nodal_u_old_neighbor; }
287  const VariableValue & nodalSlnOlderNeighbor() { return _nodal_u_older_neighbor; }
289  {
290  _need_nodal_u_previous_nl_neighbor = true;
291  return _nodal_u_previous_nl_neighbor;
292  }
293  const VariableValue & nodalSlnDotNeighbor() { return _nodal_u_dot_neighbor; }
294  const VariableValue & nodalSlnDuDotDuNeighbor() { return _nodal_du_dot_du_neighbor; }
295  const DenseVector<Number> & solutionDoFs()
296  {
297  _need_solution_dofs = true;
298  return _solution_dofs;
299  }
300  const DenseVector<Number> & solutionDoFsOld()
301  {
302  _need_solution_dofs_old = true;
303  return _solution_dofs_old;
304  }
305  const DenseVector<Number> & solutionDoFsOlder()
306  {
307  _need_solution_dofs_older = true;
308  return _solution_dofs_older;
309  }
310  const DenseVector<Number> & solutionDoFsNeighbor()
311  {
312  _need_solution_dofs_neighbor = true;
313  return _solution_dofs_neighbor;
314  }
315  const DenseVector<Number> & solutionDoFsOldNeighbor()
316  {
317  _need_solution_dofs_old_neighbor = true;
318  return _solution_dofs_old_neighbor;
319  }
320  const DenseVector<Number> & solutionDoFsOlderNeighbor()
321  {
322  _need_solution_dofs_older_neighbor = true;
323  return _solution_dofs_older_neighbor;
324  }
325 
333  void computePerturbedElemValues(unsigned i, Real scale, Real & h);
334 
339  void restoreUnperturbedElemValues();
340 
344  virtual void computeElemValues();
348  virtual void computeElemValuesFace();
352  virtual void computeNeighborValuesFace();
356  virtual void computeNeighborValues();
360  void computeNodalValues();
364  void computeNodalNeighborValues();
368  void setNodalValue(Number value, unsigned int idx = 0);
372  void setNodalValue(const DenseVector<Number> & value);
373 
377  void setNodalValueNeighbor(Number value);
381  void setNodalValueNeighbor(const DenseVector<Number> & value);
382 
386  void computeIncrementAtQps(const NumericVector<Number> & increment_vec);
387 
391  void computeIncrementAtNode(const NumericVector<Number> & increment_vec);
392 
397  std::vector<dof_id_type> & dofIndicesNeighbor() { return _dof_indices_neighbor; }
398 
399  unsigned int numberOfDofsNeighbor() { return _dof_indices_neighbor.size(); }
400 
401  void insert(NumericVector<Number> & residual);
402  void add(NumericVector<Number> & residual);
403 
407  Number getNodalValue(const Node & node);
411  Number getNodalValueOld(const Node & node);
415  Number getNodalValueOlder(const Node & node);
416 
423  Real getValue(const Elem * elem, const std::vector<std::vector<Real>> & phi) const;
424  RealGradient getGradient(const Elem * elem,
425  const std::vector<std::vector<RealGradient>> & phi) const;
426 
432  Number getElementalValue(const Elem * elem, unsigned int idx = 0) const;
433 
439  bool usesPhi() { return true; }
440 
446  bool usesGradPhi() { return true; }
447 
451  bool usesSecondPhi() { return _need_second || _need_second_old || _need_second_older; }
452 
453 protected:
459  void getDofIndices(const Elem * elem, std::vector<dof_id_type> & dof_indices);
460 
461 protected:
464 
466  QBase *& _qrule;
468  QBase *& _qrule_face;
470  QBase *& _qrule_neighbor;
471 
473  const Elem *& _elem;
475  unsigned int & _current_side;
476 
478  const Elem *& _neighbor;
479 
481  std::vector<dof_id_type> _dof_indices_neighbor;
482 
486 
490 
495 
499 
503 
508 
514 
520 
527 
528  // Shape function values, gradients. second derivatives
532 
533  // Values, gradients and second derivatives of shape function on faces
537 
538  // Values, gradients and second derivatives of shape function
542 
543  // Values, gradients and second derivatives of shape function on faces
547 
550 
563 
576 
577  // time derivatives
578 
581  VariableValue _u_dot_neighbor, _u_dot_bak_neighbor;
582 
585  VariableValue _du_dot_du_neighbor, _du_dot_du_bak_neighbor;
586 
587  // nodal stuff
588 
594  const Node *& _node;
595  dof_id_type _nodal_dof_index;
600 
605 
608  const Node *& _node_neighbor;
616 
618  DenseVector<Number> _solution_dofs;
619  DenseVector<Number> _solution_dofs_old;
620  DenseVector<Number> _solution_dofs_older;
621  DenseVector<Number> _solution_dofs_neighbor;
622  DenseVector<Number> _solution_dofs_old_neighbor;
623  DenseVector<Number> _solution_dofs_older_neighbor;
624 
626  bool _is_nodal;
627 
628  // damping
630 
631  friend class NodeFaceConstraint;
632  friend class ValueThresholdMarker;
633  friend class ValueRangeMarker;
634 };
635 
636 #endif /* MOOSEVARIABLE_H */
const VariableGradient & gradSlnNeighbor()
bool _need_grad_older_neighbor
VariableValue _increment
const VariablePhiGradient & _grad_phi
const VariableSecond & secondSlnOldNeighbor()
bool _need_grad_previous_nl_neighbor
RealVectorValue RealGradient
Definition: Assembly.h:43
bool usesPhi()
Whether or not this variable is actually using the shape function value.
const VariableSecond & secondSlnPreviousNL()
bool _need_second_old_neighbor
bool _need_u_old_neighbor
VariableSecond _second_u_previous_nl_neighbor
bool _need_solution_dofs_older_neighbor
const VariableGradient & gradSlnOld()
const VariablePhiSecond * _second_phi_neighbor
const DenseVector< Number > & solutionDoFsOld()
bool _need_second_older
VariableValue _nodal_u_old
bool _need_nodal_u_older
dof_id_type _nodal_dof_index_neighbor
DenseVector< Number > _solution_dofs_older
VariableValue _u_older_bak
Keeps track of stuff related to assembling.
Definition: Assembly.h:63
bool _need_nodal_u_previous_nl
VariableGradient _grad_u_old_bak
Class for stuff related to variables.
Definition: MooseVariable.h:43
const VariableValue & slnOlder()
const VariablePhiSecond * _second_phi
subdomain_id_type SubdomainID
Definition: MooseTypes.h:77
const DenseVector< Number > & solutionDoFsNeighbor()
bool _need_grad_previous_nl
VariableGradient _grad_u_neighbor
bool _need_u_previous_nl_neighbor
VariableGradient _grad_u_old_neighbor
const VariableSecond & secondSln()
const VariableGradient & gradSlnOldNeighbor()
bool _is_defined_neighbor
If the variable is defined at the neighbor node (used in compute nodal values)
const VariableValue & nodalSlnPreviousNL()
DenseVector< Number > _solution_dofs_neighbor
const VariableValue & uDotNeighbor()
bool _is_defined
If the variable is defined at the node (used in compute nodal values)
const VariableValue & slnPreviousNLNeighbor()
bool _need_grad_old_neighbor
QBase *& _qrule_face
Quadrature rule for the face.
const VariablePhiValue & _phi_neighbor
bool _has_nodal_value
If true, the nodal value gets inserted on calling insert()
const VariableValue & nodalSlnNeighbor()
const Node *& _node
VariableValue _nodal_u_older_neighbor
bool _need_solution_dofs_old
const VariablePhiGradient & _grad_phi_face
VariableValue _u_old_neighbor
VariableValue _u_previous_nl_neighbor
unsigned int & _current_side
the side of the current element (valid when doing face assembly)
dof_id_type & nodalDofIndexNeighbor()
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
bool _need_second_previous_nl
const VariableGradient & gradSlnPreviousNLNeighbor()
const VariableSecond & secondSlnOld()
bool usesSecondPhi()
Whether or not this variable is actually using the shape function second derivative.
Base class for a system (of equations)
Definition: SystemBase.h:91
VariableValue _u_older_neighbor
const VariableSecond & secondSlnOlder()
const VariableValue & nodalSlnDuDotDuNeighbor()
const DenseVector< Number > & solutionDoFs()
const VariablePhiSecond * _second_phi_face_neighbor
bool _need_nodal_u_older_neighbor
const VariableValue & nodalSlnDot()
const Elem *& currentElem()
Current element this variable is evaluated at.
Definition: MooseVariable.h:91
VariableValue _nodal_u
VariableValue _nodal_u_neighbor
dof_id_type & nodalDofIndex()
VariableValue _nodal_du_dot_du
nodal values of derivative of u_dot wrt u
const DenseVector< Number > & solutionDoFsOlder()
VariableSecond _second_u_old_neighbor
const VariableValue & nodalSlnOldNeighbor()
const VariableValue & nodalSlnPreviousNLNeighbor()
QBase *& _qrule
Quadrature rule for interior.
const VariableValue & slnOld()
VariableValue _u_dot_neighbor
const VariablePhiSecond * _second_phi_face
const VariablePhiGradient & _grad_phi_neighbor
const MooseArray< Point > & normals()
const Node *& _node_neighbor
bool _need_nodal_u_previous_nl_neighbor
A NodeFaceConstraint is used when you need to create constraints between two surfaces in a mesh...
dof_id_type _nodal_dof_index
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:156
const Node *& node()
const MooseArray< Point > & _normals
Normals at QPs on faces.
std::vector< dof_id_type > & dofIndicesNeighbor()
Get DOF indices for currently selected element.
const VariableValue & nodalSlnOlderNeighbor()
bool _need_nodal_u_old_neighbor
const VariableGradient & gradSlnPreviousNL()
bool _need_nodal_u_old
VariableValue _nodal_u_previous_nl
const VariableValue & nodalSlnOld()
VariableGradient _grad_u_previous_nl
unsigned int & currentSide()
Current side this variable is being evaluated on.
Definition: MooseVariable.h:96
bool _has_nodal_value_neighbor
const VariablePhiValue & _phi_face_neighbor
const VariableValue & slnOlderNeighbor()
VariableValue _u_previous_nl
VariableValue _nodal_u_previous_nl_neighbor
QBase *& _qrule_neighbor
Quadrature rule for the neighbor.
const VariableValue & duDotDu()
VariableValue _u_neighbor
VariableSecond _second_u_old_bak
const VariableSecond & secondSlnOlderNeighbor()
const VariablePhiValue & _phi_face
VariableValue & increment()
bool _need_solution_dofs_neighbor
VariableValue _u_bak
DenseVector< Number > _solution_dofs_old_neighbor
const VariableValue & duDotDuNeighbor()
const VariableValue & nodalSlnOlder()
DenseVector< Number > _solution_dofs_old
bool _need_nodal_u_dot_neighbor
const VariableGradient & gradSlnOlderNeighbor()
const VariablePhiValue & _phi
const VariableValue & slnOldNeighbor()
bool _need_second_neighbor
bool _need_nodal_u_dot
bool _is_nodal
if variable is nodal
const VariableValue & nodalSln()
bool _need_solution_dofs_older
VariableValue _nodal_du_dot_du_neighbor
VariableValue _nodal_u_dot
nodal values of u_dot
const VariableValue & nodalSlnDuDotDu()
bool computingSecond()
Whether or not this variable is computing any second derivatives.
VariableValue _u_dot_bak
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
VariableSecond _second_u_older_bak
bool _need_second_previous_nl_neighbor
DenseVector< Number > _solution_dofs_older_neighbor
const Elem *& _elem
current element
VariableGradient _grad_u_bak
VariableSecond _second_u_previous_nl
VariableValue _du_dot_du_bak
VariableGradient _grad_u_previous_nl_neighbor
const VariableValue & sln()
const VariableValue & uDot()
VariableValue _nodal_u_old_neighbor
VariableSecond _second_u_older_neighbor
const DenseVector< Number > & solutionDoFsOlderNeighbor()
THREAD_ID _tid
Thread ID.
const VariableValue & nodalSlnDotNeighbor()
const VariableGradient & gradSlnOlder()
bool _need_solution_dofs
VariableValue _nodal_u_older
const Elem *& _neighbor
neighboring element
const VariablePhiGradient & _grad_phi_face_neighbor
DenseVector< Number > _solution_dofs
local elemental DoFs
const DenseVector< Number > & solutionDoFsOldNeighbor()
bool usesGradPhi()
Whether or not this variable is actually using the shape function gradient.
VariableSecond _second_u_bak
const VariableValue & slnNeighbor()
const VariableSecond & secondSlnPreviousNLNeighbor()
const Elem *& neighbor()
Current neighboring element.
const VariableSecond & secondSlnNeighbor()
std::vector< dof_id_type > _dof_indices_neighbor
DOF indices (neighbor)
const Node *& nodeNeighbor()
bool _need_solution_dofs_old_neighbor
bool _need_second_older_neighbor
bool _need_u_older_neighbor
VariableSecond _second_u_neighbor
const VariableGradient & gradSln()
VariableGradient _grad_u_older_bak
VariableValue _du_dot_du_neighbor
bool _need_u_previous_nl
const VariableValue & slnPreviousNL()
VariableGradient _grad_u_older_neighbor
VariableValue _u_old_bak
unsigned int THREAD_ID
Definition: MooseTypes.h:79
bool isNodalNeighborDefined()
unsigned int numberOfDofsNeighbor()
bool isNodalDefined()
VariableValue _nodal_u_dot_neighbor
bool _need_nodal_u_neighbor