Reference documentation for deal.II version Git 292c2606a1 2021-03-03 00:48:14 +0100
\(\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 - 2020 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  {
270  cell_diameter_over_24 = 0,
275  cell_diameter
276  };
277 
338  template <typename InputVector>
339  static void
340  estimate(
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
362  estimate(
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
392  estimate(
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
414  estimate(
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
436  estimate(
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
459  estimate(
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
481  estimate(
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
504  estimate(
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 
522  DeclExceptionMsg(ExcInvalidComponentMask,
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.");
534  ExcInvalidCoefficient,
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.");
544  DeclException3(ExcInvalidBoundaryFunction,
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.");
557  DeclException2(ExcIncompatibleNumberOfElements,
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.");
567  DeclExceptionMsg(ExcNoSolutions,
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  {
594  cell_diameter_over_24 = 0,
599  cell_diameter
600  };
601 
624  template <typename InputVector>
625  static void
626  estimate(
627  const Mapping<1, spacedim> & mapping,
628  const DoFHandler<1, spacedim> &dof,
629  const Quadrature<0> & quadrature,
630  const std::map<types::boundary_id,
632  & neumann_bc,
633  const InputVector & solution,
634  Vector<float> & error,
635  const ComponentMask & component_mask = ComponentMask(),
636  const Function<spacedim> *coefficient = nullptr,
637  const unsigned int n_threads = numbers::invalid_unsigned_int,
640  const Strategy strategy = cell_diameter_over_24);
641 
646  template <typename InputVector>
647  static void
648  estimate(
649  const DoFHandler<1, spacedim> &dof,
650  const Quadrature<0> & quadrature,
651  const std::map<types::boundary_id,
653  & neumann_bc,
654  const InputVector & solution,
655  Vector<float> & error,
656  const ComponentMask & component_mask = ComponentMask(),
657  const Function<spacedim> *coefficients = nullptr,
658  const unsigned int n_threads = numbers::invalid_unsigned_int,
661  const Strategy strategy = cell_diameter_over_24);
662 
676  template <typename InputVector>
677  static void
678  estimate(
679  const Mapping<1, spacedim> & mapping,
680  const DoFHandler<1, spacedim> &dof,
681  const Quadrature<0> & quadrature,
682  const std::map<types::boundary_id,
684  & neumann_bc,
685  const std::vector<const InputVector *> &solutions,
686  std::vector<Vector<float> *> & errors,
687  const ComponentMask & component_mask = ComponentMask(),
688  const Function<spacedim> * coefficients = nullptr,
689  const unsigned int n_threads = numbers::invalid_unsigned_int,
692  const Strategy strategy = cell_diameter_over_24);
693 
698  template <typename InputVector>
699  static void
700  estimate(
701  const DoFHandler<1, spacedim> &dof,
702  const Quadrature<0> & quadrature,
703  const std::map<types::boundary_id,
705  & neumann_bc,
706  const std::vector<const InputVector *> &solutions,
707  std::vector<Vector<float> *> & errors,
708  const ComponentMask & component_mask = ComponentMask(),
709  const Function<spacedim> * coefficients = nullptr,
710  const unsigned int n_threads = numbers::invalid_unsigned_int,
713  const Strategy strategy = cell_diameter_over_24);
714 
715 
720  template <typename InputVector>
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,
728  & neumann_bc,
729  const InputVector & solution,
730  Vector<float> & error,
731  const ComponentMask & component_mask = ComponentMask(),
732  const Function<spacedim> *coefficients = nullptr,
733  const unsigned int n_threads = numbers::invalid_unsigned_int,
736  const Strategy strategy = cell_diameter_over_24);
737 
738 
743  template <typename InputVector>
744  static void
745  estimate(
746  const DoFHandler<1, spacedim> &dof,
747  const hp::QCollection<0> & quadrature,
748  const std::map<types::boundary_id,
750  & neumann_bc,
751  const InputVector & solution,
752  Vector<float> & error,
753  const ComponentMask & component_mask = ComponentMask(),
754  const Function<spacedim> *coefficients = nullptr,
755  const unsigned int n_threads = numbers::invalid_unsigned_int,
758  const Strategy strategy = cell_diameter_over_24);
759 
760 
765  template <typename InputVector>
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,
773  & neumann_bc,
774  const std::vector<const InputVector *> &solutions,
775  std::vector<Vector<float> *> & errors,
776  const ComponentMask & component_mask = ComponentMask(),
777  const Function<spacedim> * coefficients = nullptr,
778  const unsigned int n_threads = numbers::invalid_unsigned_int,
781  const Strategy strategy = cell_diameter_over_24);
782 
783 
788  template <typename InputVector>
789  static void
790  estimate(
791  const DoFHandler<1, spacedim> &dof,
792  const hp::QCollection<0> & quadrature,
793  const std::map<types::boundary_id,
795  & neumann_bc,
796  const std::vector<const InputVector *> &solutions,
797  std::vector<Vector<float> *> & errors,
798  const ComponentMask & component_mask = ComponentMask(),
799  const Function<spacedim> * coefficients = nullptr,
800  const unsigned int n_threads = numbers::invalid_unsigned_int,
803  const Strategy strategy = cell_diameter_over_24);
804 
808  DeclExceptionMsg(ExcInvalidComponentMask,
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.");
820  ExcInvalidCoefficient,
821  "If you do specify the argument for a (possibly "
822  "spatially variable) coefficient function for this function, "
823  "then it needs to refer to a coefficient that is either "
824  "scalar (has one vector component) or has as many vector "
825  "components as there are in the finite element used by "
826  "the DoFHandler argument.");
830  DeclException3(ExcInvalidBoundaryFunction,
832  int,
833  int,
834  << "You provided a function map that for boundary indicator "
835  << arg1 << " specifies a function with " << arg2
836  << " vector components. However, the finite "
837  "element in use has "
838  << arg3
839  << " components, and these two numbers need to match.");
843  DeclException2(ExcIncompatibleNumberOfElements,
844  int,
845  int,
846  << "The number of input vectors, " << arg1
847  << " needs to be equal to the number of output vectors, "
848  << arg2
849  << ". This is not the case in your call of this function.");
853  DeclExceptionMsg(ExcNoSolutions,
854  "You need to specify at least one solution vector as "
855  "input.");
856 };
857 
858 
859 
861 
862 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:196
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:538
const types::subdomain_id invalid_subdomain_id
Definition: types.h:276
unsigned int material_id
Definition: types.h:152
unsigned int subdomain_id
Definition: types.h:43
Abstract base class for mapping classes.
Definition: mapping.h:303
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:394
Definition: hp.h:117
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:561
const types::material_id invalid_material_id
Definition: types.h:228
unsigned int boundary_id
Definition: types.h:129