https://mooseframework.inl.gov
ComputeFVFluxThread.h
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 #pragma once
11 
12 #include "ParallelUniqueId.h"
13 #include "MooseMesh.h"
14 #include "MooseTypes.h"
15 #include "MooseException.h"
16 #include "MooseVariableFieldBase.h"
17 #include "FVFluxKernel.h"
18 #include "FVFluxBC.h"
19 #include "FVInterfaceKernel.h"
20 #include "FEProblem.h"
21 #include "DisplacedProblem.h"
22 #include "SwapBackSentinel.h"
23 #include "MaterialBase.h"
24 #include "libmesh/libmesh_exceptions.h"
25 #include "libmesh/elem.h"
26 
27 #include <set>
28 
29 class MooseVariableFVBase;
30 
51 {
52  std::function<void()> _f;
53  OnScopeExit(std::function<void()> f) noexcept : _f(std::move(f)) {}
54  OnScopeExit(OnScopeExit && other) : _f(std::move(other._f)) {}
56  {
57  _f = std::move(other._f);
58  return *this;
59  }
61  {
62  if (_f)
63  _f();
64  }
65 };
66 
72 template <typename RangeType>
74 {
75 public:
76  ThreadedFaceLoop(FEProblemBase & fe_problem,
77  const unsigned int nl_system_num,
78  const std::set<TagID> & tags,
79  bool on_displaced);
80 
82 
83  virtual ~ThreadedFaceLoop();
84 
85  virtual void operator()(const RangeType & range, bool bypass_threading = false);
86 
87  void join(const ThreadedFaceLoop & y);
88 
89  virtual void onFace(const FaceInfo & fi) = 0;
92  virtual void postFace(const FaceInfo & /*fi*/) {}
94  virtual void pre() {}
96  virtual void post() {}
97 
100  virtual void onBoundary(const FaceInfo & fi, BoundaryID boundary) = 0;
101 
105 
109  {
111  }
112 
115  {
116  Threads::spin_mutex::scoped_lock lock(threaded_element_mutex);
117  std::string what(e.what());
119  }
120 
121 protected:
123  virtual void printGeneralExecutionInformation() const {}
124 
126  virtual void printBlockExecutionInformation() const {}
127 
129  virtual void printBoundaryExecutionInformation(const BoundaryID /* bnd_id */) const {}
130 
133  {
134  _blocks_exec_printed.clear();
135  _boundaries_exec_printed.clear();
136  }
137 
140  const std::set<TagID> & _tags;
142  const unsigned int _nl_system_num;
143 
145  const bool _on_displaced;
146 
149 
152 
155 
158 
161 
163  mutable std::set<std::pair<const SubdomainID, const SubdomainID>> _blocks_exec_printed;
164 
166  mutable std::set<BoundaryID> _boundaries_exec_printed;
167 
169  std::string _error_message;
170 
171 private:
173  const bool _zeroth_copy;
174 
178 };
179 
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()),
187  _tags(tags),
188  _nl_system_num(nl_system_num),
189  _on_displaced(on_displaced),
190  _subproblem(_on_displaced ? static_cast<SubProblem &>(*_fe_problem.getDisplacedProblem())
191  : static_cast<SubProblem &>(_fe_problem)),
192  _zeroth_copy(true),
193  _incoming_throw_on_error(Moose::_throw_on_error)
194 {
195  Moose::_throw_on_error = true;
196 }
197 
198 template <typename RangeType>
200  : _fe_problem(x._fe_problem),
201  _mesh(x._mesh),
202  _tags(x._tags),
203  _nl_system_num(x._nl_system_num),
204  _on_displaced(x._on_displaced),
205  _subproblem(x._subproblem),
206  _zeroth_copy(false),
207  _incoming_throw_on_error(false)
208 {
209 }
210 
211 template <typename RangeType>
213 {
214  if (_zeroth_copy)
215  {
216  Moose::_throw_on_error = _incoming_throw_on_error;
217 
218  if (!_error_message.empty())
219  mooseError(_error_message);
220  }
221 }
222 
223 template <typename RangeType>
224 void
226 {
227  if (_error_message.empty() && !y._error_message.empty())
228  _error_message = y._error_message;
229 }
230 
231 // TODO: ensure the vector<faceinfo> data structure needs to be built such
232 // that for all sides on an interface between two subdomains, the elements of
233 // the same subdomain are used consistently for all the "elem" (i.e. not
234 // "neighbor") parameters in order to avoid jumping back and forth along the
235 // boundary between using one or the other subdomains' FV kernels
236 // unpredictably.
237 template <typename RangeType>
238 void
239 ThreadedFaceLoop<RangeType>::operator()(const RangeType & range, bool bypass_threading)
240 {
241  // TODO: make this query fv flux kernel specific or somehow integrate the
242  // fv source kernels into this loop. Also this check will need to increase
243  // in generality if/when other systems and objects besides FV stuff get
244  // added to this loop.
245  std::vector<FVKernel *> kernels;
246  _fe_problem.theWarehouse()
247  .query()
248  .template condition<AttribSysNum>(_nl_system_num)
249  .template condition<AttribSystem>("FVFluxKernel")
250  .template condition<AttribDisplaced>(_on_displaced)
251  .queryInto(kernels);
252  if (kernels.size() == 0)
253  return;
254 
255  try
256  {
257  try
258  {
259  ParallelUniqueId puid;
260  _tid = bypass_threading ? 0 : puid.id;
261 
262  pre();
263  printGeneralExecutionInformation();
264 
265  _subdomain = Moose::INVALID_BLOCK_ID;
266  _neighbor_subdomain = Moose::INVALID_BLOCK_ID;
267 
268  typename RangeType::const_iterator faceinfo = range.begin();
269  for (faceinfo = range.begin(); faceinfo != range.end(); ++faceinfo)
270  {
271  const Elem & elem = (*faceinfo)->elem();
272 
273  _fe_problem.setCurrentSubdomainID(&elem, _tid);
274 
275  _old_subdomain = _subdomain;
276  _subdomain = elem.subdomain_id();
277  if (_subdomain != _old_subdomain)
278  {
279  subdomainChanged();
280  printBlockExecutionInformation();
281  }
282 
283  _old_neighbor_subdomain = _neighbor_subdomain;
284  if (const Elem * const neighbor = (*faceinfo)->neighborPtr())
285  {
286  _fe_problem.setNeighborSubdomainID(neighbor, _tid);
287  _neighbor_subdomain = neighbor->subdomain_id();
288  }
289  else
290  _neighbor_subdomain = Moose::INVALID_BLOCK_ID;
291 
292  if (_neighbor_subdomain != _old_neighbor_subdomain)
293  {
294  neighborSubdomainChanged();
295  // This is going to cause a lot more printing
296  printBlockExecutionInformation();
297  }
298 
299  onFace(**faceinfo);
300  // Cache data now because onBoundary may clear it. E.g. there was a nasty bug for two
301  // variable FV systems where if one variable was executing an FVFluxKernel on a boundary
302  // while the other was executing an FVFluxBC, the FVFluxKernel data would get lost because
303  // onBoundary would clear the residual/Jacobian data before it was cached
304  postFace(**faceinfo);
305 
306  const std::set<BoundaryID> boundary_ids = (*faceinfo)->boundaryIDs();
307  for (auto & it : boundary_ids)
308  {
309  printBoundaryExecutionInformation(it);
310  onBoundary(**faceinfo, it);
311  }
312 
313  postFace(**faceinfo);
314 
315  } // range
316  post();
317 
318  // Clear execution printing sets to start printing on every block and boundary again
319  resetExecutionPrinting();
320  }
321  catch (libMesh::LogicError & e)
322  {
323  mooseException("We caught a libMesh error: ", e.what());
324  }
325  catch (MetaPhysicL::LogicError & e)
326  {
328  }
329  }
330  catch (MooseException & e)
331  {
332  caughtMooseException(e);
333  }
334  catch (std::runtime_error & e)
335  {
336  _error_message = e.what();
337  }
338 }
339 
343 template <typename RangeType, typename AttributeTagType>
344 class ComputeFVFluxThread : public ThreadedFaceLoop<RangeType>
345 {
346 public:
347  ComputeFVFluxThread(FEProblemBase & fe_problem,
348  const unsigned int nl_system_num,
349  const std::set<TagID> & tags,
350  bool on_displaced);
351 
353 
354  virtual ~ComputeFVFluxThread();
355 
356  virtual void onFace(const FaceInfo & fi) override;
357  virtual void pre() override;
358  virtual void post() override;
359 
360  virtual void onBoundary(const FaceInfo & fi, BoundaryID boundary) override;
361 
362  virtual void subdomainChanged() override;
363  virtual void neighborSubdomainChanged() override;
364 
365 protected:
370  virtual void compute(FVFaceResidualObject & ro, const FaceInfo & fi) = 0;
371 
376  virtual void setup(SetupInterface & obj) = 0;
377 
382  virtual void addCached() = 0;
383 
384  unsigned int _num_cached = 0;
385 
396 
397 private:
398  void reinitVariables(const FaceInfo & fi);
399  void finalizeContainers();
400 
402  virtual void printGeneralExecutionInformation() const override;
403 
405  virtual void printBlockExecutionInformation() const override;
406 
408  virtual void printBoundaryExecutionInformation(const BoundaryID bnd_id) const override;
409 
411  std::pair<SubdomainName, SubdomainName> getBlockNames() const;
412 
414  std::set<MooseVariableFieldBase *> _fv_vars;
415  std::set<MooseVariableFieldBase *> _elem_sub_fv_vars;
416  std::set<MooseVariableFieldBase *> _neigh_sub_fv_vars;
417 
419  std::set<FVFluxKernel *> _fv_flux_kernels;
420  std::set<FVFluxKernel *> _elem_sub_fv_flux_kernels;
421  std::set<FVFluxKernel *> _neigh_sub_fv_flux_kernels;
422 
424  std::vector<std::shared_ptr<MaterialBase>> _elem_face_mats;
425  std::vector<std::shared_ptr<MaterialBase>> _elem_sub_elem_face_mats;
426  std::vector<std::shared_ptr<MaterialBase>> _neigh_sub_elem_face_mats;
427 
428  // Neighbor face materials
429  std::vector<std::shared_ptr<MaterialBase>> _neigh_face_mats;
430  std::vector<std::shared_ptr<MaterialBase>> _elem_sub_neigh_face_mats;
431  std::vector<std::shared_ptr<MaterialBase>> _neigh_sub_neigh_face_mats;
432 
433  const bool _scaling_jacobian;
434  const bool _scaling_residual;
435 };
436 
437 template <typename RangeType, typename AttributeTagType>
439  unsigned int nl_system_num,
440  const std::set<TagID> & tags,
441  bool on_displaced)
442  : ThreadedFaceLoop<RangeType>(fe_problem, nl_system_num, tags, on_displaced),
443  _scaling_jacobian(fe_problem.computingScalingJacobian()),
444  _scaling_residual(fe_problem.computingScalingResidual())
445 {
446 }
447 
448 template <typename RangeType, typename AttributeTagType>
451  : ThreadedFaceLoop<RangeType>(x, split),
452  _fv_vars(x._fv_vars),
453  _scaling_jacobian(x._scaling_jacobian),
454  _scaling_residual(x._scaling_residual)
455 {
456 }
457 
458 template <typename RangeType, typename AttributeTagType>
460 {
461 }
462 
463 template <typename RangeType, typename AttributeTagType>
464 void
466 {
467  // TODO: this skips necessary FE reinit. In addition to this call, we need
468  // to conditionally do some FE-specific reinit here if we have any active FE
469  // variables. However, we still want to keep/do FV-style quadrature.
470  // Figure out how to do all this some day.
471  this->_subproblem.reinitFVFace(_tid, fi);
472 
473  // TODO: for FE variables, this is handled via setting needed vars through
474  // fe problem API which passes the value on to the system class. Then
475  // reinit is called on fe problem which forwards its calls to the system
476  // function and then fe problem also calls displaced problem reinit. All
477  // the call forwarding seems silly, but it does allow the displaced problem
478  // to be easily kept in sync. However, the displaced problem has different
479  // pointers for its own face info objects, etc, we can't just pass in fe
480  // problem's face info down to the sub-problem -- centroids are different,
481  // volumes are different, etc. To support displaced meshes correctly, we
482  // need to be able to reinit the subproblems variables using its equivalent
483  // face info object. How? See https://github.com/idaholab/moose/issues/15064
484 
485  for (auto var : _fv_vars)
486  var->computeFaceValues(fi);
487 
488  _fe_problem.resizeMaterialData(Moose::MaterialDataType::FACE_MATERIAL_DATA, /*nqp=*/1, _tid);
489 
490  for (std::shared_ptr<MaterialBase> mat : _elem_face_mats)
491  {
492  mat->setFaceInfo(fi);
493  mat->computeProperties();
494  }
495 
496  if (fi.neighborPtr())
497  {
498  _fe_problem.resizeMaterialData(
500 
501  for (std::shared_ptr<MaterialBase> mat : _neigh_face_mats)
502  {
503  mat->setFaceInfo(fi);
504  mat->computeProperties();
505  }
506  }
507 }
508 
509 template <typename RangeType, typename AttributeTagType>
510 void
512 {
513  reinitVariables(fi);
514 
515  for (auto * const k : _fv_flux_kernels)
516  compute(*k, fi);
517 }
518 
519 template <typename RangeType, typename AttributeTagType>
520 void
522 {
523  if (_scaling_jacobian || _scaling_residual)
524  return;
525 
526  std::vector<FVFluxBC *> bcs;
527  _fe_problem.theWarehouse()
528  .query()
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)
535  .queryInto(bcs);
536 
537  for (auto * const bc : bcs)
538  compute(*bc, fi);
539 
540  std::vector<FVInterfaceKernel *> iks;
541  _fe_problem.theWarehouse()
542  .query()
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)
549  .queryInto(iks);
550 
551  for (auto * const ik : iks)
552  compute(*ik, fi);
553 }
554 
555 template <typename RangeType, typename AttributeTagType>
556 void
558 {
559  // make sure we add any remaining cached residuals/jacobians to add/record
560  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
561  addCached();
562 
563  _fe_problem.clearActiveElementalMooseVariables(_tid);
564  _fe_problem.clearActiveMaterialProperties(_tid);
565 }
566 
567 template <typename RangeType, typename AttributeTagType>
568 void
570 {
571  //
572  // Finalize our variables
573  //
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()));
579 
580  //
581  // Finalize our kernels
582  //
583  const bool same_kernels = _elem_sub_fv_flux_kernels == _neigh_sub_fv_flux_kernels;
584  if (same_kernels)
585  _fv_flux_kernels = _elem_sub_fv_flux_kernels;
586  else
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;
593 
594  //
595  // Finalize our element face materials
596  //
597  _elem_face_mats = _elem_sub_elem_face_mats;
598 
599  if (need_ghosting)
600  // Add any element face materials from the neighboring subdomain that do not exist on the
601  // element subdomain
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);
607 
608  //
609  // Finalize our neighbor face materials
610  //
611  _neigh_face_mats = _neigh_sub_neigh_face_mats;
612 
613  if (need_ghosting)
614  // Add any neighbor face materials from the element subdomain that do not exist on the
615  // neighbor subdomain
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);
621 }
622 
623 template <typename RangeType, typename AttributeTagType>
624 void
626 {
628 
629  // Clear variables
630  _fv_vars.clear();
631  _elem_sub_fv_vars.clear();
632 
633  // Clear kernels
634  _fv_flux_kernels.clear();
635  _elem_sub_fv_flux_kernels.clear();
636 
637  // Clear element face materials
638  _elem_face_mats.clear();
639  _elem_sub_elem_face_mats.clear();
640 
641  // Clear neighbor face materials
642  _neigh_face_mats.clear();
643  _elem_sub_neigh_face_mats.clear();
644 
645  // TODO: do this for other relevant objects - like FV BCs, FV source term
646  // kernels, etc. - but we don't need to add them for other types of objects
647  // like FE or DG kernels because those kernels don't run in this loop. Do we
648  // really want to integrate fv source kernels into this loop?
649  std::vector<FVFluxKernel *> kernels;
650  _fe_problem.theWarehouse()
651  .query()
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)
658  .queryInto(kernels);
659 
660  _elem_sub_fv_flux_kernels = std::set<FVFluxKernel *>(kernels.begin(), kernels.end());
661 
662  for (auto * k : _elem_sub_fv_flux_kernels)
663  {
664  // TODO: we need a better way to do this - especially when FE objects begin to
665  // couple to FV vars. This code shoud be refactored out into one place
666  // where it is easy for users to say initialize all materials and
667  // variables needed by these objects for me.
668  const auto & deps = k->getMooseVariableDependencies();
669  for (auto var : deps)
670  {
671  mooseAssert(var->isFV(),
672  "We do not currently support coupling of FE variables into FV objects");
673  _elem_sub_fv_vars.insert(var);
674  }
675  }
676 
677  _fe_problem.getFVMatsAndDependencies(
678  _subdomain, _elem_sub_elem_face_mats, _elem_sub_neigh_face_mats, _elem_sub_fv_vars, _tid);
679 
680  finalizeContainers();
681 }
682 
683 template <typename RangeType, typename AttributeTagType>
684 void
686 {
688 
689  // Clear variables
690  _fv_vars.clear();
691  _neigh_sub_fv_vars.clear();
692 
693  // Clear kernels
694  _fv_flux_kernels.clear();
695  _neigh_sub_fv_flux_kernels.clear();
696 
697  // Clear element face materials
698  _elem_face_mats.clear();
699  _neigh_sub_elem_face_mats.clear();
700 
701  // Clear neighbor face materials
702  _neigh_face_mats.clear();
703  _neigh_sub_neigh_face_mats.clear();
704 
705  // TODO: do this for other relevant objects - like FV BCs, FV source term
706  // kernels, etc. - but we don't need to add them for other types of objects
707  // like FE or DG kernels because those kernels don't run in this loop. Do we
708  // really want to integrate fv source kernels into this loop?
709  std::vector<FVFluxKernel *> kernels;
710  _fe_problem.theWarehouse()
711  .query()
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)
718  .queryInto(kernels);
719 
720  _neigh_sub_fv_flux_kernels = std::set<FVFluxKernel *>(kernels.begin(), kernels.end());
721 
722  for (auto * k : _neigh_sub_fv_flux_kernels)
723  {
724  // TODO: we need a better way to do this - especially when FE objects begin to
725  // couple to FV vars. This code shoud be refactored out into one place
726  // where it is easy for users to say initialize all materials and
727  // variables needed by these objects for me.
728  const auto & deps = k->getMooseVariableDependencies();
729  for (auto var : deps)
730  {
731  mooseAssert(var->isFV(),
732  "We do not currently support coupling of FE variables into FV objects");
733  _neigh_sub_fv_vars.insert(var);
734  }
735  }
736 
737  _fe_problem.getFVMatsAndDependencies(_neighbor_subdomain,
738  _neigh_sub_elem_face_mats,
739  _neigh_sub_neigh_face_mats,
740  _neigh_sub_fv_vars,
741  _tid);
742 
743  finalizeContainers();
744 }
745 
746 template <typename RangeType, typename AttributeTagType>
747 void
749 {
750  std::vector<FVFluxBC *> bcs;
751  _fe_problem.theWarehouse()
752  .query()
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)
758  .queryInto(bcs);
759 
760  std::vector<FVInterfaceKernel *> iks;
761  _fe_problem.theWarehouse()
762  .query()
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)
768  .queryInto(iks);
769 
770  std::vector<FVFluxKernel *> kernels;
771  _fe_problem.theWarehouse()
772  .query()
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)
778  .queryInto(kernels);
779 
780  for (auto * const bc : bcs)
781  setup(*bc);
782  for (auto * const ik : iks)
783  setup(*ik);
784  for (auto * const kernel : kernels)
785  setup(*kernel);
786 
787  // Clear variables
788  _fv_vars.clear();
789  _elem_sub_fv_vars.clear();
790  _neigh_sub_fv_vars.clear();
791 
792  // Clear kernels
793  _fv_flux_kernels.clear();
794  _elem_sub_fv_flux_kernels.clear();
795  _neigh_sub_fv_flux_kernels.clear();
796 
797  // Clear element face materials
798  _elem_face_mats.clear();
799  _elem_sub_elem_face_mats.clear();
800  _neigh_sub_elem_face_mats.clear();
801 
802  // Clear neighbor face materials
803  _neigh_face_mats.clear();
804  _elem_sub_neigh_face_mats.clear();
805  _neigh_sub_neigh_face_mats.clear();
806 }
807 
808 template <typename RangeType>
809 class ComputeFVFluxResidualThread : public ComputeFVFluxThread<RangeType, AttribVectorTags>
810 {
811 public:
813  const unsigned int nl_system_num,
814  const std::set<TagID> & tags,
815  bool on_displaced);
816 
818 
819 protected:
825 
826  void postFace(const FaceInfo & fi) override;
827  void compute(FVFaceResidualObject & ro, const FaceInfo & fi) override { ro.computeResidual(fi); }
828  void setup(SetupInterface & obj) override { obj.residualSetup(); }
829  void addCached() override { this->_subproblem.SubProblem::addCachedResidual(_tid); }
830 };
831 
832 template <typename RangeType>
834  FEProblemBase & fe_problem,
835  const unsigned int nl_system_num,
836  const std::set<TagID> & tags,
837  bool on_displaced)
838  : ComputeFVFluxThread<RangeType, AttribVectorTags>(fe_problem, nl_system_num, tags, on_displaced)
839 {
840 }
841 
842 template <typename RangeType>
845  : ComputeFVFluxThread<RangeType, AttribVectorTags>(x, split)
846 {
847 }
848 
849 template <typename RangeType>
850 void
852 {
853  _num_cached++;
854  // TODO: do we need both calls - or just the neighbor one? - confirm this
855  this->_subproblem.SubProblem::cacheResidual(_tid);
856  this->_subproblem.SubProblem::cacheResidualNeighbor(_tid);
857 
858  this->_subproblem.SubProblem::addCachedResidual(_tid);
859 }
860 
861 template <typename RangeType>
862 class ComputeFVFluxJacobianThread : public ComputeFVFluxThread<RangeType, AttribMatrixTags>
863 {
864 public:
866  const unsigned int nl_system_num,
867  const std::set<TagID> & tags,
868  bool on_displaced);
869 
871 
872 protected:
876 
877  void postFace(const FaceInfo & fi) override;
878  void compute(FVFaceResidualObject & ro, const FaceInfo & fi) override { ro.computeJacobian(fi); }
879  void setup(SetupInterface & obj) override { obj.jacobianSetup(); }
880  void addCached() override { this->_subproblem.SubProblem::addCachedJacobian(_tid); }
881 };
882 
883 template <typename RangeType>
885  FEProblemBase & fe_problem,
886  const unsigned int nl_system_num,
887  const std::set<TagID> & tags,
888  bool on_displaced)
889  : ComputeFVFluxThread<RangeType, AttribMatrixTags>(fe_problem, nl_system_num, tags, on_displaced)
890 {
891 }
892 
893 template <typename RangeType>
896  : ComputeFVFluxThread<RangeType, AttribMatrixTags>(x, split)
897 {
898 }
899 
900 template <typename RangeType>
901 void
903 {
904  _num_cached++;
905  // FV objects do not store their Jacobian data in TaggingInterface data structures; instead they
906  // add into the cache directly. So we do not need the same cache calls for the Jacobian that we
907  // need for the residual
908  if (_num_cached % 20 == 0)
909  {
910  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
911  this->_subproblem.SubProblem::addCachedJacobian(_tid);
912  }
913 }
914 
915 template <typename RangeType>
916 class ComputeFVFluxRJThread : public ComputeFVFluxThread<RangeType, AttribVectorTags>
917 {
918 public:
920  const unsigned int nl_system_num,
921  const std::set<TagID> & vector_tags,
922  const std::set<TagID> & /*matrix_tags*/,
923  bool on_displaced);
924 
926 
927 protected:
931 
932  void postFace(const FaceInfo & fi) override;
933  void compute(FVFaceResidualObject & ro, const FaceInfo & fi) override
934  {
936  }
937  void setup(SetupInterface & obj) override { obj.residualSetup(); }
938  void addCached() override
939  {
940  this->_subproblem.SubProblem::addCachedResidual(_tid);
941  this->_subproblem.SubProblem::addCachedJacobian(_tid);
942  }
943 };
944 
945 template <typename RangeType>
947  const unsigned int nl_system_num,
948  const std::set<TagID> & vector_tags,
949  const std::set<TagID> & /*matrix_tags*/,
950  bool on_displaced)
952  fe_problem, nl_system_num, vector_tags, on_displaced)
953 {
954 }
955 
956 template <typename RangeType>
959  : ComputeFVFluxThread<RangeType, AttribVectorTags>(x, split)
960 {
961 }
962 
963 template <typename RangeType>
964 void
966 {
967  _num_cached++;
968  // TODO: do we need both calls - or just the neighbor one? - confirm this
969  this->_subproblem.SubProblem::cacheResidual(_tid);
970  this->_subproblem.SubProblem::cacheResidualNeighbor(_tid);
971  // FV objects do not store their Jacobian data in TaggingInterface data structures; instead they
972  // add into the cache directly. So we do not need the same cache calls for the Jacobian that we
973  // need for the residual
974  if (_num_cached % 20 == 0)
975  {
976  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
977  this->_subproblem.SubProblem::addCachedResidual(_tid);
978  this->_subproblem.SubProblem::addCachedJacobian(_tid);
979  }
980 }
981 
982 template <typename RangeType, typename AttributeTagType>
983 void
985 {
986  if (!_fe_problem.shouldPrintExecution(_tid))
987  return;
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: "
992  << std::endl;
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;);
996 }
997 
998 template <typename RangeType, typename AttributeTagType>
999 void
1001 {
1002  if (!_fe_problem.shouldPrintExecution(_tid) || !_fv_flux_kernels.size())
1003  return;
1004 
1005  // Print the location of the execution
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))
1009  return;
1010  auto & console = _fe_problem.console();
1011  console << "[DBG] Flux kernels on block " << block_pair_names.first;
1012  if (_neighbor_subdomain != Moose::INVALID_BLOCK_ID)
1013  console << " and neighbor " << block_pair_names.second << std::endl;
1014  else
1015  console << " with no neighbor block" << std::endl;
1016 
1017  // Print the list of objects
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));
1022  "[DBG]")
1023  << std::endl;
1024  _blocks_exec_printed.insert(block_pair);
1025 }
1026 
1027 template <typename RangeType, typename AttributeTagType>
1028 void
1030  const BoundaryID bnd_id) const
1031 {
1032  if (!_fe_problem.shouldPrintExecution(_tid))
1033  return;
1034  if (_boundaries_exec_printed.count(bnd_id))
1035  return;
1036  std::vector<MooseObject *> bcs;
1037  _fe_problem.theWarehouse()
1038  .query()
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)
1044  .queryInto(bcs);
1045 
1046  std::vector<MooseObject *> iks;
1047  _fe_problem.theWarehouse()
1048  .query()
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)
1054  .queryInto(iks);
1055 
1056  const auto block_pair_names = this->getBlockNames();
1057  if (bcs.size())
1058  {
1059  auto & console = _fe_problem.console();
1060  console << "[DBG] FVBCs on boundary " << bnd_id << " between subdomain "
1061  << block_pair_names.first;
1062  if (_neighbor_subdomain != Moose::INVALID_BLOCK_ID)
1063  console << " and neighbor " << block_pair_names.second << std::endl;
1064  else
1065  console << " and the exterior of the mesh " << std::endl;
1066  const std::string fv_bcs = ConsoleUtils::mooseObjectVectorToString(bcs);
1067  console << ConsoleUtils::formatString(fv_bcs, "[DBG]") << std::endl;
1068  }
1069  if (iks.size())
1070  {
1071  auto & console = _fe_problem.console();
1072  console << "[DBG] FVIKs on boundary " << bnd_id << " between subdomain "
1073  << block_pair_names.first;
1074  if (_neighbor_subdomain != Moose::INVALID_BLOCK_ID)
1075  console << " and neighbor " << block_pair_names.second << std::endl;
1076  else
1077  console << " and the exterior of the mesh " << std::endl;
1078  const std::string fv_iks = ConsoleUtils::mooseObjectVectorToString(iks);
1079  console << ConsoleUtils::formatString(fv_iks, "[DBG]") << std::endl;
1080  }
1081  _boundaries_exec_printed.insert(bnd_id);
1082 }
1083 
1084 template <typename RangeType, typename AttributeTagType>
1085 std::pair<SubdomainName, SubdomainName>
1087 {
1088  auto block_names = std::make_pair(_mesh.getSubdomainName(_subdomain),
1089  _mesh.getSubdomainName(_neighbor_subdomain));
1090  if (block_names.first == "")
1091  block_names.first = Moose::stringify(_subdomain);
1092  if (block_names.second == "")
1093  block_names.second = Moose::stringify(_neighbor_subdomain);
1094  return block_names;
1095 }
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 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.
Definition: ConsoleUtils.C:598
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...
Definition: MooseError.h:302
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
Definition: MooseError.C:112
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
MeshBase & mesh
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
Definition: MooseTypes.C:20
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 ...
Definition: FaceInfo.h:36
Base class for assembly-like calculations.
const Elem * neighborPtr() const
Definition: FaceInfo.h:84
ComputeFVFluxResidualThread(FEProblemBase &fe_problem, const unsigned int nl_system_num, const std::set< TagID > &tags, bool on_displaced)
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:582
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...
Definition: MooseMesh.h:88
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
Definition: Conversion.h:64
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&#39;s execution.
Definition: MathFVUtils.h:1169
tbb::split split
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.
Definition: SubProblem.h:78
void caughtMooseException(MooseException &e)
Called if a MooseException is caught anywhere during the computation.
virtual void subdomainSetup(SubdomainID subdomain, const THREAD_ID tid)
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...
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...
Definition: Moose.C:762
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.
unsigned int THREAD_ID
Definition: MooseTypes.h:209
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...