libMesh
Classes | Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Friends | List of all members
libMesh::PatchRecoveryErrorEstimator Class Reference

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

#include <patch_recovery_error_estimator.h>

Inheritance diagram for libMesh::PatchRecoveryErrorEstimator:
[legend]

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 >, ErrorVector * > ErrorMap
 When calculating many error vectors at once, we need a data structure to hold them all. More...
 

Public Member Functions

 PatchRecoveryErrorEstimator ()
 Constructor. More...
 
 ~PatchRecoveryErrorEstimator ()
 Destructor. More...
 
virtual void estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=libmesh_nullptr, bool estimate_parent_error=false) libmesh_override
 This function uses the Patch Recovery error estimate to estimate the error on each cell. More...
 
void set_patch_reuse (bool)
 
virtual ErrorEstimatorType type () const libmesh_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=libmesh_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. More...
 
virtual void estimate_errors (const EquationSystems &equation_systems, ErrorMap &errors_per_cell, const std::map< const System *, const NumericVector< Number > * > *solution_vectors=libmesh_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. More...
 

Public Attributes

unsigned int target_patch_size
 The PatchErrorEstimator will build patches of at least this many elements to perform estimates. More...
 
Patch::PMF patch_growth_strategy
 The PatchErrorEstimator will use this pointer to a Patch member function when growing patches. More...
 
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. More...
 

Protected Member Functions

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

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
Date
2004

Definition at line 48 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 113 of file error_estimator.h.

Constructor & Destructor Documentation

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 57 of file patch_recovery_error_estimator.h.

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

57  :
61  patch_reuse(true)
62  { error_norm = H1_SEMINORM; }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
Patch::PMF patch_growth_strategy
The PatchErrorEstimator will use this pointer to a Patch member function when growing patches...
ErrorEstimator()
Constructor.
unsigned int target_patch_size
The PatchErrorEstimator will build patches of at least this many elements to perform estimates...
void add_local_face_neighbors()
This function finds all elements on the current processor which touch the current patch at a face...
Definition: patch.C:85
libMesh::PatchRecoveryErrorEstimator::~PatchRecoveryErrorEstimator ( )

Destructor.

Definition at line 67 of file patch_recovery_error_estimator.h.

References estimate_error(), and libmesh_nullptr.

67 {}

Member Function Documentation

void libMesh::PatchRecoveryErrorEstimator::estimate_error ( const System system,
ErrorVector error_per_cell,
const NumericVector< Number > *  solution_vector = libmesh_nullptr,
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 131 of file patch_recovery_error_estimator.C.

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

Referenced by main(), and ~PatchRecoveryErrorEstimator().

135 {
136  LOG_SCOPE("estimate_error()", "PatchRecoveryErrorEstimator");
137 
138  // The current mesh
139  const MeshBase & mesh = system.get_mesh();
140 
141  // Resize the error_per_cell vector to be
142  // the number of elements, initialize it to 0.
143  error_per_cell.resize (mesh.max_elem_id());
144  std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
145 
146  // Prepare current_local_solution to localize a non-standard
147  // solution vector if necessary
148  if (solution_vector && solution_vector != system.solution.get())
149  {
150  NumericVector<Number> * newsol =
151  const_cast<NumericVector<Number> *>(solution_vector);
152  System & sys = const_cast<System &>(system);
153  newsol->swap(*sys.solution);
154  sys.update();
155  }
156 
157  //------------------------------------------------------------
158  // Iterate over all the active elements in the mesh
159  // that live on this processor.
160  Threads::parallel_for (ConstElemRange(mesh.active_local_elements_begin(),
161  mesh.active_local_elements_end(),
162  200),
163  EstimateError(system,
164  *this,
165  error_per_cell)
166  );
167 
168  // Each processor has now computed the error contributions
169  // for its local elements, and error_per_cell contains 0 for all the
170  // non-local elements. Summing the vector will provide the true
171  // value for each element, local or remote
172  this->reduce_error(error_per_cell, system.comm());
173 
174  // If we used a non-standard solution before, now is the time to fix
175  // the current_local_solution
176  if (solution_vector && solution_vector != system.solution.get())
177  {
178  NumericVector<Number> * newsol =
179  const_cast<NumericVector<Number> *>(solution_vector);
180  System & sys = const_cast<System &>(system);
181  newsol->swap(*sys.solution);
182  sys.update();
183  }
184 }
void parallel_for(const Range &range, const Body &body)
Execute the provided function object in parallel on the specified range.
Definition: threads_none.h:73
ImplicitSystem & sys
MeshBase & mesh
StoredRange< MeshBase::const_element_iterator, const Elem * > ConstElemRange
Definition: elem_range.h:34
void reduce_error(std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm LIBMESH_CAN_DEFAULT_TO_COMMWORLD) const
This method takes the local error contributions in error_per_cell from each processor and combines th...
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition: system.C:425
virtual void get(const std::vector< numeric_index_type > &index, T *values) const
Access multiple components at once.
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 = libmesh_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 48 of file error_estimator.C.

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

Referenced by libMesh::ErrorEstimator::~ErrorEstimator().

53 {
54  SystemNorm old_error_norm = this->error_norm;
55 
56  // Sum the error values from each system
57  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
58  {
59  ErrorVector system_error_per_cell;
60  const System & sys = equation_systems.get_system(s);
61  if (error_norms.find(&sys) == error_norms.end())
62  this->error_norm = old_error_norm;
63  else
64  this->error_norm = error_norms.find(&sys)->second;
65 
66  const NumericVector<Number> * solution_vector = libmesh_nullptr;
67  if (solution_vectors &&
68  solution_vectors->find(&sys) != solution_vectors->end())
69  solution_vector = solution_vectors->find(&sys)->second;
70 
71  this->estimate_error(sys, system_error_per_cell,
72  solution_vector, estimate_parent_error);
73 
74  if (s)
75  {
76  libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
77  for (std::size_t i=0; i != error_per_cell.size(); ++i)
78  error_per_cell[i] += system_error_per_cell[i];
79  }
80  else
81  error_per_cell = system_error_per_cell;
82  }
83 
84  // Restore our old state before returning
85  this->error_norm = old_error_norm;
86 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=libmesh_nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems equation_systems,
ErrorMap errors_per_cell,
const std::map< const System *, const NumericVector< Number > * > *  solution_vectors = libmesh_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.

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

98 {
99  SystemNorm old_error_norm = this->error_norm;
100 
101  // Find the requested error values from each system
102  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
103  {
104  const System & sys = equation_systems.get_system(s);
105 
106  unsigned int n_vars = sys.n_vars();
107 
108  for (unsigned int v = 0; v != n_vars; ++v)
109  {
110  // Only fill in ErrorVectors the user asks for
111  if (errors_per_cell.find(std::make_pair(&sys, v)) ==
112  errors_per_cell.end())
113  continue;
114 
115  // Calculate error in only one variable
116  std::vector<Real> weights(n_vars, 0.0);
117  weights[v] = 1.0;
118  this->error_norm =
119  SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
120  weights);
121 
122  const NumericVector<Number> * solution_vector = libmesh_nullptr;
123  if (solution_vectors &&
124  solution_vectors->find(&sys) != solution_vectors->end())
125  solution_vector = solution_vectors->find(&sys)->second;
126 
127  this->estimate_error
128  (sys, *errors_per_cell[std::make_pair(&sys, v)],
129  solution_vector, estimate_parent_error);
130  }
131  }
132 
133  // Restore our old state before returning
134  this->error_norm = old_error_norm;
135 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
const unsigned int n_vars
Definition: tecplot_io.C:68
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=libmesh_nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
void libMesh::ErrorEstimator::reduce_error ( std::vector< ErrorVectorReal > &  error_per_cell,
const Parallel::Communicator &comm  LIBMESH_CAN_DEFAULT_TO_COMMWORLD 
) 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 33 of file error_estimator.C.

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

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

35 {
36  // This function must be run on all processors at once
37  // parallel_object_only();
38 
39  // Each processor has now computed the error contributions
40  // for its local elements. We may need to sum the vector to
41  // recover the error for each element.
42 
43  comm.sum(error_per_cell);
44 }
void libMesh::PatchRecoveryErrorEstimator::set_patch_reuse ( bool  patch_reuse_flag)

Definition at line 47 of file patch_recovery_error_estimator.C.

References patch_reuse.

Referenced by build_error_estimator().

48 {
49  patch_reuse = patch_reuse_flag;
50 }
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 54 of file patch_recovery_error_estimator.C.

References libMesh::Real, and libMesh::x.

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

58 {
59  std::vector<Real> psi;
60  psi.reserve(matsize);
61  unsigned int npows = order+1;
62  std::vector<Real> xpow(npows,1.), ypow, zpow;
63  {
64  Real x = p(0);
65  for (unsigned int i=1; i != npows; ++i)
66  xpow[i] = xpow[i-1] * x;
67  }
68  if (dim > 1)
69  {
70  Real y = p(1);
71  ypow.resize(npows,1.);
72  for (unsigned int i=1; i != npows; ++i)
73  ypow[i] = ypow[i-1] * y;
74  }
75  if (dim > 2)
76  {
77  Real z = p(2);
78  zpow.resize(npows,1.);
79  for (unsigned int i=1; i != npows; ++i)
80  zpow[i] = zpow[i-1] * z;
81  }
82 
83  // builds psi vector of form 1 x y z x^2 xy xz y^2 yz z^2 etc..
84  // I haven't added 1D support here
85  for (unsigned int poly_deg=0; poly_deg <= static_cast<unsigned int>(order) ; poly_deg++)
86  { // loop over all polynomials of total degree = poly_deg
87 
88  switch (dim)
89  {
90  // 3D spectral polynomial basis functions
91  case 3:
92  {
93  for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
94  for (int yexp=poly_deg-xexp; yexp >= 0; yexp--)
95  {
96  int zexp = poly_deg - xexp - yexp;
97  psi.push_back(xpow[xexp]*ypow[yexp]*zpow[zexp]);
98  }
99  break;
100  }
101 
102  // 2D spectral polynomial basis functions
103  case 2:
104  {
105  for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
106  {
107  int yexp = poly_deg - xexp;
108  psi.push_back(xpow[xexp]*ypow[yexp]);
109  }
110  break;
111  }
112 
113  // 1D spectral polynomial basis functions
114  case 1:
115  {
116  int xexp = poly_deg;
117  psi.push_back(xpow[xexp]);
118  break;
119  }
120 
121  default:
122  libmesh_error_msg("Invalid dimension dim " << dim);
123  }
124  }
125 
126  return psi;
127 }
unsigned int dim
PetscErrorCode Vec x
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual ErrorEstimatorType libMesh::PatchRecoveryErrorEstimator::type ( ) const
virtual
Returns
The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Reimplemented in libMesh::WeightedPatchRecoveryErrorEstimator.

Definition at line 97 of file patch_recovery_error_estimator.h.

References dim, libMesh::PATCH_RECOVERY, and specpoly().

Friends And Related Function Documentation

friend class EstimateError
friend

Definition at line 137 of file patch_recovery_error_estimator.h.

Referenced by estimate_error().

Member Data Documentation

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 150 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::JumpErrorEstimator::estimate_error(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ExactErrorEstimator::ExactErrorEstimator(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::LaplacianErrorEstimator::init_context(), libMesh::DiscontinuityMeasure::init_context(), libMesh::KellyErrorEstimator::init_context(), libMesh::LaplacianErrorEstimator::internal_side_integration(), libMesh::DiscontinuityMeasure::internal_side_integration(), libMesh::KellyErrorEstimator::internal_side_integration(), libMesh::KellyErrorEstimator::KellyErrorEstimator(), libMesh::LaplacianErrorEstimator::LaplacianErrorEstimator(), main(), PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::UniformRefinementEstimator::UniformRefinementEstimator().

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 93 of file patch_recovery_error_estimator.h.

bool libMesh::PatchRecoveryErrorEstimator::patch_reuse
protected

Definition at line 110 of file patch_recovery_error_estimator.h.

Referenced by set_patch_reuse().

unsigned int libMesh::PatchRecoveryErrorEstimator::target_patch_size

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

Definition at line 85 of file patch_recovery_error_estimator.h.


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