24 #include "libmesh/libmesh_exceptions.h" 25 #include "libmesh/elem.h" 29 class MooseVariableFVBase;
52 std::function<void()>
_f;
57 _f = std::move(other._f);
72 template <
typename RangeType>
77 const unsigned int nl_system_num,
78 const std::set<TagID> & tags,
85 virtual void operator()(
const RangeType & range,
bool bypass_threading =
false);
117 std::string what(e.
what());
180 template <
typename RangeType>
182 const unsigned int nl_system_num,
183 const std::set<TagID> & tags,
184 const bool on_displaced)
185 : _fe_problem(fe_problem),
186 _mesh(fe_problem.
mesh()),
188 _nl_system_num(nl_system_num),
189 _on_displaced(on_displaced),
190 _subproblem(_on_displaced ? static_cast<
SubProblem &>(*_fe_problem.getDisplacedProblem())
198 template <
typename RangeType>
200 : _fe_problem(x._fe_problem),
203 _nl_system_num(x._nl_system_num),
204 _on_displaced(x._on_displaced),
205 _subproblem(x._subproblem),
207 _incoming_throw_on_error(false)
211 template <
typename RangeType>
218 if (!_error_message.empty())
223 template <
typename RangeType>
237 template <
typename RangeType>
245 std::vector<FVKernel *> kernels;
246 _fe_problem.theWarehouse()
248 .template condition<AttribSysNum>(_nl_system_num)
249 .
template condition<AttribSystem>(
"FVFluxKernel")
250 .template condition<AttribDisplaced>(_on_displaced)
252 if (kernels.size() == 0)
260 _tid = bypass_threading ? 0 : puid.
id;
263 printGeneralExecutionInformation();
268 typename RangeType::const_iterator faceinfo = range.begin();
269 for (faceinfo = range.begin(); faceinfo != range.end(); ++faceinfo)
271 const Elem & elem = (*faceinfo)->elem();
273 _fe_problem.setCurrentSubdomainID(&elem, _tid);
275 _old_subdomain = _subdomain;
276 _subdomain = elem.subdomain_id();
277 if (_subdomain != _old_subdomain)
280 printBlockExecutionInformation();
283 _old_neighbor_subdomain = _neighbor_subdomain;
284 if (
const Elem *
const neighbor = (*faceinfo)->neighborPtr())
286 _fe_problem.setNeighborSubdomainID(neighbor, _tid);
287 _neighbor_subdomain = neighbor->subdomain_id();
292 if (_neighbor_subdomain != _old_neighbor_subdomain)
294 neighborSubdomainChanged();
296 printBlockExecutionInformation();
304 postFace(**faceinfo);
306 const std::set<BoundaryID> boundary_ids = (*faceinfo)->boundaryIDs();
307 for (
auto & it : boundary_ids)
309 printBoundaryExecutionInformation(it);
313 postFace(**faceinfo);
319 resetExecutionPrinting();
323 mooseException(
"We caught a libMesh error: ", e.what());
325 catch (MetaPhysicL::LogicError & e)
332 caughtMooseException(e);
334 catch (std::runtime_error & e)
336 _error_message = e.what();
343 template <
typename RangeType,
typename AttributeTagType>
348 const unsigned int nl_system_num,
349 const std::set<TagID> & tags,
357 virtual void pre()
override;
358 virtual void post()
override;
411 std::pair<SubdomainName, SubdomainName>
getBlockNames()
const;
437 template <
typename RangeType,
typename AttributeTagType>
439 unsigned int nl_system_num,
440 const std::set<TagID> & tags,
442 :
ThreadedFaceLoop<RangeType>(fe_problem, nl_system_num, tags, on_displaced),
443 _scaling_jacobian(fe_problem.computingScalingJacobian()),
444 _scaling_residual(fe_problem.computingScalingResidual())
448 template <
typename RangeType,
typename AttributeTagType>
452 _fv_vars(x._fv_vars),
453 _scaling_jacobian(x._scaling_jacobian),
454 _scaling_residual(x._scaling_residual)
458 template <
typename RangeType,
typename AttributeTagType>
463 template <
typename RangeType,
typename AttributeTagType>
471 this->_subproblem.reinitFVFace(_tid, fi);
485 for (
auto var : _fv_vars)
486 var->computeFaceValues(fi);
490 for (std::shared_ptr<MaterialBase> mat : _elem_face_mats)
492 mat->setFaceInfo(fi);
493 mat->computeProperties();
498 _fe_problem.resizeMaterialData(
501 for (std::shared_ptr<MaterialBase> mat : _neigh_face_mats)
503 mat->setFaceInfo(fi);
504 mat->computeProperties();
509 template <
typename RangeType,
typename AttributeTagType>
515 for (
auto *
const k : _fv_flux_kernels)
519 template <
typename RangeType,
typename AttributeTagType>
523 if (_scaling_jacobian || _scaling_residual)
526 std::vector<FVFluxBC *> bcs;
527 _fe_problem.theWarehouse()
529 .template condition<AttribSysNum>(_nl_system_num)
530 .
template condition<AttribSystem>(
"FVFluxBC")
531 .template condition<AttribDisplaced>(_on_displaced)
532 .
template condition<AttribThread>(_tid)
533 .template condition<AttributeTagType>(_tags)
534 .
template condition<AttribBoundaries>(bnd_id)
537 for (
auto *
const bc : bcs)
540 std::vector<FVInterfaceKernel *> iks;
541 _fe_problem.theWarehouse()
543 .template condition<AttribSysNum>(_nl_system_num)
544 .
template condition<AttribSystem>(
"FVInterfaceKernel")
545 .template condition<AttribDisplaced>(_on_displaced)
546 .
template condition<AttribThread>(_tid)
547 .template condition<AttributeTagType>(_tags)
548 .
template condition<AttribBoundaries>(bnd_id)
551 for (
auto *
const ik : iks)
555 template <
typename RangeType,
typename AttributeTagType>
560 Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
563 _fe_problem.clearActiveElementalMooseVariables(_tid);
564 _fe_problem.clearActiveMaterialProperties(_tid);
567 template <
typename RangeType,
typename AttributeTagType>
574 std::set_union(_elem_sub_fv_vars.begin(),
575 _elem_sub_fv_vars.end(),
576 _neigh_sub_fv_vars.begin(),
577 _neigh_sub_fv_vars.end(),
578 std::inserter(_fv_vars, _fv_vars.begin()));
583 const bool same_kernels = _elem_sub_fv_flux_kernels == _neigh_sub_fv_flux_kernels;
585 _fv_flux_kernels = _elem_sub_fv_flux_kernels;
587 std::set_union(_elem_sub_fv_flux_kernels.begin(),
588 _elem_sub_fv_flux_kernels.end(),
589 _neigh_sub_fv_flux_kernels.begin(),
590 _neigh_sub_fv_flux_kernels.end(),
591 std::inserter(_fv_flux_kernels, _fv_flux_kernels.begin()));
592 const bool need_ghosting = !same_kernels;
597 _elem_face_mats = _elem_sub_elem_face_mats;
602 for (std::shared_ptr<MaterialBase> neigh_sub_elem_face_mat : _neigh_sub_elem_face_mats)
603 if (std::find(_elem_sub_elem_face_mats.begin(),
604 _elem_sub_elem_face_mats.end(),
605 neigh_sub_elem_face_mat) == _elem_sub_elem_face_mats.end())
606 _elem_face_mats.push_back(neigh_sub_elem_face_mat);
611 _neigh_face_mats = _neigh_sub_neigh_face_mats;
616 for (std::shared_ptr<MaterialBase> elem_sub_neigh_face_mat : _elem_sub_neigh_face_mats)
617 if (std::find(_neigh_sub_neigh_face_mats.begin(),
618 _neigh_sub_neigh_face_mats.end(),
619 elem_sub_neigh_face_mat) == _neigh_sub_neigh_face_mats.end())
620 _neigh_face_mats.push_back(elem_sub_neigh_face_mat);
623 template <
typename RangeType,
typename AttributeTagType>
631 _elem_sub_fv_vars.clear();
634 _fv_flux_kernels.clear();
635 _elem_sub_fv_flux_kernels.clear();
638 _elem_face_mats.clear();
639 _elem_sub_elem_face_mats.clear();
642 _neigh_face_mats.clear();
643 _elem_sub_neigh_face_mats.clear();
649 std::vector<FVFluxKernel *> kernels;
650 _fe_problem.theWarehouse()
652 .template condition<AttribSysNum>(_nl_system_num)
653 .
template condition<AttribSystem>(
"FVFluxKernel")
654 .template condition<AttribDisplaced>(_on_displaced)
655 .
template condition<AttribSubdomains>(_subdomain)
656 .template condition<AttribThread>(_tid)
657 .
template condition<AttributeTagType>(_tags)
660 _elem_sub_fv_flux_kernels = std::set<FVFluxKernel *>(kernels.begin(), kernels.end());
662 for (
auto * k : _elem_sub_fv_flux_kernels)
668 const auto & deps = k->getMooseVariableDependencies();
669 for (
auto var : deps)
671 mooseAssert(var->isFV(),
672 "We do not currently support coupling of FE variables into FV objects");
673 _elem_sub_fv_vars.insert(var);
677 _fe_problem.getFVMatsAndDependencies(
678 _subdomain, _elem_sub_elem_face_mats, _elem_sub_neigh_face_mats, _elem_sub_fv_vars, _tid);
680 finalizeContainers();
683 template <
typename RangeType,
typename AttributeTagType>
691 _neigh_sub_fv_vars.clear();
694 _fv_flux_kernels.clear();
695 _neigh_sub_fv_flux_kernels.clear();
698 _elem_face_mats.clear();
699 _neigh_sub_elem_face_mats.clear();
702 _neigh_face_mats.clear();
703 _neigh_sub_neigh_face_mats.clear();
709 std::vector<FVFluxKernel *> kernels;
710 _fe_problem.theWarehouse()
712 .template condition<AttribSysNum>(_nl_system_num)
713 .
template condition<AttribSystem>(
"FVFluxKernel")
714 .template condition<AttribDisplaced>(_on_displaced)
715 .
template condition<AttribSubdomains>(_neighbor_subdomain)
716 .template condition<AttribThread>(_tid)
717 .
template condition<AttributeTagType>(_tags)
720 _neigh_sub_fv_flux_kernels = std::set<FVFluxKernel *>(kernels.begin(), kernels.end());
722 for (
auto * k : _neigh_sub_fv_flux_kernels)
728 const auto & deps = k->getMooseVariableDependencies();
729 for (
auto var : deps)
731 mooseAssert(var->isFV(),
732 "We do not currently support coupling of FE variables into FV objects");
733 _neigh_sub_fv_vars.insert(var);
737 _fe_problem.getFVMatsAndDependencies(_neighbor_subdomain,
738 _neigh_sub_elem_face_mats,
739 _neigh_sub_neigh_face_mats,
743 finalizeContainers();
746 template <
typename RangeType,
typename AttributeTagType>
750 std::vector<FVFluxBC *> bcs;
751 _fe_problem.theWarehouse()
753 .template condition<AttribSysNum>(_nl_system_num)
754 .
template condition<AttribSystem>(
"FVFluxBC")
755 .template condition<AttribDisplaced>(_on_displaced)
756 .
template condition<AttribThread>(_tid)
757 .template condition<AttributeTagType>(_tags)
760 std::vector<FVInterfaceKernel *> iks;
761 _fe_problem.theWarehouse()
763 .template condition<AttribSysNum>(_nl_system_num)
764 .
template condition<AttribSystem>(
"FVInterfaceKernel")
765 .template condition<AttribDisplaced>(_on_displaced)
766 .
template condition<AttribThread>(_tid)
767 .template condition<AttributeTagType>(_tags)
770 std::vector<FVFluxKernel *> kernels;
771 _fe_problem.theWarehouse()
773 .template condition<AttribSysNum>(_nl_system_num)
774 .
template condition<AttribSystem>(
"FVFluxKernel")
775 .template condition<AttribDisplaced>(_on_displaced)
776 .
template condition<AttribThread>(_tid)
777 .template condition<AttributeTagType>(_tags)
780 for (
auto *
const bc : bcs)
782 for (
auto *
const ik : iks)
784 for (
auto *
const kernel : kernels)
789 _elem_sub_fv_vars.clear();
790 _neigh_sub_fv_vars.clear();
793 _fv_flux_kernels.clear();
794 _elem_sub_fv_flux_kernels.clear();
795 _neigh_sub_fv_flux_kernels.clear();
798 _elem_face_mats.clear();
799 _elem_sub_elem_face_mats.clear();
800 _neigh_sub_elem_face_mats.clear();
803 _neigh_face_mats.clear();
804 _elem_sub_neigh_face_mats.clear();
805 _neigh_sub_neigh_face_mats.clear();
808 template <
typename RangeType>
813 const unsigned int nl_system_num,
814 const std::set<TagID> & tags,
832 template <
typename RangeType>
835 const unsigned int nl_system_num,
836 const std::set<TagID> & tags,
842 template <
typename RangeType>
849 template <
typename RangeType>
855 this->_subproblem.SubProblem::cacheResidual(_tid);
856 this->_subproblem.SubProblem::cacheResidualNeighbor(_tid);
858 this->_subproblem.SubProblem::addCachedResidual(_tid);
861 template <
typename RangeType>
866 const unsigned int nl_system_num,
867 const std::set<TagID> & tags,
883 template <
typename RangeType>
886 const unsigned int nl_system_num,
887 const std::set<TagID> & tags,
893 template <
typename RangeType>
900 template <
typename RangeType>
908 if (_num_cached % 20 == 0)
910 Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
911 this->_subproblem.SubProblem::addCachedJacobian(_tid);
915 template <
typename RangeType>
920 const unsigned int nl_system_num,
921 const std::set<TagID> & vector_tags,
922 const std::set<TagID> & ,
945 template <
typename RangeType>
947 const unsigned int nl_system_num,
948 const std::set<TagID> & vector_tags,
949 const std::set<TagID> & ,
952 fe_problem, nl_system_num, vector_tags, on_displaced)
956 template <
typename RangeType>
963 template <
typename RangeType>
969 this->_subproblem.SubProblem::cacheResidual(_tid);
970 this->_subproblem.SubProblem::cacheResidualNeighbor(_tid);
974 if (_num_cached % 20 == 0)
976 Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
977 this->_subproblem.SubProblem::addCachedResidual(_tid);
978 this->_subproblem.SubProblem::addCachedJacobian(_tid);
982 template <
typename RangeType,
typename AttributeTagType>
986 if (!_fe_problem.shouldPrintExecution(_tid))
988 auto & console = _fe_problem.console();
989 auto execute_on = _fe_problem.getCurrentExecuteOnFlag();
990 console <<
"[DBG] Beginning finite volume flux objects loop on " << execute_on << std::endl;
991 mooseDoOnce(console <<
"[DBG] Loop on faces (FaceInfo), objects ordered on each face: " 993 console <<
"[DBG] - (finite volume) flux kernels" << std::endl;
994 console <<
"[DBG] - (finite volume) flux boundary conditions" << std::endl;
995 console <<
"[DBG] - (finite volume) interface kernels" << std::endl;);
998 template <
typename RangeType,
typename AttributeTagType>
1002 if (!_fe_problem.shouldPrintExecution(_tid) || !_fv_flux_kernels.size())
1006 const auto block_pair = std::make_pair(_subdomain, _neighbor_subdomain);
1007 const auto block_pair_names = this->getBlockNames();
1008 if (_blocks_exec_printed.count(block_pair))
1010 auto & console = _fe_problem.console();
1011 console <<
"[DBG] Flux kernels on block " << block_pair_names.first;
1013 console <<
" and neighbor " << block_pair_names.second << std::endl;
1015 console <<
" with no neighbor block" << std::endl;
1018 std::vector<MooseObject *> fv_flux_kernels;
1019 for (
const auto & fv_kernel : _fv_flux_kernels)
1020 fv_flux_kernels.push_back(dynamic_cast<MooseObject *>(fv_kernel));
1024 _blocks_exec_printed.insert(block_pair);
1027 template <
typename RangeType,
typename AttributeTagType>
1032 if (!_fe_problem.shouldPrintExecution(_tid))
1034 if (_boundaries_exec_printed.count(bnd_id))
1036 std::vector<MooseObject *> bcs;
1037 _fe_problem.theWarehouse()
1039 .template condition<AttribSystem>(
"FVFluxBC")
1040 .
template condition<AttribDisplaced>(_on_displaced)
1041 .template condition<AttribThread>(_tid)
1042 .
template condition<AttributeTagType>(_tags)
1043 .template condition<AttribBoundaries>(bnd_id)
1046 std::vector<MooseObject *> iks;
1047 _fe_problem.theWarehouse()
1049 .template condition<AttribSystem>(
"FVInterfaceKernel")
1050 .
template condition<AttribDisplaced>(_on_displaced)
1051 .template condition<AttribThread>(_tid)
1052 .
template condition<AttributeTagType>(_tags)
1053 .template condition<AttribBoundaries>(bnd_id)
1056 const auto block_pair_names = this->getBlockNames();
1059 auto & console = _fe_problem.console();
1060 console <<
"[DBG] FVBCs on boundary " << bnd_id <<
" between subdomain " 1061 << block_pair_names.first;
1063 console <<
" and neighbor " << block_pair_names.second << std::endl;
1065 console <<
" and the exterior of the mesh " << std::endl;
1071 auto & console = _fe_problem.console();
1072 console <<
"[DBG] FVIKs on boundary " << bnd_id <<
" between subdomain " 1073 << block_pair_names.first;
1075 console <<
" and neighbor " << block_pair_names.second << std::endl;
1077 console <<
" and the exterior of the mesh " << std::endl;
1081 _boundaries_exec_printed.insert(bnd_id);
1084 template <
typename RangeType,
typename AttributeTagType>
1085 std::pair<SubdomainName, SubdomainName>
1088 auto block_names = std::make_pair(_mesh.getSubdomainName(_subdomain),
1089 _mesh.getSubdomainName(_neighbor_subdomain));
1090 if (block_names.first ==
"")
1092 if (block_names.second ==
"")
void setup(SetupInterface &obj) override
call either residualSetup or jacobianSetup depending on what derived class of this class is instantia...
void compute(FVFaceResidualObject &ro, const FaceInfo &fi) override
call either computeResidual, computeJacobian, or computeResidualAndJacobian on the provided residual ...
virtual ~ComputeFVFluxThread()
virtual void post() override
This is called once after all face-looping is finished.
void postFace(const FaceInfo &fi) override
This is called once for each face after all face and boundary callbacks have been finished for that f...
std::vector< std::shared_ptr< MaterialBase > > _neigh_face_mats
const bool _incoming_throw_on_error
The value of Moose::_throw_on_error at the time of construction.
std::string mooseObjectVectorToString(const std::vector< MooseObject *> &objs, const std::string &sep=" ")
Routine to output the name of MooseObjects in a string.
virtual void neighborSubdomainChanged() override
Called every time the neighbor subdomain changes (i.e.
virtual void computeResidualAndJacobian(const FaceInfo &fi)=0
Compute the residual and Jacobian on the supplied face.
void setup(SetupInterface &obj) override
call either residualSetup or jacobianSetup depending on what derived class of this class is instantia...
const unsigned int _nl_system_num
std::set< FVFluxKernel * > _neigh_sub_fv_flux_kernels
ComputeFVFluxRJThread(FEProblemBase &fe_problem, const unsigned int nl_system_num, const std::set< TagID > &vector_tags, const std::set< TagID > &, bool on_displaced)
virtual const char * what() const
Get out the error message.
std::vector< std::shared_ptr< MaterialBase > > _neigh_sub_neigh_face_mats
std::string _error_message
Holds caught runtime error messages.
std::set< BoundaryID > _boundaries_exec_printed
Set to keep track of boundaries for which we have printed the execution pattern.
virtual void computeResidual(const FaceInfo &fi)=0
Compute the residual on the supplied face.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
std::vector< std::shared_ptr< MaterialBase > > _elem_sub_neigh_face_mats
virtual void computeJacobian(const FaceInfo &fi)=0
Compute the jacobian on the supplied face.
const bool _zeroth_copy
Whether this is the zeroth threaded copy of this body.
void compute(FVFaceResidualObject &ro, const FaceInfo &fi) override
call either computeResidual, computeJacobian, or computeResidualAndJacobian on the provided residual ...
FEProblemBase & _fe_problem
void translateMetaPhysicLError(const MetaPhysicL::LogicError &)
emit a relatively clear error message when we catch a MetaPhysicL logic error
SubdomainID _subdomain
The subdomain for the current element.
virtual void postFace(const FaceInfo &)
This is called once for each face after all face and boundary callbacks have been finished for that f...
virtual void printBoundaryExecutionInformation(const BoundaryID bnd_id) const override
Print ordering of objects exected on each boundary.
void resetExecutionPrinting()
Reset lists of blocks and boundaries for which execution printing has been done.
virtual void setException(const std::string &message)
Set an exception, which is stored at this point by toggling a member variable in this class...
std::vector< std::shared_ptr< MaterialBase > > _neigh_sub_elem_face_mats
std::set< MooseVariableFieldBase * > _neigh_sub_fv_vars
Interface class for a finite volume residual object whose residuals are based on faces.
void postFace(const FaceInfo &fi) override
This is called once for each face after all face and boundary callbacks have been finished for that f...
virtual void onBoundary(const FaceInfo &fi, BoundaryID boundary) override
This is called once for every face that is on a boundary after onFace is called for the face...
void compute(FVFaceResidualObject &ro, const FaceInfo &fi) override
call either computeResidual, computeJacobian, or computeResidualAndJacobian on the provided residual ...
ThreadedFaceLoop(FEProblemBase &fe_problem, const unsigned int nl_system_num, const std::set< TagID > &tags, bool on_displaced)
ComputeFVFluxThread(FEProblemBase &fe_problem, const unsigned int nl_system_num, const std::set< TagID > &tags, bool on_displaced)
virtual void post()
This is called once after all face-looping is finished.
virtual void printBoundaryExecutionInformation(const BoundaryID) const
Print ordering of objects exected on each boundary.
std::function< void()> _f
std::pair< SubdomainName, SubdomainName > getBlockNames() const
Utility to get the subdomain names from the ids.
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
virtual void compute(FVFaceResidualObject &ro, const FaceInfo &fi)=0
call either computeResidual, computeJacobian, or computeResidualAndJacobian on the provided residual ...
const std::set< TagID > & _tags
SubdomainID _neighbor_subdomain
The subdomain for the current neighbor.
const SubdomainID INVALID_BLOCK_ID
virtual void neighborSubdomainChanged()
Called every time the neighbor subdomain changes (i.e.
virtual void onBoundary(const FaceInfo &fi, BoundaryID boundary)=0
This is called once for every face that is on a boundary after onFace is called for the face...
This data structure is used to store geometric and variable related metadata about each cell face in ...
Base class for assembly-like calculations.
const Elem * neighborPtr() const
ComputeFVFluxResidualThread(FEProblemBase &fe_problem, const unsigned int nl_system_num, const std::set< TagID > &tags, bool on_displaced)
virtual ~ThreadedFaceLoop()
std::string formatString(std::string message, const std::string &prefix)
Add new lines and prefixes to a string for pretty display in output NOTE: This makes a copy of the st...
boundary_id_type BoundaryID
This works like the Sentinel helper classes in MOOSE, except it is more flexible and concise to use...
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
ComputeFVFluxJacobianThread(FEProblemBase &fe_problem, const unsigned int nl_system_num, const std::set< TagID > &tags, bool on_displaced)
OnScopeExit(OnScopeExit &&other)
std::vector< std::shared_ptr< MaterialBase > > _elem_face_mats
Element face materials.
SubProblem & _subproblem
FEProblemBase or DisplacedProblem depending on _on_displaced.
virtual void printGeneralExecutionInformation() const
Print list of object types executed and in which order.
std::vector< std::shared_ptr< MaterialBase > > _elem_sub_elem_face_mats
std::string stringify(const T &t)
conversion to string
virtual void pre()
This is called once before all face-looping.
virtual void onFace(const FaceInfo &fi)=0
virtual void printGeneralExecutionInformation() const override
Print list of object types executed and in which order.
std::set< FVFluxKernel * > _fv_flux_kernels
FVFluxKernels.
std::set< MooseVariableFieldBase * > _elem_sub_fv_vars
virtual void operator()(const RangeType &range, bool bypass_threading=false)
std::set< MooseVariableFieldBase * > _fv_vars
Variables.
bool onBoundary(const SubdomainRestrictable &obj, const FaceInfo &fi)
Return whether the supplied face is on a boundary of the object's execution.
std::set< FVFluxKernel * > _elem_sub_fv_flux_kernels
void postFace(const FaceInfo &fi) override
This is called once for each face after all face and boundary callbacks have been finished for that f...
This loops over a set of mesh faces (i.e.
void addCached() override
call either addCachedJacobian or addCachedResidual or both depending on what derived class of this cl...
Provides a way for users to bail out of the current solve.
virtual void setup(SetupInterface &obj)=0
call either residualSetup or jacobianSetup depending on what derived class of this class is instantia...
void setup(SetupInterface &obj) override
call either residualSetup or jacobianSetup depending on what derived class of this class is instantia...
Generic class for solving transient nonlinear problems.
void caughtMooseException(MooseException &e)
Called if a MooseException is caught anywhere during the computation.
virtual void subdomainSetup(SubdomainID subdomain, const THREAD_ID tid)
const bool _scaling_residual
void finalizeContainers()
OnScopeExit & operator=(OnScopeExit &&other)
void setup(EquationSystems &systems, Mesh &mesh, GetPot &args)
virtual void subdomainChanged()
Called every time the current subdomain changes (i.e.
const bool _on_displaced
Whether this loop is operating on the displaced mesh.
static Threads::spin_mutex threaded_element_mutex
This mutex is used by all derived classes of the ThreadedElementLoop.
virtual void subdomainChanged() override
Called every time the current subdomain changes (i.e.
virtual void onFace(const FaceInfo &fi) override
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void join(const ThreadedFaceLoop &y)
void addCached() override
call either addCachedJacobian or addCachedResidual or both depending on what derived class of this cl...
const bool _scaling_jacobian
void reinitVariables(const FaceInfo &fi)
virtual void printBlockExecutionInformation() const
Print ordering of objects executed on each block.
virtual void addCached()=0
call either addCachedJacobian or addCachedResidual or both depending on what derived class of this cl...
virtual void jacobianSetup()
Gets called just before the Jacobian is computed and before this object is asked to do its job...
virtual void pre() override
This is called once before all face-looping.
OnScopeExit(std::function< void()> f) noexcept
virtual void residualSetup()
Gets called just before the residual is computed and before this object is asked to do its job...
virtual void printBlockExecutionInformation() const override
Print ordering of objects executed on each block.
bool _throw_on_error
Variable to turn on exceptions during mooseError(), should only be used within MOOSE unit tests or wh...
std::set< std::pair< const SubdomainID, const SubdomainID > > _blocks_exec_printed
Set to keep track of blocks for which we have printed the execution pattern.
SubdomainID _old_subdomain
The subdomain for the last element.
SubdomainID _old_neighbor_subdomain
The subdomain for the last neighbor.
virtual void neighborSubdomainSetup(SubdomainID subdomain, const THREAD_ID tid)
void addCached() override
call either addCachedJacobian or addCachedResidual or both depending on what derived class of this cl...