www.mooseframework.org
DisplacedProblem.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 // MOOSE includes
16 #include "DisplacedProblem.h"
17 
18 #include "Assembly.h"
19 #include "AuxiliarySystem.h"
20 #include "FEProblem.h"
21 #include "MooseApp.h"
22 #include "MooseMesh.h"
23 #include "NonlinearSystem.h"
24 #include "Problem.h"
26 #include "SubProblem.h"
28 
29 #include "libmesh/numeric_vector.h"
30 
31 template <>
34 {
36  params.addPrivateParam<std::vector<std::string>>("displacements");
37  return params;
38 }
39 
41  : SubProblem(parameters),
42  _mproblem(parameters.have_parameter<FEProblemBase *>("_fe_problem_base")
43  ? *getParam<FEProblemBase *>("_fe_problem_base")
44  : *getParam<FEProblem *>("_fe_problem")),
45  _mesh(*getParam<MooseMesh *>("mesh")),
46  _eq(_mesh),
47  _ref_mesh(_mproblem.mesh()),
48  _displacements(getParam<std::vector<std::string>>("displacements")),
49  _displaced_nl(*this,
50  _mproblem.getNonlinearSystemBase(),
51  _mproblem.getNonlinearSystemBase().name() + "_displaced",
53  _displaced_aux(*this,
54  _mproblem.getAuxiliarySystem(),
55  _mproblem.getAuxiliarySystem().name() + "_displaced",
57  _geometric_search_data(_mproblem, _mesh)
58 {
59  // TODO: Move newAssemblyArray further up to SubProblem so that we can use it here
60  unsigned int n_threads = libMesh::n_threads();
61 
62  _assembly.reserve(n_threads);
63  for (unsigned int i = 0; i < n_threads; ++i)
64  _assembly.emplace_back(libmesh_make_unique<Assembly>(_displaced_nl, i));
65 }
66 
67 bool
69 {
70  return _mproblem.isTransient();
71 }
72 
75 {
76  return _mproblem.getCoordSystem(sid);
77 }
78 
79 std::set<dof_id_type> &
81 {
82  return _mproblem.ghostedElems();
83 }
84 
85 void
87  Order order,
88  Order volume_order,
89  Order face_order)
90 {
91  for (unsigned int tid = 0; tid < libMesh::n_threads(); ++tid)
92  _assembly[tid]->createQRules(type, order, volume_order, face_order);
93 }
94 
95 void
97 {
98  unsigned int n_threads = libMesh::n_threads();
99 
100  for (unsigned int i = 0; i < n_threads; ++i)
101  _assembly[i]->useFECache(
102  fe_cache); // fe caching is turned off for now for the displaced system.
103 }
104 
105 void
107 {
108  for (THREAD_ID tid = 0; tid < libMesh::n_threads(); ++tid)
110 
111  _displaced_nl.dofMap().attach_extra_send_list_function(&extraSendList, &_displaced_nl);
112  _displaced_aux.dofMap().attach_extra_send_list_function(&extraSendList, &_displaced_aux);
113 
116 
117  Moose::perf_log.push("DisplacedProblem::init::eq.init()", "Setup");
118  _eq.init();
119  Moose::perf_log.pop("DisplacedProblem::init::eq.init()", "Setup");
120 
121  Moose::perf_log.push("DisplacedProblem::init::meshChanged()", "Setup");
122  _mesh.meshChanged();
123  Moose::perf_log.pop("DisplacedProblem::init::meshChanged()", "Setup");
124 }
125 
126 void
128 {
129 }
130 
131 void
133 {
136 }
137 
138 void
140 {
143 }
144 
145 void
147 {
152 }
153 
154 void
155 DisplacedProblem::syncSolutions(const NumericVector<Number> & soln,
156  const NumericVector<Number> & aux_soln)
157 {
158  (*_displaced_nl.sys().solution) = soln;
159  (*_displaced_aux.sys().solution) = aux_soln;
162 }
163 
164 void
166 {
167  Moose::perf_log.push("updateDisplacedMesh()", "Execution");
168 
169  unsigned int n_threads = libMesh::n_threads();
170 
171  syncSolutions();
172 
173  for (unsigned int i = 0; i < n_threads; ++i)
174  _assembly[i]->invalidateCache();
175 
178 
179  // If the displaced mesh has been serialized to one processor (as
180  // may have occurred if it was used for Exodus output), then we need
181  // the reference mesh to be also. For that matter, did anyone
182  // somehow serialize the whole mesh? Hopefully not but let's avoid
183  // causing errors if so.
184  if (_mesh.getMesh().is_serial() && !this->refMesh().getMesh().is_serial())
185  this->refMesh().getMesh().allgather();
186 
187  if (_mesh.getMesh().is_serial_on_zero() && !this->refMesh().getMesh().is_serial_on_zero())
188  this->refMesh().getMesh().gather_to_zero();
189 
191 
192  // We displace all nodes, not just semilocal nodes, because
193  // parallel-inconsistent mesh geometry makes libMesh cry.
194  NodeRange node_range(_mesh.getMesh().nodes_begin(),
195  _mesh.getMesh().nodes_end(),
196  /*grainsize=*/1);
197 
198  Threads::parallel_reduce(node_range, udmt);
199 
200  // Update the geometric searches that depend on the displaced mesh
202 
203  // Since the Mesh changed, update the PointLocator object used by DiracKernels.
205 
206  Moose::perf_log.pop("updateDisplacedMesh()", "Execution");
207 }
208 
209 void
210 DisplacedProblem::updateMesh(const NumericVector<Number> & soln,
211  const NumericVector<Number> & aux_soln)
212 {
213  Moose::perf_log.push("updateDisplacedMesh()", "Execution");
214 
215  unsigned int n_threads = libMesh::n_threads();
216 
217  syncSolutions(soln, aux_soln);
218 
219  for (unsigned int i = 0; i < n_threads; ++i)
220  _assembly[i]->invalidateCache();
221 
222  _nl_solution = &soln;
223  _aux_solution = &aux_soln;
224 
226 
227  // We displace all nodes, not just semilocal nodes, because
228  // parallel-inconsistent mesh geometry makes libMesh cry.
229  NodeRange node_range(_mesh.getMesh().nodes_begin(),
230  _mesh.getMesh().nodes_end(),
231  /*grainsize=*/1);
232 
233  Threads::parallel_reduce(node_range, udmt);
234 
235  // Update the geometric searches that depend on the displaced mesh
237 
238  // Since the Mesh changed, update the PointLocator object used by DiracKernels.
240 
241  Moose::perf_log.pop("updateDisplacedMesh()", "Execution");
242 }
243 
244 bool
245 DisplacedProblem::hasVariable(const std::string & var_name)
246 {
247  if (_displaced_nl.hasVariable(var_name))
248  return true;
249  else if (_displaced_aux.hasVariable(var_name))
250  return true;
251  else
252  return false;
253 }
254 
256 DisplacedProblem::getVariable(THREAD_ID tid, const std::string & var_name)
257 {
258  if (_displaced_nl.hasVariable(var_name))
259  return _displaced_nl.getVariable(tid, var_name);
260  else if (!_displaced_aux.hasVariable(var_name))
261  mooseError("No variable with name '" + var_name + "'");
262 
263  return _displaced_aux.getVariable(tid, var_name);
264 }
265 
266 bool
267 DisplacedProblem::hasScalarVariable(const std::string & var_name)
268 {
269  if (_displaced_nl.hasScalarVariable(var_name))
270  return true;
271  else if (_displaced_aux.hasScalarVariable(var_name))
272  return true;
273  else
274  return false;
275 }
276 
278 DisplacedProblem::getScalarVariable(THREAD_ID tid, const std::string & var_name)
279 {
280  if (_displaced_nl.hasScalarVariable(var_name))
281  return _displaced_nl.getScalarVariable(tid, var_name);
282  else if (_displaced_aux.hasScalarVariable(var_name))
283  return _displaced_aux.getScalarVariable(tid, var_name);
284  else
285  mooseError("No variable with name '" + var_name + "'");
286 }
287 
288 System &
289 DisplacedProblem::getSystem(const std::string & var_name)
290 {
291  if (_displaced_nl.hasVariable(var_name))
292  return _displaced_nl.system();
293  else if (_displaced_aux.hasVariable(var_name))
294  return _displaced_aux.system();
295  else
296  mooseError("Unable to find a system containing the variable " + var_name);
297 }
298 
299 void
300 DisplacedProblem::addVariable(const std::string & var_name,
301  const FEType & type,
302  Real scale_factor,
303  const std::set<SubdomainID> * const active_subdomains)
304 {
305  _displaced_nl.addVariable(var_name, type, scale_factor, active_subdomains);
306 }
307 
308 void
309 DisplacedProblem::addAuxVariable(const std::string & var_name,
310  const FEType & type,
311  const std::set<SubdomainID> * const active_subdomains)
312 {
313  _displaced_aux.addVariable(var_name, type, 1.0, active_subdomains);
314 }
315 
316 void
317 DisplacedProblem::addScalarVariable(const std::string & var_name,
318  Order order,
319  Real scale_factor,
320  const std::set<SubdomainID> * const active_subdomains)
321 {
322  _displaced_nl.addScalarVariable(var_name, order, scale_factor, active_subdomains);
323 }
324 
325 void
326 DisplacedProblem::addAuxScalarVariable(const std::string & var_name,
327  Order order,
328  Real scale_factor,
329  const std::set<SubdomainID> * const active_subdomains)
330 {
331  _displaced_aux.addScalarVariable(var_name, order, scale_factor, active_subdomains);
332 }
333 
334 void
335 DisplacedProblem::prepare(const Elem * elem, THREAD_ID tid)
336 {
337  _assembly[tid]->reinit(elem);
338 
339  _displaced_nl.prepare(tid);
340  _displaced_aux.prepare(tid);
341  _assembly[tid]->prepare();
342 }
343 
344 void
346 {
347  _assembly[tid]->prepareNonlocal();
348 }
349 
350 void
351 DisplacedProblem::prepareFace(const Elem * /*elem*/, THREAD_ID tid)
352 {
353  _displaced_nl.prepareFace(tid, true);
354  _displaced_aux.prepareFace(tid, false);
355 }
356 
357 void
358 DisplacedProblem::prepare(const Elem * elem,
359  unsigned int ivar,
360  unsigned int jvar,
361  const std::vector<dof_id_type> & dof_indices,
362  THREAD_ID tid)
363 {
364  _assembly[tid]->reinit(elem);
365 
366  _displaced_nl.prepare(tid);
367  _displaced_aux.prepare(tid);
368  _assembly[tid]->prepareBlock(ivar, jvar, dof_indices);
369 }
370 
371 void
373 {
374  SubdomainID did = elem->subdomain_id();
375  _assembly[tid]->setCurrentSubdomainID(did);
376 }
377 
378 void
379 DisplacedProblem::setNeighborSubdomainID(const Elem * elem, unsigned int side, THREAD_ID tid)
380 {
381  SubdomainID did = elem->neighbor_ptr(side)->subdomain_id();
382  _assembly[tid]->setCurrentNeighborSubdomainID(did);
383 }
384 
385 void
387  unsigned int jvar,
388  const std::vector<dof_id_type> & idof_indices,
389  const std::vector<dof_id_type> & jdof_indices,
390  THREAD_ID tid)
391 {
392  _assembly[tid]->prepareBlockNonlocal(ivar, jvar, idof_indices, jdof_indices);
393 }
394 
395 void
397 {
398  _assembly[tid]->prepare();
399 }
400 
401 void
403 {
404  _assembly[tid]->prepareNeighbor();
405 }
406 
407 bool
409 {
410  std::vector<Point> & points = _dirac_kernel_info.getPoints()[elem].first;
411 
412  unsigned int n_points = points.size();
413 
414  if (n_points)
415  {
416  _assembly[tid]->reinitAtPhysical(elem, points);
417 
418  _displaced_nl.prepare(tid);
419  _displaced_aux.prepare(tid);
420 
421  reinitElem(elem, tid);
422  }
423 
424  _assembly[tid]->prepare();
425 
426  return n_points > 0;
427 }
428 
429 void
431 {
432  _displaced_nl.reinitElem(elem, tid);
433  _displaced_aux.reinitElem(elem, tid);
434 }
435 
436 void
438  std::vector<Point> phys_points_in_elem,
439  THREAD_ID tid)
440 {
441  std::vector<Point> points(phys_points_in_elem.size());
442  std::copy(phys_points_in_elem.begin(), phys_points_in_elem.end(), points.begin());
443 
444  _assembly[tid]->reinitAtPhysical(elem, points);
445 
446  _displaced_nl.prepare(tid);
447  _displaced_aux.prepare(tid);
448  _assembly[tid]->prepare();
449 
450  reinitElem(elem, tid);
451 }
452 
453 void
455  unsigned int side,
456  BoundaryID bnd_id,
457  THREAD_ID tid)
458 {
459  _assembly[tid]->reinit(elem, side);
460  _displaced_nl.reinitElemFace(elem, side, bnd_id, tid);
461  _displaced_aux.reinitElemFace(elem, side, bnd_id, tid);
462 }
463 
464 void
466 {
467  _assembly[tid]->reinit(node);
468  _displaced_nl.reinitNode(node, tid);
469  _displaced_aux.reinitNode(node, tid);
470 }
471 
472 void
474 {
475  _assembly[tid]->reinit(node);
476  _displaced_nl.reinitNodeFace(node, bnd_id, tid);
477  _displaced_aux.reinitNodeFace(node, bnd_id, tid);
478 }
479 
480 void
481 DisplacedProblem::reinitNodes(const std::vector<dof_id_type> & nodes, THREAD_ID tid)
482 {
483  _displaced_nl.reinitNodes(nodes, tid);
484  _displaced_aux.reinitNodes(nodes, tid);
485 }
486 
487 void
488 DisplacedProblem::reinitNodesNeighbor(const std::vector<dof_id_type> & nodes, THREAD_ID tid)
489 {
492 }
493 
494 void
495 DisplacedProblem::reinitNeighbor(const Elem * elem, unsigned int side, THREAD_ID tid)
496 {
497  const Elem * neighbor = elem->neighbor_ptr(side);
498  unsigned int neighbor_side = neighbor->which_neighbor_am_i(elem);
499 
500  _assembly[tid]->reinitElemAndNeighbor(elem, side, neighbor, neighbor_side);
501 
504 
505  _assembly[tid]->prepareNeighbor();
506 
507  BoundaryID bnd_id = 0; // some dummy number (it is not really used for anything, right now)
508  _displaced_nl.reinitElemFace(elem, side, bnd_id, tid);
509  _displaced_aux.reinitElemFace(elem, side, bnd_id, tid);
510 
511  _displaced_nl.reinitNeighborFace(neighbor, neighbor_side, bnd_id, tid);
512  _displaced_aux.reinitNeighborFace(neighbor, neighbor_side, bnd_id, tid);
513 }
514 
515 void
517  unsigned int neighbor_side,
518  const std::vector<Point> & physical_points,
519  THREAD_ID tid)
520 {
521  // Reinit shape functions
522  _assembly[tid]->reinitNeighborAtPhysical(neighbor, neighbor_side, physical_points);
523 
524  // Set the neighbor dof indices
527 
529 
530  // Compute values at the points
531  _displaced_nl.reinitNeighborFace(neighbor, neighbor_side, 0, tid);
532  _displaced_aux.reinitNeighborFace(neighbor, neighbor_side, 0, tid);
533 }
534 
535 void
537  const std::vector<Point> & physical_points,
538  THREAD_ID tid)
539 {
540  // Reinit shape functions
541  _assembly[tid]->reinitNeighborAtPhysical(neighbor, physical_points);
542 
543  // Set the neighbor dof indices
546 
548 
549  // Compute values at the points
550  _displaced_nl.reinitNeighbor(neighbor, tid);
551  _displaced_aux.reinitNeighbor(neighbor, tid);
552 }
553 
554 void
556 {
557  _assembly[tid]->reinitNodeNeighbor(node);
560 }
561 
562 void
564 {
567 }
568 
569 void
571 {
572  _assembly[tid]->prepareOffDiagScalar();
573 }
574 
575 void
576 DisplacedProblem::getDiracElements(std::set<const Elem *> & elems)
577 {
579 }
580 
581 void
583 {
585 }
586 
587 void
589 {
594 }
595 
596 void
598 {
604 }
605 
606 void
608 {
609  _assembly[tid]->cacheResidual();
610 }
611 
612 void
614 {
615  _assembly[tid]->cacheResidualNeighbor();
616 }
617 
618 void
620 {
621  _assembly[tid]->addCachedResiduals();
622 }
623 
624 void
625 DisplacedProblem::addCachedResidualDirectly(NumericVector<Number> & residual, THREAD_ID tid)
626 {
627  _assembly[tid]->addCachedResidual(residual, Moose::KT_TIME);
628  _assembly[tid]->addCachedResidual(residual, Moose::KT_NONTIME);
629 }
630 
631 void
632 DisplacedProblem::setResidual(NumericVector<Number> & residual, THREAD_ID tid)
633 {
634  _assembly[tid]->setResidual(residual);
635 }
636 
637 void
638 DisplacedProblem::setResidualNeighbor(NumericVector<Number> & residual, THREAD_ID tid)
639 {
640  _assembly[tid]->setResidualNeighbor(residual);
641 }
642 
643 void
644 DisplacedProblem::addJacobian(SparseMatrix<Number> & jacobian, THREAD_ID tid)
645 {
646  _assembly[tid]->addJacobian(jacobian);
647 }
648 
649 void
650 DisplacedProblem::addJacobianNonlocal(SparseMatrix<Number> & jacobian, THREAD_ID tid)
651 {
652  _assembly[tid]->addJacobianNonlocal(jacobian);
653 }
654 
655 void
656 DisplacedProblem::addJacobianNeighbor(SparseMatrix<Number> & jacobian, THREAD_ID tid)
657 {
658  _assembly[tid]->addJacobianNeighbor(jacobian);
659 }
660 
661 void
663 {
664  _assembly[tid]->cacheJacobian();
665 }
666 
667 void
669 {
670  _assembly[tid]->cacheJacobianNonlocal();
671 }
672 
673 void
675 {
676  _assembly[tid]->cacheJacobianNeighbor();
677 }
678 
679 void
680 DisplacedProblem::addCachedJacobian(SparseMatrix<Number> & jacobian, THREAD_ID tid)
681 {
682  _assembly[tid]->addCachedJacobian(jacobian);
683 }
684 
685 void
686 DisplacedProblem::addJacobianBlock(SparseMatrix<Number> & jacobian,
687  unsigned int ivar,
688  unsigned int jvar,
689  const DofMap & dof_map,
690  std::vector<dof_id_type> & dof_indices,
691  THREAD_ID tid)
692 {
693  _assembly[tid]->addJacobianBlock(jacobian, ivar, jvar, dof_map, dof_indices);
694 }
695 
696 void
697 DisplacedProblem::addJacobianBlockNonlocal(SparseMatrix<Number> & jacobian,
698  unsigned int ivar,
699  unsigned int jvar,
700  const DofMap & dof_map,
701  const std::vector<dof_id_type> & idof_indices,
702  const std::vector<dof_id_type> & jdof_indices,
703  THREAD_ID tid)
704 {
705  _assembly[tid]->addJacobianBlockNonlocal(
706  jacobian, ivar, jvar, dof_map, idof_indices, jdof_indices);
707 }
708 
709 void
710 DisplacedProblem::addJacobianNeighbor(SparseMatrix<Number> & jacobian,
711  unsigned int ivar,
712  unsigned int jvar,
713  const DofMap & dof_map,
714  std::vector<dof_id_type> & dof_indices,
715  std::vector<dof_id_type> & neighbor_dof_indices,
716  THREAD_ID tid)
717 {
718  _assembly[tid]->addJacobianNeighbor(
719  jacobian, ivar, jvar, dof_map, dof_indices, neighbor_dof_indices);
720 }
721 
722 void
724 {
725  _assembly[tid]->copyShapes(var);
726 }
727 
728 void
730 {
731  _assembly[tid]->copyFaceShapes(var);
732 }
733 
734 void
736 {
737  _assembly[tid]->copyNeighborShapes(var);
738 }
739 
740 void
742 {
744 }
745 
746 void
748 {
749  // mesh changed
750  _eq.reinit();
751  _mesh.meshChanged();
752 
753  // Since the Mesh changed, update the PointLocator object used by DiracKernels.
755 
756  unsigned int n_threads = libMesh::n_threads();
757 
758  for (unsigned int i = 0; i < n_threads; ++i)
759  _assembly[i]->invalidateCache();
761 }
762 
763 void
765 {
766  _mproblem.addGhostedElem(elem_id);
767 }
768 
769 void
771 {
772  _mproblem.addGhostedBoundary(boundary_id);
773 }
774 
775 void
777 {
779 }
780 
781 MooseMesh &
783 {
784  return _ref_mesh;
785 }
786 
787 void
789 {
790 }
791 
792 bool
794 {
795  return _mproblem.converged();
796 }
797 
798 bool
800 {
802 }
803 
804 void
806 {
807 }
808 
809 void
811 {
812 }
813 
814 void
816 {
817  // If undisplaceMesh() is called during initial adaptivity, it is
818  // not valid to call _mesh.getActiveSemiLocalNodeRange() since it is
819  // not set up yet. So we are creating the Range by hand.
820  //
821  // We must undisplace *all* our nodes to the _ref_mesh
822  // configuration, not just the local ones, since the partitioners
823  // require this. We are using the GRAIN_SIZE=1 from MooseMesh.C,
824  // not sure how this value was decided upon.
825  //
826  // (DRG: The grainsize parameter is ultimately passed to TBB to help
827  // it choose how to split up the range. A grainsize of 1 says "split
828  // it as much as you want". Years ago I experimentally found that it
829  // didn't matter much and that using 1 was fine.)
830  //
831  // Note: we don't have to invalidate/update as much stuff as
832  // DisplacedProblem::updateMesh() does, since this will be handled
833  // by a later call to updateMesh().
834  NodeRange node_range(_mesh.getMesh().nodes_begin(),
835  _mesh.getMesh().nodes_end(),
836  /*grainsize=*/1);
837 
838  ResetDisplacedMeshThread rdmt(_mproblem, *this);
839 
840  // Undisplace the mesh using threads.
841  Threads::parallel_reduce(node_range, rdmt);
842 }
virtual void cacheJacobian(THREAD_ID tid) override
virtual const NumericVector< Number > *& currentSolution() override
The solution vector that is currently being operated on.
GeometricSearchData _geometric_search_data
virtual void cacheResidualNeighbor(THREAD_ID tid) override
virtual void reinitNode(const Node *node, THREAD_ID tid)
Reinit nodal assembly info.
Definition: SystemBase.C:329
virtual void addGhostedElem(dof_id_type elem_id) override
Will make sure that all dofs connected to elem_id are ghosted to this processor.
virtual bool converged() override
virtual void getDiracElements(std::set< const Elem * > &elems) override
Fills "elems" with the elements that should be looped over for Dirac Kernels.
virtual void setResidual(NumericVector< Number > &residual, THREAD_ID tid) override
DisplacedProblem(const InputParameters &parameters)
virtual System & getSystem(const std::string &var_name) override
Returns the equation system containing the variable provided.
const CouplingMatrix * couplingMatrix()
virtual void reinitScalars(THREAD_ID tid) override
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name) override
Returns the scalar variable reference from whichever system contains it.
virtual void addGhostedElem(dof_id_type elem_id) override
Will make sure that all dofs connected to elem_id are ghosted to this processor.
virtual void saveOldSolutions()
Allocate vectors and save old solutions into them.
virtual void reinitScalars(THREAD_ID tid)
Reinit scalar varaibles.
Definition: SystemBase.C:393
virtual void reinitNeighborFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid)
Compute the values of the variables at all the current points.
Definition: SystemBase.C:310
virtual void addJacobianNonlocal(SparseMatrix< Number > &jacobian, THREAD_ID tid)
void undisplaceMesh()
Resets the displaced mesh to the reference mesh.
virtual void useFECache(bool fe_cache) override
Whether or not this problem should utilize FE shape function caching.
virtual Moose::CoordinateSystemType getCoordSystem(SubdomainID sid) override
Class for stuff related to variables.
Definition: MooseVariable.h:43
NonlinearSystemBase & getNonlinearSystemBase()
virtual TransientExplicitSystem & sys()
virtual void reinitElemFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid)
Reinit assembly info for a side of an element.
Definition: SystemBase.C:299
void addPrivateParam(const std::string &name, const T &value)
These method add a parameter to the InputParameters object which can be retrieved like any other para...
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
Definition: FEProblem.h:30
virtual void initAdaptivity()
subdomain_id_type SubdomainID
Definition: MooseTypes.h:77
virtual void prepareFace(const Elem *elem, THREAD_ID tid) override
virtual bool hasScalarVariable(const std::string &var_name)
Definition: SystemBase.C:565
virtual void prepareAssembly(THREAD_ID tid) override
virtual void addVariable(const std::string &var_name, const FEType &type, Real scale_factor, const std::set< SubdomainID > *const active_subdomains=NULL)
Adds a variable to the system.
Definition: SystemBase.C:500
virtual bool hasScalarVariable(const std::string &var_name) override
Returns a Boolean indicating whether any system contains a variable with the name provided...
virtual void reinitElem(const Elem *elem, THREAD_ID tid)
Reinit an element assembly info.
Definition: SystemBase.C:279
FEProblemBase & _mproblem
virtual void reinitNodeFace(const Node *node, BoundaryID bnd_id, THREAD_ID tid) override
std::vector< std::unique_ptr< Assembly > > _assembly
virtual void addAuxScalarVariable(const std::string &var_name, Order order, Real scale_factor=1., const std::set< SubdomainID > *const active_subdomains=NULL)
DisplacedSystem _displaced_aux
virtual void setResidualNeighbor(NumericVector< Number > &residual, THREAD_ID tid) override
virtual void addGhostedBoundary(BoundaryID boundary_id) override
Will make sure that all necessary elements from boundary_id are ghosted to this processor.
DisplacedSystem _displaced_nl
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void cacheResidual(THREAD_ID tid) override
virtual void reinitNodesNeighbor(const std::vector< dof_id_type > &nodes, THREAD_ID tid) override
virtual void addCachedJacobian(SparseMatrix< Number > &jacobian, THREAD_ID tid) override
virtual void prepareFace(THREAD_ID tid, bool resize_data)
Prepare the system for use on sides.
Definition: SystemBase.C:237
virtual void createQRules(QuadratureType type, Order order, Order volume_order, Order face_order)
virtual const NumericVector< Number > *& currentSolution() override
The solution vector that is currently being operated on.
virtual void reinitNodes(const std::vector< dof_id_type > &nodes, THREAD_ID tid)
Reinit variables at a set of nodes.
Definition: SystemBase.C:371
virtual void addJacobianBlockNonlocal(SparseMatrix< Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices, THREAD_ID tid)
virtual void reinitNodesNeighbor(const std::vector< dof_id_type > &nodes, THREAD_ID tid)
Reinit variables at a set of neighbor nodes.
Definition: SystemBase.C:382
virtual void reinitOffDiagScalars(THREAD_ID tid) override
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
virtual void onTimestepBegin() override
virtual void syncSolutions()
Copy the solutions on the undisplaced systems to the displaced systems.
virtual std::set< dof_id_type > & ghostedElems() override
Return the list of elements that should have their DoFs ghosted to this processor.
virtual void reinitNeighbor(const Elem *elem, unsigned int side, THREAD_ID tid) override
virtual void solve() override
InputParameters validParams< SubProblem >()
Definition: SubProblem.C:26
virtual void meshChanged() override
virtual void update()
Update the system (doing libMesh magic)
Definition: SystemBase.C:665
EquationSystems _eq
virtual bool converged() override
virtual DofMap & dofMap()
Gets the dof map.
Definition: SystemBase.C:610
virtual void addResidual(THREAD_ID tid) override
virtual void cacheJacobianNonlocal(THREAD_ID tid)
virtual void prepareNeighborShapes(unsigned int var, THREAD_ID tid) override
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:103
virtual void setCurrentSubdomainID(const Elem *elem, THREAD_ID tid) override
MultiPointMap & getPoints()
Returns a writeable reference to the _points container.
virtual void reinitNodeNeighbor(const Node *node, THREAD_ID tid)
Reinit nodal assembly info for neighbor node.
Definition: SystemBase.C:357
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2408
MooseMesh & _ref_mesh
reference mesh
void extraSendList(std::vector< dof_id_type > &send_list, void *context)
///< Type of coordinate system
Definition: SystemBase.C:33
virtual void prepareAssemblyNeighbor(THREAD_ID tid)
virtual void prepare(const Elem *elem, THREAD_ID tid) override
virtual bool isTransient() const override
virtual void addJacobianBlock(SparseMatrix< Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, std::vector< dof_id_type > &dof_indices, THREAD_ID tid) override
virtual void clearDiracInfo() override
Gets called before Dirac Kernels are asked to add the points they are supposed to be evaluated in...
virtual void restoreOldSolutions()
Restore old solutions from the backup vectors and deallocate them.
virtual void prepareNeighbor(THREAD_ID tid)
Prepare the system for use.
Definition: SystemBase.C:271
DofMap & dof_map
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name) override
Returns the variable reference for requested variable which may be in any system. ...
PerfLog perf_log
Perflog to be used by applications.
virtual void prepareBlockNonlocal(unsigned int ivar, unsigned int jvar, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices, THREAD_ID tid)
virtual bool isTransient() const override
void reinit()
Completely redo all geometric search objects.
virtual void reinitNeighborPhys(const Elem *neighbor, unsigned int neighbor_side, const std::vector< Point > &physical_points, THREAD_ID tid) override
virtual std::set< dof_id_type > & ghostedElems()
Return the list of elements that should have their DoFs ghosted to this processor.
Definition: SubProblem.h:400
void clearPoints()
Remove all of the current points and elements.
virtual bool hasVariable(const std::string &var_name)
Query a system for a variable.
Definition: SystemBase.C:556
const NumericVector< Number > * _aux_solution
virtual void addResidualNeighbor(THREAD_ID tid) override
AuxiliarySystem & getAuxiliarySystem()
virtual void reinitNeighbor(const Elem *elem, THREAD_ID tid)
Compute the values of the variables at all the current points.
Definition: SystemBase.C:321
virtual void updateMesh()
Copy the solutions on the undisplaced systems to the displaced systems and reinitialize the geometry ...
GeometricSearchType
Used to select groups of geometric search objects to update.
virtual void reinitNodeNeighbor(const Node *node, THREAD_ID tid) override
CoordinateSystemType
Definition: MooseTypes.h:212
virtual void onTimestepEnd() override
virtual void reinitElemPhys(const Elem *elem, std::vector< Point > phys_points_in_elem, THREAD_ID tid) override
virtual void addJacobianNeighbor(SparseMatrix< Number > &jacobian, THREAD_ID tid) override
virtual void cacheJacobianNeighbor(THREAD_ID tid) override
virtual void prepareFaceShapes(unsigned int var, THREAD_ID tid) override
virtual void saveOldSolutions()
Save the old and older solutions.
Definition: SystemBase.C:459
MatType type
virtual void prepareShapes(unsigned int var, THREAD_ID tid) override
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
InputParameters validParams< DisplacedProblem >()
virtual void reinitNodeFace(const Node *node, BoundaryID bnd_id, THREAD_ID tid)
Reinit nodal assembly info on a face.
Definition: SystemBase.C:343
virtual void addVariable(const std::string &var_name, const FEType &type, Real scale_factor, const std::set< SubdomainID > *const active_subdomains=NULL)
virtual System & system() override
Get the reference to the libMesh system.
virtual Moose::CoordinateSystemType getCoordSystem(SubdomainID sid) override
void update(GeometricSearchType type=ALL)
Update all of the search objects.
Class for scalar variables (they are different).
virtual void addCachedResidualDirectly(NumericVector< Number > &residual, THREAD_ID tid)
virtual void addJacobian(SparseMatrix< Number > &jacobian, THREAD_ID tid) override
virtual void prepare(THREAD_ID tid)
Prepare the system for use.
Definition: SystemBase.C:214
virtual void reinitNodes(const std::vector< dof_id_type > &nodes, THREAD_ID tid) override
virtual void addCachedResidual(THREAD_ID tid) override
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
virtual void updateGeomSearch(GeometricSearchData::GeometricSearchType type=GeometricSearchData::ALL) override
virtual bool computingInitialResidual() override
Returns true if the problem is in the process of computing it&#39;s initial residual. ...
const NumericVector< Number > * _nl_solution
virtual void setNeighborSubdomainID(const Elem *elem, unsigned int side, THREAD_ID tid) override
std::set< const Elem * > & getElements()
Returns a writeable reference to the _elements container.
virtual bool hasVariable(const std::string &var_name) override
MooseMesh & refMesh()
Definition: Moose.h:84
virtual void init() override
Initialize the system.
virtual bool hasResidualVector(Moose::KernelType type) const override
virtual void reinitNode(const Node *node, THREAD_ID tid) override
DiracKernelInfo _dirac_kernel_info
Definition: SubProblem.h:425
virtual void reinitElemFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid) override
virtual bool reinitDirac(const Elem *elem, THREAD_ID tid) override
Returns true if the Problem has Dirac kernels it needs to compute on elem.
virtual NonlinearSystem & getNonlinearSystem()
virtual void addGhostedBoundary(BoundaryID boundary_id) override
Will make sure that all necessary elements from boundary_id are ghosted to this processor.
virtual void reinitElem(const Elem *elem, THREAD_ID tid) override
virtual void prepareNonlocal(THREAD_ID tid)
void updatePointLocator(const MooseMesh &mesh)
Called during FEProblemBase::meshChanged() to update the PointLocator object used by the DiracKernels...
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:121
virtual void restoreOldSolutions()
Restore the old and older solutions when the saved solutions present.
Definition: SystemBase.C:473
virtual bool computingInitialResidual() override
Returns true if the problem is in the process of computing it&#39;s initial residual. ...
virtual void addScalarVariable(const std::string &var_name, Order order, Real scale_factor, const std::set< SubdomainID > *const active_subdomains=NULL)
Adds a scalar variable.
Definition: SystemBase.C:530
virtual NumericVector< Number > & residualVector(Moose::KernelType type)
void meshChanged()
Declares that the MooseMesh has changed, invalidates cached data and rebuilds caches.
Definition: MooseMesh.C:491
boundary_id_type BoundaryID
Definition: MooseTypes.h:75
unsigned int THREAD_ID
Definition: MooseTypes.h:79
virtual void ghostGhostedBoundaries() override
Causes the boundaries added using addGhostedBoundary to actually be ghosted.
virtual void init() override
virtual void ghostGhostedBoundaries() override
Causes the boundaries added using addGhostedBoundary to actually be ghosted.
virtual void addAuxVariable(const std::string &var_name, const FEType &type, const std::set< SubdomainID > *const active_subdomains=NULL)
virtual void addScalarVariable(const std::string &var_name, Order order, Real scale_factor=1., const std::set< SubdomainID > *const active_subdomains=NULL)