Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
parsed_convergence_table.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2019 - 2023 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
25
27
28#include <deal.II/fe/mapping.h>
29
30#include <deal.II/grid/tria.h>
31
32#include <deal.II/lac/vector.h>
33
35
37
38#ifndef DOXYGEN
40#endif
41
137{
138public:
176 const std::vector<std::string> & component_names = {"u"},
177 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms = {
179
245 const std::vector<std::string> & component_names,
246 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
247 const double exponent,
248 const std::set<std::string> & extra_columns,
249 const std::string & rate_key,
250 const std::string & rate_mode,
251 const std::string & error_file_name,
252 const unsigned int precision,
253 const bool compute_error);
254
260 void
262
274 template <int dim, int spacedim, typename VectorType>
275 void
277 const VectorType & solution,
278 const Function<spacedim> & exact,
279 const Function<spacedim> * weight = nullptr);
280
284 template <int dim, int spacedim, typename VectorType>
285 void
287 const DoFHandler<dim, spacedim> &vspace,
288 const VectorType & solution,
289 const Function<spacedim> & exact,
290 const Function<spacedim> * weight = nullptr);
291
358 void
359 add_extra_column(const std::string & column_name,
360 const std::function<double()> &custom_function,
361 const bool compute_rate = true);
362
366 template <int dim, int spacedim, typename VectorType>
367 void
369 const VectorType &,
370 const VectorType &,
371 const Function<spacedim> *weight = nullptr);
372
376 template <int dim, int spacedim, typename VectorType>
377 void
380 const VectorType &,
381 const VectorType &,
382 const Function<spacedim> *weight = nullptr);
383
389 void
390 output_table(std::ostream &out);
391
398 void
399 output_table();
400
401private:
405 void
407
411 const std::vector<std::string> component_names;
412
416 const std::vector<std::string> unique_component_names;
417
421 const std::vector<ComponentMask> unique_component_masks;
422
426 std::map<std::string, std::pair<std::function<double()>, bool>>
428
432 std::vector<std::set<VectorTools::NormType>> norms_per_unique_component;
433
437 double exponent;
438
443
447 std::set<std::string> extra_columns;
448
452 std::string rate_key;
453
457 std::string rate_mode;
458
462 unsigned int precision;
463
467 std::string error_file_name;
468
474};
475
476
477
478#ifndef DOXYGEN
479// ============================================================
480// Template functions
481// ============================================================
482template <int dim, int spacedim, typename VectorType>
483void
485 const VectorType & solution1,
486 const VectorType & solution2,
487 const Function<spacedim> * weight)
488{
489 AssertThrow(solution1.size() == solution2.size(),
490 ExcDimensionMismatch(solution1.size(), solution2.size()));
491 VectorType solution(solution1);
492 solution -= solution2;
494 dh,
495 solution,
497 weight);
498}
499
500
501
502template <int dim, int spacedim, typename VectorType>
503void
506 const VectorType & solution1,
507 const VectorType & solution2,
508 const Function<spacedim> * weight)
509{
510 AssertThrow(solution1.size() == solution2.size(),
511 ExcDimensionMismatch(solution1.size(), solution2.size()));
512 VectorType solution(solution1);
513 solution -= solution2;
515 mapping,
516 dh,
517 solution,
519 weight);
520}
521
522
523
524template <int dim, int spacedim, typename VectorType>
525void
527 const VectorType & solution,
528 const Function<spacedim> &exact,
529 const Function<spacedim> *weight)
530{
532 dh,
533 solution,
534 exact,
535 weight);
536}
537
538
539
540template <int dim, int spacedim, typename VectorType>
541void
544 const VectorType & solution,
545 const Function<spacedim> &exact,
546 const Function<spacedim> *weight)
547{
548 const auto n_components = component_names.size();
549
550 if (compute_error)
551 {
552 AssertDimension(exact.n_components, n_components);
553 AssertDimension(dh.get_fe().n_components(), n_components);
554
557 const unsigned int n_dofs = dh.n_dofs();
558
559 for (const auto &col : extra_columns)
560 if (col == "cells")
561 {
562 table.add_value("cells", n_active_cells);
563 table.set_tex_caption("cells", "\\# cells");
564 table.set_tex_format("cells", "r");
565 }
566 else if (col == "dofs")
567 {
568 table.add_value("dofs", n_dofs);
569 table.set_tex_caption("dofs", "\\# dofs");
570 table.set_tex_format("dofs", "r");
571 }
572
573 // A vector of zero std::functions with n_components components
574 const std::vector<std::function<double(const Point<spacedim> &)>>
575 zero_components(n_components,
576 [](const Point<spacedim> &) { return 0.0; });
577
578 // The default weight function, with n_components components
579 std::vector<std::function<double(const Point<spacedim> &)>>
580 weight_components(n_components,
581 [](const Point<spacedim> &) { return 1.0; });
582
583 if (weight != nullptr)
584 {
585 if (weight->n_components == 1)
586 {
587 for (auto &f : weight_components)
588 f = [&](const Point<spacedim> &p) { return weight->value(p); };
589 }
590 else
591 {
592 AssertDimension(weight->n_components, n_components);
593 for (unsigned int i = 0; i < n_components; ++i)
594 weight_components[i] = [&](const Point<spacedim> &p) {
595 return weight->value(p, i);
596 };
597 }
598 }
599
600 for (unsigned int i = 0; i < norms_per_unique_component.size(); ++i)
601 {
602 std::map<VectorTools::NormType, double> errors;
603
604 const auto &norms = norms_per_unique_component[i];
605 const auto &mask = unique_component_masks[i];
606
607 // Simple case first
608 if (norms.empty())
609 continue;
610
611 auto components_expr = zero_components;
612 for (unsigned int j = 0; j < n_components; ++j)
613 if (mask[j] == true)
614 components_expr[j] = weight_components[j];
615
617 components_expr);
618
619 Vector<float> difference_per_cell(
621
622 QGauss<dim> q_gauss((dh.get_fe().degree + 1) * 2);
623
624 for (const auto &norm : norms)
625 {
626 difference_per_cell = 0;
628 dh,
629 solution,
630 exact,
631 difference_per_cell,
632 q_gauss,
633 norm,
634 &select_component,
635 exponent);
636
638 dh.get_triangulation(), difference_per_cell, norm, exponent);
639
640 std::string name = unique_component_names[i] + "_" +
642 std::string latex_name = "@f$\\| " + unique_component_names[i] +
643 " - " + unique_component_names[i] +
644 "_h \\|_{" +
645 Patterns::Tools::to_string(norm) + "}@f$";
646
647 table.add_value(name, errors[norm]);
649 table.set_scientific(name, true);
650 table.set_tex_caption(name, latex_name);
651 }
652 }
653
654 for (const auto &extra_col : extra_column_functions)
655 {
656 const double custom_error = extra_col.second.first();
657
658 std::string name = extra_col.first;
659 table.add_value(name, custom_error);
661 table.set_scientific(name, true);
662 }
663 }
664}
665
666#endif
667
669
670#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:453
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:317
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:112
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#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:285
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &Du)
Definition divergence.h:472
std::string to_string(const T &t)
Definition patterns.h:2392
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 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.)
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition tria.cc:13833