ThreadedElementLoopBase.h

Go to the documentation of this file.

//* This file is part of the MOOSE framework
//* https://mooseframework.inl.gov
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html
 
#pragma once
 
#include "ParallelUniqueId.h"
#include "MooseMesh.h"
#include "MooseTypes.h"
#include "MooseException.h"
#include "libmesh/libmesh_exceptions.h"
#include "libmesh/elem.h"
 
// C++
#include <cstring> // for "Jacobian" exception test
 
template <typename RangeType>
class ThreadedElementLoopBase
{
public:
  ThreadedElementLoopBase(MooseMesh & mesh);
 
  ThreadedElementLoopBase(ThreadedElementLoopBase & x, Threads::split split);
 
  virtual ~ThreadedElementLoopBase();
 
  virtual void operator()(const RangeType & range, bool bypass_threading = false);
 
  virtual void pre();
 
  virtual void post();
 
  virtual void onElement(const Elem * elem);
 
  virtual void preElement(const Elem * elem);
 
  virtual void postElement(const Elem * elem);
 
  virtual void preBoundary(const Elem * elem,
                           unsigned int side,
                           BoundaryID bnd_id,
                           const Elem * lower_d_elem = nullptr);
 
  virtual void onBoundary(const Elem * elem,
                          unsigned int side,
                          BoundaryID bnd_id,
                          const Elem * lower_d_elem = nullptr);
 
  virtual void preInternalSide(const Elem * elem, unsigned int side);
 
  virtual void postInternalSide(const Elem * elem, unsigned int side);
 
  virtual void onInternalSide(const Elem * elem, unsigned int side);
 
  virtual void onExternalSide(const Elem * elem, unsigned int side);
 
  virtual void onInterface(const Elem * elem, unsigned int side, BoundaryID bnd_id);
 
  virtual void subdomainChanged();
 
  virtual void neighborSubdomainChanged();
 
  virtual void caughtMooseException(MooseException &) {}
 
  virtual bool keepGoing() { return true; }
 
protected:
  MooseMesh & _mesh;
  THREAD_ID _tid;
 
  SubdomainID _subdomain;
 
  SubdomainID _old_subdomain;
 
  SubdomainID _neighbor_subdomain;
 
  SubdomainID _old_neighbor_subdomain;
 
  virtual void printGeneralExecutionInformation() const {}
 
  virtual void printBlockExecutionInformation() const {}
 
  virtual void printBoundaryExecutionInformation(const unsigned int /*bid*/) const {}
 
  mutable std::set<SubdomainID> _blocks_exec_printed;
 
  mutable std::set<BoundaryID> _boundaries_exec_printed;
 
  void resetExecPrintedSets() const;
 
  virtual bool shouldComputeInternalSide(const Elem & elem, const Elem & neighbor) const;
};
 
template <typename RangeType>
ThreadedElementLoopBase<RangeType>::ThreadedElementLoopBase(MooseMesh & mesh) : _mesh(mesh)
{
}
 
template <typename RangeType>
ThreadedElementLoopBase<RangeType>::ThreadedElementLoopBase(ThreadedElementLoopBase & x,
                                                            Threads::split /*split*/)
  : _mesh(x._mesh)
{
}
 
