Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3087-ga35b476a78 2025-04-19 20:30:01+00:00
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
parsed_convergence_table.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2019 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_base_parsed_convergence_table_h
16#define dealii_base_parsed_convergence_table_h
17
18#include <deal.II/base/config.h>
19
24
26
27#include <deal.II/fe/mapping.h>
28
29#include <deal.II/grid/tria.h>
30
31#include <deal.II/lac/vector.h>
32
34
35#include <set>
36#include <string>
37#include <vector>
38
39
41
42#ifndef DOXYGEN
44#endif
45
129{
130public:
168 const std::vector<std::string> &component_names = {"u"},
169 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms = {
171
237 const std::vector<std::string> &component_names,
238 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
239 const double exponent,
240 const std::set<std::string> &extra_columns,
241 const std::string &rate_key,
242 const std::string &rate_mode,
243 const std::string &error_file_name,
244 const unsigned int precision,
245 const bool compute_error);
246
252 void
254
266 template <int dim, int spacedim, typename VectorType>
267 void
269 const VectorType &solution,
270 const Function<spacedim> &exact,
271 const Function<spacedim> *weight = nullptr);
272
276 template <int dim, int spacedim, typename VectorType>
277 void
279 const DoFHandler<dim, spacedim> &vspace,
280 const VectorType &solution,
281 const Function<spacedim> &exact,
282 const Function<spacedim> *weight = nullptr);
283
350 void
351 add_extra_column(const std::string &column_name,
352 const std::function<double()> &custom_function,
353 const bool compute_rate = true);
354
358 template <int dim, int spacedim, typename VectorType>
359 void
361 const VectorType &,
362 const VectorType &,
363 const Function<spacedim> *weight = nullptr);
364
368 template <int dim, int spacedim, typename VectorType>
369 void
372 const VectorType &,
373 const VectorType &,
374 const Function<spacedim> *weight = nullptr);
375
381 void
382 output_table(std::ostream &out);
383
390 void
391 output_table();
392
393private:
397 void
399
403 const std::vector<std::string> component_names;
404
408 const std::vector<std::string> unique_component_names;
409
413 const std::vector<ComponentMask> unique_component_masks;
414
418 std::map<std::string, std::pair<std::function<double()>, bool>>
420
424 std::vector<std::set<VectorTools::NormType>> norms_per_unique_component;
425
429 double exponent;
430
435
439 std::set<std::string> extra_columns;
440
444 std::string rate_key;
445
449 std::string rate_mode;
450
454 unsigned int precision;
455
459 std::string error_file_name;
460
466};
467
468
469
470#ifndef DOXYGEN
471// ============================================================
472// Template functions
473// ============================================================
474template <int dim, int spacedim, typename VectorType>
475void
477 const VectorType &solution1,
478 const VectorType &solution2,
479 const Function<spacedim> *weight)
480{
481 AssertThrow(solution1.size() == solution2.size(),
482 ExcDimensionMismatch(solution1.size(), solution2.size()));
483 VectorType solution(solution1);
484 solution -= solution2;
486 dh,
487 solution,
489 weight);
490}
491
492
493
494template <int dim, int spacedim, typename VectorType>
495void
498 const VectorType &solution1,
499 const VectorType &solution2,
500 const Function<spacedim> *weight)
501{
502 AssertThrow(solution1.size() == solution2.size(),
503 ExcDimensionMismatch(solution1.size(), solution2.size()));
504 VectorType solution(solution1);
505 solution -= solution2;
507 mapping,
508 dh,
509 solution,
511 weight);
512}
513
514
515
516template <int dim, int spacedim, typename VectorType>
517void
519 const VectorType &solution,
520 const Function<spacedim> &exact,
521 const Function<spacedim> *weight)
522{
524 dh,
525 solution,
526 exact,
527 weight);
528}
529
530
531
532template <int dim, int spacedim, typename VectorType>
533void
536 const VectorType &solution,
537 const Function<spacedim> &exact,
538 const Function<spacedim> *weight)
539{
540 const auto n_components = component_names.size();
541
542 if (compute_error)
543 {
544 AssertDimension(exact.n_components, n_components);
545 AssertDimension(dh.get_fe().n_components(), n_components);
546
549 const unsigned int n_dofs = dh.n_dofs();
550
551 for (const auto &col : extra_columns)
552 if (col == "cells")
553 {
554 table.add_value("cells", n_active_cells);
555 table.set_tex_caption("cells", "\\# cells");
556 table.set_tex_format("cells", "r");
557 }
558 else if (col == "dofs")
559 {
560 table.add_value("dofs", n_dofs);
561 table.set_tex_caption("dofs", "\\# dofs");
562 table.set_tex_format("dofs", "r");
563 }
564
565 // A vector of zero std::functions with n_components components
566 const std::vector<std::function<double(const Point<spacedim> &)>>
567 zero_components(n_components,
568 [](const Point<spacedim> &) { return 0.0; });
569
570 // The default weight function, with n_components components
571 std::vector<std::function<double(const Point<spacedim> &)>>
572 weight_components(n_components,
573 [](const Point<spacedim> &) { return 1.0; });
574
575 if (weight != nullptr)
576 {
577 if (weight->n_components == 1)
578 {
579 for (auto &f : weight_components)
580 f = [&](const Point<spacedim> &p) { return weight->value(p); };
581 }
582 else
583 {
584 AssertDimension(weight->n_components, n_components);
585 for (unsigned int i = 0; i < n_components; ++i)
586 weight_components[i] = [&](const Point<spacedim> &p) {
587 return weight->value(p, i);
588 };
589 }
590 }
591
592 for (unsigned int i = 0; i < norms_per_unique_component.size(); ++i)
593 {
594 std::map<VectorTools::NormType, double> errors;
595
596 const auto &norms = norms_per_unique_component[i];
597 const auto &mask = unique_component_masks[i];
598
599 // Simple case first
600 if (norms.empty())
601 continue;
602
603 auto components_expr = zero_components;
604 for (unsigned int j = 0; j < n_components; ++j)
605 if (mask[j] == true)
606 components_expr[j] = weight_components[j];
607
609 components_expr);
610
611 Vector<float> difference_per_cell(
613
614 QGauss<dim> q_gauss((dh.get_fe().degree + 1) * 2);
615
616 for (const auto &norm : norms)
617 {
618 difference_per_cell = 0;
620 dh,
621 solution,
622 exact,
623 difference_per_cell,
624 q_gauss,
625 norm,
626 &select_component,
627 exponent);
628
630 dh.get_triangulation(), difference_per_cell, norm, exponent);
631
632 std::string name = unique_component_names[i] + "_" +
634 std::string latex_name = "@f$\\| " + unique_component_names[i] +
635 " - " + unique_component_names[i] +
636 "_h \\|_{" +
637 Patterns::Tools::to_string(norm) + "}@f$";
638
639 table.add_value(name, errors[norm]);
641 table.set_scientific(name, true);
642 table.set_tex_caption(name, latex_name);
643 }
644 }
645
646 for (const auto &extra_col : extra_column_functions)
647 {
648 const double custom_error = extra_col.second.first();
649
650 std::string name = extra_col.first;
651 table.add_value(name, custom_error);
653 table.set_scientific(name, true);
654 }
655 }
656}
657
658#endif
659
661
662#endif
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
Definition fe_data.h:452
unsigned int n_components() const
const unsigned int n_components
Definition function.h:164
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
Abstract base class for mapping classes.
Definition mapping.h:320
The ParsedConvergenceTable class.
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)
const std::vector< std::string > component_names
std::set< std::string > extra_columns
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
Definition point.h:113
void set_tex_format(const std::string &key, const std::string &format="c")
void add_value(const std::string &key, const T value)
void set_tex_caption(const std::string &key, const std::string &tex_caption)
void set_scientific(const std::string &key, const bool scientific)
void set_precision(const std::string &key, const unsigned int precision)
virtual types::global_cell_index n_global_active_cells() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
#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)
Definition mapping.cc:316
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &Du)
Definition divergence.h:471
std::string to_string(const T &t)
Definition patterns.h:2452
double compute_global_error(const Triangulation< dim, spacedim > &tria, const InVector &cellwise_error, const NormType &norm, const double exponent=2.)
void integrate_difference(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const ReadVector< Number > &fe_function, const Function< spacedim, Number > &exact_solution, OutVector &difference, const Quadrature< dim > &q, const NormType &norm, const Function< spacedim, double > *weight=nullptr, const double exponent=2.)
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition tria.cc:14911