Poster of Linux kernelThe best gift for a Linux geek
PatchRecoveryErrorEstimator

PatchRecoveryErrorEstimator

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

PatchRecoveryErrorEstimator -  

SYNOPSIS


#include <patch_recovery_error_estimator.h>

Inherits ErrorEstimator.  

Classes


class EstimateError
 

Public Types


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

Public Member Functions


PatchRecoveryErrorEstimator ()

~PatchRecoveryErrorEstimator ()

virtual void estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=NULL, bool estimate_parent_error=false)

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=NULL, bool estimate_parent_error=false)

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

Public Attributes


unsigned int target_patch_size

Patch::PMF patch_growth_strategy

SystemNorm error_norm
 

Protected Member Functions


void reduce_error (std::vector< float > &error_per_cell) const
 

Static Private Member Functions


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

Friends


class EsitmateError
 

Detailed Description

This class implements the Patch Recovery error indicator.

Author:

Varis Carey, Benjamin S. Kirk, 2004.

Definition at line 45 of file patch_recovery_error_estimator.h.  

Member Typedef Documentation

 

typedef std::map<std::pair<const System*, unsigned int>, ErrorVector*> ErrorEstimator::ErrorMap [inherited]When calculating many error vectors at once, we need a data structure to hold them all

Definition at line 105 of file error_estimator.h.  

Constructor & Destructor Documentation

 

PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator () [inline]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 54 of file patch_recovery_error_estimator.h.

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

                                :
    target_patch_size(20),
    patch_growth_strategy(&Patch::add_local_face_neighbors) 
  { error_norm = H1_SEMINORM; }
 

PatchRecoveryErrorEstimator::~PatchRecoveryErrorEstimator () [inline]Destructor.

Definition at line 62 of file patch_recovery_error_estimator.h.

{}
 

Member Function Documentation

 

void PatchRecoveryErrorEstimator::estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector = NULL, boolestimate_parent_error = false) [virtual]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 ErrorEstimator.

Definition at line 111 of file patch_recovery_error_estimator.C.

References MeshBase::active_local_elements_begin(), MeshBase::active_local_elements_end(), AutoPtr< Tp >::get(), System::get_mesh(), MeshBase::max_elem_id(), Threads::parallel_for(), ErrorEstimator::reduce_error(), System::solution, and NumericVector< T >::swap().

{
  START_LOG('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 contribuions
  // for its local elements, and error_per_cell contains 0 for all the
  // non-local elements.  Summing the vector will provide the L_oo value
  // for each element, local or remote
  this->reduce_error(error_per_cell);
  
  // 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();
    }

  STOP_LOG('estimate_error()', 'PatchRecoveryErrorEstimator');
}
 

void ErrorEstimator::estimate_errors (const EquationSystems &equation_systems, ErrorMap &errors_per_cell, const std::map< const System *, const NumericVector< Number > * > *solution_vectors = NULL, boolestimate_parent_error = false) [virtual, inherited]This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors.

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

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

Reimplemented in UniformRefinementEstimator.

Definition at line 92 of file error_estimator.C.

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

{
  SystemNorm old_error_norm = this->error_norm;

  // Find the requested error values from each system
  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
    {
      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.find(std::make_pair(&sys, v)) ==
              errors_per_cell.end())
            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(0)),
                       weights);

          const NumericVector<Number>* solution_vector = NULL;
          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;
}
 

void 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 = NULL, boolestimate_parent_error = false) [virtual, inherited]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 UniformRefinementEstimator.

Definition at line 46 of file error_estimator.C.

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

