libMesh
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | List of all members
libMesh::UniformRefinementEstimator Class Reference

This class implements a `‘brute force’' error estimator which integrates differences between the current solution and the solution on a uniformly refined (in h and/or p, for an arbitrary number of levels) grid. More...

#include <uniform_refinement_estimator.h>

Inheritance diagram for libMesh::UniformRefinementEstimator:
[legend]

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

 UniformRefinementEstimator ()
 Constructor. More...
 
 UniformRefinementEstimator (const UniformRefinementEstimator &)=default
 Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this simple class. More...
 
 UniformRefinementEstimator (UniformRefinementEstimator &&)=default
 
UniformRefinementEstimatoroperator= (const UniformRefinementEstimator &)=default
 
UniformRefinementEstimatoroperator= (UniformRefinementEstimator &&)=default
 
virtual ~UniformRefinementEstimator ()=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 does uniform refinements and a solve to get an improved solution on each cell, then estimates the error by integrating differences between the coarse and fine solutions. More...
 
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) override
 Currently this function ignores the error_norm member variable, and uses the function argument error_norms instead. 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) override
 Currently this function ignores the component_scale member variable, because it calculates each error individually, unscaled. More...
 
virtual ErrorEstimatorType type () const override
 

Public Attributes

unsigned char number_h_refinements
 How many h refinements to perform to get the fine grid. More...
 
unsigned char number_p_refinements
 How many p refinements to perform to get the fine grid. 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

virtual void _estimate_error (const EquationSystems *equation_systems, const System *system, ErrorVector *error_per_cell, std::map< std::pair< const System *, unsigned int >, ErrorVector *> *errors_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)
 The code for estimate_error and both estimate_errors versions is very similar, so we use the same function for all three. More...
 
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...
 

Detailed Description

This class implements a `‘brute force’' error estimator which integrates differences between the current solution and the solution on a uniformly refined (in h and/or p, for an arbitrary number of levels) grid.

Author
Roy H. Stogner
Date
2006

Definition at line 45 of file uniform_refinement_estimator.h.

Member Typedef Documentation

◆ ErrorMap

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

Constructor & Destructor Documentation

◆ UniformRefinementEstimator() [1/3]

libMesh::UniformRefinementEstimator::UniformRefinementEstimator ( )

Constructor.

Sets the most common default parameter values.

◆ UniformRefinementEstimator() [2/3]

libMesh::UniformRefinementEstimator::UniformRefinementEstimator ( const UniformRefinementEstimator )
default

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

◆ UniformRefinementEstimator() [3/3]

libMesh::UniformRefinementEstimator::UniformRefinementEstimator ( UniformRefinementEstimator &&  )
default

◆ ~UniformRefinementEstimator()

virtual libMesh::UniformRefinementEstimator::~UniformRefinementEstimator ( )
virtualdefault

Member Function Documentation

◆ _estimate_error()

virtual void libMesh::UniformRefinementEstimator::_estimate_error ( const EquationSystems equation_systems,
const System system,
ErrorVector error_per_cell,
std::map< std::pair< const System *, unsigned int >, ErrorVector *> *  errors_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 
)
protectedvirtual

The code for estimate_error and both estimate_errors versions is very similar, so we use the same function for all three.

◆ estimate_error()

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

This function does uniform refinements and a solve to get an improved solution on each cell, then estimates the error by integrating differences between the coarse and fine solutions.

system.solve() must be called, and so should have no side effects.

Only the provided system is solved on the refined mesh; for problems decoupled into multiple systems, use of estimate_errors() should be more reliable.

The estimated error is output in the vector error_per_cell

Implements libMesh::ErrorEstimator.

◆ estimate_errors() [1/2]

virtual void libMesh::UniformRefinementEstimator::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 
)
overridevirtual

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 from libMesh::ErrorEstimator.

◆ estimate_errors() [2/2]

virtual void libMesh::UniformRefinementEstimator::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 
)
overridevirtual

Currently this function ignores the component_scale member variable, because it calculates each 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 from libMesh::ErrorEstimator.

◆ operator=() [1/2]

UniformRefinementEstimator& libMesh::UniformRefinementEstimator::operator= ( const UniformRefinementEstimator )
default

◆ operator=() [2/2]

UniformRefinementEstimator& libMesh::UniformRefinementEstimator::operator= ( UniformRefinementEstimator &&  )
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.

◆ type()

virtual ErrorEstimatorType libMesh::UniformRefinementEstimator::type ( ) const
overridevirtual
Returns
The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Member Data Documentation

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

◆ number_h_refinements

unsigned char libMesh::UniformRefinementEstimator::number_h_refinements

How many h refinements to perform to get the fine grid.

Definition at line 116 of file uniform_refinement_estimator.h.

◆ number_p_refinements

unsigned char libMesh::UniformRefinementEstimator::number_p_refinements

How many p refinements to perform to get the fine grid.

Definition at line 121 of file uniform_refinement_estimator.h.


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