libMesh::PatchRecoveryErrorEstimator Class Reference

This class implements the Patch Recovery error indicator. More...

#include <patch_recovery_error_estimator.h>

Inheritance diagram for libMesh::PatchRecoveryErrorEstimator:

Classes
class	EstimateError
	Class to compute the error contribution for a range of elements. More...

Public Types
typedef std::map< std::pair< const System *, unsigned int >, std::unique_ptr< ErrorVector > >	ErrorMap
	When calculating many error vectors at once, we need a data structure to hold them all.

Public Member Functions
	PatchRecoveryErrorEstimator ()
	Constructor.
	PatchRecoveryErrorEstimator (const PatchRecoveryErrorEstimator &)=default
	Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this class.
	PatchRecoveryErrorEstimator (PatchRecoveryErrorEstimator &&)=default
PatchRecoveryErrorEstimator &	operator= (const PatchRecoveryErrorEstimator &)=default
PatchRecoveryErrorEstimator &	operator= (PatchRecoveryErrorEstimator &&)=default
virtual	~PatchRecoveryErrorEstimator ()=default
virtual void	estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false) override
	This function uses the Patch Recovery error estimate to estimate the error on each cell.
void	set_patch_reuse (bool)
void	extra_quadrature_order (const int extraorder)
	Increases or decreases the order of the quadrature rule used for numerical integration.
virtual ErrorEstimatorType	type () const override
virtual void	estimate_errors (const EquationSystems &equation_systems, ErrorVector &error_per_cell, const std::map< const System *, SystemNorm > &error_norms, const std::map< const System *, const NumericVector< Number > * > *solution_vectors=nullptr, bool estimate_parent_error=false)
	This virtual function can be redefined in derived classes, but by default computes the sum of the error_per_cell for each system in the equation_systems.
virtual void	estimate_errors (const EquationSystems &equation_systems, ErrorMap &errors_per_cell, const std::map< const System *, const NumericVector< Number > * > *solution_vectors=nullptr, bool estimate_parent_error=false)
	This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors.

Public Attributes
unsigned int	target_patch_size
	The PatchErrorEstimator will build patches of at least this many elements to perform estimates.
Patch::PMF	patch_growth_strategy
	The PatchErrorEstimator will use this pointer to a Patch member function when growing patches.
SystemNorm	error_norm
	When estimating the error in a single system, the error_norm is used to control the scaling and norm choice for each variable.

Protected Member Functions
void	reduce_error (std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm) const
	This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector.

Static Protected Member Functions
static std::vector< Real >	specpoly (const unsigned int dim, const Order order, const Point p, const unsigned int matsize)

Protected Attributes
bool	patch_reuse
int	_extra_order
	Extra order to use for quadrature rule.

Friends
class	EstimateError

Detailed Description

This class implements the Patch Recovery error indicator.

Author

Varis Carey

Benjamin S. Kirk

Date

2004

Definition at line 47 of file patch_recovery_error_estimator.h.

Member Typedef Documentation

ErrorMap

typedef std::map<std::pair<const System *, unsigned int>, std::unique_ptr<ErrorVector> > libMesh::ErrorEstimator::ErrorMap

inherited

When calculating many error vectors at once, we need a data structure to hold them all.

Definition at line 121 of file error_estimator.h.

Constructor & Destructor Documentation

PatchRecoveryErrorEstimator() [1/3]

libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator

(

)

Constructor.

Defaults to H1 seminorm. All Hilbert norms and seminorms should be supported now. W1,p and W2,p norms would be natural to support if any contributors make the effort.

Definition at line 65 of file patch_recovery_error_estimator.C.

                                                         :
    ErrorEstimator(),
    target_patch_size(20),
    patch_growth_strategy(&Patch::add_local_face_neighbors),
    patch_reuse(true),
    _extra_order(1)
{
  error_norm = H1_SEMINORM;
}

References libMesh::ErrorEstimator::error_norm, and libMesh::H1_SEMINORM.

PatchRecoveryErrorEstimator() [2/3]

libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator ( const PatchRecoveryErrorEstimator &  )

default

Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this class.

PatchRecoveryErrorEstimator() [3/3]

libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator ( PatchRecoveryErrorEstimator &&  )

default

~PatchRecoveryErrorEstimator()

virtual libMesh::PatchRecoveryErrorEstimator::~PatchRecoveryErrorEstimator ( )

virtualdefault

Member Function Documentation

estimate_error()

