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 >, std::unique_ptr< 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 (const PatchRecoveryErrorEstimator &)=default
 Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this class. More...
 
 PatchRecoveryErrorEstimator (PatchRecoveryErrorEstimator &&)=default
 
PatchRecoveryErrorEstimatoroperator= (const PatchRecoveryErrorEstimator &)=default
 
PatchRecoveryErrorEstimatoroperator= (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. More...
 
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. More...
 
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. More...
 
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. 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) 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
 
int _extra_order
 Extra order to use for quadrature rule. More...
 

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.

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

65  :
69  patch_reuse(true),
70  _extra_order(1)
71 {
73 }
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...
unsigned int target_patch_size
The PatchErrorEstimator will build patches of at least this many elements to perform estimates...
ErrorEstimator()=default
Constructor.
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:75
int _extra_order
Extra order to use for quadrature rule.

◆ 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.

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

Referenced by main().

158 {
159  LOG_SCOPE("estimate_error()", "PatchRecoveryErrorEstimator");
160 
161  // The current mesh
162  const MeshBase & mesh = system.get_mesh();
163 
164  // Resize the error_per_cell vector to be
165  // the number of elements, initialize it to 0.
166  error_per_cell.resize (mesh.max_elem_id());
167  std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
168 
169  // Prepare current_local_solution to localize a non-standard
170  // solution vector if necessary
171  if (solution_vector && solution_vector != system.solution.get())
172  {
173  NumericVector<Number> * newsol =
174  const_cast<NumericVector<Number> *>(solution_vector);
175  System & sys = const_cast<System &>(system);
176  newsol->swap(*sys.solution);
177  sys.update();
178  }
179 
180  //------------------------------------------------------------
181  // Iterate over all the active elements in the mesh
182  // that live on this processor.
183  Threads::parallel_for (ConstElemRange(mesh.active_local_elements_begin(),
184  mesh.active_local_elements_end(),
185  200),
186  EstimateError(system,
187  *this,
188  error_per_cell)
189  );
190 
191  // Each processor has now computed the error contributions
192  // for its local elements, and error_per_cell contains 0 for all the
193  // non-local elements. Summing the vector will provide the true
194  // value for each element, local or remote
195  this->reduce_error(error_per_cell, system.comm());
196 
197  // If we used a non-standard solution before, now is the time to fix
198  // the current_local_solution
199  if (solution_vector && solution_vector != system.solution.get())
200  {
201  NumericVector<Number> * newsol =
202  const_cast<NumericVector<Number> *>(solution_vector);
203  System & sys = const_cast<System &>(system);
204  newsol->swap(*sys.solution);
205  sys.update();
206  }
207 }
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
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) const
This method takes the local error contributions in error_per_cell from each processor and combines th...
template class LIBMESH_EXPORT NumericVector< Number >

