libMesh
|
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>
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 | |
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 | |
UniformRefinementEstimator & | operator= (const UniformRefinementEstimator &)=default |
UniformRefinementEstimator & | operator= (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... | |
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 |
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, ErrorMap *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... | |
Protected Attributes | |
int | _extra_order |
Extra order to use for quadrature rule. More... | |
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.
Definition at line 45 of file uniform_refinement_estimator.h.
|
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.
libMesh::UniformRefinementEstimator::UniformRefinementEstimator | ( | ) |
Constructor.
Sets the most common default parameter values.
Definition at line 54 of file uniform_refinement_estimator.C.
References libMesh::ErrorEstimator::error_norm, and libMesh::H1.
|
default |
Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this simple class.
|
default |
|
virtualdefault |
|
protectedvirtual |
The code for estimate_error and both estimate_errors versions is very similar, so we use the same function for all three.
Definition at line 119 of file uniform_refinement_estimator.C.
References _extra_order, libMesh::EquationSystems::adjoint_solve(), libMesh::System::adjoint_solve(), libMesh::FEGenericBase< OutputType >::build(), libMesh::NumericVector< T >::build(), libMesh::NumericVector< T >::clear(), libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, libMesh::System::current_solution(), libMesh::FEType::default_quadrature_rule(), dim, libMesh::System::disable_cache(), libMesh::DofMap::dof_indices(), libMesh::ErrorEstimator::error_norm, libMesh::ErrorVectorReal, libMesh::NumericVector< T >::get(), libMesh::System::get_adjoint_solution(), libMesh::System::get_dof_map(), libMesh::System::get_equation_systems(), libMesh::EquationSystems::get_mesh(), libMesh::DofMap::get_send_list(), libMesh::EquationSystems::get_system(), libMesh::System::get_vector(), libMesh::GHOSTED, libMesh::H1, libMesh::H1_SEMINORM, libMesh::H2, libMesh::H2_SEMINORM, libMesh::DofObject::id(), libMesh::index_range(), libMesh::invalid_uint, libMesh::L2, libMesh::libmesh_assert(), libMesh::make_range(), mesh, libMesh::System::n_qois(), libMesh::EquationSystems::n_systems(), n_vars, libMesh::System::n_vars(), libMesh::TensorTools::norm_sq(), libMesh::TypeVector< T >::norm_sq(), libMesh::TypeTensor< T >::norm_sq(), number_h_refinements, number_p_refinements, libMesh::Elem::parent(), libMesh::System::project_solution_on_reinit(), libMesh::Real, libMesh::ErrorEstimator::reduce_error(), libMesh::EquationSystems::reinit(), libMesh::System::request_vector(), libMesh::System::solution, libMesh::EquationSystems::solve(), libMesh::System::solve(), std::sqrt(), libMesh::NumericVector< T >::swap(), libMesh::SystemNorm::type(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::MeshRefinement::uniformly_p_coarsen(), libMesh::MeshRefinement::uniformly_p_refine(), libMesh::MeshRefinement::uniformly_refine(), libMesh::System::update(), libMesh::DofMap::variable_type(), libMesh::System::vectors_begin(), libMesh::System::vectors_end(), and libMesh::SystemNorm::weight_sq().
Referenced by estimate_error(), and estimate_errors().
|
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.
Definition at line 71 of file uniform_refinement_estimator.C.
References _estimate_error().
Referenced by main().
|
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.
Definition at line 88 of file uniform_refinement_estimator.C.
References _estimate_error().
|
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.
Definition at line 104 of file uniform_refinement_estimator.C.
References _estimate_error().
|
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 118 of file uniform_refinement_estimator.h.
References _extra_order.
|
default |
|
default |
|
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 _estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), and libMesh::ExactErrorEstimator::estimate_error().
|
overridevirtual |
Implements libMesh::ErrorEstimator.
Definition at line 65 of file uniform_refinement_estimator.C.
References libMesh::UNIFORM_REFINEMENT.
|
protected |
Extra order to use for quadrature rule.
Definition at line 138 of file uniform_refinement_estimator.h.
Referenced by _estimate_error(), and extra_quadrature_order().
|
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 _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()(), libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and UniformRefinementEstimator().
unsigned char libMesh::UniformRefinementEstimator::number_h_refinements |
How many h refinements to perform to get the fine grid.
Definition at line 126 of file uniform_refinement_estimator.h.
Referenced by _estimate_error().
unsigned char libMesh::UniformRefinementEstimator::number_p_refinements |
How many p refinements to perform to get the fine grid.
Definition at line 131 of file uniform_refinement_estimator.h.
Referenced by _estimate_error().