16 #ifndef dealii_base_parsed_convergence_table_h 17 #define dealii_base_parsed_convergence_table_h 19 #include <deal.II/base/config.h> 21 #include <deal.II/base/convergence_table.h> 22 #include <deal.II/base/function.h> 23 #include <deal.II/base/function_parser.h> 24 #include <deal.II/base/parameter_handler.h> 25 #include <deal.II/base/quadrature_lib.h> 26 #include <deal.II/base/utilities.h> 28 #include <deal.II/dofs/dof_handler.h> 30 #include <deal.II/fe/mapping_q.h> 32 #include <deal.II/grid/grid_tools.h> 33 #include <deal.II/grid/tria.h> 35 #include <deal.II/lac/vector.h> 37 #include <deal.II/numerics/vector_tools.h> 39 DEAL_II_NAMESPACE_OPEN
178 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms = {
247 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
275 template <
typename DoFHandlerType,
typename VectorType>
278 const DoFHandlerType & vspace,
279 const VectorType & solution,
286 template <
typename DoFHandlerType,
typename VectorType>
291 const DoFHandlerType & vspace,
292 const VectorType & solution,
364 const std::function<
double()> &custom_function,
365 const bool & compute_rate =
true);
370 template <
typename DoFHandlerType,
typename VectorType>
380 template <
typename DoFHandlerType,
typename VectorType>
383 DoFHandlerType::space_dimension> &mapping,
384 const DoFHandlerType &,
431 std::map<std::string, std::pair<std::function<double()>,
bool>>
487 template <
typename DoFHandlerType,
typename VectorType>
490 const DoFHandlerType & dh,
491 const VectorType & solution1,
492 const VectorType & solution2,
497 VectorType solution(solution1);
498 solution -= solution2;
508 template <
typename DoFHandlerType,
typename VectorType>
513 const DoFHandlerType & dh,
514 const VectorType & solution1,
515 const VectorType & solution2,
520 VectorType solution(solution1);
521 solution -= solution2;
532 template <
typename DoFHandlerType,
typename VectorType>
535 const DoFHandlerType & dh,
536 const VectorType & solution,
541 DoFHandlerType::space_dimension>::mapping,
550 template <
typename DoFHandlerType,
typename VectorType>
555 const DoFHandlerType & dh,
556 const VectorType & solution,
560 const int dim = DoFHandlerType::dimension;
561 const int spacedim = DoFHandlerType::space_dimension;
569 const unsigned int n_active_cells =
570 dh.get_triangulation().n_global_active_cells();
571 const unsigned int n_dofs = dh.n_dofs();
580 else if (col ==
"dofs")
588 const std::vector<std::function<double(const Point<spacedim> &)>>
589 zero_components(n_components,
593 std::vector<std::function<double(const Point<spacedim> &)>>
594 weight_components(n_components,
597 if (weight !=
nullptr)
601 for (
auto &f : weight_components)
609 return weight->
value(p, i);
616 std::map<VectorTools::NormType, double> errors;
625 auto components_expr = zero_components;
628 components_expr[i] = weight_components[i];
633 Vector<float> difference_per_cell(
634 dh.get_triangulation().n_global_active_cells());
638 for (
const auto &norm : norms)
640 difference_per_cell = 0;
652 dh.get_triangulation(), difference_per_cell,
norm,
exponent);
670 const double custom_error = extra_col.second.first();
672 std::string name = extra_col.first;
682 DEAL_II_NAMESPACE_CLOSE
void set_precision(const std::string &key, const unsigned int precision)
std::set< std::string > extra_columns
#define AssertDimension(dim1, dim2)
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 unsigned int n_components
void add_parameters(ParameterHandler &prm)
void add_value(const std::string &key, const T value)
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim >>> &Du)
void set_tex_caption(const std::string &key, const std::string &tex_caption)
#define AssertThrow(cond, exc)
void set_tex_format(const std::string &key, const std::string &format="c")
const std::vector< std::string > component_names
void difference(const DoFHandlerType &, const VectorType &, const VectorType &, const Function< DoFHandlerType::space_dimension > *weight=nullptr)
void error_from_exact(const DoFHandlerType &vspace, const VectorType &solution, const Function< DoFHandlerType::space_dimension > &exact, const Function< DoFHandlerType::space_dimension > *weight=nullptr)
std::vector< std::set< VectorTools::NormType > > norms_per_unique_component
void prepare_table_for_output()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
std::string error_file_name
const std::vector< ComponentMask > unique_component_masks
void set_scientific(const std::string &key, const bool scientific)
const std::vector< std::string > unique_component_names
The ParsedConvergenceTable class.
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
void add_extra_column(const std::string &column_name, const std::function< double()> &custom_function, const bool &compute_rate=true)
std::map< std::string, std::pair< std::function< double()>, bool > > extra_column_functions