template <typename RangeType>
ThreadedElementLoopBase<RangeType>::~ThreadedElementLoopBase()
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::operator()(const RangeType & range, bool bypass_threading)
{
  try
  {
    try
    {
      ParallelUniqueId puid;
      _tid = bypass_threading ? 0 : puid.id;
 
      pre();
      printGeneralExecutionInformation();
 
      _subdomain = Moose::INVALID_BLOCK_ID;
      _neighbor_subdomain = Moose::INVALID_BLOCK_ID;
      typename RangeType::const_iterator el = range.begin();
      for (el = range.begin(); el != range.end(); ++el)
      {
        if (!keepGoing())
          break;
 
        const Elem * elem = *el;
 
        preElement(elem);
 
        _old_subdomain = _subdomain;
        _subdomain = elem->subdomain_id();
        if (_subdomain != _old_subdomain)
        {
          subdomainChanged();
          printBlockExecutionInformation();
        }
 
        onElement(elem);
 
        if (_mesh.interiorLowerDBlocks().count(elem->subdomain_id()) > 0 ||
            _mesh.boundaryLowerDBlocks().count(elem->subdomain_id()) > 0)
        {
          postElement(elem);
          continue;
        }
 
        const auto elem_boundary_ids = _mesh.getBoundaryIDs(elem);
        for (unsigned int side = 0; side < elem->n_sides(); side++)
        {
          const auto & boundary_ids = elem_boundary_ids[side];
          const Elem * lower_d_elem = _mesh.getLowerDElem(elem, side);
 
          for (const auto bnd_id : boundary_ids)
          {
            preBoundary(elem, side, bnd_id, lower_d_elem);
            printBoundaryExecutionInformation(bnd_id);
            onBoundary(elem, side, bnd_id, lower_d_elem);
          }
 
          const Elem * neighbor = elem->neighbor_ptr(side);
          if (neighbor)
          {
            preInternalSide(elem, side);
 
            _old_neighbor_subdomain = _neighbor_subdomain;
            _neighbor_subdomain = neighbor->subdomain_id();
            if (_neighbor_subdomain != _old_neighbor_subdomain)
              neighborSubdomainChanged();
 
            if (shouldComputeInternalSide(*elem, *neighbor))
              onInternalSide(elem, side);
 
            for (const auto bnd_id : boundary_ids)
              onInterface(elem, side, bnd_id);
 
            postInternalSide(elem, side);
          }
          else
            onExternalSide(elem, side);
        } // sides
 
        postElement(elem);
      } // range
 
      post();
      resetExecPrintedSets();
    }
    catch (MetaPhysicL::LogicError & e)
    {
      moose::translateMetaPhysicLError(e);
    }
    catch (std::exception & e)
    {
      // Continue if we find a libMesh degenerate map exception, but
      // just re-throw for any real error
      if (!strstr(e.what(), "Jacobian") && !strstr(e.what(), "singular") &&
          !strstr(e.what(), "det != 0"))
        throw; // not "throw e;" - that destroys type info!
 
      mooseException("We caught a libMesh degeneracy exception in ThreadedElementLoopBase:\n",
                     e.what());
    }
  }
  catch (MooseException & e)
  {
    caughtMooseException(e);
  }
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::pre()
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::post()
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::onElement(const Elem * /*elem*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::preElement(const Elem * /*elem*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::postElement(const Elem * /*elem*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::preBoundary(const Elem * /*elem*/,
                                                unsigned int /*side*/,
                                                BoundaryID /*bnd_id*/,
                                                const Elem * /*lower_d_elem = nullptr*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::onBoundary(const Elem * /*elem*/,
                                               unsigned int /*side*/,
                                               BoundaryID /*bnd_id*/,
                                               const Elem * /*lower_d_elem = nullptr*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::preInternalSide(const Elem * /*elem*/, unsigned int /*side*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::postInternalSide(const Elem * /*elem*/, unsigned int /*side*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::onInternalSide(const Elem * /*elem*/, unsigned int /*side*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::onExternalSide(const Elem * /*elem*/, unsigned int /*side*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::onInterface(const Elem * /*elem*/,
                                                unsigned int /*side*/,
                                                BoundaryID /*bnd_id*/)
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::subdomainChanged()
{
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::neighborSubdomainChanged()
{
}
 
template <typename RangeType>
bool
ThreadedElementLoopBase<RangeType>::shouldComputeInternalSide(const Elem & elem,
                                                              const Elem & neighbor) const
{
  // If we're going to compute the internal side with this elem-neighbor pair, then they must both
  // be active. Note that if elem is an active coarse element at an interface with finer elements,
  // then its neighbor will be an equal-level inactive element, hence the following
  // neighbor.active() check. In that case we'll catch current 'elem' when we come around and
  // examine one of the finer elements as 'elem'
  mooseAssert(elem.active(), "This method should never be called with an inactive element");
  if (!neighbor.active())
    return false;
 
  //
  // Define an ordering: first prefer finer by h (higher level()), then prefer finer by p (higher
  // p_level()), then prefer smaller id()"
  //
 
  const auto elem_level = elem.level(), neighbor_level = neighbor.level();
  if (elem_level != neighbor_level)
    return elem_level > neighbor_level;
 
  const auto elem_p_level = elem.p_level(), neighbor_p_level = neighbor.p_level();
  if (elem_p_level != neighbor_p_level)
    return elem_p_level > neighbor_p_level;
 
  return elem.id() < neighbor.id();
}
 
template <typename RangeType>
void
ThreadedElementLoopBase<RangeType>::resetExecPrintedSets() const
{
  _blocks_exec_printed.clear();
  _boundaries_exec_printed.clear();
}