Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
error_estimator.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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 
22 #include <deal.II/base/exceptions.h>
23 #include <deal.II/base/function.h>
24 
25 #include <deal.II/dofs/deprecated_function_map.h>
26 
27 #include <deal.II/fe/component_mask.h>
28 
29 #include <map>
30 
31 DEAL_II_NAMESPACE_OPEN
32 
33 
34 template <int, int>
35 class Mapping;
36 template <int>
37 class Quadrature;
38 
39 namespace hp
40 {
41  template <int>
42  class QCollection;
43 }
44 
45 
46 
260 template <int dim, int spacedim = dim>
262 {
263 public:
268  enum Strategy
269  {
277  };
278 
339  template <typename InputVector, typename DoFHandlerType>
340  static void
341  estimate(
342  const Mapping<dim, spacedim> &mapping,
343  const DoFHandlerType & 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, typename DoFHandlerType>
362  static void
363  estimate(
364  const DoFHandlerType & 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, typename DoFHandlerType>
392  static void
393  estimate(
394  const Mapping<dim, spacedim> &mapping,
395  const DoFHandlerType & 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, typename DoFHandlerType>
414  static void
415  estimate(
416  const DoFHandlerType & 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, typename DoFHandlerType>
436  static void
437  estimate(
438  const Mapping<dim, spacedim> & mapping,
439  const DoFHandlerType & 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, typename DoFHandlerType>
459  static void
460  estimate(
461  const DoFHandlerType & 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, typename DoFHandlerType>
481  static void
482  estimate(
483  const Mapping<dim, spacedim> & mapping,
484  const DoFHandlerType & 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, typename DoFHandlerType>
504  static void
505  estimate(
506  const DoFHandlerType & 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  << arg2
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 
586 template <int spacedim>
587 class KellyErrorEstimator<1, spacedim>
588 {
589 public:
594  enum Strategy
595  {
603  };
604 
627  template <typename InputVector, typename DoFHandlerType>
628  static void
629  estimate(
630  const Mapping<1, spacedim> &mapping,
631  const DoFHandlerType & 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 
649  template <typename InputVector, typename DoFHandlerType>
650  static void
651  estimate(
652  const DoFHandlerType &dof,
653  const Quadrature<0> & quadrature,
654  const std::map<types::boundary_id,
656  & neumann_bc,
657  const InputVector & solution,
658  Vector<float> & error,
659  const ComponentMask & component_mask = ComponentMask(),
660  const Function<spacedim> *coefficients = nullptr,
661  const unsigned int n_threads = numbers::invalid_unsigned_int,
664  const Strategy strategy = cell_diameter_over_24);
665 
679  template <typename InputVector, typename DoFHandlerType>
680  static void
681  estimate(
682  const Mapping<1, spacedim> &mapping,
683  const DoFHandlerType & dof,
684  const Quadrature<0> & quadrature,
685  const std::map<types::boundary_id,
687  & neumann_bc,
688  const std::vector<const InputVector *> &solutions,
689  std::vector<Vector<float> *> & errors,
690  const ComponentMask & component_mask = ComponentMask(),
691  const Function<spacedim> * coefficients = nullptr,
692  const unsigned int n_threads = numbers::invalid_unsigned_int,
695  const Strategy strategy = cell_diameter_over_24);
696 
701  template <typename InputVector, typename DoFHandlerType>
702  static void
703  estimate(
704  const DoFHandlerType &dof,
705  const Quadrature<0> & quadrature,
706  const std::map<types::boundary_id,
708  & neumann_bc,
709  const std::vector<const InputVector *> &solutions,
710  std::vector<Vector<float> *> & errors,
711  const ComponentMask & component_mask = ComponentMask(),
712  const Function<spacedim> * coefficients = nullptr,
713  const unsigned int n_threads = numbers::invalid_unsigned_int,
716  const Strategy strategy = cell_diameter_over_24);
717 
718 
723  template <typename InputVector, typename DoFHandlerType>
724  static void
725  estimate(
726  const Mapping<1, spacedim> &mapping,
727  const DoFHandlerType & dof,
728  const hp::QCollection<0> & quadrature,
729  const std::map<types::boundary_id,
731  & neumann_bc,
732  const InputVector & solution,
733  Vector<float> & error,
734  const ComponentMask & component_mask = ComponentMask(),
735  const Function<spacedim> *coefficients = nullptr,
736  const unsigned int n_threads = numbers::invalid_unsigned_int,
739  const Strategy strategy = cell_diameter_over_24);
740 
741 
746  template <typename InputVector, typename DoFHandlerType>
747  static void
748  estimate(
749  const DoFHandlerType & dof,
750  const hp::QCollection<0> &quadrature,
751  const std::map<types::boundary_id,
753  & neumann_bc,
754  const InputVector & solution,
755  Vector<float> & error,
756  const ComponentMask & component_mask = ComponentMask(),
757  const Function<spacedim> *coefficients = nullptr,
758  const unsigned int n_threads = numbers::invalid_unsigned_int,
761  const Strategy strategy = cell_diameter_over_24);
762 
763 
768  template <typename InputVector, typename DoFHandlerType>
769  static void
770  estimate(
771  const Mapping<1, spacedim> &mapping,
772  const DoFHandlerType & dof,
773  const hp::QCollection<0> & quadrature,
774  const std::map<types::boundary_id,
776  & neumann_bc,
777  const std::vector<const InputVector *> &solutions,
778  std::vector<Vector<float> *> & errors,
779  const ComponentMask & component_mask = ComponentMask(),
780  const Function<spacedim> * coefficients = nullptr,
781  const unsigned int n_threads = numbers::invalid_unsigned_int,
784  const Strategy strategy = cell_diameter_over_24);
785 
786 
791  template <typename InputVector, typename DoFHandlerType>
792  static void
793  estimate(
794  const DoFHandlerType & dof,
795  const hp::QCollection<0> &quadrature,
796  const std::map<types::boundary_id,
798  & neumann_bc,
799  const std::vector<const InputVector *> &solutions,
800  std::vector<Vector<float> *> & errors,
801  const ComponentMask & component_mask = ComponentMask(),
802  const Function<spacedim> * coefficients = nullptr,
803  const unsigned int n_threads = numbers::invalid_unsigned_int,
806  const Strategy strategy = cell_diameter_over_24);
807 
812  "You provided a ComponentMask argument that is invalid. "
813  "Component masks need to be either default constructed "
814  "(in which case they indicate that every component is "
815  "selected) or need to have a length equal to the number "
816  "of vector components of the finite element in use "
817  "by the DoFHandler object. In the latter case, at "
818  "least one component needs to be selected.");
824  "If you do specify the argument for a (possibly "
825  "spatially variable) coefficient function for this function, "
826  "then it needs to refer to a coefficient that is either "
827  "scalar (has one vector component) or has as many vector "
828  "components as there are in the finite element used by "
829  "the DoFHandler argument.");
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.");
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.");
857  "You need to specify at least one solution vector as "
858  "input.");
859 };
860 
861 
862 
863 DEAL_II_NAMESPACE_CLOSE
864 
865 #endif
static ::ExceptionBase & ExcInvalidBoundaryFunction(types::boundary_id arg1, int arg2, int arg3)
static const unsigned int invalid_unsigned_int
Definition: types.h:173
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
static ::ExceptionBase & ExcInvalidCoefficient()
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandlerType &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)
const types::subdomain_id invalid_subdomain_id
Definition: types.h:258
static ::ExceptionBase & ExcInvalidComponentMask()
unsigned int material_id
Definition: types.h:134
unsigned int subdomain_id
Definition: types.h:43
Abstract base class for mapping classes.
Definition: dof_tools.h:57
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
static ::ExceptionBase & ExcIncompatibleNumberOfElements(int arg1, int arg2)
Definition: hp.h:117
Kelly error estimator with the factor .
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:564
Kelly error estimator with the factor .
const types::material_id invalid_material_id
Definition: types.h:196
unsigned int boundary_id
Definition: types.h:111
static ::ExceptionBase & ExcNoSolutions()