libMesh::PatchRecoveryErrorEstimator Class Reference

#include <patch_recovery_error_estimator.h>

Inheritance diagram for libMesh::PatchRecoveryErrorEstimator:

List of all members.

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)
void set_patch_reuse (bool)
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 Protected Member Functions

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

Protected Attributes

bool patch_reuse

Friends

class EstimateError

Detailed Description

This class implements the Patch Recovery error indicator.

Author:
Varis Carey, Benjamin S. Kirk, 2004.

Definition at line 47 of file patch_recovery_error_estimator.h.


Member Typedef Documentation

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

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

Definition at line 107 of file error_estimator.h.


Constructor & Destructor Documentation

libMesh::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 56 of file patch_recovery_error_estimator.h.

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

libMesh::PatchRecoveryErrorEstimator::~PatchRecoveryErrorEstimator (  )  [inline]

Destructor.

Definition at line 65 of file patch_recovery_error_estimator.h.

00065 {}


Member Function Documentation

void libMesh::PatchRecoveryErrorEstimator::estimate_error ( const System system,
ErrorVector error_per_cell,
const NumericVector< Number > *  solution_vector = NULL,
bool  estimate_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 libMesh::ErrorEstimator.

Reimplemented in libMesh::WeightedPatchRecoveryErrorEstimator.

Definition at line 129 of file patch_recovery_error_estimator.C.

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

00133 {
00134   START_LOG("estimate_error()", "PatchRecoveryErrorEstimator");
00135 
00136   // The current mesh
00137   const MeshBase& mesh = system.get_mesh();
00138 
00139   // Resize the error_per_cell vector to be
00140   // the number of elements, initialize it to 0.
00141   error_per_cell.resize (mesh.max_elem_id());
00142   std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
00143 
00144   // Prepare current_local_solution to localize a non-standard
00145   // solution vector if necessary
00146   if (solution_vector && solution_vector != system.solution.get())
00147     {
00148       NumericVector<Number>* newsol =
00149         const_cast<NumericVector<Number>*>(solution_vector);
00150       System &sys = const_cast<System&>(system);
00151       newsol->swap(*sys.solution);
00152       sys.update();
00153     }
00154 
00155   //------------------------------------------------------------
00156   // Iterate over all the active elements in the mesh
00157   // that live on this processor.
00158   Threads::parallel_for (ConstElemRange(mesh.active_local_elements_begin(),
00159                                         mesh.active_local_elements_end(),
00160                                         200),
00161                          EstimateError(system,
00162                                        *this,
00163                                        error_per_cell)
00164                          );
00165 
00166   // Each processor has now computed the error contributions
00167   // for its local elements, and error_per_cell contains 0 for all the
00168   // non-local elements.  Summing the vector will provide the true
00169   // value for each element, local or remote
00170   this->reduce_error(error_per_cell);
00171 
00172   // If we used a non-standard solution before, now is the time to fix
00173   // the current_local_solution
00174   if (solution_vector && solution_vector != system.solution.get())
00175     {
00176       NumericVector<Number>* newsol =
00177         const_cast<NumericVector<Number>*>(solution_vector);
00178       System &sys = const_cast<System&>(system);
00179       newsol->swap(*sys.solution);
00180       sys.update();
00181     }
00182 
00183   STOP_LOG("estimate_error()", "PatchRecoveryErrorEstimator");
00184 }

void libMesh::ErrorEstimator::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 
) [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 libMesh::UniformRefinementEstimator.

Definition at line 93 of file error_estimator.C.

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

00097 {
00098   SystemNorm old_error_norm = this->error_norm;
00099 
00100   // Find the requested error values from each system
00101   for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
00102     {
00103       const System &sys = equation_systems.get_system(s);
00104 
00105       unsigned int n_vars = sys.n_vars();
00106 
00107       for (unsigned int v = 0; v != n_vars; ++v)
00108         {
00109           // Only fill in ErrorVectors the user asks for
00110           if (errors_per_cell.find(std::make_pair(&sys, v)) ==
00111               errors_per_cell.end())
00112             continue;
00113 
00114           // Calculate error in only one variable
00115           std::vector<Real> weights(n_vars, 0.0);
00116           weights[v] = 1.0;
00117           this->error_norm =
00118             SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
00119                        weights);
00120 
00121           const NumericVector<Number>* solution_vector = NULL;
00122           if (solution_vectors &&
00123               solution_vectors->find(&sys) != solution_vectors->end())
00124             solution_vector = solution_vectors->find(&sys)->second;
00125 
00126           this->estimate_error
00127             (sys, *errors_per_cell[std::make_pair(&sys, v)],
00128              solution_vector, estimate_parent_error);
00129         }
00130     }
00131 
00132   // Restore our old state before returning
00133   this->error_norm = old_error_norm;
00134 }

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 = NULL,
bool  estimate_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 libMesh::UniformRefinementEstimator.

