Reference documentation for deal.II version GIT 9042b9283b 2023-12-02 14:50: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\}}\)
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 
23 #include <deal.II/base/function.h>
24 
26 
28 
29 #include <map>
30 
32 
33 // Forward declarations
34 #ifndef DOXYGEN
35 template <int dim, int spacedim>
36 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
37 class DoFHandler;
38 template <int, int>
39 class Mapping;
40 template <int>
41 class Quadrature;
42 
43 namespace hp
44 {
45  template <int>
46  class QCollection;
47 }
48 #endif
49 
50 
262 template <int dim, int spacedim = dim>
264 {
265 public:
270  enum Strategy
271  {
279  };
280 
341  template <typename Number>
342  static void
344  const Mapping<dim, spacedim> &mapping,
345  const DoFHandler<dim, spacedim> &dof,
346  const Quadrature<dim - 1> &quadrature,
347  const std::map<types::boundary_id, const Function<spacedim, Number> *>
348  &neumann_bc,
349  const ReadVector<Number> &solution,
350  Vector<float> &error,
351  const ComponentMask &component_mask = {},
352  const Function<spacedim> *coefficients = nullptr,
353  const unsigned int n_threads = numbers::invalid_unsigned_int,
356  const Strategy strategy = cell_diameter_over_24);
357 
362  template <typename Number>
363  static void
365  const DoFHandler<dim, spacedim> &dof,
366  const Quadrature<dim - 1> &quadrature,
367  const std::map<types::boundary_id, const Function<spacedim, Number> *>
368  &neumann_bc,
369  const ReadVector<Number> &solution,
370  Vector<float> &error,
371  const ComponentMask &component_mask = {},
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 Number>
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, const Function<spacedim, Number> *>
398  &neumann_bc,
399  const ArrayView<const ReadVector<Number> *> &solutions,
400  ArrayView<Vector<float> *> &errors,
401  const ComponentMask &component_mask = {},
402  const Function<spacedim> *coefficients = nullptr,
403  const unsigned int n_threads = numbers::invalid_unsigned_int,
406  const Strategy strategy = cell_diameter_over_24);
407 
412  template <typename Number>
413  static void
415  const DoFHandler<dim, spacedim> &dof,
416  const Quadrature<dim - 1> &quadrature,
417  const std::map<types::boundary_id, const Function<spacedim, Number> *>
418  &neumann_bc,
419  const ArrayView<const ReadVector<Number> *> &solutions,
420  ArrayView<Vector<float> *> &errors,
421  const ComponentMask &component_mask = {},
422  const Function<spacedim> *coefficients = nullptr,
423  const unsigned int n_threads = numbers::invalid_unsigned_int,
426  const Strategy strategy = cell_diameter_over_24);
427 
428 
433  template <typename Number>
434  static void
436  const Mapping<dim, spacedim> &mapping,
437  const DoFHandler<dim, spacedim> &dof,
438  const hp::QCollection<dim - 1> &quadrature,
439  const std::map<types::boundary_id, const Function<spacedim, Number> *>
440  &neumann_bc,
441  const ReadVector<Number> &solution,
442  Vector<float> &error,
443  const ComponentMask &component_mask = {},
444  const Function<spacedim> *coefficients = nullptr,
445  const unsigned int n_threads = numbers::invalid_unsigned_int,
448  const Strategy strategy = cell_diameter_over_24);
449 
450 
455  template <typename Number>
456  static void
458  const DoFHandler<dim, spacedim> &dof,
459  const hp::QCollection<dim - 1> &quadrature,
460  const std::map<types::boundary_id, const Function<spacedim, Number> *>
461  &neumann_bc,
462  const ReadVector<Number> &solution,
463  Vector<float> &error,
464  const ComponentMask &component_mask = {},
465  const Function<spacedim> *coefficients = nullptr,
466  const unsigned int n_threads = numbers::invalid_unsigned_int,
469  const Strategy strategy = cell_diameter_over_24);
470 
471 
476  template <typename Number>
477  static void
479  const Mapping<dim, spacedim> &mapping,
480  const DoFHandler<dim, spacedim> &dof,
481  const hp::QCollection<dim - 1> &quadrature,
482  const std::map<types::boundary_id, const Function<spacedim, Number> *>
483  &neumann_bc,
484  const ArrayView<const ReadVector<Number> *> &solutions,
485  ArrayView<Vector<float> *> &errors,
486  const ComponentMask &component_mask = {},
487  const Function<spacedim> *coefficients = nullptr,
488  const unsigned int n_threads = numbers::invalid_unsigned_int,
491  const Strategy strategy = cell_diameter_over_24);
492 
493 
498  template <typename Number>
499  static void
501  const DoFHandler<dim, spacedim> &dof,
502  const hp::QCollection<dim - 1> &quadrature,
503  const std::map<types::boundary_id, const Function<spacedim, Number> *>
504  &neumann_bc,
505  const ArrayView<const ReadVector<Number> *> &solutions,
506  ArrayView<Vector<float> *> &errors,
507  const ComponentMask &component_mask = {},
508  const Function<spacedim> *coefficients = nullptr,
509  const unsigned int n_threads = numbers::invalid_unsigned_int,
512  const Strategy strategy = cell_diameter_over_24);
513 
518  "You provided a ComponentMask argument that is invalid. "
519  "Component masks need to be either default constructed "
520  "(in which case they indicate that every component is "
521  "selected) or need to have a length equal to the number "
522  "of vector components of the finite element in use "
523  "by the DoFHandler object. In the latter case, at "
524  "least one component needs to be selected.");
530  "If you do specify the argument for a (possibly "
531  "spatially variable) coefficient function for this function, "
532  "then it needs to refer to a coefficient that is either "
533  "scalar (has one vector component) or has as many vector "
534  "components as there are in the finite element used by "
535  "the DoFHandler argument.");
541  int,
542  int,
543  << "You provided a function map that for boundary indicator "
544  << arg1 << " specifies a function with " << arg2
545  << " vector components. However, the finite "
546  "element in use has "
547  << arg3
548  << " components, and these two numbers need to match.");
553  int,
554  int,
555  << "The number of input vectors, " << arg1
556  << " needs to be equal to the number of output vectors, "
557  << arg2
558  << ". This is not the case in your call of this function.");
563  "You need to specify at least one solution vector as "
564  "input.");
565 };
566 
567 
568 
578 template <int spacedim>
579 class KellyErrorEstimator<1, spacedim>
580 {
581 public:
586  enum Strategy
587  {
595  };
596 
597 
598 
621  template <typename Number>
622  static void
623  estimate(
624  const Mapping<1, spacedim> &mapping,
625  const DoFHandler<1, spacedim> &dof,
626  const Quadrature<0> &quadrature,
627  const std::map<types::boundary_id, const Function<spacedim, Number> *>
628  &neumann_bc,
629  const ReadVector<Number> &solution,
630  Vector<float> &error,
631  const ComponentMask &component_mask = {},
632  const Function<spacedim> *coefficient = nullptr,
633  const unsigned int n_threads = numbers::invalid_unsigned_int,
636  const Strategy strategy = cell_diameter_over_24);
637 
638 
639 
644  template <typename Number>
645  static void
646  estimate(
647  const DoFHandler<1, spacedim> &dof,
648  const Quadrature<0> &quadrature,
649  const std::map<types::boundary_id, const Function<spacedim, Number> *>
650  &neumann_bc,
651  const ReadVector<Number> &solution,
652  Vector<float> &error,
653  const ComponentMask &component_mask = {},
654  const Function<spacedim> *coefficients = nullptr,
655  const unsigned int n_threads = numbers::invalid_unsigned_int,
658  const Strategy strategy = cell_diameter_over_24);
659 
660 
661 
675  template <typename Number>
676  static void
677  estimate(
678  const Mapping<1, spacedim> &mapping,
679  const DoFHandler<1, spacedim> &dof,
680  const Quadrature<0> &quadrature,
681  const std::map<types::boundary_id, const Function<spacedim, Number> *>
682  &neumann_bc,
683  const ArrayView<const ReadVector<Number> *> &solutions,
684  ArrayView<Vector<float> *> &errors,
685  const ComponentMask &component_mask = {},
686  const Function<spacedim> *coefficients = nullptr,
687  const unsigned int n_threads = numbers::invalid_unsigned_int,
690  const Strategy strategy = cell_diameter_over_24);
691 
692 
693 
698  template <typename Number>
699  static void
700  estimate(
701  const DoFHandler<1, spacedim> &dof,
702  const Quadrature<0> &quadrature,
703  const std::map<types::boundary_id, const Function<spacedim, Number> *>
704  &neumann_bc,
705  const ArrayView<const ReadVector<Number> *> &solutions,
706  ArrayView<Vector<float> *> &errors,
707  const ComponentMask &component_mask = {},
708  const Function<spacedim> *coefficients = nullptr,
709  const unsigned int n_threads = numbers::invalid_unsigned_int,
712  const Strategy strategy = cell_diameter_over_24);
713 
714 
715 
720  template <typename Number>
721  static void
722  estimate(
723  const Mapping<1, spacedim> &mapping,
724  const DoFHandler<1, spacedim> &dof,
725  const hp::QCollection<0> &quadrature,
726  const std::map<types::boundary_id, const Function<spacedim, Number> *>
727  &neumann_bc,
728  const ReadVector<Number> &solution,
729  Vector<float> &error,
730  const ComponentMask &component_mask = {},
731  const Function<spacedim> *coefficients = nullptr,
732  const unsigned int n_threads = numbers::invalid_unsigned_int,
735  const Strategy strategy = cell_diameter_over_24);
736 
737 
738 
743  template <typename Number>
744  static void
745  estimate(
746  const DoFHandler<1, spacedim> &dof,
747  const hp::QCollection<0> &quadrature,
748  const std::map<types::boundary_id, const Function<spacedim, Number> *>
749  &neumann_bc,
750  const ReadVector<Number> &solution,
751  Vector<float> &error,
752  const ComponentMask &component_mask = {},
753  const Function<spacedim> *coefficients = nullptr,
754  const unsigned int n_threads = numbers::invalid_unsigned_int,
757  const Strategy strategy = cell_diameter_over_24);
758 
759 
760 
765  template <typename Number>
766  static void
767  estimate(
768  const Mapping<1, spacedim> &mapping,
769  const DoFHandler<1, spacedim> &dof,
770  const hp::QCollection<0> &quadrature,
771  const std::map<types::boundary_id, const Function<spacedim, Number> *>
772  &neumann_bc,
773  const ArrayView<const ReadVector<Number> *> &solutions,
774  ArrayView<Vector<float> *> &errors,
775  const ComponentMask &component_mask = {},
776  const Function<spacedim> *coefficients = nullptr,
777  const unsigned int n_threads = numbers::invalid_unsigned_int,
780  const Strategy strategy = cell_diameter_over_24);
781 
782 
783 
788  template <typename Number>
789  static void
790  estimate(
791  const DoFHandler<1, spacedim> &dof,
792  const hp::QCollection<0> &quadrature,
793  const std::map<types::boundary_id, const Function<spacedim, Number> *>
794  &neumann_bc,
795  const ArrayView<const ReadVector<Number> *> &solutions,
796  ArrayView<Vector<float> *> &errors,
797  const ComponentMask &component_mask = {},
798  const Function<spacedim> *coefficients = nullptr,
799  const unsigned int n_threads = numbers::invalid_unsigned_int,
802  const Strategy strategy = cell_diameter_over_24);
803 
804 
805 
810  "You provided a ComponentMask argument that is invalid. "
811  "Component masks need to be either default constructed "
812  "(in which case they indicate that every component is "
813  "selected) or need to have a length equal to the number "
814  "of vector components of the finite element in use "
815  "by the DoFHandler object. In the latter case, at "
816  "least one component needs to be selected.");
817 
823  "If you do specify the argument for a (possibly "
824  "spatially variable) coefficient function for this function, "
825  "then it needs to refer to a coefficient that is either "
826  "scalar (has one vector component) or has as many vector "
827  "components as there are in the finite element used by "
828  "the DoFHandler argument.");
829 
835  int,
836  int,
837  << "You provided a function map that for boundary indicator "
838  << arg1 << " specifies a function with " << arg2
839  << " vector components. However, the finite "
840  "element in use has "
841  << arg3
842  << " components, and these two numbers need to match.");
843 
848  int,
849  int,
850  << "The number of input vectors, " << arg1
851  << " needs to be equal to the number of output vectors, "
852  << arg2
853  << ". This is not the case in your call of this function.");
854 
859  "You need to specify at least one solution vector as "
860  "input.");
861 };
862 
863 
865 
866 #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:317
Definition: vector.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_CXX20_REQUIRES(condition)
Definition: config.h:166
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcIncompatibleNumberOfElements(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidCoefficient()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:540
static ::ExceptionBase & ExcNoSolutions()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:495
static ::ExceptionBase & ExcInvalidBoundaryFunction(types::boundary_id arg1, int arg2, int arg3)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:563
static ::ExceptionBase & ExcInvalidComponentMask()
Definition: hp.h:118
const types::material_id invalid_material_id
Definition: types.h:278
const types::subdomain_id invalid_subdomain_id
Definition: types.h:342
static const unsigned int invalid_unsigned_int
Definition: types.h:221
unsigned int subdomain_id
Definition: types.h:44
unsigned int material_id
Definition: types.h:168
unsigned int boundary_id
Definition: types.h:145