https://mooseframework.inl.gov
NonlinearThread.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "NonlinearThread.h"
11 #include "NonlinearSystem.h"
12 #include "Problem.h"
13 #include "FEProblem.h"
14 #include "KernelBase.h"
15 #include "IntegratedBCBase.h"
16 #include "DGKernelBase.h"
17 #include "InterfaceKernelBase.h"
18 #include "Material.h"
19 #include "TimeKernel.h"
20 #include "SwapBackSentinel.h"
21 #include "FVTimeKernel.h"
22 #include "ComputeJacobianThread.h"
23 
24 #include "libmesh/threads.h"
25 
27  : ThreadedElementLoop<ConstElemRange>(fe_problem),
28  _nl(fe_problem.currentNonlinearSystem()),
29  _num_cached(0),
30  _integrated_bcs(_nl.getIntegratedBCWarehouse()),
31  _dg_kernels(_nl.getDGKernelWarehouse()),
32  _interface_kernels(_nl.getInterfaceKernelWarehouse()),
33  _kernels(_nl.getKernelWarehouse()),
34  _hdg_kernels(_nl.getHDGKernelWarehouse()),
35  _has_active_objects(_integrated_bcs.hasActiveObjects() || _dg_kernels.hasActiveObjects() ||
36  _interface_kernels.hasActiveObjects() || _kernels.hasActiveObjects() ||
37  _fe_problem.haveFV()),
38  _should_execute_dg(false)
39 {
40 }
41 
42 // Splitting Constructor
45  _nl(x._nl),
46  _num_cached(x._num_cached),
47  _integrated_bcs(x._integrated_bcs),
48  _dg_kernels(x._dg_kernels),
49  _interface_kernels(x._interface_kernels),
50  _kernels(x._kernels),
51  _tag_kernels(x._tag_kernels),
52  _hdg_kernels(x._hdg_kernels),
53  _has_active_objects(x._has_active_objects),
54  _should_execute_dg(x._should_execute_dg)
55 {
56 }
57 
59 
60 void
61 NonlinearThread::operator()(const ConstElemRange & range, bool bypass_threading)
62 {
64  ThreadedElementLoop<ConstElemRange>::operator()(range, bypass_threading);
65 }
66 
67 void
69 {
70  // This should come first to setup the residual objects before we do dependency determination of
71  // material properties and variables
73 
75 
76  // Update variable Dependencies
77  std::set<MooseVariableFEBase *> needed_moose_vars;
82 
83  // Update FE variable coupleable vector tags
84  std::set<TagID> needed_fe_var_vector_tags;
86  _subdomain, needed_fe_var_vector_tags, _tid);
88  _subdomain, needed_fe_var_vector_tags, _tid);
90  _subdomain, needed_fe_var_vector_tags, _tid);
91 
92  // Update material dependencies
93  std::unordered_set<unsigned int> needed_mat_props;
98 
99  if (_fe_problem.haveFV())
100  {
101  // Re-query the finite volume elemental kernels
102  _fv_kernels.clear();
104  .query()
105  .template condition<AttribSysNum>(_nl.number())
106  .template condition<AttribSystem>("FVElementalKernel")
107  .template condition<AttribSubdomains>(_subdomain)
108  .template condition<AttribThread>(_tid)
110  for (const auto fv_kernel : _fv_kernels)
111  {
112  const auto & fv_mv_deps = fv_kernel->getMooseVariableDependencies();
113  needed_moose_vars.insert(fv_mv_deps.begin(), fv_mv_deps.end());
114  const auto & fv_mp_deps = fv_kernel->getMatPropDependencies();
115  needed_mat_props.insert(fv_mp_deps.begin(), fv_mp_deps.end());
116  }
117  }
118 
119  // Cache these to avoid computing them on every side
122 
124  _fe_problem.setActiveFEVariableCoupleableVectorTags(needed_fe_var_vector_tags, _tid);
125  _fe_problem.prepareMaterials(needed_mat_props, _subdomain, _tid);
126 }
127 
128 void
129 NonlinearThread::onElement(const Elem * const elem)
130 {
131  // Set up Sentinel class so that, even if reinitMaterials() throws in prepareElement, we
132  // still remember to swap back during stack unwinding.
134 
135  prepareElement(elem);
136 
137  if (dynamic_cast<ComputeJacobianThread *>(this))
138  if (_nl.getScalarVariables(_tid).size() > 0)
140 
142 }
143 
144 void
146 {
148  {
149  const auto & kernels = _tag_kernels->getActiveBlockObjects(_subdomain, _tid);
150  for (const auto & kernel : kernels)
151  compute(*kernel);
152  }
153 
154  if (_fe_problem.haveFV())
155  for (auto kernel : _fv_kernels)
156  compute(*kernel);
157 }
158 
159 void
160 NonlinearThread::onBoundary(const Elem * const elem,
161  const unsigned int side,
162  const BoundaryID bnd_id,
163  const Elem * const lower_d_elem /*=nullptr*/)
164 {
166  {
167  // Set up Sentinel class so that, after we swap in reinitMaterialsFace in prepareFace, even if
168  // one of our callees throws we remember to swap back during stack unwinding. We put our
169  // sentinel here as opposed to in prepareFace because we certainly don't want our materials
170  // swapped back before we proceed to residual/Jacobian computation
172 
173  prepareFace(elem, side, bnd_id, lower_d_elem);
174  computeOnBoundary(bnd_id, lower_d_elem);
175 
176  if (lower_d_elem)
177  {
178  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
179  accumulateLower();
180  }
181  }
182 }
183 
184 void
185 NonlinearThread::computeOnBoundary(BoundaryID bnd_id, const Elem * /*lower_d_elem*/)
186 {
187  const auto & bcs = _ibc_warehouse->getActiveBoundaryObjects(bnd_id, _tid);
188  for (const auto & bc : bcs)
189  if (bc->shouldApply())
190  compute(*bc);
191 }
192 
193 void
194 NonlinearThread::onInterface(const Elem * elem, unsigned int side, BoundaryID bnd_id)
195 {
197  {
198 
199  // Pointer to the neighbor we are currently working on.
200  const Elem * neighbor = elem->neighbor_ptr(side);
201 
202  if (neighbor->active())
203  {
204  _fe_problem.reinitNeighbor(elem, side, _tid);
205 
206  // Set up Sentinels so that, even if one of the reinitMaterialsXXX() calls throws, we
207  // still remember to swap back during stack unwinding. Note that face, boundary, and interface
208  // all operate with the same MaterialData object
210  _fe_problem.reinitMaterialsFaceOnBoundary(bnd_id, elem->subdomain_id(), _tid);
212 
214  _fe_problem.reinitMaterialsNeighborOnBoundary(bnd_id, neighbor->subdomain_id(), _tid);
215 
216  // Has to happen after face and neighbor properties have been computed. Note that we don't use
217  // a sentinel here because FEProblem::swapBackMaterialsFace is going to handle face materials,
218  // boundary materials, and interface materials (e.g. it queries the boundary material data
219  // with the current element and side
221 
222  computeOnInterface(bnd_id);
223 
224  {
225  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
227  }
228  }
229  }
230 }
231 
232 void
234 {
235  const auto & int_ks = _ik_warehouse->getActiveBoundaryObjects(bnd_id, _tid);
236  for (const auto & interface_kernel : int_ks)
237  compute(*interface_kernel);
238 }
239 
240 void
241 NonlinearThread::onInternalSide(const Elem * elem, unsigned int side)
242 {
243  if (_should_execute_dg)
244  {
245  // Pointer to the neighbor we are currently working on.
246  const Elem * neighbor = elem->neighbor_ptr(side);
247 
249 
250  // Set up Sentinels so that, even if one of the reinitMaterialsXXX() calls throws, we
251  // still remember to swap back during stack unwinding.
253  _fe_problem.reinitMaterialsFace(elem->subdomain_id(), _tid);
254 
256  _fe_problem.reinitMaterialsNeighbor(neighbor->subdomain_id(), _tid);
257 
258  computeOnInternalFace(neighbor);
259 
260  {
261  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
263  }
264  }
265  if (_subdomain_has_hdg)
266  {
267  // Set up Sentinel class so that, after we swap in reinitMaterialsFace in prepareFace, even if
268  // one of our callees throws we remember to swap back during stack unwinding. We put our
269  // sentinel here as opposed to in prepareFace because we certainly don't want our materials
270  // swapped back before we proceed to residual/Jacobian computation
272 
273  prepareFace(elem, side, Moose::INVALID_BOUNDARY_ID, nullptr);
275  }
276 }
277 
278 void
280 {
281  const auto & dgks = _dg_warehouse->getActiveBlockObjects(_subdomain, _tid);
282  for (const auto & dg_kernel : dgks)
283  if (dg_kernel->hasBlocks(neighbor->subdomain_id()))
284  compute(*dg_kernel, neighbor);
285 }
286 
287 void
289 {
290  compute(static_cast<ResidualObject &>(kernel));
291 }
292 
293 void
295 {
296  compute(static_cast<ResidualObject &>(kernel));
297 }
298 
299 void
301 {
302  compute(static_cast<ResidualObject &>(bc));
303 }
304 
305 void
306 NonlinearThread::compute(DGKernelBase & dg, const Elem * /*neighbor*/)
307 {
308  compute(static_cast<ResidualObject &>(dg));
309 }
310 
311 void
313 {
314  compute(static_cast<ResidualObject &>(ik));
315 }
316 
317 void
318 NonlinearThread::postElement(const Elem * /*elem*/)
319 {
320  accumulate();
321 }
322 
323 void
325 {
327 }
328 
329 void
331 {
333  {
334  const auto & console = _fe_problem.console();
335  const auto execute_on = _fe_problem.getCurrentExecuteOnFlag();
336  console << "[DBG] Beginning elemental loop to compute " + objectType() + " on " << execute_on
337  << std::endl;
338  mooseDoOnce(
339  console << "[DBG] Execution order on each element:" << std::endl;
340  console << "[DBG] - kernels on element quadrature points" << std::endl;
341  console << "[DBG] - finite volume elemental kernels on element" << std::endl;
342  console << "[DBG] - integrated boundary conditions on element side quadrature points"
343  << std::endl;
344  console << "[DBG] - DG kernels on element side quadrature points" << std::endl;
345  console << "[DBG] - interface kernels on element side quadrature points" << std::endl;);
346  }
347 }
348 
349 void
351 {
352  // Number of objects executing is approximated by size of warehouses
353  const int num_objects = _kernels.size() + _fv_kernels.size() + _integrated_bcs.size() +
355  const auto & console = _fe_problem.console();
356  const auto block_name = _mesh.getSubdomainName(_subdomain);
357 
358  if (_fe_problem.shouldPrintExecution(_tid) && num_objects > 0)
359  {
360  if (_blocks_exec_printed.count(_subdomain))
361  return;
362  console << "[DBG] Ordering of " + objectType() + " Objects on block " << block_name << " ("
363  << _subdomain << ")" << std::endl;
365  {
366  console << "[DBG] Ordering of kernels:" << std::endl;
367  console << _kernels.activeObjectsToFormattedString() << std::endl;
368  }
369  if (_fv_kernels.size())
370  {
371  console << "[DBG] Ordering of FV elemental kernels:" << std::endl;
372  std::string fvkernels =
373  std::accumulate(_fv_kernels.begin() + 1,
374  _fv_kernels.end(),
375  _fv_kernels[0]->name(),
376  [](const std::string & str_out, FVElementalKernel * kernel)
377  { return str_out + " " + kernel->name(); });
378  console << ConsoleUtils::formatString(fvkernels, "[DBG]") << std::endl;
379  }
381  {
382  console << "[DBG] Ordering of DG kernels:" << std::endl;
383  console << _dg_kernels.activeObjectsToFormattedString() << std::endl;
384  }
385  }
386  else if (_fe_problem.shouldPrintExecution(_tid) && num_objects == 0 &&
388  console << "[DBG] No Active " + objectType() + " Objects on block " << block_name << " ("
389  << _subdomain << ")" << std::endl;
390 
392 }
393 
394 void
396 {
400  return;
401 
402  const auto & console = _fe_problem.console();
403  const auto b_name = _mesh.getBoundaryName(bid);
404  console << "[DBG] Ordering of " + objectType() + " Objects on boundary " << b_name << " (" << bid
405  << ")" << std::endl;
406 
408  {
409  console << "[DBG] Ordering of integrated boundary conditions:" << std::endl;
410  console << _integrated_bcs.activeObjectsToFormattedString() << std::endl;
411  }
412 
413  // We have not checked if we have a neighbor. This could be premature for saying we are executing
414  // interface kernels. However, we should assume the execution will happen on another side of the
415  // same boundary
417  {
418  console << "[DBG] Ordering of interface kernels:" << std::endl;
419  console << _interface_kernels.activeObjectsToFormattedString() << std::endl;
420  }
421 
422  _boundaries_exec_printed.insert(bid);
423 }
424 
425 void
426 NonlinearThread::prepareFace(const Elem * const elem,
427  const unsigned int side,
428  const BoundaryID bnd_id,
429  const Elem * const lower_d_elem)
430 {
431  _fe_problem.reinitElemFace(elem, side, _tid);
432 
433  // Needed to use lower-dimensional variables on Materials
434  if (lower_d_elem)
435  _fe_problem.reinitLowerDElem(lower_d_elem, _tid);
436 
437  if (bnd_id != Moose::INVALID_BOUNDARY_ID)
438  {
439  _fe_problem.reinitMaterialsFaceOnBoundary(bnd_id, elem->subdomain_id(), _tid);
441  }
442  // Currently only used by HDG
443  else
444  _fe_problem.reinitMaterialsFace(elem->subdomain_id(), _tid);
445 }
446 
447 bool
448 NonlinearThread::shouldComputeInternalSide(const Elem & elem, const Elem & neighbor) const
449 {
450  // ThreadedElementLoop<ConstElemRange>::shouldComputeInternalSide gets expensive on high
451  // h-refinement cases so we avoid it if possible
452  _should_execute_dg = false;
453  if (_subdomain_has_dg)
457 }
std::string activeObjectsToFormattedString(THREAD_ID tid=0, const std::string &prefix="[DBG]") const
Output the active content of the warehouse to a string, meant to be output to the console...
unsigned int size(THREAD_ID tid=0) const
Return how many kernels we store in the current warehouse.
virtual bool shouldComputeInternalSide(const Elem &elem, const Elem &neighbor) const
Whether to compute the internal side for the provided element-neighbor pair.
Base class for assembly-like calculations.
virtual void accumulate()=0
Add element residual/Jacobian into assembly global data.
virtual void accumulateLower()=0
Add lower-d residual/Jacobian into assembly global data.
const std::vector< MooseVariableScalar * > & getScalarVariables(THREAD_ID tid)
Definition: SystemBase.h:756
bool hasActiveBlockObjects(THREAD_ID tid=0) const
virtual void compute(ResidualObject &ro)=0
Will dispatch to computeResidual/computeJacobian/computeResidualAndJacobian based on the derived clas...
const std::map< SubdomainID, std::vector< std::shared_ptr< T > > > & getActiveBlockObjects(THREAD_ID tid=0) const
MooseObjectWarehouse< InterfaceKernelBase > * _ik_warehouse
virtual void computeOnBoundary(BoundaryID bnd_id, const Elem *lower_d_elem)
virtual bool haveFV() const override
returns true if this problem includes/needs finite volume functionality.
const BoundaryID INVALID_BOUNDARY_ID
Definition: MooseTypes.C:22
MooseObjectTagWarehouse< IntegratedBCBase > & _integrated_bcs
Reference to BC storage structures.
void reinitMaterialsFaceOnBoundary(const BoundaryID boundary_id, const SubdomainID blk_id, const THREAD_ID tid, const bool swap_stateful=true, const std::deque< MaterialBase *> *const reinit_mats=nullptr)
reinit materials on element faces on a boundary (internal or external) This specific routine helps us...
const ExecFlagType & getCurrentExecuteOnFlag() const
Return/set the current execution flag.
std::vector< T * > & queryInto(std::vector< T *> &results, Args &&... args)
queryInto executes the query and stores the results in the given vector.
Definition: TheWarehouse.h:311
const std::string & getBoundaryName(BoundaryID boundary_id)
Return the name of the boundary given the id.
Definition: MooseMesh.C:1830
void reinitMaterialsBoundary(BoundaryID boundary_id, const THREAD_ID tid, bool swap_stateful=true, const std::deque< MaterialBase *> *reinit_mats=nullptr)
reinit materials on a boundary
const MaterialWarehouse & getMaterialWarehouse() const
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
void prepareFace(const Elem *elem, unsigned int side, BoundaryID bnd_id=Moose::INVALID_BOUNDARY_ID, const Elem *lower_d_elem=nullptr)
Reinitialize variables and materials on a face.
void updateBlockVariableDependency(SubdomainID id, std::set< MooseVariableFieldBase *> &needed_moose_vars, THREAD_ID tid=0) const
const std::string & getSubdomainName(SubdomainID subdomain_id) const
Return the name of a block given an id.
Definition: MooseMesh.C:1801
void reinitMaterialsFace(SubdomainID blk_id, const THREAD_ID tid, bool swap_stateful=true, const std::deque< MaterialBase *> *reinit_mats=nullptr)
reinit materials on element faces
bool shouldComputeInternalSide(const Elem &elem, const Elem &neighbor) const override
Whether to compute the internal side for the provided element-neighbor pair.
Serves as a base class for DGKernel and ADDGKernel.
Definition: DGKernelBase.h:32
virtual void reinitElemNeighborAndLowerD(const Elem *elem, unsigned int side, const THREAD_ID tid) override
virtual void onBoundary(const Elem *elem, unsigned int side, BoundaryID bnd_id, const Elem *lower_d_elem=nullptr) override
Called when doing boundary assembling.
virtual void accumulateNeighbor()=0
Add neighbor residual/Jacobian into assembly global data.
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
virtual void swapBackMaterialsFace(const THREAD_ID tid)
FVElemental is used for calculating residual contributions from volume integral terms of a PDE where ...
virtual void setActiveElementalMooseVariables(const std::set< MooseVariableFEBase *> &moose_vars, const THREAD_ID tid) override
Set the MOOSE variables to be reinited on each element.
const bool _has_active_objects
Whether there are any active residual objects; otherwise we will do an early return.
bool shouldPrintExecution(const THREAD_ID tid) const
Check whether the problem should output execution orders at this time.
MooseObjectWarehouse< HDGKernel > * _hdg_warehouse
This is the common base class for the three main kernel types implemented in MOOSE, Kernel, VectorKernel and ArrayKernel.
Definition: KernelBase.h:23
TheWarehouse & theWarehouse() const
void printBlockExecutionInformation() const override
Print list of specific objects executed on each block and in which order.
void reinitMaterialsNeighbor(SubdomainID blk_id, const THREAD_ID tid, bool swap_stateful=true, const std::deque< MaterialBase *> *reinit_mats=nullptr)
reinit materials on the neighboring element face
MooseObjectTagWarehouse< KernelBase > & _kernels
void reinitMaterialsNeighborOnBoundary(const BoundaryID boundary_id, const SubdomainID blk_id, const THREAD_ID tid, const bool swap_stateful=true, const std::deque< MaterialBase *> *const reinit_mats=nullptr)
reinit materials on neighbor element (usually faces) on a boundary (internal or external) This specif...
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...
Definition: ConsoleUtils.C:581
boundary_id_type BoundaryID
virtual void computeOnElement()
bool _subdomain_has_hdg
Whether the subdomain has HDGKernels.
virtual std::string objectType() const
Return what the loops is meant to compute.
NonlinearThread(FEProblemBase &fe_problem)
void printGeneralExecutionInformation() const override
Print information about the loop, mostly order of execution of objects.
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1157
virtual void computeOnInterface(BoundaryID bnd_id)
NonlinearSystemBase & _nl
virtual void postElement(const Elem *) override
Called after the element assembly is done (including surface assembling)
virtual void operator()(const RangeType &range, bool bypass_threading=false)
virtual void operator()(const ConstElemRange &range, bool bypass_threading=false) override
void updateBoundaryMatPropDependency(std::unordered_set< unsigned int > &needed_mat_props, THREAD_ID tid=0) const
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getActiveBoundaryObjects(THREAD_ID tid=0) const
virtual void swapBackMaterialsNeighbor(const THREAD_ID tid)
virtual void onInternalSide(const Elem *elem, unsigned int side) override
Called when doing internal edge assembling.
virtual void computeOnInternalFace()=0
std::set< SubdomainID > _blocks_exec_printed
Keep track of which blocks were visited.
tbb::split split
MooseObjectWarehouse< IntegratedBCBase > * _ibc_warehouse
virtual void reinitLowerDElem(const Elem *lower_d_elem, const THREAD_ID tid, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr) override
const ConsoleStream & console() const
Return console handle.
Definition: Problem.h:48
virtual void swapBackMaterials(const THREAD_ID tid)
void updateBlockMatPropDependency(SubdomainID id, std::unordered_set< unsigned int > &needed_mat_props, THREAD_ID tid=0) const
void updateBlockFEVariableCoupledVectorTagDependency(SubdomainID id, std::set< TagID > &needed_fe_var_vector_tags, THREAD_ID tid=0) const
Update FE variable coupleable vector tag vector.
std::set< BoundaryID > _boundaries_exec_printed
Keep track of which boundaries were visited.
virtual void subdomainSetup(SubdomainID subdomain, const THREAD_ID tid)
bool _should_execute_dg
Whether DG kernels should be executed for a given elem-neighbor pairing.
Query query()
query creates and returns an initialized a query object for querying objects from the warehouse...
Definition: TheWarehouse.h:466
virtual void accumulateNeighborLower()=0
Add neighbor and lower residual/Jacobian into assembly global data.
MooseObjectTagWarehouse< InterfaceKernelBase > & _interface_kernels
Reference to interface kernel storage structure.
Base class for deriving any boundary condition of a integrated type.
void printBoundaryExecutionInformation(const unsigned int bid) const override
Print list of specific objects executed on each boundary and in which order.
virtual void onInterface(const Elem *elem, unsigned int side, BoundaryID bnd_id) override
Called when doing interface assembling.
void reinitElemFace(const Elem *elem, unsigned int side, BoundaryID, const THREAD_ID tid)
InterfaceKernelBase is the base class for all InterfaceKernel type classes.
virtual void determineObjectWarehouses()=0
Determine the objects we will actually compute based on vector/matrix tag information.
SubdomainID _subdomain
The subdomain for the current element.
MooseObjectWarehouse< KernelBase > * _tag_kernels
MooseObjectWarehouse< DGKernelBase > * _dg_warehouse
virtual void reinitNeighbor(const Elem *elem, unsigned int side, const THREAD_ID tid) override
virtual void setActiveFEVariableCoupleableVectorTags(std::set< TagID > &vtags, const THREAD_ID tid) override
virtual void onElement(const Elem *elem) override
Assembly of the element (not including surface assembly)
void reinitMaterialsInterface(BoundaryID boundary_id, const THREAD_ID tid, bool swap_stateful=true)
void updateBoundaryVariableDependency(std::set< MooseVariableFieldBase *> &needed_moose_vars, THREAD_ID tid=0) const
std::vector< FVElementalKernel * > _fv_kernels
Current subdomain FVElementalKernels.
void prepareMaterials(const std::unordered_set< unsigned int > &consumer_needed_mat_props, const SubdomainID blk_id, const THREAD_ID tid)
Add the MooseVariables and the material properties that the current materials depend on to the depend...
virtual void reinitOffDiagScalars(const THREAD_ID tid) override
virtual void post() override
Called after the element range loop.
MooseObjectTagWarehouse< DGKernelBase > & _dg_kernels
Reference to DGKernel storage structure.
virtual ~NonlinearThread()
bool _subdomain_has_dg
Whether the subdomain has DGKernels.
virtual void subdomainChanged() override
Called every time the current subdomain changes (i.e.
The "SwapBackSentinel" class&#39;s destructor guarantees that FEProblemBase::swapBackMaterials{Face,Neighbor}() is called even when an exception is thrown from FEProblemBase::reinitMaterials{Face,Neighbor}.