www.mooseframework.org
InitialCondition.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 "InitialCondition.h"
16 #include "FEProblem.h"
17 #include "SystemBase.h"
18 #include "Assembly.h"
19 #include "MooseVariable.h"
20 
21 #include "libmesh/fe_interface.h"
22 #include "libmesh/quadrature.h"
23 
24 template <>
27 {
31 
32  params.addRequiredParam<VariableName>("variable",
33  "The variable this initial condition is "
34  "supposed to provide values for.");
35  params.addParam<bool>("ignore_uo_dependency",
36  false,
37  "When set to true, a UserObject retrieved "
38  "by this IC will not be executed before the "
39  "this IC");
40 
41  params.addParamNamesToGroup("ignore_uo_dependency", "Advanced");
42 
43  params.registerBase("InitialCondition");
44 
45  return params;
46 }
47 
49  : MooseObject(parameters),
50  BlockRestrictable(this),
51  Coupleable(this,
52  getParam<SystemBase *>("_sys")
53  ->getVariable(parameters.get<THREAD_ID>("_tid"),
54  parameters.get<VariableName>("variable"))
55  .isNodal()),
56  FunctionInterface(this),
57  UserObjectInterface(this),
58  BoundaryRestrictable(this, _c_nodal),
60  Restartable(parameters, "InitialConditions"),
61  ZeroInterface(parameters),
62  _fe_problem(*parameters.getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
63  _sys(*parameters.getCheckedPointerParam<SystemBase *>("_sys")),
64  _tid(getParam<THREAD_ID>("_tid")),
65  _assembly(_fe_problem.assembly(_tid)),
66  _t(_fe_problem.time()),
67  _coord_sys(_assembly.coordSystem()),
68  _var(_sys.getVariable(_tid, getParam<VariableName>("variable"))),
69 
70  _current_elem(_var.currentElem()),
71  _current_node(NULL),
72  _qp(0),
73  _ignore_uo_dependency(getParam<bool>("ignore_uo_dependency"))
74 {
75  _supplied_vars.insert(getParam<VariableName>("variable"));
76 
77  std::map<std::string, std::vector<MooseVariable *>> coupled_vars = getCoupledVars();
78  for (const auto & it : coupled_vars)
79  for (const auto & var : it.second)
80  _depend_vars.insert(var->name());
81 }
82 
84 
85 const std::set<std::string> &
87 {
88  return _depend_vars;
89 }
90 
91 const std::set<std::string> &
93 {
94  return _supplied_vars;
95 }
96 
97 const UserObject &
99 {
101  _depend_uo.insert(_pars.get<UserObjectName>(name));
103 }
104 
105 void
107 {
108  // -- NOTE ----
109  // The following code is a copy from libMesh project_vector.C plus it adds some features, so we
110  // can couple variable values
111  // and we also do not call any callbacks, but we use our initial condition system directly.
112  // ------------
113 
114  // The element matrix and RHS for projections.
115  // Note that Ke is always real-valued, whereas Fe may be complex valued if complex number support
116  // is enabled
117  DenseMatrix<Real> Ke;
118  DenseVector<Number> Fe;
119  // The new element coefficients
120  DenseVector<Number> Ue;
121 
122  const FEType & fe_type = _var.feType();
123 
124  // The dimension of the current element
125  const unsigned int dim = _current_elem->dim();
126  // The element type
127  const ElemType elem_type = _current_elem->type();
128  // The number of nodes on the new element
129  const unsigned int n_nodes = _current_elem->n_nodes();
130  // The global DOF indices
131  std::vector<dof_id_type> dof_indices;
132  // Side/edge DOF indices
133  std::vector<unsigned int> side_dofs;
134 
135  // Get FE objects of the appropriate type
136  // We cannot use the FE object in Assembly, since the following code is messing with the
137  // quadrature rules
138  // for projections and would screw it up. However, if we implement projections from one mesh to
139  // another,
140  // this code should use that implementation.
141  std::unique_ptr<FEBase> fe(FEBase::build(dim, fe_type));
142 
143  // Prepare variables for projection
144  std::unique_ptr<QBase> qrule(fe_type.default_quadrature_rule(dim));
145  std::unique_ptr<QBase> qedgerule(fe_type.default_quadrature_rule(1));
146  std::unique_ptr<QBase> qsiderule(fe_type.default_quadrature_rule(dim - 1));
147 
148  // The values of the shape functions at the quadrature points
149  const std::vector<std::vector<Real>> & phi = fe->get_phi();
150 
151  // The gradients of the shape functions at the quadrature points on the child element.
152  const std::vector<std::vector<RealGradient>> * dphi = NULL;
153 
154  const FEContinuity cont = fe->get_continuity();
155 
156  if (cont == C_ONE)
157  {
158  const std::vector<std::vector<RealGradient>> & ref_dphi = fe->get_dphi();
159  dphi = &ref_dphi;
160  }
161 
162  // The Jacobian * quadrature weight at the quadrature points
163  const std::vector<Real> & JxW = fe->get_JxW();
164  // The XYZ locations of the quadrature points
165  const std::vector<Point> & xyz_values = fe->get_xyz();
166 
167  // Update the DOF indices for this element based on the current mesh
168  _var.prepareIC();
169  dof_indices = _var.dofIndices();
170 
171  // The number of DOFs on the element
172  const unsigned int n_dofs = dof_indices.size();
173  if (n_dofs == 0)
174  return;
175 
176  // Fixed vs. free DoFs on edge/face projections
177  std::vector<char> dof_is_fixed(n_dofs, false); // bools
178  std::vector<int> free_dof(n_dofs, 0);
179 
180  // Zero the interpolated values
181  Ue.resize(n_dofs);
182  Ue.zero();
183 
184  // In general, we need a series of
185  // projections to ensure a unique and continuous
186  // solution. We start by interpolating nodes, then
187  // hold those fixed and project edges, then
188  // hold those fixed and project faces, then
189  // hold those fixed and project interiors
190 
191  // Interpolate node values first
192  unsigned int current_dof = 0;
193  for (unsigned int n = 0; n != n_nodes; ++n)
194  {
195  // FIXME: this should go through the DofMap,
196  // not duplicate dof_indices code badly!
197  const unsigned int nc = FEInterface::n_dofs_at_node(dim, fe_type, elem_type, n);
198  if (!_current_elem->is_vertex(n))
199  {
200  current_dof += nc;
201  continue;
202  }
203  if (cont == DISCONTINUOUS)
204  {
205  libmesh_assert(nc == 0);
206  }
207  // Assume that C_ZERO elements have a single nodal
208  // value shape function
209  else if (cont == C_ZERO)
210  {
211  libmesh_assert(nc == 1);
212  _qp = n;
213  _current_node = _current_elem->node_ptr(n);
214  Ue(current_dof) = value(*_current_node);
215  dof_is_fixed[current_dof] = true;
216  current_dof++;
217  }
218  // The hermite element vertex shape functions are weird
219  else if (fe_type.family == HERMITE)
220  {
221  _qp = n;
222  _current_node = _current_elem->node_ptr(n);
223  Ue(current_dof) = value(*_current_node);
224  dof_is_fixed[current_dof] = true;
225  current_dof++;
226  Gradient grad = gradient(*_current_node);
227  // x derivative
228  Ue(current_dof) = grad(0);
229  dof_is_fixed[current_dof] = true;
230  current_dof++;
231  if (dim > 1)
232  {
233  // We'll finite difference mixed derivatives
234  Point nxminus = _current_elem->point(n), nxplus = _current_elem->point(n);
235  nxminus(0) -= TOLERANCE;
236  nxplus(0) += TOLERANCE;
237  Gradient gxminus = gradient(nxminus);
238  Gradient gxplus = gradient(nxplus);
239  // y derivative
240  Ue(current_dof) = grad(1);
241  dof_is_fixed[current_dof] = true;
242  current_dof++;
243  // xy derivative
244  Ue(current_dof) = (gxplus(1) - gxminus(1)) / 2. / TOLERANCE;
245  dof_is_fixed[current_dof] = true;
246  current_dof++;
247 
248  if (dim > 2)
249  {
250  // z derivative
251  Ue(current_dof) = grad(2);
252  dof_is_fixed[current_dof] = true;
253  current_dof++;
254  // xz derivative
255  Ue(current_dof) = (gxplus(2) - gxminus(2)) / 2. / TOLERANCE;
256  dof_is_fixed[current_dof] = true;
257  current_dof++;
258  // We need new points for yz
259  Point nyminus = _current_elem->point(n), nyplus = _current_elem->point(n);
260  nyminus(1) -= TOLERANCE;
261  nyplus(1) += TOLERANCE;
262  Gradient gyminus = gradient(nyminus);
263  Gradient gyplus = gradient(nyplus);
264  // xz derivative
265  Ue(current_dof) = (gyplus(2) - gyminus(2)) / 2. / TOLERANCE;
266  dof_is_fixed[current_dof] = true;
267  current_dof++;
268  // Getting a 2nd order xyz is more tedious
269  Point nxmym = _current_elem->point(n), nxmyp = _current_elem->point(n),
270  nxpym = _current_elem->point(n), nxpyp = _current_elem->point(n);
271  nxmym(0) -= TOLERANCE;
272  nxmym(1) -= TOLERANCE;
273  nxmyp(0) -= TOLERANCE;
274  nxmyp(1) += TOLERANCE;
275  nxpym(0) += TOLERANCE;
276  nxpym(1) -= TOLERANCE;
277  nxpyp(0) += TOLERANCE;
278  nxpyp(1) += TOLERANCE;
279  Gradient gxmym = gradient(nxmym);
280  Gradient gxmyp = gradient(nxmyp);
281  Gradient gxpym = gradient(nxpym);
282  Gradient gxpyp = gradient(nxpyp);
283  Number gxzplus = (gxpyp(2) - gxmyp(2)) / 2. / TOLERANCE;
284  Number gxzminus = (gxpym(2) - gxmym(2)) / 2. / TOLERANCE;
285  // xyz derivative
286  Ue(current_dof) = (gxzplus - gxzminus) / 2. / TOLERANCE;
287  dof_is_fixed[current_dof] = true;
288  current_dof++;
289  }
290  }
291  }
292  // Assume that other C_ONE elements have a single nodal
293  // value shape function and nodal gradient component
294  // shape functions
295  else if (cont == C_ONE)
296  {
297  libmesh_assert(nc == 1 + dim);
298  _current_node = _current_elem->node_ptr(n);
299  Ue(current_dof) = value(*_current_node);
300  dof_is_fixed[current_dof] = true;
301  current_dof++;
302  Gradient grad = gradient(*_current_node);
303  for (unsigned int i = 0; i != dim; ++i)
304  {
305  Ue(current_dof) = grad(i);
306  dof_is_fixed[current_dof] = true;
307  current_dof++;
308  }
309  }
310  else
311  libmesh_error();
312  } // loop over nodes
313 
314  // From here on out we won't be sampling at nodes anymore
315  _current_node = NULL;
316 
317  // In 3D, project any edge values next
318  if (dim > 2 && cont != DISCONTINUOUS)
319  for (unsigned int e = 0; e != _current_elem->n_edges(); ++e)
320  {
321  FEInterface::dofs_on_edge(_current_elem, dim, fe_type, e, side_dofs);
322 
323  // Some edge dofs are on nodes and already
324  // fixed, others are free to calculate
325  unsigned int free_dofs = 0;
326  for (unsigned int i = 0; i != side_dofs.size(); ++i)
327  if (!dof_is_fixed[side_dofs[i]])
328  free_dof[free_dofs++] = i;
329 
330  // There may be nothing to project
331  if (!free_dofs)
332  continue;
333 
334  Ke.resize(free_dofs, free_dofs);
335  Ke.zero();
336  Fe.resize(free_dofs);
337  Fe.zero();
338  // The new edge coefficients
339  DenseVector<Number> Uedge(free_dofs);
340 
341  // Initialize FE data on the edge
342  fe->attach_quadrature_rule(qedgerule.get());
343  fe->edge_reinit(_current_elem, e);
344  const unsigned int n_qp = qedgerule->n_points();
345 
346  // Loop over the quadrature points
347  for (unsigned int qp = 0; qp < n_qp; qp++)
348  {
349  // solution at the quadrature point
350  Number fineval = value(xyz_values[qp]);
351  // solution grad at the quadrature point
352  Gradient finegrad;
353  if (cont == C_ONE)
354  finegrad = gradient(xyz_values[qp]);
355 
356  // Form edge projection matrix
357  for (unsigned int sidei = 0, freei = 0; sidei != side_dofs.size(); ++sidei)
358  {
359  unsigned int i = side_dofs[sidei];
360  // fixed DoFs aren't test functions
361  if (dof_is_fixed[i])
362  continue;
363  for (unsigned int sidej = 0, freej = 0; sidej != side_dofs.size(); ++sidej)
364  {
365  unsigned int j = side_dofs[sidej];
366  if (dof_is_fixed[j])
367  Fe(freei) -= phi[i][qp] * phi[j][qp] * JxW[qp] * Ue(j);
368  else
369  Ke(freei, freej) += phi[i][qp] * phi[j][qp] * JxW[qp];
370  if (cont == C_ONE)
371  {
372  if (dof_is_fixed[j])
373  Fe(freei) -= ((*dphi)[i][qp] * (*dphi)[j][qp]) * JxW[qp] * Ue(j);
374  else
375  Ke(freei, freej) += ((*dphi)[i][qp] * (*dphi)[j][qp]) * JxW[qp];
376  }
377  if (!dof_is_fixed[j])
378  freej++;
379  }
380  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
381  if (cont == C_ONE)
382  Fe(freei) += (finegrad * (*dphi)[i][qp]) * JxW[qp];
383  freei++;
384  }
385  }
386 
387  Ke.cholesky_solve(Fe, Uedge);
388 
389  // Transfer new edge solutions to element
390  for (unsigned int i = 0; i != free_dofs; ++i)
391  {
392  Number & ui = Ue(side_dofs[free_dof[i]]);
393  libmesh_assert(std::abs(ui) < TOLERANCE || std::abs(ui - Uedge(i)) < TOLERANCE);
394  ui = Uedge(i);
395  dof_is_fixed[side_dofs[free_dof[i]]] = true;
396  }
397  }
398 
399  // Project any side values (edges in 2D, faces in 3D)
400  if (dim > 1 && cont != DISCONTINUOUS)
401  for (unsigned int s = 0; s != _current_elem->n_sides(); ++s)
402  {
403  FEInterface::dofs_on_side(_current_elem, dim, fe_type, s, side_dofs);
404 
405  // Some side dofs are on nodes/edges and already
406  // fixed, others are free to calculate
407  unsigned int free_dofs = 0;
408  for (unsigned int i = 0; i != side_dofs.size(); ++i)
409  if (!dof_is_fixed[side_dofs[i]])
410  free_dof[free_dofs++] = i;
411 
412  // There may be nothing to project
413  if (!free_dofs)
414  continue;
415 
416  Ke.resize(free_dofs, free_dofs);
417  Ke.zero();
418  Fe.resize(free_dofs);
419  Fe.zero();
420  // The new side coefficients
421  DenseVector<Number> Uside(free_dofs);
422 
423  // Initialize FE data on the side
424  fe->attach_quadrature_rule(qsiderule.get());
425  fe->reinit(_current_elem, s);
426  const unsigned int n_qp = qsiderule->n_points();
427 
428  // Loop over the quadrature points
429  for (unsigned int qp = 0; qp < n_qp; qp++)
430  {
431  // solution at the quadrature point
432  Number fineval = value(xyz_values[qp]);
433  // solution grad at the quadrature point
434  Gradient finegrad;
435  if (cont == C_ONE)
436  finegrad = gradient(xyz_values[qp]);
437 
438  // Form side projection matrix
439  for (unsigned int sidei = 0, freei = 0; sidei != side_dofs.size(); ++sidei)
440  {
441  unsigned int i = side_dofs[sidei];
442  // fixed DoFs aren't test functions
443  if (dof_is_fixed[i])
444  continue;
445  for (unsigned int sidej = 0, freej = 0; sidej != side_dofs.size(); ++sidej)
446  {
447  unsigned int j = side_dofs[sidej];
448  if (dof_is_fixed[j])
449  Fe(freei) -= phi[i][qp] * phi[j][qp] * JxW[qp] * Ue(j);
450  else
451  Ke(freei, freej) += phi[i][qp] * phi[j][qp] * JxW[qp];
452  if (cont == C_ONE)
453  {
454  if (dof_is_fixed[j])
455  Fe(freei) -= ((*dphi)[i][qp] * (*dphi)[j][qp]) * JxW[qp] * Ue(j);
456  else
457  Ke(freei, freej) += ((*dphi)[i][qp] * (*dphi)[j][qp]) * JxW[qp];
458  }
459  if (!dof_is_fixed[j])
460  freej++;
461  }
462  Fe(freei) += (fineval * phi[i][qp]) * JxW[qp];
463  if (cont == C_ONE)
464  Fe(freei) += (finegrad * (*dphi)[i][qp]) * JxW[qp];
465  freei++;
466  }
467  }
468 
469  Ke.cholesky_solve(Fe, Uside);
470 
471  // Transfer new side solutions to element
472  for (unsigned int i = 0; i != free_dofs; ++i)
473  {
474  Number & ui = Ue(side_dofs[free_dof[i]]);
475  libmesh_assert(std::abs(ui) < TOLERANCE || std::abs(ui - Uside(i)) < TOLERANCE);
476  ui = Uside(i);
477  dof_is_fixed[side_dofs[free_dof[i]]] = true;
478  }
479  }
480 
481  // Project the interior values, finally
482 
483  // Some interior dofs are on nodes/edges/sides and
484  // already fixed, others are free to calculate
485  unsigned int free_dofs = 0;
486  for (unsigned int i = 0; i != n_dofs; ++i)
487  if (!dof_is_fixed[i])
488  free_dof[free_dofs++] = i;
489 
490  // There may be nothing to project
491  if (free_dofs)
492  {
493  Ke.resize(free_dofs, free_dofs);
494  Ke.zero();
495  Fe.resize(free_dofs);
496  Fe.zero();
497  // The new interior coefficients
498  DenseVector<Number> Uint(free_dofs);
499 
500  // Initialize FE data
501  fe->attach_quadrature_rule(qrule.get());
502  fe->reinit(_current_elem);
503  const unsigned int n_qp = qrule->n_points();
504 
505  // Loop over the quadrature points
506  for (unsigned int qp = 0; qp < n_qp; qp++)
507  {
508  // solution at the quadrature point
509  Number fineval = value(xyz_values[qp]);
510  // solution grad at the quadrature point
511  Gradient finegrad;
512  if (cont == C_ONE)
513  finegrad = gradient(xyz_values[qp]);
514 
515  // Form interior projection matrix
516  for (unsigned int i = 0, freei = 0; i != n_dofs; ++i)
517  {
518  // fixed DoFs aren't test functions
519  if (dof_is_fixed[i])
520  continue;
521  for (unsigned int j = 0, freej = 0; j != n_dofs; ++j)
522  {
523  if (dof_is_fixed[j])
524  Fe(freei) -= phi[i][qp] * phi[j][qp] * JxW[qp] * Ue(j);
525  else
526  Ke(freei, freej) += phi[i][qp] * phi[j][qp] * JxW[qp];
527  if (cont == C_ONE)
528  {
529  if (dof_is_fixed[j])
530  Fe(freei) -= ((*dphi)[i][qp] * (*dphi)[j][qp]) * JxW[qp] * Ue(j);
531  else
532  Ke(freei, freej) += ((*dphi)[i][qp] * (*dphi)[j][qp]) * JxW[qp];
533  }
534  if (!dof_is_fixed[j])
535  freej++;
536  }
537  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
538  if (cont == C_ONE)
539  Fe(freei) += (finegrad * (*dphi)[i][qp]) * JxW[qp];
540  freei++;
541  }
542  }
543  Ke.cholesky_solve(Fe, Uint);
544 
545  // Transfer new interior solutions to element
546  for (unsigned int i = 0; i != free_dofs; ++i)
547  {
548  Number & ui = Ue(free_dof[i]);
549  libmesh_assert(std::abs(ui) < TOLERANCE || std::abs(ui - Uint(i)) < TOLERANCE);
550  ui = Uint(i);
551  dof_is_fixed[free_dof[i]] = true;
552  }
553  } // if there are free interior dofs
554 
555  // Make sure every DoF got reached!
556  for (unsigned int i = 0; i != n_dofs; ++i)
557  libmesh_assert(dof_is_fixed[i]);
558 
559  NumericVector<Number> & solution = _var.sys().solution();
560 
561  // 'first' and 'last' are no longer used, see note about subdomain-restricted variables below
562  // const dof_id_type
563  // first = solution.first_local_index(),
564  // last = solution.last_local_index();
565 
566  // Lock the new_vector since it is shared among threads.
567  {
568  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
569 
570  for (unsigned int i = 0; i < n_dofs; i++)
571  // We may be projecting a new zero value onto
572  // an old nonzero approximation - RHS
573  // if (Ue(i) != 0.)
574 
575  // This is commented out because of subdomain restricted variables.
576  // It can be the case that if a subdomain restricted variable's boundary
577  // aligns perfectly with a processor boundary that the variable will get
578  // no value. To counteract this we're going to let every processor set a
579  // value at every node and then let PETSc figure it out.
580  // Later we can choose to do something different / better.
581  // if ((dof_indices[i] >= first) && (dof_indices[i] < last))
582  {
583  solution.set(dof_indices[i], Ue(i));
584  }
585  _var.setNodalValue(Ue);
586  }
587 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
A class for creating restricted objects.
Definition: Restartable.h:31
const FEType & feType() const
Get the type of finite element object.
InputParameters validParams< BlockRestrictable >()
virtual ~InitialCondition()
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
/class BoundaryRestrictable /brief Provides functionality for limiting the object to certain boundary...
InputParameters validParams< InitialCondition >()
Base class for a system (of equations)
Definition: SystemBase.h:91
std::set< std::string > _depend_uo
Depend UserObjects.
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
const Elem *& _current_elem
The current element we are on will retrieving values at specific points in the domain.
const std::map< std::string, std::vector< MooseVariable * > > & getCoupledVars()
Get the list of coupled variables.
Definition: Coupleable.h:54
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:36
virtual const std::set< std::string > & getRequestedItems() override
Return a set containing the names of items requested by the object.
void setNodalValue(Number value, unsigned int idx=0)
Set the nodal value for this variable to keep everything up to date.
unsigned int _qp
The current quadrature point, contains the "nth" node number when visiting nodes. ...
MooseVariable & _var
The variable that this initial condition is acting upon.
Interface for objects that need to use UserObjects.
const InputParameters & _pars
Parameters of this object, references the InputParameters stored in the InputParametersWarehouse.
Definition: MooseObject.h:111
InputParameters validParams< MooseObject >()
Definition: MooseObject.C:22
virtual const std::set< std::string > & getSuppliedItems() override
Return a set containing the names of items owned by the object.
const UserObject & getUserObjectBase(const std::string &name)
Interface for objects that needs coupling capabilities.
Definition: Coupleable.h:35
std::vector< dof_id_type > & dofIndices()
std::set< std::string > _supplied_vars
const bool _ignore_uo_dependency
If set, UOs retrieved by this IC will not be executed before this IC.
Interface to bring zero values inside objects.
Definition: ZeroInterface.h:35
PetscInt n
const Node * _current_node
The current node if the point we are evaluating at also happens to be a node.
InputParameters validParams< BoundaryRestrictable >()
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
virtual NumericVector< Number > & solution()=0
Interface for sorting dependent vectors of objects.
virtual void compute()
virtual Real value(const Point &p)=0
The value of the variable at a point.
std::set< std::string > _depend_vars
InitialCondition(const InputParameters &parameters)
Constructor.
SystemBase & sys()
Get the system this variable is part of.
const UserObject & getUserObjectBase(const std::string &name)
Get an user object with a given parameter name.
Interface for objects that need to use functions.
Base class for user-specific data.
Definition: UserObject.h:42
virtual RealGradient gradient(const Point &)
The gradient of the variable at a point.
unsigned int THREAD_ID
Definition: MooseTypes.h:79