Reference documentation for deal.II version GIT relicensing-245-g36f19064f7 2024-03-29 07:20:02+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\}}\)
Loading...
Searching...
No Matches
error_estimator.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1998 - 2023 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_error_estimator_h
16#define dealii_error_estimator_h
17
18
19#include <deal.II/base/config.h>
20
23
25
27
28#include <map>
29
31
32// Forward declarations
33#ifndef DOXYGEN
34template <int dim, int spacedim>
36class DoFHandler;
37template <int, int>
38class Mapping;
39template <int>
40class Quadrature;
41
42namespace hp
43{
44 template <int>
45 class QCollection;
46}
47#endif
48
49
261template <int dim, int spacedim = dim>
263{
264public:
279
340 template <typename Number>
341 static void
343 const Mapping<dim, spacedim> &mapping,
344 const DoFHandler<dim, spacedim> &dof,
345 const Quadrature<dim - 1> &quadrature,
346 const std::map<types::boundary_id, const Function<spacedim, Number> *>
347 &neumann_bc,
348 const ReadVector<Number> &solution,
349 Vector<float> &error,
350 const ComponentMask &component_mask = {},
351 const Function<spacedim> *coefficients = nullptr,
352 const unsigned int n_threads = numbers::invalid_unsigned_int,
355 const Strategy strategy = cell_diameter_over_24);
356
361 template <typename Number>
362 static void
364 const DoFHandler<dim, spacedim> &dof,
365 const Quadrature<dim - 1> &quadrature,
366 const std::map<types::boundary_id, const Function<spacedim, Number> *>
367 &neumann_bc,
368 const ReadVector<Number> &solution,
369 Vector<float> &error,
370 const ComponentMask &component_mask = {},
371 const Function<spacedim> *coefficients = nullptr,
372 const unsigned int n_threads = numbers::invalid_unsigned_int,
375 const Strategy strategy = cell_diameter_over_24);
376
390 template <typename Number>
391 static void
393 const Mapping<dim, spacedim> &mapping,
394 const DoFHandler<dim, spacedim> &dof,
395 const Quadrature<dim - 1> &quadrature,
396 const std::map<types::boundary_id, const Function<spacedim, Number> *>
397 &neumann_bc,
398 const ArrayView<const ReadVector<Number> *> &solutions,
399 ArrayView<Vector<float> *> &errors,
400 const ComponentMask &component_mask = {},
401 const Function<spacedim> *coefficients = nullptr,
402 const unsigned int n_threads = numbers::invalid_unsigned_int,
405 const Strategy strategy = cell_diameter_over_24);
406
411 template <typename Number>
412 static void
414 const DoFHandler<dim, spacedim> &dof,
415 const Quadrature<dim - 1> &quadrature,
416 const std::map<types::boundary_id, const Function<spacedim, Number> *>
417 &neumann_bc,
418 const ArrayView<const ReadVector<Number> *> &solutions,
419 ArrayView<Vector<float> *> &errors,
420 const ComponentMask &component_mask = {},
421 const Function<spacedim> *coefficients = nullptr,
422 const unsigned int n_threads = numbers::invalid_unsigned_int,
425 const Strategy strategy = cell_diameter_over_24);
426
427
432 template <typename Number>
433 static void
435 const Mapping<dim, spacedim> &mapping,
436 const DoFHandler<dim, spacedim> &dof,
437 const hp::QCollection<dim - 1> &quadrature,
438 const std::map<types::boundary_id, const Function<spacedim, Number> *>
439 &neumann_bc,
440 const ReadVector<Number> &solution,
441 Vector<float> &error,
442 const ComponentMask &component_mask = {},
443 const Function<spacedim> *coefficients = nullptr,
444 const unsigned int n_threads = numbers::invalid_unsigned_int,
447 const Strategy strategy = cell_diameter_over_24);
448
449
454 template <typename Number>
455 static void
457 const DoFHandler<dim, spacedim> &dof,
458 const hp::QCollection<dim - 1> &quadrature,
459 const std::map<types::boundary_id, const Function<spacedim, Number> *>
460 &neumann_bc,
461 const ReadVector<Number> &solution,
462 Vector<float> &error,
463 const ComponentMask &component_mask = {},
464 const Function<spacedim> *coefficients = nullptr,
465 const unsigned int n_threads = numbers::invalid_unsigned_int,
468 const Strategy strategy = cell_diameter_over_24);
469
470
475 template <typename Number>
476 static void
478 const Mapping<dim, spacedim> &mapping,
479 const DoFHandler<dim, spacedim> &dof,
480 const hp::QCollection<dim - 1> &quadrature,
481 const std::map<types::boundary_id, const Function<spacedim, Number> *>
482 &neumann_bc,
483 const ArrayView<const ReadVector<Number> *> &solutions,
484 ArrayView<Vector<float> *> &errors,
485 const ComponentMask &component_mask = {},
486 const Function<spacedim> *coefficients = nullptr,
487 const unsigned int n_threads = numbers::invalid_unsigned_int,
490 const Strategy strategy = cell_diameter_over_24);
491
492
497 template <typename Number>
498 static void
500 const DoFHandler<dim, spacedim> &dof,
501 const hp::QCollection<dim - 1> &quadrature,
502 const std::map<types::boundary_id, const Function<spacedim, Number> *>
503 &neumann_bc,
504 const ArrayView<const ReadVector<Number> *> &solutions,
505 ArrayView<Vector<float> *> &errors,
506 const ComponentMask &component_mask = {},
507 const Function<spacedim> *coefficients = nullptr,
508 const unsigned int n_threads = numbers::invalid_unsigned_int,
511 const Strategy strategy = cell_diameter_over_24);
512
517 "You provided a ComponentMask argument that is invalid. "
518 "Component masks need to be either default constructed "
519 "(in which case they indicate that every component is "
520 "selected) or need to have a length equal to the number "
521 "of vector components of the finite element in use "
522 "by the DoFHandler object. In the latter case, at "
523 "least one component needs to be selected.");
529 "If you do specify the argument for a (possibly "
530 "spatially variable) coefficient function for this function, "
531 "then it needs to refer to a coefficient that is either "
532 "scalar (has one vector component) or has as many vector "
533 "components as there are in the finite element used by "
534 "the DoFHandler argument.");
540 int,
541 int,
542 << "You provided a function map that for boundary indicator "
543 << arg1 << " specifies a function with " << arg2
544 << " vector components. However, the finite "
545 "element in use has "
546 << arg3
547 << " components, and these two numbers need to match.");
552 int,
553 int,
554 << "The number of input vectors, " << arg1
555 << " needs to be equal to the number of output vectors, "
556 << arg2
557 << ". This is not the case in your call of this function.");
562 "You need to specify at least one solution vector as "
563 "input.");
564};
565
566
567
577template <int spacedim>
578class KellyErrorEstimator<1, spacedim>
579{
580public:
595
596
597
620 template <typename Number>
621 static void
622 estimate(
623 const Mapping<1, spacedim> &mapping,
624 const DoFHandler<1, spacedim> &dof,
625 const Quadrature<0> &quadrature,
626 const std::map<types::boundary_id, const Function<spacedim, Number> *>
627 &neumann_bc,
628 const ReadVector<Number> &solution,
629 Vector<float> &error,
630 const ComponentMask &component_mask = {},
631 const Function<spacedim> *coefficient = nullptr,
632 const unsigned int n_threads = numbers::invalid_unsigned_int,
635 const Strategy strategy = cell_diameter_over_24);
636
637
638
643 template <typename Number>
644 static void
645 estimate(
646 const DoFHandler<1, spacedim> &dof,
647 const Quadrature<0> &quadrature,
648 const std::map<types::boundary_id, const Function<spacedim, Number> *>
649 &neumann_bc,
650 const ReadVector<Number> &solution,
651 Vector<float> &error,
652 const ComponentMask &component_mask = {},
653 const Function<spacedim> *coefficients = nullptr,
654 const unsigned int n_threads = numbers::invalid_unsigned_int,
657 const Strategy strategy = cell_diameter_over_24);
658
659
660
674 template <typename Number>
675 static void
676 estimate(
677 const Mapping<1, spacedim> &mapping,
678 const DoFHandler<1, spacedim> &dof,
679 const Quadrature<0> &quadrature,
680 const std::map<types::boundary_id, const Function<spacedim, Number> *>
681 &neumann_bc,
682 const ArrayView<const ReadVector<Number> *> &solutions,
683 ArrayView<Vector<float> *> &errors,
684 const ComponentMask &component_mask = {},
685 const Function<spacedim> *coefficients = nullptr,
686 const unsigned int n_threads = numbers::invalid_unsigned_int,
689 const Strategy strategy = cell_diameter_over_24);
690
691
692
697 template <typename Number>
698 static void
699 estimate(
700 const DoFHandler<1, spacedim> &dof,
701 const Quadrature<0> &quadrature,
702 const std::map<types::boundary_id, const Function<spacedim, Number> *>
703 &neumann_bc,
704 const ArrayView<const ReadVector<Number> *> &solutions,
705 ArrayView<Vector<float> *> &errors,
706 const ComponentMask &component_mask = {},
707 const Function<spacedim> *coefficients = nullptr,
708 const unsigned int n_threads = numbers::invalid_unsigned_int,
711 const Strategy strategy = cell_diameter_over_24);
712
713
714
719 template <typename Number>
720 static void
721 estimate(
722 const Mapping<1, spacedim> &mapping,
723 const DoFHandler<1, spacedim> &dof,
724 const hp::QCollection<0> &quadrature,
725 const std::map<types::boundary_id, const Function<spacedim, Number> *>
726 &neumann_bc,
727 const ReadVector<Number> &solution,
728 Vector<float> &error,
729 const ComponentMask &component_mask = {},
730 const Function<spacedim> *coefficients = nullptr,
731 const unsigned int n_threads = numbers::invalid_unsigned_int,
734 const Strategy strategy = cell_diameter_over_24);
735
736
737
742 template <typename Number>
743 static void
744 estimate(
745 const DoFHandler<1, spacedim> &dof,
746 const hp::QCollection<0> &quadrature,
747 const std::map<types::boundary_id, const Function<spacedim, Number> *>
748 &neumann_bc,
749 const ReadVector<Number> &solution,
750 Vector<float> &error,
751 const ComponentMask &component_mask = {},
752 const Function<spacedim> *coefficients = nullptr,
753 const unsigned int n_threads = numbers::invalid_unsigned_int,
756 const Strategy strategy = cell_diameter_over_24);
757
758
759
764 template <typename Number>
765 static void
766 estimate(
767 const Mapping<1, spacedim> &mapping,
768 const DoFHandler<1, spacedim> &dof,
769 const hp::QCollection<0> &quadrature,
770 const std::map<types::boundary_id, const Function<spacedim, Number> *>
771 &neumann_bc,
772 const ArrayView<const ReadVector<Number> *> &solutions,
773 ArrayView<Vector<float> *> &errors,
774 const ComponentMask &component_mask = {},
775 const Function<spacedim> *coefficients = nullptr,
776 const unsigned int n_threads = numbers::invalid_unsigned_int,
779 const Strategy strategy = cell_diameter_over_24);
780
781
782
787 template <typename Number>
788 static void
789 estimate(
790 const DoFHandler<1, spacedim> &dof,
791 const hp::QCollection<0> &quadrature,
792 const std::map<types::boundary_id, const Function<spacedim, Number> *>
793 &neumann_bc,
794 const ArrayView<const ReadVector<Number> *> &solutions,
795 ArrayView<Vector<float> *> &errors,
796 const ComponentMask &component_mask = {},
797 const Function<spacedim> *coefficients = nullptr,
798 const unsigned int n_threads = numbers::invalid_unsigned_int,
801 const Strategy strategy = cell_diameter_over_24);
802
803
804
809 "You provided a ComponentMask argument that is invalid. "
810 "Component masks need to be either default constructed "
811 "(in which case they indicate that every component is "
812 "selected) or need to have a length equal to the number "
813 "of vector components of the finite element in use "
814 "by the DoFHandler object. In the latter case, at "
815 "least one component needs to be selected.");
816
822 "If you do specify the argument for a (possibly "
823 "spatially variable) coefficient function for this function, "
824 "then it needs to refer to a coefficient that is either "
825 "scalar (has one vector component) or has as many vector "
826 "components as there are in the finite element used by "
827 "the DoFHandler argument.");
828
834 int,
835 int,
836 << "You provided a function map that for boundary indicator "
837 << arg1 << " specifies a function with " << arg2
838 << " vector components. However, the finite "
839 "element in use has "
840 << arg3
841 << " components, and these two numbers need to match.");
842
847 int,
848 int,
849 << "The number of input vectors, " << arg1
850 << " needs to be equal to the number of output vectors, "
851 << arg2
852 << ". This is not the case in your call of this function.");
853
858 "You need to specify at least one solution vector as "
859 "input.");
860};
861
862
864
865#endif
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ArrayView< const ReadVector< Number > * > &solutions, ArrayView< Vector< float > * > &errors, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
@ cell_diameter_over_24
Kelly error estimator with the factor .
@ cell_diameter
Kelly error estimator with the factor .
static void estimate(const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
static void estimate(const DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, Number > * > &neumann_bc, const ReadVector< Number > &solution, Vector< float > &error, const ComponentMask &component_mask={}, const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
Abstract base class for mapping classes.
Definition mapping.h:316
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:177
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcInvalidBoundaryFunction(types::boundary_id arg1, int arg2, int arg3)
static ::ExceptionBase & ExcIncompatibleNumberOfElements(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:539
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:494
static ::ExceptionBase & ExcInvalidComponentMask()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:562
static ::ExceptionBase & ExcInvalidCoefficient()
static ::ExceptionBase & ExcNoSolutions()
Definition hp.h:117
const types::material_id invalid_material_id
Definition types.h:277
const types::subdomain_id invalid_subdomain_id
Definition types.h:341
static const unsigned int invalid_unsigned_int
Definition types.h:220
unsigned int subdomain_id
Definition types.h:43
unsigned int material_id
Definition types.h:167