Reference documentation for deal.II version GIT c415589cf0 2022-08-14 18: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 - 2021 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 
27 #include <map>
28 
30 
31 // Forward declarations
32 #ifndef DOXYGEN
33 template <int, int>
34 class DoFHandler;
35 template <int, int>
36 class Mapping;
37 template <int>
38 class Quadrature;
39 
40 namespace hp
41 {
42  template <int>
43  class QCollection;
44 }
45 #endif
46 
47 
259 template <int dim, int spacedim = dim>
261 {
262 public:
267  enum Strategy
268  {
276  };
277 
338  template <typename InputVector>
339  static void
341  const Mapping<dim, spacedim> & mapping,
342  const DoFHandler<dim, spacedim> &dof,
343  const Quadrature<dim - 1> & quadrature,
344  const std::map<types::boundary_id,
346  & neumann_bc,
347  const InputVector & solution,
348  Vector<float> & error,
349  const ComponentMask & component_mask = ComponentMask(),
350  const Function<spacedim> *coefficients = nullptr,
351  const unsigned int n_threads = numbers::invalid_unsigned_int,
354  const Strategy strategy = cell_diameter_over_24);
355 
360  template <typename InputVector>
361  static void
363  const DoFHandler<dim, spacedim> &dof,
364  const Quadrature<dim - 1> & quadrature,
365  const std::map<types::boundary_id,
367  & neumann_bc,
368  const InputVector & solution,
369  Vector<float> & error,
370  const ComponentMask & component_mask = ComponentMask(),
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 InputVector>
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,
398  & neumann_bc,
399  const std::vector<const InputVector *> &solutions,
400  std::vector<Vector<float> *> & errors,
401  const ComponentMask & component_mask = ComponentMask(),
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 InputVector>
413  static void
415  const DoFHandler<dim, spacedim> &dof,
416  const Quadrature<dim - 1> & quadrature,
417  const std::map<types::boundary_id,
419  & neumann_bc,
420  const std::vector<const InputVector *> &solutions,
421  std::vector<Vector<float> *> & errors,
422  const ComponentMask & component_mask = ComponentMask(),
423  const Function<spacedim> * coefficients = nullptr,
424  const unsigned int n_threads = numbers::invalid_unsigned_int,
427  const Strategy strategy = cell_diameter_over_24);
428 
429 
434  template <typename InputVector>
435  static void
437  const Mapping<dim, spacedim> & mapping,
438  const DoFHandler<dim, spacedim> &dof,
439  const hp::QCollection<dim - 1> & quadrature,
440  const std::map<types::boundary_id,
442  & neumann_bc,
443  const InputVector & solution,
444  Vector<float> & error,
445  const ComponentMask & component_mask = ComponentMask(),
446  const Function<spacedim> *coefficients = nullptr,
447  const unsigned int n_threads = numbers::invalid_unsigned_int,
450  const Strategy strategy = cell_diameter_over_24);
451 
452 
457  template <typename InputVector>
458  static void
460  const DoFHandler<dim, spacedim> &dof,
461  const hp::QCollection<dim - 1> & quadrature,
462  const std::map<types::boundary_id,
464  & neumann_bc,
465  const InputVector & solution,
466  Vector<float> & error,
467  const ComponentMask & component_mask = ComponentMask(),
468  const Function<spacedim> *coefficients = nullptr,
469  const unsigned int n_threads = numbers::invalid_unsigned_int,
472  const Strategy strategy = cell_diameter_over_24);
473 
474 
479  template <typename InputVector>
480  static void
482  const Mapping<dim, spacedim> & mapping,
483  const DoFHandler<dim, spacedim> &dof,
484  const hp::QCollection<dim - 1> & quadrature,
485  const std::map<types::boundary_id,
487  & neumann_bc,
488  const std::vector<const InputVector *> &solutions,
489  std::vector<Vector<float> *> & errors,
490  const ComponentMask & component_mask = ComponentMask(),
491  const Function<spacedim> * coefficients = nullptr,
492  const unsigned int n_threads = numbers::invalid_unsigned_int,
495  const Strategy strategy = cell_diameter_over_24);
496 
497 
502  template <typename InputVector>
503  static void
505  const DoFHandler<dim, spacedim> &dof,
506  const hp::QCollection<dim - 1> & quadrature,
507  const std::map<types::boundary_id,
509  & neumann_bc,
510  const std::vector<const InputVector *> &solutions,
511  std::vector<Vector<float> *> & errors,
512  const ComponentMask & component_mask = ComponentMask(),
513  const Function<spacedim> * coefficients = nullptr,
514  const unsigned int n_threads = numbers::invalid_unsigned_int,
517  const Strategy strategy = cell_diameter_over_24);
518 
523  "You provided a ComponentMask argument that is invalid. "
524  "Component masks need to be either default constructed "
525  "(in which case they indicate that every component is "
526  "selected) or need to have a length equal to the number "
527  "of vector components of the finite element in use "
528  "by the DoFHandler object. In the latter case, at "
529  "least one component needs to be selected.");
535  "If you do specify the argument for a (possibly "
536  "spatially variable) coefficient function for this function, "
537  "then it needs to refer to a coefficient that is either "
538  "scalar (has one vector component) or has as many vector "
539  "components as there are in the finite element used by "
540  "the DoFHandler argument.");
546  int,
547  int,
548  << "You provided a function map that for boundary indicator "
549  << arg1 << " specifies a function with " << arg2
550  << " vector components. However, the finite "
551  "element in use has "
552  << arg3
553  << " components, and these two numbers need to match.");
558  int,
559  int,
560  << "The number of input vectors, " << arg1
561  << " needs to be equal to the number of output vectors, "
562  << arg2
563  << ". This is not the case in your call of this function.");
568  "You need to specify at least one solution vector as "
569  "input.");
570 };
571 
572 
573 
583 template <int spacedim>
584 class KellyErrorEstimator<1, spacedim>
585 {
586 public:
591  enum Strategy
592  {
600  };
601 
602 
603 
626  template <typename InputVector>
627  static void
628  estimate(
629  const Mapping<1, spacedim> & mapping,
630  const DoFHandler<1, spacedim> &dof,
631  const Quadrature<0> & quadrature,
632  const std::map<types::boundary_id,
634  & neumann_bc,
635  const InputVector & solution,
636  Vector<float> & error,
637  const ComponentMask & component_mask = ComponentMask(),
638  const Function<spacedim> *coefficient = nullptr,
639  const unsigned int n_threads = numbers::invalid_unsigned_int,
642  const Strategy strategy = cell_diameter_over_24);
643 
644 
645 
650  template <typename InputVector>
651  static void
652  estimate(
653  const DoFHandler<1, spacedim> &dof,
654  const Quadrature<0> & quadrature,
655  const std::map<types::boundary_id,
657  & neumann_bc,
658  const InputVector & solution,
659  Vector<float> & error,
660  const ComponentMask & component_mask = ComponentMask(),
661  const Function<spacedim> *coefficients = nullptr,
662  const unsigned int n_threads = numbers::invalid_unsigned_int,
665  const Strategy strategy = cell_diameter_over_24);
666 
667 
668 
682  template <typename InputVector>
683  static void
684  estimate(
685  const Mapping<1, spacedim> & mapping,
686  const DoFHandler<1, spacedim> &dof,
687  const Quadrature<0> & quadrature,
688  const std::map<types::boundary_id,
690  & neumann_bc,
691  const std::vector<const InputVector *> &solutions,
692  std::vector<Vector<float> *> & errors,
693  const ComponentMask & component_mask = ComponentMask(),
694  const Function<spacedim> * coefficients = nullptr,
695  const unsigned int n_threads = numbers::invalid_unsigned_int,
698  const Strategy strategy = cell_diameter_over_24);
699 
700 
701 
706  template <typename InputVector>
707  static void
708  estimate(
709  const DoFHandler<1, spacedim> &dof,
710  const Quadrature<0> & quadrature,
711  const std::map<types::boundary_id,
713  & neumann_bc,
714  const std::vector<const InputVector *> &solutions,
715  std::vector<Vector<float> *> & errors,
716  const ComponentMask & component_mask = ComponentMask(),
717  const Function<spacedim> * coefficients = nullptr,
718  const unsigned int n_threads = numbers::invalid_unsigned_int,
721  const Strategy strategy = cell_diameter_over_24);
722 
723 
724 
729  template <typename InputVector>
730  static void
731  estimate(
732  const Mapping<1, spacedim> & mapping,
733  const DoFHandler<1, spacedim> &dof,
734  const hp::QCollection<0> & quadrature,
735  const std::map<types::boundary_id,
737  & neumann_bc,
738  const InputVector & solution,
739  Vector<float> & error,
740  const ComponentMask & component_mask = ComponentMask(),
741  const Function<spacedim> *coefficients = nullptr,
742  const unsigned int n_threads = numbers::invalid_unsigned_int,
745  const Strategy strategy = cell_diameter_over_24);
746 
747 
748 
753  template <typename InputVector>
754  static void
755  estimate(
756  const DoFHandler<1, spacedim> &dof,
757  const hp::QCollection<0> & quadrature,
758  const std::map<types::boundary_id,
760  & neumann_bc,
761  const InputVector & solution,
762  Vector<float> & error,
763  const ComponentMask & component_mask = ComponentMask(),
764  const Function<spacedim> *coefficients = nullptr,
765  const unsigned int n_threads = numbers::invalid_unsigned_int,
768  const Strategy strategy = cell_diameter_over_24);
769 
770 
771 
776  template <typename InputVector>
777  static void
778  estimate(
779  const Mapping<1, spacedim> & mapping,
780  const DoFHandler<1, spacedim> &dof,
781  const hp::QCollection<0> & quadrature,
782  const std::map<types::boundary_id,
784  & neumann_bc,
785  const std::vector<const InputVector *> &solutions,
786  std::vector<Vector<float> *> & errors,
787  const ComponentMask & component_mask = ComponentMask(),
788  const Function<spacedim> * coefficients = nullptr,
789  const unsigned int n_threads = numbers::invalid_unsigned_int,
792  const Strategy strategy = cell_diameter_over_24);
793 
794 
795 
800  template <typename InputVector>
801  static void
802  estimate(
803  const DoFHandler<1, spacedim> &dof,
804  const hp::QCollection<0> & quadrature,
805  const std::map<types::boundary_id,
807  & neumann_bc,
808  const std::vector<const InputVector *> &solutions,
809  std::vector<Vector<float> *> & errors,
810  const ComponentMask & component_mask = ComponentMask(),
811  const Function<spacedim> * coefficients = nullptr,
812  const unsigned int n_threads = numbers::invalid_unsigned_int,
815  const Strategy strategy = cell_diameter_over_24);
816 
817 
818 
823  "You provided a ComponentMask argument that is invalid. "
824  "Component masks need to be either default constructed "
825  "(in which case they indicate that every component is "
826  "selected) or need to have a length equal to the number "
827  "of vector components of the finite element in use "
828  "by the DoFHandler object. In the latter case, at "
829  "least one component needs to be selected.");
830 
836  "If you do specify the argument for a (possibly "
837  "spatially variable) coefficient function for this function, "
838  "then it needs to refer to a coefficient that is either "
839  "scalar (has one vector component) or has as many vector "
840  "components as there are in the finite element used by "
841  "the DoFHandler argument.");
842 
848  int,
849  int,
850  << "You provided a function map that for boundary indicator "
851  << arg1 << " specifies a function with " << arg2
852  << " vector components. However, the finite "
853  "element in use has "
854  << arg3
855  << " components, and these two numbers need to match.");
856 
861  int,
862  int,
863  << "The number of input vectors, " << arg1
864  << " needs to be equal to the number of output vectors, "
865  << arg2
866  << ". This is not the case in your call of this function.");
867 
872  "You need to specify at least one solution vector as "
873  "input.");
874 };
875 
876 
878 
879 #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:311
Definition: vector.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
static ::ExceptionBase & ExcIncompatibleNumberOfElements(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidCoefficient()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
static ::ExceptionBase & ExcNoSolutions()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcInvalidBoundaryFunction(types::boundary_id arg1, int arg2, int arg3)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:555
static ::ExceptionBase & ExcInvalidComponentMask()
Definition: hp.h:118
const types::material_id invalid_material_id
Definition: types.h:233
const types::subdomain_id invalid_subdomain_id
Definition: types.h:281
static const unsigned int invalid_unsigned_int
Definition: types.h:201
unsigned int subdomain_id
Definition: types.h:43
unsigned int material_id
Definition: types.h:152
unsigned int boundary_id
Definition: types.h:129