void libMesh::PatchRecoveryErrorEstimator::estimate_error ( const System &  system, ErrorVector &  error_per_cell, const NumericVector< Number > *  solution_vector = nullptr, bool  estimate_parent_error = false  )

overridevirtual

This function uses the Patch Recovery error estimate to estimate the error on each cell.

The estimated error is output in the vector error_per_cell

Implements libMesh::ErrorEstimator.

Reimplemented in libMesh::WeightedPatchRecoveryErrorEstimator.

Definition at line 154 of file patch_recovery_error_estimator.C.

{
  LOG_SCOPE("estimate_error()", "PatchRecoveryErrorEstimator");
 
  // The current mesh
  const MeshBase & mesh = system.get_mesh();
 
  // Resize the error_per_cell vector to be
  // the number of elements, initialize it to 0.
  error_per_cell.resize (mesh.max_elem_id());
  std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
 
  // Prepare current_local_solution to localize a non-standard
  // solution vector if necessary
  if (solution_vector && solution_vector != system.solution.get())
    {
      NumericVector<Number> * newsol =
        const_cast<NumericVector<Number> *>(solution_vector);
      System & sys = const_cast<System &>(system);
      newsol->swap(*sys.solution);
      sys.update();
    }
 
  //------------------------------------------------------------
  // Iterate over all the active elements in the mesh
  // that live on this processor.
  Threads::parallel_for (ConstElemRange(mesh.active_local_elements_begin(),
                                        mesh.active_local_elements_end(),
                                        200),
                         EstimateError(system,
                                       *this,
                                       error_per_cell)
                         );
 
  // Each processor has now computed the error contributions
  // for its local elements, and error_per_cell contains 0 for all the
  // non-local elements.  Summing the vector will provide the true
  // value for each element, local or remote
  this->reduce_error(error_per_cell, system.comm());
 
  // If we used a non-standard solution before, now is the time to fix
  // the current_local_solution
  if (solution_vector && solution_vector != system.solution.get())
    {
      NumericVector<Number> * newsol =
        const_cast<NumericVector<Number> *>(solution_vector);
      System & sys = const_cast<System &>(system);
      newsol->swap(*sys.solution);
      sys.update();
    }
}

References libMesh::ParallelObject::comm(), libMesh::System::get_mesh(), mesh, libMesh::Threads::parallel_for(), libMesh::ErrorEstimator::reduce_error(), libMesh::System::solution, and libMesh::NumericVector< T >::swap().

Referenced by main().

estimate_errors() [1/2]

void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems &  equation_systems, ErrorMap &  errors_per_cell, const std::map< const System *, const NumericVector< Number > * > *  solution_vectors = nullptr, bool  estimate_parent_error = false  )

virtualinherited

This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors.

FIXME: This is a default implementation - derived classes should reimplement it for efficiency.

Currently this function ignores the error_norm.weight() values because it calculates each variable's error individually, unscaled.

The user selects which errors get computed by filling a map with error vectors: If errors_per_cell[&system][v] exists, it will be filled with the error values in variable v of system

Reimplemented in libMesh::UniformRefinementEstimator.

Definition at line 94 of file error_estimator.C.

{
  SystemNorm old_error_norm = this->error_norm;
 
  // Find the requested error values from each system
  for (auto s : make_range(equation_systems.n_systems()))
    {
      const System & sys = equation_systems.get_system(s);
 
      unsigned int n_vars = sys.n_vars();
 
      for (unsigned int v = 0; v != n_vars; ++v)
        {
          // Only fill in ErrorVectors the user asks for
          if (!errors_per_cell.count(std::make_pair(&sys, v)))
            continue;
 
          // Calculate error in only one variable
          std::vector<Real> weights(n_vars, 0.0);
          weights[v] = 1.0;
          this->error_norm =
            SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
                       weights);
 
          const NumericVector<Number> * solution_vector = nullptr;
          if (solution_vectors &&
              solution_vectors->find(&sys) != solution_vectors->end())
            solution_vector = solution_vectors->find(&sys)->second;
 
          this->estimate_error
            (sys, *errors_per_cell[std::make_pair(&sys, v)],
             solution_vector, estimate_parent_error);
        }
    }
 
  // Restore our old state before returning
  this->error_norm = old_error_norm;
}

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libMesh::make_range(), libMesh::EquationSystems::n_systems(), n_vars, libMesh::System::n_vars(), and libMesh::SystemNorm::type().

estimate_errors() [2/2]

void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems &  equation_systems, ErrorVector &  error_per_cell, const std::map< const System *, SystemNorm > &  error_norms, const std::map< const System *, const NumericVector< Number > * > *  solution_vectors = nullptr, bool  estimate_parent_error = false  )

