Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
parsed_convergence_table.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2019 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_base_parsed_convergence_table_h
17 #define dealii_base_parsed_convergence_table_h
18 
19 #include <deal.II/base/config.h>
20 
22 #include <deal.II/base/function.h>
26 #include <deal.II/base/utilities.h>
27 
29 
30 #include <deal.II/fe/mapping_q.h>
31 
33 #include <deal.II/grid/tria.h>
34 
35 #include <deal.II/lac/vector.h>
36 
38 
40 
138 {
139 public:
177  const std::vector<std::string> & component_names = {"u"},
178  const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms = {
180 
246  const std::vector<std::string> & component_names,
247  const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
248  const double exponent,
249  const std::set<std::string> & extra_columns,
250  const std::string & rate_key,
251  const std::string & rate_mode,
252  const std::string & error_file_name,
253  const unsigned int precision,
254  const bool compute_error);
255 
261  void
263 
275  template <typename DoFHandlerType, typename VectorType>
276  void
278  const DoFHandlerType & vspace,
279  const VectorType & solution,
281  const Function<DoFHandlerType::space_dimension> *weight = nullptr);
282 
286  template <typename DoFHandlerType, typename VectorType>
287  void
290  & mapping,
291  const DoFHandlerType & vspace,
292  const VectorType & solution,
294  const Function<DoFHandlerType::space_dimension> *weight = nullptr);
295 
362  void
363  add_extra_column(const std::string & column_name,
364  const std::function<double()> &custom_function,
365  const bool compute_rate = true);
366 
370  template <typename DoFHandlerType, typename VectorType>
371  void
372  difference(const DoFHandlerType &,
373  const VectorType &,
374  const VectorType &,
375  const Function<DoFHandlerType::space_dimension> *weight = nullptr);
376 
380  template <typename DoFHandlerType, typename VectorType>
381  void
382  difference(const Mapping<DoFHandlerType::dimension,
383  DoFHandlerType::space_dimension> &mapping,
384  const DoFHandlerType &,
385  const VectorType &,
386  const VectorType &,
387  const Function<DoFHandlerType::space_dimension> *weight = nullptr);
388 
394  void
395  output_table(std::ostream &out);
396 
403  void
404  output_table();
405 
406 private:
410  void
412 
416  const std::vector<std::string> component_names;
417 
421  const std::vector<std::string> unique_component_names;
422 
426  const std::vector<ComponentMask> unique_component_masks;
427 
431  std::map<std::string, std::pair<std::function<double()>, bool>>
433 
437  std::vector<std::set<VectorTools::NormType>> norms_per_unique_component;
438 
442  double exponent;
443 
448 
452  std::set<std::string> extra_columns;
453 
457  std::string rate_key;
458 
462  std::string rate_mode;
463 
467  unsigned int precision;
468 
472  std::string error_file_name;
473 
479 };
480 
481 
482 
483 #ifndef DOXYGEN
484 // ============================================================
485 // Template functions
486 // ============================================================
487 template <typename DoFHandlerType, typename VectorType>
488 void
490  const DoFHandlerType & dh,
491  const VectorType & solution1,
492  const VectorType & solution2,
494 {
495  AssertThrow(solution1.size() == solution2.size(),
496  ExcDimensionMismatch(solution1.size(), solution2.size()));
497  VectorType solution(solution1);
498  solution -= solution2;
499  error_from_exact(dh,
500  solution,
502  0, component_names.size()),
503  weight);
504 }
505 
506 
507 
508 template <typename DoFHandlerType, typename VectorType>
509 void
512  & mapping,
513  const DoFHandlerType & dh,
514  const VectorType & solution1,
515  const VectorType & solution2,
517 {
518  AssertThrow(solution1.size() == solution2.size(),
519  ExcDimensionMismatch(solution1.size(), solution2.size()));
520  VectorType solution(solution1);
521  solution -= solution2;
522  error_from_exact(mapping,
523  dh,
524  solution,
526  0, component_names.size()),
527  weight);
528 }
529 
530 
531 
532 template <typename DoFHandlerType, typename VectorType>
533 void
535  const DoFHandlerType & dh,
536  const VectorType & solution,
539 {
540  error_from_exact(StaticMappingQ1<DoFHandlerType::dimension,
541  DoFHandlerType::space_dimension>::mapping,
542  dh,
543  solution,
544  exact,
545  weight);
546 }
547 
548 
549 
550 template <typename DoFHandlerType, typename VectorType>
551 void
554  & mapping,
555  const DoFHandlerType & dh,
556  const VectorType & solution,
559 {
560  const int dim = DoFHandlerType::dimension;
561  const int spacedim = DoFHandlerType::space_dimension;
562  const auto n_components = component_names.size();
563 
564  if (compute_error)
565  {
567  AssertDimension(dh.get_fe().n_components(), n_components);
568 
570  dh.get_triangulation().n_global_active_cells();
571  const unsigned int n_dofs = dh.n_dofs();
572 
573  for (const auto &col : extra_columns)
574  if (col == "cells")
575  {
576  table.add_value("cells", n_active_cells);
577  table.set_tex_caption("cells", "\\# cells");
578  table.set_tex_format("cells", "r");
579  }
580  else if (col == "dofs")
581  {
582  table.add_value("dofs", n_dofs);
583  table.set_tex_caption("dofs", "\\# dofs");
584  table.set_tex_format("dofs", "r");
585  }
586 
587  // A vector of zero std::functions with n_components components
588  const std::vector<std::function<double(const Point<spacedim> &)>>
589  zero_components(n_components,
590  [](const Point<spacedim> &) { return 0.0; });
591 
592  // The default weight function, with n_components components
593  std::vector<std::function<double(const Point<spacedim> &)>>
594  weight_components(n_components,
595  [](const Point<spacedim> &) { return 1.0; });
596 
597  if (weight != nullptr)
598  {
599  if (weight->n_components == 1)
600  {
601  for (auto &f : weight_components)
602  f = [&](const Point<spacedim> &p) { return weight->value(p); };
603  }
604  else
605  {
607  for (unsigned int i = 0; i < n_components; ++i)
608  weight_components[i] = [&](const Point<spacedim> &p) {
609  return weight->value(p, i);
610  };
611  }
612  }
613 
614  for (unsigned int i = 0; i < norms_per_unique_component.size(); ++i)
615  {
616  std::map<VectorTools::NormType, double> errors;
617 
618  const auto &norms = norms_per_unique_component[i];
619  const auto &mask = unique_component_masks[i];
620 
621  // Simple case first
622  if (norms.empty())
623  continue;
624 
625  auto components_expr = zero_components;
626  for (unsigned int i = 0; i < n_components; ++i)
627  if (mask[i] == true)
628  components_expr[i] = weight_components[i];
629 
630  FunctionFromFunctionObjects<spacedim> select_component(
631  components_expr);
632 
633  Vector<float> difference_per_cell(
634  dh.get_triangulation().n_global_active_cells());
635 
636  QGauss<dim> q_gauss((dh.get_fe().degree + 1) * 2);
637 
638  for (const auto &norm : norms)
639  {
640  difference_per_cell = 0;
642  dh,
643  solution,
644  exact,
645  difference_per_cell,
646  q_gauss,
647  norm,
648  &select_component,
649  exponent);
650 
652  dh.get_triangulation(), difference_per_cell, norm, exponent);
653 
654  std::string name = unique_component_names[i] + "_" +
656  std::string latex_name = "@f$\\| " + unique_component_names[i] +
657  " - " + unique_component_names[i] +
658  "_h \\|_{" +
660 
661  table.add_value(name, errors[norm]);
663  table.set_scientific(name, true);
664  table.set_tex_caption(name, latex_name);
665  }
666  }
667 
668  for (const auto &extra_col : extra_column_functions)
669  {
670  const double custom_error = extra_col.second.first();
671 
672  std::string name = extra_col.first;
673  table.add_value(name, custom_error);
675  table.set_scientific(name, true);
676  }
677  }
678 }
679 
680 #endif
681 
683 
684 #endif
TableHandler::add_value
void add_value(const std::string &key, const T value)
Definition: table_handler.h:937
ParsedConvergenceTable
The ParsedConvergenceTable class.
Definition: parsed_convergence_table.h:137
TableHandler::set_tex_caption
void set_tex_caption(const std::string &key, const std::string &tex_caption)
Definition: table_handler.cc:309
ParsedConvergenceTable::error_file_name
std::string error_file_name
Definition: parsed_convergence_table.h:472
mapping_q.h
VectorTools::L2_norm
@ L2_norm
Definition: vector_tools_common.h:113
vector_tools_integrate_difference.h
TableHandler::set_scientific
void set_scientific(const std::string &key, const bool scientific)
Definition: table_handler.cc:371
ParsedConvergenceTable::rate_mode
std::string rate_mode
Definition: parsed_convergence_table.h:462
StaticMappingQ1
Definition: mapping_q1.h:88
tria.h
FunctionFromFunctionObjects
Definition: function.h:885
utilities.h
convergence_table.h
VectorTools::H1_norm
@ H1_norm
Definition: vector_tools_common.h:204
Functions::ConstantFunction
Definition: function.h:412
VectorType
Function::n_components
const unsigned int n_components
Definition: function.h:165
Patterns::Tools::to_string
std::string to_string(const T &t)
Definition: patterns.h:2360
VectorTools::integrate_difference
void integrate_difference(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InVector &fe_function, const Function< spacedim, typename InVector::value_type > &exact_solution, OutVector &difference, const Quadrature< dim > &q, const NormType &norm, const Function< spacedim, double > *weight=nullptr, const double exponent=2.)
ParsedConvergenceTable::unique_component_masks
const std::vector< ComponentMask > unique_component_masks
Definition: parsed_convergence_table.h:426
VectorTools::compute_global_error
double compute_global_error(const Triangulation< dim, spacedim > &tria, const InVector &cellwise_error, const NormType &norm, const double exponent=2.)
ParsedConvergenceTable::exponent
double exponent
Definition: parsed_convergence_table.h:442
ParsedConvergenceTable::unique_component_names
const std::vector< std::string > unique_component_names
Definition: parsed_convergence_table.h:421
DoFTools::n_components
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
ParsedConvergenceTable::add_extra_column
void add_extra_column(const std::string &column_name, const std::function< double()> &custom_function, const bool compute_rate=true)
Definition: parsed_convergence_table.cc:262
quadrature_lib.h
internal::TriangulationImplementation::n_active_cells
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12625
ParsedConvergenceTable::precision
unsigned int precision
Definition: parsed_convergence_table.h:467
Mapping
Abstract base class for mapping classes.
Definition: mapping.h:302
double
grid_tools.h
LocalIntegrators::Divergence::norm
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim >>> &Du)
Definition: divergence.h:548
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
ParsedConvergenceTable::extra_column_functions
std::map< std::string, std::pair< std::function< double()>, bool > > extra_column_functions
Definition: parsed_convergence_table.h:432
ParsedConvergenceTable::compute_error
bool compute_error
Definition: parsed_convergence_table.h:478
AssertDimension
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1579
parameter_handler.h
Function::value
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
function_parser.h
ParsedConvergenceTable::extra_columns
std::set< std::string > extra_columns
Definition: parsed_convergence_table.h:452
ParsedConvergenceTable::ParsedConvergenceTable
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}})
Definition: parsed_convergence_table.cc:62
ParsedConvergenceTable::add_parameters
void add_parameters(ParameterHandler &prm)
Definition: parsed_convergence_table.cc:107
ParsedConvergenceTable::output_table
void output_table()
Definition: parsed_convergence_table.cc:228
QGauss
Definition: quadrature_lib.h:40
ParsedConvergenceTable::table
ConvergenceTable table
Definition: parsed_convergence_table.h:447
unsigned int
function.h
ParsedConvergenceTable::rate_key
std::string rate_key
Definition: parsed_convergence_table.h:457
vector.h
dof_handler.h
VectorTools::Linfty_norm
@ Linfty_norm
Definition: vector_tools_common.h:148
StandardExceptions::ExcDimensionMismatch
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Point< spacedim >
ParsedConvergenceTable::component_names
const std::vector< std::string > component_names
Definition: parsed_convergence_table.h:416
ParsedConvergenceTable::error_from_exact
void error_from_exact(const DoFHandlerType &vspace, const VectorType &solution, const Function< DoFHandlerType::space_dimension > &exact, const Function< DoFHandlerType::space_dimension > *weight=nullptr)
ParameterHandler
Definition: parameter_handler.h:845
ParsedConvergenceTable::difference
void difference(const DoFHandlerType &, const VectorType &, const VectorType &, const Function< DoFHandlerType::space_dimension > *weight=nullptr)
config.h
Function
Definition: function.h:151
TableHandler::set_tex_format
void set_tex_format(const std::string &key, const std::string &format="c")
Definition: table_handler.cc:346
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
ParsedConvergenceTable::norms_per_unique_component
std::vector< std::set< VectorTools::NormType > > norms_per_unique_component
Definition: parsed_convergence_table.h:437
ConvergenceTable
Definition: convergence_table.h:63
AssertThrow
#define AssertThrow(cond, exc)
Definition: exceptions.h:1531
ParsedConvergenceTable::prepare_table_for_output
void prepare_table_for_output()
Definition: parsed_convergence_table.cc:161
TableHandler::set_precision
void set_precision(const std::string &key, const unsigned int precision)
Definition: table_handler.cc:358