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
error_estimator.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1998 - 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_error_estimator_h
17#define dealii_error_estimator_h
18
19
20#include <deal.II/base/config.h>
21
24
26
27#include <map>
28
30
31// Forward declarations
32#ifndef DOXYGEN
33template <int dim, int spacedim>
35class DoFHandler;
36template <int, int>
37class Mapping;
38template <int>
39class Quadrature;
40
41namespace hp
42{
43 template <int>
44 class QCollection;
45}
46#endif
47
48
260template <int dim, int spacedim = dim>
262{
263public:
269 {
277 };
278
339 template <typename InputVector>
340 static void
342 const Mapping<dim, spacedim> & mapping,
343 const DoFHandler<dim, spacedim> &dof,
344 const Quadrature<dim - 1> & quadrature,
345 const std::map<types::boundary_id,
347 & neumann_bc,
348 const InputVector & solution,
349 Vector<float> & error,
350 const ComponentMask & component_mask = ComponentMask(),
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 InputVector>
362 static void
364 const DoFHandler<dim, spacedim> &dof,
365 const Quadrature<dim - 1> & quadrature,
366 const std::map<types::boundary_id,
368 & neumann_bc,
369 const InputVector & solution,
370 Vector<float> & error,
371 const ComponentMask & component_mask = ComponentMask(),
372 const Function<spacedim> *coefficients = nullptr,
373 const unsigned int n_threads = numbers::invalid_unsigned_int,
376 const Strategy strategy = cell_diameter_over_24);
377
391 template <typename InputVector>
392 static void
394 const Mapping<dim, spacedim> & mapping,
395 const DoFHandler<dim, spacedim> &dof,
396 const Quadrature<dim - 1> & quadrature,
397 const std::map<types::boundary_id,
399 & neumann_bc,
400 const std::vector<const InputVector *> &solutions,
401 std::vector<Vector<float> *> & errors,
402 const ComponentMask & component_mask = ComponentMask(),
403 const Function<spacedim> * coefficients = nullptr,
404 const unsigned int n_threads = numbers::invalid_unsigned_int,
407 const Strategy strategy = cell_diameter_over_24);
408
413 template <typename InputVector>
414 static void
416 const DoFHandler<dim, spacedim> &dof,
417 const Quadrature<dim - 1> & quadrature,
418 const std::map<types::boundary_id,
420 & neumann_bc,
421 const std::vector<const InputVector *> &solutions,
422 std::vector<Vector<float> *> & errors,
423 const ComponentMask & component_mask = ComponentMask(),
424 const Function<spacedim> * coefficients = nullptr,
425 const unsigned int n_threads = numbers::invalid_unsigned_int,
428 const Strategy strategy = cell_diameter_over_24);
429
430
435 template <typename InputVector>
436 static void
438 const Mapping<dim, spacedim> & mapping,
439 const DoFHandler<dim, spacedim> &dof,
440 const hp::QCollection<dim - 1> & quadrature,
441 const std::map<types::boundary_id,
443 & neumann_bc,
444 const InputVector & solution,
445 Vector<float> & error,
446 const ComponentMask & component_mask = ComponentMask(),
447 const Function<spacedim> *coefficients = nullptr,
448 const unsigned int n_threads = numbers::invalid_unsigned_int,
451 const Strategy strategy = cell_diameter_over_24);
452
453
458 template <typename InputVector>
459 static void
461 const DoFHandler<dim, spacedim> &dof,
462 const hp::QCollection<dim - 1> & quadrature,
463 const std::map<types::boundary_id,
465 & neumann_bc,
466 const InputVector & solution,
467 Vector<float> & error,
468 const ComponentMask & component_mask = ComponentMask(),
469 const Function<spacedim> *coefficients = nullptr,
470 const unsigned int n_threads = numbers::invalid_unsigned_int,
473 const Strategy strategy = cell_diameter_over_24);
474
475
480 template <typename InputVector>
481 static void
483 const Mapping<dim, spacedim> & mapping,
484 const DoFHandler<dim, spacedim> &dof,
485 const hp::QCollection<dim - 1> & quadrature,
486 const std::map<types::boundary_id,
488 & neumann_bc,
489 const std::vector<const InputVector *> &solutions,
490 std::vector<Vector<float> *> & errors,
491 const ComponentMask & component_mask = ComponentMask(),
492 const Function<spacedim> * coefficients = nullptr,
493 const unsigned int n_threads = numbers::invalid_unsigned_int,
496 const Strategy strategy = cell_diameter_over_24);
497
498
503 template <typename InputVector>
504 static void
506 const DoFHandler<dim, spacedim> &dof,
507 const hp::QCollection<dim - 1> & quadrature,
508 const std::map<types::boundary_id,
510 & neumann_bc,
511 const std::vector<const InputVector *> &solutions,
512 std::vector<Vector<float> *> & errors,
513 const ComponentMask & component_mask = ComponentMask(),
514 const Function<spacedim> * coefficients = nullptr,
515 const unsigned int n_threads = numbers::invalid_unsigned_int,
518 const Strategy strategy = cell_diameter_over_24);
519
524 "You provided a ComponentMask argument that is invalid. "
525 "Component masks need to be either default constructed "
526 "(in which case they indicate that every component is "
527 "selected) or need to have a length equal to the number "
528 "of vector components of the finite element in use "
529 "by the DoFHandler object. In the latter case, at "
530 "least one component needs to be selected.");
536 "If you do specify the argument for a (possibly "
537 "spatially variable) coefficient function for this function, "
538 "then it needs to refer to a coefficient that is either "
539 "scalar (has one vector component) or has as many vector "
540 "components as there are in the finite element used by "
541 "the DoFHandler argument.");
547 int,
548 int,
549 << "You provided a function map that for boundary indicator "
550 << arg1 << " specifies a function with " << arg2
551 << " vector components. However, the finite "
552 "element in use has "
553 << arg3
554 << " components, and these two numbers need to match.");
559 int,
560 int,
561 << "The number of input vectors, " << arg1
562 << " needs to be equal to the number of output vectors, "
563 << arg2
564 << ". This is not the case in your call of this function.");
569 "You need to specify at least one solution vector as "
570 "input.");
571};
572
573
574
584template <int spacedim>
585class KellyErrorEstimator<1, spacedim>
586{
587public:
593 {
601 };
602
603
604
627 template <typename InputVector>
628 static void
629 estimate(
630 const Mapping<1, spacedim> & mapping,
631 const DoFHandler<1, spacedim> &dof,
632 const Quadrature<0> & quadrature,
633 const std::map<types::boundary_id,
635 & neumann_bc,
636 const InputVector & solution,
637 Vector<float> & error,
638 const ComponentMask & component_mask = ComponentMask(),
639 const Function<spacedim> *coefficient = nullptr,
640 const unsigned int n_threads = numbers::invalid_unsigned_int,
643 const Strategy strategy = cell_diameter_over_24);
644
645
646
651 template <typename InputVector>
652 static void
653 estimate(
654 const DoFHandler<1, spacedim> &dof,
655 const Quadrature<0> & quadrature,
656 const std::map<types::boundary_id,
658 & neumann_bc,
659 const InputVector & solution,
660 Vector<float> & error,
661 const ComponentMask & component_mask = ComponentMask(),
662 const Function<spacedim> *coefficients = nullptr,
663 const unsigned int n_threads = numbers::invalid_unsigned_int,
666 const Strategy strategy = cell_diameter_over_24);
667
668
669
683 template <typename InputVector>
684 static void
685 estimate(
686 const Mapping<1, spacedim> & mapping,
687 const DoFHandler<1, spacedim> &dof,
688 const Quadrature<0> & quadrature,
689 const std::map<types::boundary_id,
691 & neumann_bc,
692 const std::vector<const InputVector *> &solutions,
693 std::vector<Vector<float> *> & errors,
694 const ComponentMask & component_mask = ComponentMask(),
695 const Function<spacedim> * coefficients = nullptr,
696 const unsigned int n_threads = numbers::invalid_unsigned_int,
699 const Strategy strategy = cell_diameter_over_24);
700
701
702
707 template <typename InputVector>
708 static void
709 estimate(
710 const DoFHandler<1, spacedim> &dof,
711 const Quadrature<0> & quadrature,
712 const std::map<types::boundary_id,
714 & neumann_bc,
715 const std::vector<const InputVector *> &solutions,
716 std::vector<Vector<float> *> & errors,
717 const ComponentMask & component_mask = ComponentMask(),
718 const Function<spacedim> * coefficients = nullptr,
719 const unsigned int n_threads = numbers::invalid_unsigned_int,
722 const Strategy strategy = cell_diameter_over_24);
723
724
725
730 template <typename InputVector>
731 static void
732 estimate(
733 const Mapping<1, spacedim> & mapping,
734 const DoFHandler<1, spacedim> &dof,
735 const hp::QCollection<0> & quadrature,
736 const std::map<types::boundary_id,
738 & neumann_bc,
739 const InputVector & solution,
740 Vector<float> & error,
741 const ComponentMask & component_mask = ComponentMask(),
742 const Function<spacedim> *coefficients = nullptr,
743 const unsigned int n_threads = numbers::invalid_unsigned_int,
746 const Strategy strategy = cell_diameter_over_24);
747
748
749
754 template <typename InputVector>
755 static void
756 estimate(
757 const DoFHandler<1, spacedim> &dof,
758 const hp::QCollection<0> & quadrature,
759 const std::map<types::boundary_id,
761 & neumann_bc,
762 const InputVector & solution,
763 Vector<float> & error,
764 const ComponentMask & component_mask = ComponentMask(),
765 const Function<spacedim> *coefficients = nullptr,
766 const unsigned int n_threads = numbers::invalid_unsigned_int,
769 const Strategy strategy = cell_diameter_over_24);
770
771
772
777 template <typename InputVector>
778 static void
779 estimate(
780 const Mapping<1, spacedim> & mapping,
781 const DoFHandler<1, spacedim> &dof,
782 const hp::QCollection<0> & quadrature,
783 const std::map<types::boundary_id,
785 & neumann_bc,
786 const std::vector<const InputVector *> &solutions,
787 std::vector<Vector<float> *> & errors,
788 const ComponentMask & component_mask = ComponentMask(),
789 const Function<spacedim> * coefficients = nullptr,
790 const unsigned int n_threads = numbers::invalid_unsigned_int,
793 const Strategy strategy = cell_diameter_over_24);
794
795
796
801 template <typename InputVector>
802 static void
803 estimate(
804 const DoFHandler<1, spacedim> &dof,
805 const hp::QCollection<0> & quadrature,
806 const std::map<types::boundary_id,
808 & neumann_bc,
809 const std::vector<const InputVector *> &solutions,
810 std::vector<Vector<float> *> & errors,
811 const ComponentMask & component_mask = ComponentMask(),
812 const Function<spacedim> * coefficients = nullptr,
813 const unsigned int n_threads = numbers::invalid_unsigned_int,
816 const Strategy strategy = cell_diameter_over_24);
817
818
819
824 "You provided a ComponentMask argument that is invalid. "
825 "Component masks need to be either default constructed "
826 "(in which case they indicate that every component is "
827 "selected) or need to have a length equal to the number "
828 "of vector components of the finite element in use "
829 "by the DoFHandler object. In the latter case, at "
830 "least one component needs to be selected.");
831
837 "If you do specify the argument for a (possibly "
838 "spatially variable) coefficient function for this function, "
839 "then it needs to refer to a coefficient that is either "
840 "scalar (has one vector component) or has as many vector "
841 "components as there are in the finite element used by "
842 "the DoFHandler argument.");
843
849 int,
850 int,
851 << "You provided a function map that for boundary indicator "
852 << arg1 << " specifies a function with " << arg2
853 << " vector components. However, the finite "
854 "element in use has "
855 << arg3
856 << " components, and these two numbers need to match.");
857
862 int,
863 int,
864 << "The number of input vectors, " << arg1
865 << " needs to be equal to the number of output vectors, "
866 << arg2
867 << ". This is not the case in your call of this function.");
868
873 "You need to specify at least one solution vector as "
874 "input.");
875};
876
877
879
880#endif
static void estimate(const DoFHandler< dim, spacedim > &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, typename InputVector::value_type > * > &neumann_bc, const std::vector< const InputVector * > &solutions, std::vector< Vector< float > * > &errors, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const std::vector< const InputVector * > &solutions, std::vector< Vector< float > * > &errors, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const InputVector &solution, Vector< float > &error, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const InputVector &solution, Vector< float > &error, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const std::vector< const InputVector * > &solutions, std::vector< Vector< float > * > &errors, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const std::vector< const InputVector * > &solutions, std::vector< Vector< float > * > &errors, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const InputVector &solution, Vector< float > &error, const ComponentMask &component_mask=ComponentMask(), 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, typename InputVector::value_type > * > &neumann_bc, const InputVector &solution, Vector< float > &error, const ComponentMask &component_mask=ComponentMask(), 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 .
Abstract base class for mapping classes.
Definition mapping.h:317
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
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:533
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
static ::ExceptionBase & ExcInvalidComponentMask()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:556
static ::ExceptionBase & ExcInvalidCoefficient()
static ::ExceptionBase & ExcNoSolutions()
Definition hp.h:118
const types::material_id invalid_material_id
Definition types.h:250
const types::subdomain_id invalid_subdomain_id
Definition types.h:298
static const unsigned int invalid_unsigned_int
Definition types.h:213