16 #ifndef dealii_error_estimator_h
17 #define dealii_error_estimator_h
35 template <
int dim,
int spacedim>
262 template <
int dim,
int spacedim = dim>
341 template <
typename Number>
362 template <
typename Number>
391 template <
typename Number>
412 template <
typename Number>
433 template <
typename Number>
455 template <
typename Number>
476 template <
typename Number>
498 template <
typename Number>
518 "You provided a ComponentMask argument that is invalid. "
519 "Component masks need to be either default constructed "
520 "(in which case they indicate that every component is "
521 "selected) or need to have a length equal to the number "
522 "of vector components of the finite element in use "
523 "by the DoFHandler object. In the latter case, at "
524 "least one component needs to be selected.");
530 "If you do specify the argument for a (possibly "
531 "spatially variable) coefficient function for this function, "
532 "then it needs to refer to a coefficient that is either "
533 "scalar (has one vector component) or has as many vector "
534 "components as there are in the finite element used by "
535 "the DoFHandler argument.");
543 <<
"You provided a function map that for boundary indicator "
544 << arg1 <<
" specifies a function with " << arg2
545 <<
" vector components. However, the finite "
546 "element in use has "
548 <<
" components, and these two numbers need to match.");
555 <<
"The number of input vectors, " << arg1
556 <<
" needs to be equal to the number of output vectors, "
558 <<
". This is not the case in your call of this function.");
563 "You need to specify at least one solution vector as "
578 template <
int spacedim>
621 template <
typename Number>
644 template <
typename Number>
675 template <
typename Number>
698 template <
typename Number>
720 template <
typename Number>
743 template <
typename Number>
765 template <
typename Number>
788 template <
typename Number>
810 "You provided a ComponentMask argument that is invalid. "
811 "Component masks need to be either default constructed "
812 "(in which case they indicate that every component is "
813 "selected) or need to have a length equal to the number "
814 "of vector components of the finite element in use "
815 "by the DoFHandler object. In the latter case, at "
816 "least one component needs to be selected.");
823 "If you do specify the argument for a (possibly "
824 "spatially variable) coefficient function for this function, "
825 "then it needs to refer to a coefficient that is either "
826 "scalar (has one vector component) or has as many vector "
827 "components as there are in the finite element used by "
828 "the DoFHandler argument.");
837 <<
"You provided a function map that for boundary indicator "
838 << arg1 <<
" specifies a function with " << arg2
839 <<
" vector components. However, the finite "
840 "element in use has "
842 <<
" components, and these two numbers need to match.");
850 <<
"The number of input vectors, " << arg1
851 <<
" needs to be equal to the number of output vectors, "
853 <<
". This is not the case in your call of this function.");
859 "You need to specify at least one solution vector as "
@ face_diameter_over_twice_max_degree
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
@ face_diameter_over_twice_max_degree
@ cell_diameter_over_24
Kelly error estimator with the factor .
@ cell_diameter
Kelly error estimator with the factor .
static void estimate(const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
Abstract base class for mapping classes.
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_CXX20_REQUIRES(condition)
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcIncompatibleNumberOfElements(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidCoefficient()
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcNoSolutions()
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInvalidBoundaryFunction(types::boundary_id arg1, int arg2, int arg3)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
static ::ExceptionBase & ExcInvalidComponentMask()
const types::material_id invalid_material_id
const types::subdomain_id invalid_subdomain_id
static const unsigned int invalid_unsigned_int
unsigned int subdomain_id