Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
error_estimator.h
1 // ---------------------------------------------------------------------
2 // @f$Id: error_estimator.h 31072 2013-10-02 17:35:40Z heister @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__error_estimator_h
18 #define __deal2__error_estimator_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/function.h>
24 #include <deal.II/base/multithread_info.h>
25 #include <deal.II/dofs/function_map.h>
26 #include <deal.II/fe/component_mask.h>
27 #include <map>
28 
30 
31 
32 template <int, int> class DoFHandler;
33 template <int, int> class Mapping;
34 template <int> class Quadrature;
35 
36 namespace hp
37 {
38  template <int, int> class DoFHandler;
39  template <int> class QCollection;
40 }
41 
42 
43 
232 template <int dim, int spacedim=dim>
234 {
235 public:
292  template <typename InputVector, class DH>
293  static void estimate (const Mapping<dim, spacedim> &mapping,
294  const DH &dof,
295  const Quadrature<dim-1> &quadrature,
296  const typename FunctionMap<spacedim>::type &neumann_bc,
297  const InputVector &solution,
298  Vector<float> &error,
299  const ComponentMask &component_mask = ComponentMask(),
300  const Function<spacedim> *coefficients = 0,
301  const unsigned int n_threads = numbers::invalid_unsigned_int,
304 
309  template <typename InputVector, class DH>
310  static void estimate (const DH &dof,
311  const Quadrature<dim-1> &quadrature,
312  const typename FunctionMap<spacedim>::type &neumann_bc,
313  const InputVector &solution,
314  Vector<float> &error,
315  const ComponentMask &component_mask = ComponentMask(),
316  const Function<spacedim> *coefficients = 0,
317  const unsigned int n_threads = numbers::invalid_unsigned_int,
320 
334  template <typename InputVector, class DH>
335  static void estimate (const Mapping<dim, spacedim> &mapping,
336  const DH &dof,
337  const Quadrature<dim-1> &quadrature,
338  const typename FunctionMap<spacedim>::type &neumann_bc,
339  const std::vector<const InputVector *> &solutions,
340  std::vector<Vector<float>*> &errors,
341  const ComponentMask &component_mask = ComponentMask(),
342  const Function<spacedim> *coefficients = 0,
343  const unsigned int n_threads = numbers::invalid_unsigned_int,
346 
351  template <typename InputVector, class DH>
352  static void estimate (const DH &dof,
353  const Quadrature<dim-1> &quadrature,
354  const typename FunctionMap<spacedim>::type &neumann_bc,
355  const std::vector<const InputVector *> &solutions,
356  std::vector<Vector<float>*> &errors,
357  const ComponentMask &component_mask = ComponentMask(),
358  const Function<spacedim> *coefficients = 0,
359  const unsigned int n_threads = numbers::invalid_unsigned_int,
362 
363 
368  template <typename InputVector, class DH>
369  static void estimate (const Mapping<dim, spacedim> &mapping,
370  const DH &dof,
371  const hp::QCollection<dim-1> &quadrature,
372  const typename FunctionMap<spacedim>::type &neumann_bc,
373  const InputVector &solution,
374  Vector<float> &error,
375  const ComponentMask &component_mask = ComponentMask(),
376  const Function<spacedim> *coefficients = 0,
377  const unsigned int n_threads = numbers::invalid_unsigned_int,
380 
381 
386  template <typename InputVector, class DH>
387  static void estimate (const DH &dof,
388  const hp::QCollection<dim-1> &quadrature,
389  const typename FunctionMap<spacedim>::type &neumann_bc,
390  const InputVector &solution,
391  Vector<float> &error,
392  const ComponentMask &component_mask = ComponentMask(),
393  const Function<spacedim> *coefficients = 0,
394  const unsigned int n_threads = numbers::invalid_unsigned_int,
397 
398 
403  template <typename InputVector, class DH>
404  static void estimate (const Mapping<dim, spacedim> &mapping,
405  const DH &dof,
406  const hp::QCollection<dim-1> &quadrature,
407  const typename FunctionMap<spacedim>::type &neumann_bc,
408  const std::vector<const InputVector *> &solutions,
409  std::vector<Vector<float>*> &errors,
410  const ComponentMask &component_mask = ComponentMask(),
411  const Function<spacedim> *coefficients = 0,
412  const unsigned int n_threads = numbers::invalid_unsigned_int,
415 
416 
421  template <typename InputVector, class DH>
422  static void estimate (const DH &dof,
423  const hp::QCollection<dim-1> &quadrature,
424  const typename FunctionMap<spacedim>::type &neumann_bc,
425  const std::vector<const InputVector *> &solutions,
426  std::vector<Vector<float>*> &errors,
427  const ComponentMask &component_mask = ComponentMask(),
428  const Function<spacedim> *coefficients = 0,
429  const unsigned int n_threads = numbers::invalid_unsigned_int,
432 
433 
437  DeclException0 (ExcInvalidBoundaryIndicator);
441  DeclException0 (ExcInvalidComponentMask);
445  DeclException0 (ExcInvalidCoefficient);
449  DeclException0 (ExcInvalidBoundaryFunction);
453  DeclException2 (ExcIncompatibleNumberOfElements,
454  int, int,
455  << "The number of elements " << arg1 << " and " << arg2
456  << " of the vectors do not match!");
460  DeclException0 (ExcInvalidSolutionVector);
464  DeclException0 (ExcNoSolutions);
465 };
466 
467 
468 
480 template <int spacedim>
481 class KellyErrorEstimator<1,spacedim>
482 {
483 public:
506  template <typename InputVector, class DH>
507  static void estimate (const Mapping<1,spacedim> &mapping,
508  const DH &dof,
509  const Quadrature<0> &quadrature,
510  const typename FunctionMap<spacedim>::type &neumann_bc,
511  const InputVector &solution,
512  Vector<float> &error,
513  const ComponentMask &component_mask = ComponentMask(),
514  const Function<spacedim> *coefficients = 0,
515  const unsigned int n_threads = numbers::invalid_unsigned_int,
518 
523  template <typename InputVector, class DH>
524  static void estimate (const DH &dof,
525  const Quadrature<0> &quadrature,
526  const typename FunctionMap<spacedim>::type &neumann_bc,
527  const InputVector &solution,
528  Vector<float> &error,
529  const ComponentMask &component_mask = ComponentMask(),
530  const Function<spacedim> *coefficients = 0,
531  const unsigned int n_threads = numbers::invalid_unsigned_int,
534 
548  template <typename InputVector, class DH>
549  static void estimate (const Mapping<1,spacedim> &mapping,
550  const DH &dof,
551  const Quadrature<0> &quadrature,
552  const typename FunctionMap<spacedim>::type &neumann_bc,
553  const std::vector<const InputVector *> &solutions,
554  std::vector<Vector<float>*> &errors,
555  const ComponentMask &component_mask = ComponentMask(),
556  const Function<spacedim> *coefficients = 0,
557  const unsigned int n_threads = numbers::invalid_unsigned_int,
560 
565  template <typename InputVector, class DH>
566  static void estimate (const DH &dof,
567  const Quadrature<0> &quadrature,
568  const typename FunctionMap<spacedim>::type &neumann_bc,
569  const std::vector<const InputVector *> &solutions,
570  std::vector<Vector<float>*> &errors,
571  const ComponentMask &component_mask = ComponentMask(),
572  const Function<spacedim> *coefficients = 0,
573  const unsigned int n_threads = numbers::invalid_unsigned_int,
576 
577 
582  template <typename InputVector, class DH>
583  static void estimate (const Mapping<1,spacedim> &mapping,
584  const DH &dof,
585  const hp::QCollection<0> &quadrature,
586  const typename FunctionMap<spacedim>::type &neumann_bc,
587  const InputVector &solution,
588  Vector<float> &error,
589  const ComponentMask &component_mask = ComponentMask(),
590  const Function<spacedim> *coefficients = 0,
591  const unsigned int n_threads = numbers::invalid_unsigned_int,
594 
595 
600  template <typename InputVector, class DH>
601  static void estimate (const DH &dof,
602  const hp::QCollection<0> &quadrature,
603  const typename FunctionMap<spacedim>::type &neumann_bc,
604  const InputVector &solution,
605  Vector<float> &error,
606  const ComponentMask &component_mask = ComponentMask(),
607  const Function<spacedim> *coefficients = 0,
608  const unsigned int n_threads = numbers::invalid_unsigned_int,
611 
612 
617  template <typename InputVector, class DH>
618  static void estimate (const Mapping<1,spacedim> &mapping,
619  const DH &dof,
620  const hp::QCollection<0> &quadrature,
621  const typename FunctionMap<spacedim>::type &neumann_bc,
622  const std::vector<const InputVector *> &solutions,
623  std::vector<Vector<float>*> &errors,
624  const ComponentMask &component_mask = ComponentMask(),
625  const Function<spacedim> *coefficients = 0,
626  const unsigned int n_threads = numbers::invalid_unsigned_int,
629 
630 
635  template <typename InputVector, class DH>
636  static void estimate (const DH &dof,
637  const hp::QCollection<0> &quadrature,
638  const typename FunctionMap<spacedim>::type &neumann_bc,
639  const std::vector<const InputVector *> &solutions,
640  std::vector<Vector<float>*> &errors,
641  const ComponentMask &component_mask = ComponentMask(),
642  const Function<spacedim> *coefficients = 0,
643  const unsigned int n_threads = numbers::invalid_unsigned_int,
646 
650  DeclException0 (ExcInvalidBoundaryIndicator);
654  DeclException0 (ExcInvalidComponentMask);
658  DeclException0 (ExcInvalidCoefficient);
662  DeclException0 (ExcInvalidBoundaryFunction);
666  DeclException2 (ExcIncompatibleNumberOfElements,
667  int, int,
668  << "The number of elements " << arg1 << " and " << arg2
669  << " of the vectors do not match!");
673  DeclException0 (ExcInvalidSolutionVector);
677  DeclException0 (ExcNoSolutions);
678 };
679 
680 
681 
682 DEAL_II_NAMESPACE_CLOSE
683 
684 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:191
const types::subdomain_id invalid_subdomain_id
Definition: types.h:267
unsigned char material_id
Definition: types.h:142
std::map< types::boundary_id, const Function< dim > * > type
Definition: function_map.h:56
DeclException0(ExcInvalidBoundaryIndicator)
unsigned int subdomain_id
Definition: types.h:43
Definition: hp.h:103
DeclException2(ExcIncompatibleNumberOfElements, int, int,<< "The number of elements "<< arg1<< " and "<< arg2<< " of the vectors do not match!")
const types::material_id invalid_material_id
Definition: types.h:226
static void estimate(const Mapping< dim, spacedim > &mapping, const DH &dof, const Quadrature< dim-1 > &quadrature, const typename FunctionMap< spacedim >::type &neumann_bc, const InputVector &solution, Vector< float > &error, const ComponentMask &component_mask=ComponentMask(), const Function< spacedim > *coefficients=0, 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)