15#ifndef dealii_base_parsed_convergence_table_h
16#define dealii_base_parsed_convergence_table_h
164 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms = {
233 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
261 template <
int dim,
int spacedim,
typename VectorType>
264 const VectorType &solution,
271 template <
int dim,
int spacedim,
typename VectorType>
275 const VectorType &solution,
347 const std::function<
double()> &custom_function,
348 const bool compute_rate =
true);
353 template <
int dim,
int spacedim,
typename VectorType>
363 template <
int dim,
int spacedim,
typename VectorType>
413 std::map<std::string, std::pair<std::function<double()>,
bool>>
469template <
int dim,
int spacedim,
typename VectorType>
472 const VectorType &solution1,
473 const VectorType &solution2,
478 VectorType solution(solution1);
479 solution -= solution2;
489template <
int dim,
int spacedim,
typename VectorType>
493 const VectorType &solution1,
494 const VectorType &solution2,
500 solution -= solution2;
511template <
int dim,
int spacedim,
typename VectorType>
514 const VectorType &solution,
527template <
int dim,
int spacedim,
typename VectorType>
531 const VectorType &solution,
544 const unsigned int n_dofs = dh.
n_dofs();
549 table.add_value(
"cells", n_active_cells);
550 table.set_tex_caption(
"cells",
"\\# cells");
551 table.set_tex_format(
"cells",
"r");
553 else if (col ==
"dofs")
555 table.add_value(
"dofs", n_dofs);
556 table.set_tex_caption(
"dofs",
"\\# dofs");
557 table.set_tex_format(
"dofs",
"r");
561 const std::vector<std::function<double(
const Point<spacedim> &)>>
562 zero_components(n_components,
563 [](
const Point<spacedim> &) {
return 0.0; });
566 std::vector<std::function<double(
const Point<spacedim> &)>>
567 weight_components(n_components,
568 [](
const Point<spacedim> &) {
return 1.0; });
570 if (weight !=
nullptr)
574 for (
auto &f : weight_components)
575 f = [&](
const Point<spacedim> &p) {
return weight->
value(p); };
580 for (
unsigned int i = 0; i < n_components; ++i)
581 weight_components[i] = [&](
const Point<spacedim> &p) {
582 return weight->
value(p, i);
589 std::map<VectorTools::NormType, double> errors;
598 auto components_expr = zero_components;
599 for (
unsigned int j = 0; j < n_components; ++j)
601 components_expr[j] = weight_components[j];
603 FunctionFromFunctionObjects<spacedim> select_component(
606 Vector<float> difference_per_cell(
611 for (
const auto &norm : norms)
613 difference_per_cell = 0;
634 table.add_value(name, errors[norm]);
636 table.set_scientific(name,
true);
637 table.set_tex_caption(name, latex_name);
643 const double custom_error = extra_col.second.first();
645 std::string name = extra_col.first;
646 table.add_value(name, custom_error);
648 table.set_scientific(name,
true);
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const Triangulation< dim, spacedim > & get_triangulation() const
types::global_dof_index n_dofs() const
const unsigned int degree
unsigned int n_components() const
const unsigned int n_components
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
Abstract base class for mapping classes.
void add_extra_column(const std::string &column_name, const std::function< double()> &custom_function, const bool compute_rate=true)
void difference(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &, const VectorType &, const VectorType &, const Function< spacedim > *weight=nullptr)
void difference(const DoFHandler< dim, spacedim > &, const VectorType &, const VectorType &, const Function< spacedim > *weight=nullptr)
void add_parameters(ParameterHandler &prm)
void error_from_exact(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &vspace, const VectorType &solution, const Function< spacedim > &exact, const Function< spacedim > *weight=nullptr)
void error_from_exact(const DoFHandler< dim, spacedim > &vspace, const VectorType &solution, const Function< spacedim > &exact, const Function< spacedim > *weight=nullptr)
std::string error_file_name
void prepare_table_for_output()
const std::vector< std::string > component_names
std::set< std::string > extra_columns
ParsedConvergenceTable(const std::vector< std::string > &component_names={"u"}, const std::vector< std::set< VectorTools::NormType > > &list_of_error_norms={ {VectorTools::H1_norm, VectorTools::L2_norm, VectorTools::Linfty_norm}})
const std::vector< std::string > unique_component_names
std::vector< std::set< VectorTools::NormType > > norms_per_unique_component
const std::vector< ComponentMask > unique_component_masks
std::map< std::string, std::pair< std::function< double()>, bool > > extra_column_functions
virtual types::global_cell_index n_global_active_cells() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define AssertThrow(cond, exc)
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
MappingQ< dim, spacedim > StaticMappingQ1< dim, spacedim >::mapping
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &Du)
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
unsigned int global_cell_index