{
  SystemNorm old_error_norm = this->error_norm;

  // Sum the error values from each system
  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
    {
      ErrorVector system_error_per_cell;
      const System &sys = equation_systems.get_system(s);
      if (error_norms.find(&sys) == error_norms.end())
        this->error_norm = old_error_norm;
      else
        this->error_norm = error_norms.find(&sys)->second;

      const NumericVector<Number>* solution_vector = NULL;
      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(error_per_cell.size() == system_error_per_cell.size());
          for (unsigned int i=0; i != error_per_cell.size(); ++i)
            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;
}
 

void ErrorEstimator::reduce_error (std::vector< float > &error_per_cell) const [protected, inherited]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.

Referenced by UniformRefinementEstimator::_estimate_error(), estimate_error(), and JumpErrorEstimator::estimate_error().

{
  // This function must be run on all processors at once
  parallel_only();

  // Each processor has now computed the error contribuions
  // for its local elements.  We may need to sum the vector to
  // recover the error for each element.
  
  Parallel::sum(error_per_cell);
}
 

std::vector< Real > PatchRecoveryErrorEstimator::specpoly (const unsigned intdim, const Orderorder, const Pointp, const unsigned intmatsize) [static, private]Returns the spectral polynomial basis function values at a point x,y,z

Definition at line 47 of file patch_recovery_error_estimator.C.

References Utility::pow().

Referenced by PatchRecoveryErrorEstimator::EstimateError::operator()().

{
  std::vector<Real> psi;
  psi.reserve(matsize);
  Real x = p(0), y=0., z=0.;
  if (dim > 1)
    y = p(1);
  if (dim > 2)
    z = p(2);
    
  // 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 degreee = 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(std::pow(x,static_cast<Real>(xexp))*
                                std::pow(y,static_cast<Real>(yexp))*
                                std::pow(z,static_cast<Real>(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(std::pow(x,static_cast<Real>(xexp))*
                              std::pow(y,static_cast<Real>(yexp)));
              }
            break;
          }

          // 1D spectral polynomial basis functions
        case 1:
          {
            int xexp = poly_deg;
            psi.push_back(std::pow(x,static_cast<Real>(xexp)));
            break;
          }
          
        default:
          libmesh_error();
        }
    }

  return psi;
}
 

Friends And Related Function Documentation

 

friend class EsitmateError [friend]

Definition at line 125 of file patch_recovery_error_estimator.h.  

Member Data Documentation

 

SystemNorm 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 137 of file error_estimator.h.

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

Patch::PMF PatchRecoveryErrorEstimator::patch_growth_strategyThe 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 88 of file patch_recovery_error_estimator.h.

Referenced by PatchRecoveryErrorEstimator::EstimateError::operator()().  

unsigned int PatchRecoveryErrorEstimator::target_patch_sizeThe PatchErrorEstimator will build patches of at least this many elements to perform estimates

Definition at line 80 of file patch_recovery_error_estimator.h.

Referenced by PatchRecoveryErrorEstimator::EstimateError::operator()().

 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Classes
Public Types
Public Member Functions
Public Attributes
Protected Member Functions
Static Private Member Functions
Friends
Detailed Description
Member Typedef Documentation
typedef std::map<std::pair<const System*, unsigned int>, ErrorVector*> ErrorEstimator::ErrorMap [inherited]When calculating many error vectors at once, we need a data structure to hold them all
Constructor & Destructor Documentation
PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator () [inline]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.
PatchRecoveryErrorEstimator::~PatchRecoveryErrorEstimator () [inline]Destructor.
Member Function Documentation
void PatchRecoveryErrorEstimator::estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector = NULL, boolestimate_parent_error = false) [virtual]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
void ErrorEstimator::estimate_errors (const EquationSystems &equation_systems, ErrorMap &errors_per_cell, const std::map< const System *, const NumericVector< Number > * > *solution_vectors = NULL, boolestimate_parent_error = false) [virtual, inherited]This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors.
void 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 = NULL, boolestimate_parent_error = false) [virtual, inherited]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.
void ErrorEstimator::reduce_error (std::vector< float > &error_per_cell) const [protected, inherited]This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector.
std::vector< Real > PatchRecoveryErrorEstimator::specpoly (const unsigned intdim, const Orderorder, const Pointp, const unsigned intmatsize) [static, private]Returns the spectral polynomial basis function values at a point x,y,z
Friends And Related Function Documentation
friend class EsitmateError [friend]
Member Data Documentation
SystemNorm 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.
Patch::PMF PatchRecoveryErrorEstimator::patch_growth_strategyThe 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.
unsigned int PatchRecoveryErrorEstimator::target_patch_sizeThe PatchErrorEstimator will build patches of at least this many elements to perform estimates
Author

This document was created by man2html, using the manual pages.
Time: 21:52:04 GMT, April 16, 2011