virtualinherited

This virtual function can be redefined in derived classes, but by default computes the sum of the error_per_cell for each system in the equation_systems.

Currently this function ignores the error_norm member variable, and uses the function argument error_norms instead.

This function is named estimate_errors instead of estimate_error because otherwise C++ can get confused.

Reimplemented in libMesh::UniformRefinementEstimator.

Definition at line 47 of file error_estimator.C.

{
  SystemNorm old_error_norm = this->error_norm;
 
  // Sum the error values from each system
  for (auto s : make_range(equation_systems.n_systems()))
    {
      ErrorVector system_error_per_cell;
      const System & sys = equation_systems.get_system(s);
      if (const auto it = error_norms.find(&sys);
          it == error_norms.end())
        this->error_norm = old_error_norm;
      else
        this->error_norm = it->second;
 
      const NumericVector<Number> * solution_vector = nullptr;
      if (solution_vectors &&
          solution_vectors->find(&sys) != solution_vectors->end())
        solution_vector = solution_vectors->find(&sys)->second;
 
      this->estimate_error(sys, system_error_per_cell,
                           solution_vector, estimate_parent_error);
 
      if (s)
        {
          libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
          for (auto i : index_range(error_per_cell))
            error_per_cell[i] += system_error_per_cell[i];
        }
      else
        error_per_cell = system_error_per_cell;
    }
 
  // Restore our old state before returning
  this->error_norm = old_error_norm;
}

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libMesh::index_range(), libMesh::make_range(), and libMesh::EquationSystems::n_systems().

extra_quadrature_order()

void libMesh::PatchRecoveryErrorEstimator::extra_quadrature_order ( const int  extraorder )

inline

Increases or decreases the order of the quadrature rule used for numerical integration.

The default extraorder is 1, because properly integrating L2 error requires integrating the squares of terms with order p+1, and 2p+2 is 1 higher than what we default to using for reasonable mass matrix integration.

Definition at line 102 of file patch_recovery_error_estimator.h.

  { _extra_order = extraorder; }

References _extra_order.

operator=() [1/2]

PatchRecoveryErrorEstimator & libMesh::PatchRecoveryErrorEstimator::operator= ( const PatchRecoveryErrorEstimator &  )

default

operator=() [2/2]

PatchRecoveryErrorEstimator & libMesh::PatchRecoveryErrorEstimator::operator= ( PatchRecoveryErrorEstimator &&  )

default

reduce_error()

void libMesh::ErrorEstimator::reduce_error ( std::vector< ErrorVectorReal > &  error_per_cell, const Parallel::Communicator &  comm  ) const

protectedinherited

This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector.

Definition at line 32 of file error_estimator.C.

{
  // This function must be run on all processors at once
  // parallel_object_only();
 
  // Each processor has now computed the error contributions
  // for its local elements.  We may need to sum the vector to
  // recover the error for each element.
 
  comm.sum(error_per_cell);
}

References libMesh::Parallel::Communicator::sum().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), estimate_error(), and libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error().

set_patch_reuse()

void libMesh::PatchRecoveryErrorEstimator::set_patch_reuse

(

bool

patch_reuse_flag

)

Definition at line 58 of file patch_recovery_error_estimator.C.

{
  patch_reuse = patch_reuse_flag;
}

References patch_reuse.

specpoly()

std::vector< Real > libMesh::PatchRecoveryErrorEstimator::specpoly ( const unsigned int  dim, const Order  order, const Point  p, const unsigned int  matsize  )

staticprotected

Returns

The spectral polynomial basis function values at a point (x,y,z).

Definition at line 77 of file patch_recovery_error_estimator.C.

