www.mooseframework.org
MooseVariableConstMonomial.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 
16 #include "SubProblem.h"
17 #include "SystemBase.h"
18 #include "Assembly.h"
19 
20 #include "libmesh/quadrature.h"
21 
23  const FEType & fe_type,
24  SystemBase & sys,
25  Assembly & assembly,
26  Moose::VarKindType var_kind)
27  : MooseVariable(var_num, fe_type, sys, assembly, var_kind)
28 {
29 }
30 
31 void
32 MooseVariableConstMonomial::computeElemValuesHelper(const unsigned & nqp, const Real & phi)
33 {
34  bool is_transient = _subproblem.isTransient();
35 
36  _u.resize(nqp);
37  _grad_u.resize(nqp);
38 
39  if (_need_second)
40  _second_u.resize(nqp);
41 
44 
47 
50 
51  if (is_transient)
52  {
53  _u_dot.resize(nqp);
54  _du_dot_du.resize(nqp);
55 
56  if (_need_u_old)
57  _u_old.resize(nqp);
58 
59  if (_need_u_older)
60  _u_older.resize(nqp);
61 
62  if (_need_grad_old)
63  _grad_u_old.resize(nqp);
64 
65  if (_need_grad_older)
66  _grad_u_older.resize(nqp);
67 
68  if (_need_second_old)
69  _second_u_old.resize(nqp);
70 
73  }
74 
75  if (_need_nodal_u)
76  _nodal_u.resize(1);
77 
80 
81  if (is_transient)
82  {
89  }
90 
92  _solution_dofs.resize(1);
93 
95  _solution_dofs_old.resize(1);
96 
98  _solution_dofs_older.resize(1);
99 
100  const dof_id_type & idx = _dof_indices[0];
101  const Real & soln = (*_sys.currentSolution())(idx);
102  Real soln_old;
103  Real soln_older;
104  Real soln_previous_nl;
105  Real u_dot;
106  const Real & du_dot_du = _sys.duDotDu();
107 
108  if (_need_nodal_u)
109  _nodal_u[0] = soln;
110 
113  soln_previous_nl = (*_sys.solutionPreviousNewton())(idx);
114 
116  _nodal_u_previous_nl[0] = soln_previous_nl;
117 
119  _solution_dofs(0) = soln;
120 
121  if (is_transient)
122  {
124  soln_old = _sys.solutionOld()(idx);
125 
127  soln_older = _sys.solutionOlder()(idx);
128 
129  if (_need_nodal_u_old)
130  _nodal_u_old[0] = soln_old;
131 
133  _nodal_u_older[0] = soln_older;
134 
135  u_dot = _sys.solutionUDot()(idx);
136 
137  if (_need_nodal_u_dot)
138  _nodal_u_dot[0] = u_dot;
139 
141  _solution_dofs_old(0) = soln_old;
142 
144  _solution_dofs_older(0) = soln_older;
145  }
146 
147  _u[0] = phi * soln;
148 
150  _u_previous_nl[0] = phi * soln_previous_nl;
151 
152  if (is_transient)
153  {
154  _u_dot[0] = phi * u_dot;
155  _du_dot_du[0] = du_dot_du;
156 
157  if (_need_u_old)
158  _u_old[0] = phi * soln_old;
159 
160  if (_need_u_older)
161  _u_older[0] = phi * soln_older;
162  }
163 
164  for (unsigned qp = 1; qp < nqp; ++qp)
165  {
166  _u[qp] = _u[0];
167 
170 
171  if (is_transient)
172  {
173  _u_dot[qp] = _u_dot[0];
174  _du_dot_du[qp] = _du_dot_du[0];
175 
176  if (_need_u_old)
177  _u_old[qp] = _u_old[0];
178 
179  if (_need_u_older)
180  _u_older[qp] = _u_older[qp];
181  }
182  }
183 }
184 
185 void
187 {
188  bool is_transient = _subproblem.isTransient();
189 
190  _u_neighbor.resize(nqp);
192 
195 
196  if (is_transient)
197  {
199  _u_old_neighbor.resize(nqp);
200 
203 
206 
209 
212 
215  }
216 
219  if (is_transient)
220  {
227  }
228 
230  _solution_dofs_neighbor.resize(1);
231 
233  _solution_dofs_old_neighbor.resize(1);
234 
237 
238  const dof_id_type & idx = _dof_indices_neighbor[0];
239  const Real & soln = (*_sys.currentSolution())(idx);
240  Real soln_old;
241  Real soln_older;
242  Real u_dot;
243 
245  _nodal_u_neighbor[0] = soln;
246 
248  _solution_dofs_neighbor(0) = soln;
249 
250  if (is_transient)
251  {
253  soln_old = _sys.solutionOld()(idx);
254 
256  soln_older = _sys.solutionOlder()(idx);
257 
259  _nodal_u_old_neighbor[0] = soln_old;
261  _nodal_u_older_neighbor[0] = soln_older;
262 
263  u_dot = _sys.solutionUDot()(idx);
264 
266  _nodal_u_dot_neighbor[0] = u_dot;
267 
269  _solution_dofs_old_neighbor(0) = soln_old;
270 
272  _solution_dofs_older_neighbor(0) = soln_older;
273  }
274 
275  _u_neighbor[0] = phi * soln;
276 
277  if (is_transient)
278  {
280  _u_old_neighbor[0] = phi * soln_old;
281 
283  _u_older_neighbor[0] = phi * soln_older;
284  }
285 
286  for (unsigned qp = 1; qp < nqp; ++qp)
287  {
288  _u_neighbor[qp] = _u_neighbor[0];
289 
290  if (is_transient)
291  {
294 
297  }
298  }
299 }
300 
301 void
303 {
304  if (_dof_indices.size() == 0)
305  return;
306 
307  computeElemValuesHelper(_qrule->n_points(), _phi[0][0]);
308 }
309 
310 void
312 {
313  if (_dof_indices.size() == 0)
314  return;
315 
316  computeElemValuesHelper(_qrule_face->n_points(), _phi_face[0][0]);
317 }
318 
319 void
321 {
322  if (_dof_indices_neighbor.size() == 0)
323  return;
324 
326 }
327 
328 void
330 {
331  if (_dof_indices_neighbor.size() == 0)
332  return;
333 
335 }
bool _need_grad_older_neighbor
bool _need_second_old_neighbor
bool _need_u_old_neighbor
bool _need_solution_dofs_older_neighbor
bool _need_second_older
bool _need_nodal_u_older
VariableValue _nodal_u_old
DenseVector< Number > _solution_dofs_older
Keeps track of stuff related to assembling.
Definition: Assembly.h:63
bool _need_nodal_u_previous_nl
Class for stuff related to variables.
Definition: MooseVariable.h:43
VariableSecond _second_u_older
virtual void computeElemValues() override
Compute values at interior quadrature points.
void resize(const unsigned int size)
Change the number of elements the array can store.
Definition: MooseArray.h:205
VariableSecond _second_u_old
virtual NumericVector< Number > & solutionOld()=0
bool _need_grad_previous_nl
VariableGradient _grad_u_neighbor
VariableGradient _grad_u_old_neighbor
DenseVector< Number > _solution_dofs_neighbor
QBase *& _qrule_face
Quadrature rule for the face.
bool _need_grad_old_neighbor
virtual bool isTransient() const =0
const VariablePhiValue & _phi_neighbor
VariableValue _nodal_u_older_neighbor
bool _need_solution_dofs_old
VariableValue _u_old_neighbor
VariableValue _du_dot_du
derivative of u_dot wrt u
VariableGradient _grad_u_older
bool _need_second_previous_nl
Base class for a system (of equations)
Definition: SystemBase.h:91
VariableValue _u_older_neighbor
VariableGradient _grad_u_old
virtual NumericVector< Number > * solutionPreviousNewton()=0
bool _need_nodal_u_older_neighbor
std::vector< dof_id_type > _dof_indices
DOF indices.
virtual void computeNeighborValuesFace() override
Compute values at facial quadrature points for the neighbor.
VariableValue _nodal_u
VariableSecond _second_u
VariableValue _nodal_u_neighbor
virtual Number & duDotDu()
Definition: SystemBase.h:156
VariableSecond _second_u_old_neighbor
MooseVariableConstMonomial(unsigned int var_num, const FEType &fe_type, SystemBase &sys, Assembly &assembly, Moose::VarKindType var_kind)
VariableValue _u
QBase *& _qrule
Quadrature rule for interior.
VariableGradient _grad_u
SubProblem & _subproblem
Problem this variable is part of.
virtual NumericVector< Number > & solutionUDot()
Definition: SystemBase.h:157
const VariablePhiValue & phi()
SystemBase & _sys
System this variable is part of.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:156
bool _need_nodal_u_old_neighbor
bool _need_nodal_u_old
VariableValue _nodal_u_previous_nl
VariableGradient _grad_u_previous_nl
const VariablePhiValue & _phi_face_neighbor
VariableValue _u_previous_nl
VariableValue _u_old
QBase *& _qrule_neighbor
Quadrature rule for the neighbor.
VariableValue _u_neighbor
const VariablePhiValue & _phi_face
bool _need_solution_dofs_neighbor
DenseVector< Number > _solution_dofs_old_neighbor
DenseVector< Number > _solution_dofs_old
bool _need_nodal_u_dot_neighbor
virtual void computeNeighborValues() override
Compute values at quadrature points for the neighbor.
const VariablePhiValue & _phi
virtual void computeElemValuesHelper(const unsigned &nqp, const Real &phi)
bool _need_second_neighbor
bool _need_nodal_u_dot
bool _need_solution_dofs_older
VariableValue _nodal_u_dot
nodal values of u_dot
DenseVector< Number > _solution_dofs_older_neighbor
VariableSecond _second_u_previous_nl
VariableValue _nodal_u_old_neighbor
VariableSecond _second_u_older_neighbor
virtual NumericVector< Number > & solutionOlder()=0
bool _need_solution_dofs
VariableValue _nodal_u_older
DenseVector< Number > _solution_dofs
local elemental DoFs
virtual void computeNeighborValuesHelper(const unsigned &nqp, const Real &phi)
virtual void computeElemValuesFace() override
Compute values at facial quadrature points.
VariableValue _u_dot
u_dot (time derivative)
std::vector< dof_id_type > _dof_indices_neighbor
DOF indices (neighbor)
bool _need_solution_dofs_old_neighbor
bool _need_second_older_neighbor
bool _need_u_older_neighbor
VariableSecond _second_u_neighbor
bool _need_u_previous_nl
VariableGradient _grad_u_older_neighbor
VariableValue _u_older
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