Reference documentation for deal.II version 9.3.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 - 2021 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 
136 {
137 public:
175  const std::vector<std::string> & component_names = {"u"},
176  const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms = {
178 
244  const std::vector<std::string> & component_names,
245  const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
246  const double exponent,
247  const std::set<std::string> & extra_columns,
248  const std::string & rate_key,
249  const std::string & rate_mode,
250  const std::string & error_file_name,
251  const unsigned int precision,
252  const bool compute_error);
253 
259  void
261 
273  template <int dim, int spacedim, typename VectorType>
274  void
276  const VectorType & solution,
277  const Function<spacedim> & exact,
278  const Function<spacedim> * weight = nullptr);
279 
283  template <int dim, int spacedim, typename VectorType>
284  void
286  const DoFHandler<dim, spacedim> &vspace,
287  const VectorType & solution,
288  const Function<spacedim> & exact,
289  const Function<spacedim> * weight = nullptr);
290 
357  void
358  add_extra_column(const std::string & column_name,
359  const std::function<double()> &custom_function,
360  const bool compute_rate = true);
361 
365  template <int dim, int spacedim, typename VectorType>
366  void
368  const VectorType &,
369  const VectorType &,
370  const Function<spacedim> *weight = nullptr);
371 
375  template <int dim, int spacedim, typename VectorType>
376  void
377  difference(const Mapping<dim, spacedim> &mapping,
379  const VectorType &,
380  const VectorType &,
381  const Function<spacedim> *weight = nullptr);
382 
388  void
389  output_table(std::ostream &out);
390 
397  void
398  output_table();
399 
400 private:
404  void
406 
410  const std::vector<std::string> component_names;
411 
415  const std::vector<std::string> unique_component_names;
416 
420  const std::vector<ComponentMask> unique_component_masks;
421 
425  std::map<std::string, std::pair<std::function<double()>, bool>>
427 
431  std::vector<std::set<VectorTools::NormType>> norms_per_unique_component;
432 
436  double exponent;
437 
442 
446  std::set<std::string> extra_columns;
447 
451  std::string rate_key;
452 
456  std::string rate_mode;
457 
461  unsigned int precision;
462 
466  std::string error_file_name;
467 
473 };
474 
475 
476 
477 #ifndef DOXYGEN
478 // ============================================================
479 // Template functions
480 // ============================================================
481 template <int dim, int spacedim, typename VectorType>
482 void
484  const VectorType & solution1,
485  const VectorType & solution2,
486  const Function<spacedim> * weight)
487 {
488  AssertThrow(solution1.size() == solution2.size(),
489  ExcDimensionMismatch(solution1.size(), solution2.size()));
490  VectorType solution(solution1);
491  solution -= solution2;
493  dh,
494  solution,
496  weight);
497 }
498 
499 
500 
501 template <int dim, int spacedim, typename VectorType>
502 void
504  const DoFHandler<dim, spacedim> &dh,
505  const VectorType & solution1,
506  const VectorType & solution2,
507  const Function<spacedim> * weight)
508 {
509  AssertThrow(solution1.size() == solution2.size(),
510  ExcDimensionMismatch(solution1.size(), solution2.size()));
511  VectorType solution(solution1);
512  solution -= solution2;
514  mapping,
515  dh,
516  solution,
518  weight);
519 }
520 
521 
522 
523 template <int dim, int spacedim, typename VectorType>
524 void
526  const VectorType & solution,
527  const Function<spacedim> &exact,
528  const Function<spacedim> *weight)
529 {
531  dh,
532  solution,
533  exact,
534  weight);
535 }
536 
537 
538 
539 template <int dim, int spacedim, typename VectorType>
540 void
542  const DoFHandler<dim, spacedim> &dh,
543  const VectorType & solution,
544  const Function<spacedim> &exact,
545  const Function<spacedim> *weight)
546 {
547  const auto n_components = component_names.size();
548 
549  if (compute_error)
550  {
551  AssertDimension(exact.n_components, n_components);
552  AssertDimension(dh.get_fe().n_components(), n_components);
553 
555  dh.get_triangulation().n_global_active_cells();
556  const unsigned int n_dofs = dh.n_dofs();
557 
558  for (const auto &col : extra_columns)
559  if (col == "cells")
560  {
561  table.add_value("cells", n_active_cells);
562  table.set_tex_caption("cells", "\\# cells");
563  table.set_tex_format("cells", "r");
564  }
565  else if (col == "dofs")
566  {
567  table.add_value("dofs", n_dofs);
568  table.set_tex_caption("dofs", "\\# dofs");
569  table.set_tex_format("dofs", "r");
570  }
571 
572  // A vector of zero std::functions with n_components components
573  const std::vector<std::function<double(const Point<spacedim> &)>>
574  zero_components(n_components,
575  [](const Point<spacedim> &) { return 0.0; });
576 
577  // The default weight function, with n_components components
578  std::vector<std::function<double(const Point<spacedim> &)>>
579  weight_components(n_components,
580  [](const Point<spacedim> &) { return 1.0; });
581 
582  if (weight != nullptr)
583  {
584  if (weight->n_components == 1)
585  {
586  for (auto &f : weight_components)
587  f = [&](const Point<spacedim> &p) { return weight->value(p); };
588  }
589  else
590  {
591  AssertDimension(weight->n_components, n_components);
592  for (unsigned int i = 0; i < n_components; ++i)
593  weight_components[i] = [&](const Point<spacedim> &p) {
594  return weight->value(p, i);
595  };
596  }
597  }
598 
599  for (unsigned int i = 0; i < norms_per_unique_component.size(); ++i)
600  {
601  std::map<VectorTools::NormType, double> errors;
602 
603  const auto &norms = norms_per_unique_component[i];
604  const auto &mask = unique_component_masks[i];
605 
606  // Simple case first
607  if (norms.empty())
608  continue;
609 
610  auto components_expr = zero_components;
611  for (unsigned int i = 0; i < n_components; ++i)
612  if (mask[i] == true)
613  components_expr[i] = weight_components[i];
614 
615  FunctionFromFunctionObjects<spacedim> select_component(
616  components_expr);
617 
618  Vector<float> difference_per_cell(
619  dh.get_triangulation().n_global_active_cells());
620 
621  QGauss<dim> q_gauss((dh.get_fe().degree + 1) * 2);
622 
623  for (const auto &norm : norms)
624  {
625  difference_per_cell = 0;
627  dh,
628  solution,
629  exact,
630  difference_per_cell,
631  q_gauss,
632  norm,
633  &select_component,
634  exponent);
635 
637  dh.get_triangulation(), difference_per_cell, norm, exponent);
638 
639  std::string name = unique_component_names[i] + "_" +
641  std::string latex_name = "@f$\\| " + unique_component_names[i] +
642  " - " + unique_component_names[i] +
643  "_h \\|_{" +
644  Patterns::Tools::to_string(norm) + "}@f$";
645 
646  table.add_value(name, errors[norm]);
648  table.set_scientific(name, true);
649  table.set_tex_caption(name, latex_name);
650  }
651  }
652 
653  for (const auto &extra_col : extra_column_functions)
654  {
655  const double custom_error = extra_col.second.first();
656 
657  std::string name = extra_col.first;
658  table.add_value(name, custom_error);
660  table.set_scientific(name, true);
661  }
662  }
663 }
664 
665 #endif
666 
668 
669 #endif
void set_precision(const std::string &key, const unsigned int precision)
std::set< std::string > extra_columns
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
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
Definition: function.h:164
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
Definition: mapping.cc:240
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)
Definition: divergence.h:472
void set_tex_caption(const std::string &key, const std::string &tex_caption)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1575
void set_tex_format(const std::string &key, const std::string &format="c")
void difference(const DoFHandler< dim, spacedim > &, const VectorType &, const VectorType &, const Function< spacedim > *weight=nullptr)
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const std::vector< std::string > component_names
std::vector< std::set< VectorTools::NormType > > norms_per_unique_component
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
types::global_dof_index n_dofs() const
Abstract base class for mapping classes.
Definition: mapping.h:303
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
std::string to_string(const T &t)
Definition: patterns.h:2329
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12594
const std::vector< ComponentMask > unique_component_masks
void set_scientific(const std::string &key, const bool scientific)
const Triangulation< dim, spacedim > & get_triangulation() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
const std::vector< std::string > unique_component_names
void add_extra_column(const std::string &column_name, const std::function< double()> &custom_function, const bool compute_rate=true)
The ParsedConvergenceTable class.
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
void error_from_exact(const DoFHandler< dim, spacedim > &vspace, const VectorType &solution, const Function< spacedim > &exact, const Function< spacedim > *weight=nullptr)
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.)
double compute_global_error(const Triangulation< dim, spacedim > &tria, const InVector &cellwise_error, const NormType &norm, const double exponent=2.)
std::map< std::string, std::pair< std::function< double()>, bool > > extra_column_functions