{
  std::vector<Real> psi;
  psi.reserve(matsize);
  int npows = order+1;
  std::vector<Real> xpow(npows,1.), ypow, zpow;
  {
    Real x = p(0);
    for (auto i : make_range(1, npows))
      xpow[i] = xpow[i-1] * x;
  }
  if (dim > 1)
    {
      Real y = p(1);
      ypow.resize(npows,1.);
      for (auto i : make_range(1, npows))
        ypow[i] = ypow[i-1] * y;
    }
  if (dim > 2)
    {
      Real z = p(2);
      zpow.resize(npows,1.);
      for (auto i : make_range(1, npows))
        zpow[i] = zpow[i-1] * z;
    }
 
  // builds psi vector of form 1 x y z x^2 xy xz y^2 yz z^2 etc..
  // I haven't added 1D support here
  for (unsigned int poly_deg=0; poly_deg <= static_cast<unsigned int>(order) ; poly_deg++)
    { // loop over all polynomials of total degree = poly_deg
 
      switch (dim)
        {
          // 3D spectral polynomial basis functions
        case 3:
          {
            for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
              for (int yexp=poly_deg-xexp; yexp >= 0; yexp--)
                {
                  int zexp = poly_deg - xexp - yexp;
                  psi.push_back(xpow[xexp]*ypow[yexp]*zpow[zexp]);
                }
            break;
          }
 
          // 2D spectral polynomial basis functions
        case 2:
          {
            for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
              {
                int yexp = poly_deg - xexp;
                psi.push_back(xpow[xexp]*ypow[yexp]);
              }
            break;
          }
 
          // 1D spectral polynomial basis functions
        case 1:
          {
            int xexp = poly_deg;
            psi.push_back(xpow[xexp]);
            break;
          }
 
        default:
          libmesh_error_msg("Invalid dimension dim " << dim);
        }
    }
 
  return psi;
}

References dim, libMesh::make_range(), and libMesh::Real.

Referenced by libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), and libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()().

type()

ErrorEstimatorType libMesh::PatchRecoveryErrorEstimator::type ( ) const

overridevirtual

Returns

The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Reimplemented in libMesh::WeightedPatchRecoveryErrorEstimator.

Definition at line 50 of file patch_recovery_error_estimator.C.

{
  return PATCH_RECOVERY;
}

References libMesh::PATCH_RECOVERY.

Friends And Related Symbol Documentation

EstimateError

friend class EstimateError

friend

Definition at line 149 of file patch_recovery_error_estimator.h.

Member Data Documentation

_extra_order

int libMesh::PatchRecoveryErrorEstimator::_extra_order

protected

Extra order to use for quadrature rule.

Definition at line 122 of file patch_recovery_error_estimator.h.

Referenced by extra_quadrature_order(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), and libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()().

error_norm

SystemNorm libMesh::ErrorEstimator::error_norm

inherited

When estimating the error in a single system, the error_norm is used to control the scaling and norm choice for each variable.

Not all estimators will support all norm choices. The default scaling is for all variables to be weighted equally. The default norm choice depends on the error estimator.

Part of this functionality was supported via component_scale and sobolev_order in older libMesh versions, and a small part was supported via component_mask in even older versions. Hopefully the encapsulation here will allow us to avoid changing this API again.

Definition at line 158 of file error_estimator.h.

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::AdjointRefinementEstimator::AdjointRefinementEstimator(), libMesh::DiscontinuityMeasure::boundary_side_integration(), libMesh::KellyErrorEstimator::boundary_side_integration(), libMesh::DiscontinuityMeasure::DiscontinuityMeasure(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ExactErrorEstimator::ExactErrorEstimator(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::DiscontinuityMeasure::init_context(), libMesh::LaplacianErrorEstimator::init_context(), libMesh::KellyErrorEstimator::init_context(), libMesh::DiscontinuityMeasure::internal_side_integration(), libMesh::LaplacianErrorEstimator::internal_side_integration(), libMesh::KellyErrorEstimator::internal_side_integration(), libMesh::KellyErrorEstimator::KellyErrorEstimator(), libMesh::LaplacianErrorEstimator::LaplacianErrorEstimator(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::UniformRefinementEstimator::UniformRefinementEstimator().

patch_growth_strategy

Patch::PMF libMesh::PatchRecoveryErrorEstimator::patch_growth_strategy

The PatchErrorEstimator will use this pointer to a Patch member function when growing patches.

The default strategy used is Patch::add_local_face_neighbors. Patch::add_local_point_neighbors may be more reliable but slower.

Definition at line 91 of file patch_recovery_error_estimator.h.

Referenced by libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), and libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()().

patch_reuse

bool libMesh::PatchRecoveryErrorEstimator::patch_reuse

protected

Definition at line 117 of file patch_recovery_error_estimator.h.

Referenced by libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), and set_patch_reuse().

target_patch_size

unsigned int libMesh::PatchRecoveryErrorEstimator::target_patch_size

The PatchErrorEstimator will build patches of at least this many elements to perform estimates.

Definition at line 83 of file patch_recovery_error_estimator.h.

Referenced by libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), and libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()().

---

The documentation for this class was generated from the following files:

• include/error_estimation/patch_recovery_error_estimator.h
• src/error_estimation/patch_recovery_error_estimator.C