Definition at line 47 of file error_estimator.C.

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

00052 {
00053   SystemNorm old_error_norm = this->error_norm;
00054 
00055   // Sum the error values from each system
00056   for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
00057     {
00058       ErrorVector system_error_per_cell;
00059       const System &sys = equation_systems.get_system(s);
00060       if (error_norms.find(&sys) == error_norms.end())
00061         this->error_norm = old_error_norm;
00062       else
00063         this->error_norm = error_norms.find(&sys)->second;
00064 
00065       const NumericVector<Number>* solution_vector = NULL;
00066       if (solution_vectors &&
00067           solution_vectors->find(&sys) != solution_vectors->end())
00068         solution_vector = solution_vectors->find(&sys)->second;
00069 
00070       this->estimate_error(sys, system_error_per_cell,
00071                            solution_vector, estimate_parent_error);
00072 
00073       if (s)
00074         {
00075           libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
00076           for (unsigned int i=0; i != error_per_cell.size(); ++i)
00077             error_per_cell[i] += system_error_per_cell[i];
00078         }
00079       else
00080         error_per_cell = system_error_per_cell;
00081     }
00082 
00083   // Restore our old state before returning
00084   this->error_norm = old_error_norm;
00085 }

void libMesh::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.

Referenced by libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), estimate_error(), and libMesh::JumpErrorEstimator::estimate_error().

void libMesh::PatchRecoveryErrorEstimator::set_patch_reuse ( bool  patch_reuse_flag  ) 

Definition at line 47 of file patch_recovery_error_estimator.C.

References patch_reuse.

00048   {
00049     patch_reuse = patch_reuse_flag;
00050   }

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

Returns the spectral polynomial basis function values at a point x,y,z

Definition at line 54 of file patch_recovery_error_estimator.C.

References libMesh::Real.

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

00058 {
00059   std::vector<Real> psi;
00060   psi.reserve(matsize);
00061   Real x = p(0), y=0., z=0.;
00062   unsigned int npows = order+1;
00063   std::vector<Real> xpow(npows,1.), ypow, zpow;
00064   for (unsigned int i=1; i != npows; ++i)
00065     xpow[i] = xpow[i-1] * x;
00066   if (dim > 1)
00067     {
00068       y = p(1);
00069       ypow.resize(npows,1.);
00070       for (unsigned int i=1; i != npows; ++i)
00071         ypow[i] = ypow[i-1] * y;
00072     }
00073   if (dim > 2)
00074     {
00075       z = p(2);
00076       zpow.resize(npows,1.);
00077       for (unsigned int i=1; i != npows; ++i)
00078         zpow[i] = zpow[i-1] * z;
00079     }
00080 
00081   // builds psi vector of form 1 x y z x^2 xy xz y^2 yz z^2 etc..
00082   // I haven't added 1D support here
00083   for (unsigned int poly_deg=0; poly_deg <= static_cast<unsigned int>(order) ; poly_deg++)
00084     { // loop over all polynomials of total degreee = poly_deg
00085 
00086       switch (dim)
00087         {
00088           // 3D spectral polynomial basis functions
00089         case 3:
00090           {
00091             for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
00092               for (int yexp=poly_deg-xexp; yexp >= 0; yexp--)
00093                 {
00094                   int zexp = poly_deg - xexp - yexp;
00095                   psi.push_back(xpow[xexp]*ypow[yexp]*zpow[zexp]);
00096                 }
00097             break;
00098           }
00099 
00100           // 2D spectral polynomial basis functions
00101         case 2:
00102           {
00103             for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
00104               {
00105                 int yexp = poly_deg - xexp;
00106                 psi.push_back(xpow[xexp]*ypow[yexp]);
00107               }
00108             break;
00109           }
00110 
00111           // 1D spectral polynomial basis functions
00112         case 1:
00113           {
00114             int xexp = poly_deg;
00115             psi.push_back(xpow[xexp]);
00116             break;
00117           }
00118 
00119         default:
00120           libmesh_error();
00121         }
00122     }
00123 
00124   return psi;
00125 }


Friends And Related Function Documentation

friend class EstimateError [friend]

Reimplemented in libMesh::WeightedPatchRecoveryErrorEstimator.

Definition at line 139 of file patch_recovery_error_estimator.h.

Referenced by estimate_error().


Member Data Documentation

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

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

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::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), and libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()().

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::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), and libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()().


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

Site Created By: libMesh Developers
Last modified: February 05 2013 19:55:32 UTC

Hosted By:
SourceForge.net Logo