◆ estimate_errors() [1/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.

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

52 {
53  SystemNorm old_error_norm = this->error_norm;
54 
55  // Sum the error values from each system
56  for (auto s : make_range(equation_systems.n_systems()))
57  {
58  ErrorVector system_error_per_cell;
59  const System & sys = equation_systems.get_system(s);
60  if (error_norms.find(&sys) == error_norms.end())
61  this->error_norm = old_error_norm;
62  else
63  this->error_norm = error_norms.find(&sys)->second;
64 
65  const NumericVector<Number> * solution_vector = nullptr;
66  if (solution_vectors &&
67  solution_vectors->find(&sys) != solution_vectors->end())
68  solution_vector = solution_vectors->find(&sys)->second;
69 
70  this->estimate_error(sys, system_error_per_cell,
71  solution_vector, estimate_parent_error);
72 
73  if (s)
74  {
75  libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
76  for (auto i : index_range(error_per_cell))
77  error_per_cell[i] += system_error_per_cell[i];
78  }
79  else
80  error_per_cell = system_error_per_cell;
81  }
82 
83  // Restore our old state before returning
84  this->error_norm = old_error_norm;
85 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
template class LIBMESH_EXPORT NumericVector< Number >
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ estimate_errors() [2/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 93 of file error_estimator.C.

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().

97 {
98  SystemNorm old_error_norm = this->error_norm;
99 
100  // Find the requested error values from each system
101  for (auto s : make_range(equation_systems.n_systems()))
102  {
103  const System & sys = equation_systems.get_system(s);
104 
105  unsigned int n_vars = sys.n_vars();
106 
107  for (unsigned int v = 0; v != n_vars; ++v)
108  {
109  // Only fill in ErrorVectors the user asks for
110  if (!errors_per_cell.count(std::make_pair(&sys, v)))
111  continue;
112 
113  // Calculate error in only one variable
114  std::vector<Real> weights(n_vars, 0.0);
115  weights[v] = 1.0;
116  this->error_norm =
117  SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
118  weights);
119 
120  const NumericVector<Number> * solution_vector = nullptr;
121  if (solution_vectors &&
122  solution_vectors->find(&sys) != solution_vectors->end())
123  solution_vector = solution_vectors->find(&sys)->second;
124 
125  this->estimate_error
126  (sys, *errors_per_cell[std::make_pair(&sys, v)],
127  solution_vector, estimate_parent_error);
128  }
129  }
130 
131  // Restore our old state before returning
132  this->error_norm = old_error_norm;
133 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
unsigned int n_vars
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
template class LIBMESH_EXPORT NumericVector< Number >

◆ 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.

References _extra_order.

103  { _extra_order = extraorder; }
int _extra_order
Extra order to use for quadrature rule.

◆ 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.

References TIMPI::Communicator::sum().

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

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

◆ set_patch_reuse()

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

Definition at line 58 of file patch_recovery_error_estimator.C.

References patch_reuse.

59 {
60  patch_reuse = patch_reuse_flag;
61 }

◆ 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.

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

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

81 {
82  std::vector<Real> psi;
83  psi.reserve(matsize);
84  int npows = order+1;
85  std::vector<Real> xpow(npows,1.), ypow, zpow;
86  {
87  Real x = p(0);
88  for (auto i : make_range(1, npows))
89  xpow[i] = xpow[i-1] * x;
90  }
91  if (dim > 1)
92  {
93  Real y = p(1);
94  ypow.resize(npows,1.);
95  for (auto i : make_range(1, npows))
96  ypow[i] = ypow[i-1] * y;
97  }
98  if (dim > 2)
99  {
100  Real z = p(2);
101  zpow.resize(npows,1.);
102  for (auto i : make_range(1, npows))
103  zpow[i] = zpow[i-1] * z;
104  }
105 
106  // builds psi vector of form 1 x y z x^2 xy xz y^2 yz z^2 etc..
107  // I haven't added 1D support here
108  for (unsigned int poly_deg=0; poly_deg <= static_cast<unsigned int>(order) ; poly_deg++)
109  { // loop over all polynomials of total degree = poly_deg
110 
111  switch (dim)
112  {
113  // 3D spectral polynomial basis functions
114  case 3:
115  {
116  for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
117  for (int yexp=poly_deg-xexp; yexp >= 0; yexp--)
118  {
119  int zexp = poly_deg - xexp - yexp;
120  psi.push_back(xpow[xexp]*ypow[yexp]*zpow[zexp]);
121  }
122  break;
123  }
124 
125  // 2D spectral polynomial basis functions
126  case 2:
127  {
128  for (int xexp=poly_deg; xexp >= 0; xexp--) // use an int for xexp since we -- it
129  {
130  int yexp = poly_deg - xexp;
131  psi.push_back(xpow[xexp]*ypow[yexp]);
132  }
133  break;
134  }
135 
136  // 1D spectral polynomial basis functions
137  case 1:
138  {
139  int xexp = poly_deg;
140  psi.push_back(xpow[xexp]);
141  break;
142  }
143 
144  default:
145  libmesh_error_msg("Invalid dimension dim " << dim);
146  }
147  }
148 
149  return psi;
150 }
unsigned int dim
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134

◆ 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.

References libMesh::PATCH_RECOVERY.

Friends And Related Function Documentation

◆ EstimateError

friend class EstimateError
friend

Definition at line 149 of file patch_recovery_error_estimator.h.

Referenced by estimate_error().

Member Data Documentation

◆ _extra_order

int libMesh::PatchRecoveryErrorEstimator::_extra_order
protected

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

◆ patch_reuse

bool libMesh::PatchRecoveryErrorEstimator::patch_reuse
protected

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


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