include/deal.II/numerics/error_estimator.h

00001 //---------------------------------------------------------------------------
00002 //    @f$Id: error_estimator.h 25505 2012-05-10 16:26:18Z bangerth @f$
00003 //
00004 //    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
00005 //
00006 //    This file is subject to QPL and may not be  distributed
00007 //    without copyright and license information. Please refer
00008 //    to the file deal.II/doc/license.html for the  text  and
00009 //    further information on this license.
00010 //
00011 //---------------------------------------------------------------------------
00012 #ifndef __deal2__error_estimator_h
00013 #define __deal2__error_estimator_h
00014 
00015 
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/exceptions.h>
00018 #include <deal.II/base/function.h>
00019 #include <deal.II/base/multithread_info.h>
00020 #include <deal.II/dofs/function_map.h>
00021 #include <map>
00022 
00023 DEAL_II_NAMESPACE_OPEN
00024 
00025 
00026 template <int, int> class DoFHandler;
00027 template <int, int> class Mapping;
00028 template <int> class Quadrature;
00029 
00030 namespace hp
00031 {
00032   template <int, int> class DoFHandler;
00033   template <int> class QCollection;
00034 }
00035 
00036 
00037 
00226 template <int dim, int spacedim=dim>
00227 class KellyErrorEstimator
00228 {
00229   public:
00339     template <typename InputVector, class DH>
00340     static void estimate (const Mapping<dim, spacedim>      &mapping,
00341                           const DH                &dof,
00342                           const Quadrature<dim-1> &quadrature,
00343                           const typename FunctionMap<spacedim>::type &neumann_bc,
00344                           const InputVector       &solution,
00345                           Vector<float>           &error,
00346                           const std::vector<bool> &component_mask = std::vector<bool>(),
00347                           const Function<spacedim>     *coefficients   = 0,
00348                           const unsigned int       n_threads = numbers::invalid_unsigned_int,
00349                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00350                           const types::material_id_t       material_id = types::invalid_material_id);
00351 
00357     template <typename InputVector, class DH>
00358     static void estimate (const DH                &dof,
00359                           const Quadrature<dim-1> &quadrature,
00360                           const typename FunctionMap<spacedim>::type &neumann_bc,
00361                           const InputVector       &solution,
00362                           Vector<float>           &error,
00363                           const std::vector<bool> &component_mask = std::vector<bool>(),
00364                           const Function<spacedim>     *coefficients   = 0,
00365                           const unsigned int       n_threads = multithread_info.n_default_threads,
00366                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00367                           const types::material_id_t       material_id = types::invalid_material_id);
00368 
00396     template <typename InputVector, class DH>
00397     static void estimate (const Mapping<dim, spacedim>          &mapping,
00398                           const DH                    &dof,
00399                           const Quadrature<dim-1>     &quadrature,
00400                           const typename FunctionMap<spacedim>::type &neumann_bc,
00401                           const std::vector<const InputVector *> &solutions,
00402                           std::vector<Vector<float>*> &errors,
00403                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00404                           const Function<spacedim>         *coefficients   = 0,
00405                           const unsigned int           n_threads = multithread_info.n_default_threads,
00406                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00407                           const types::material_id_t           material_id = types::invalid_material_id);
00408 
00414     template <typename InputVector, class DH>
00415     static void estimate (const DH                    &dof,
00416                           const Quadrature<dim-1>     &quadrature,
00417                           const typename FunctionMap<spacedim>::type &neumann_bc,
00418                           const std::vector<const InputVector *> &solutions,
00419                           std::vector<Vector<float>*> &errors,
00420                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00421                           const Function<spacedim>         *coefficients   = 0,
00422                           const unsigned int           n_threads = multithread_info.n_default_threads,
00423                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00424                           const types::material_id_t          material_id = types::invalid_material_id);
00425 
00426 
00433     template <typename InputVector, class DH>
00434     static void estimate (const Mapping<dim, spacedim>      &mapping,
00435                           const DH                &dof,
00436                           const hp::QCollection<dim-1> &quadrature,
00437                           const typename FunctionMap<spacedim>::type &neumann_bc,
00438                           const InputVector       &solution,
00439                           Vector<float>           &error,
00440                           const std::vector<bool> &component_mask = std::vector<bool>(),
00441                           const Function<spacedim>     *coefficients   = 0,
00442                           const unsigned int       n_threads = multithread_info.n_default_threads,
00443                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00444                           const types::material_id_t       material_id = types::invalid_material_id);
00445 
00446 
00453     template <typename InputVector, class DH>
00454     static void estimate (const DH                &dof,
00455                           const hp::QCollection<dim-1> &quadrature,
00456                           const typename FunctionMap<spacedim>::type &neumann_bc,
00457                           const InputVector       &solution,
00458                           Vector<float>           &error,
00459                           const std::vector<bool> &component_mask = std::vector<bool>(),
00460                           const Function<spacedim>     *coefficients   = 0,
00461                           const unsigned int       n_threads = multithread_info.n_default_threads,
00462                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00463                           const types::material_id_t       material_id = types::invalid_material_id);
00464 
00465 
00472     template <typename InputVector, class DH>
00473     static void estimate (const Mapping<dim, spacedim>          &mapping,
00474                           const DH                    &dof,
00475                           const hp::QCollection<dim-1> &quadrature,
00476                           const typename FunctionMap<spacedim>::type &neumann_bc,
00477                           const std::vector<const InputVector *> &solutions,
00478                           std::vector<Vector<float>*> &errors,
00479                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00480                           const Function<spacedim>         *coefficients   = 0,
00481                           const unsigned int           n_threads = multithread_info.n_default_threads,
00482                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00483                           const types::material_id_t           material_id = types::invalid_material_id);
00484 
00485 
00492     template <typename InputVector, class DH>
00493     static void estimate (const DH                    &dof,
00494                           const hp::QCollection<dim-1> &quadrature,
00495                           const typename FunctionMap<spacedim>::type &neumann_bc,
00496                           const std::vector<const InputVector *> &solutions,
00497                           std::vector<Vector<float>*> &errors,
00498                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00499                           const Function<spacedim>*    coefficients   = 0,
00500                           const unsigned int           n_threads = multithread_info.n_default_threads,
00501                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00502                           const types::material_id_t           material_id = types::invalid_material_id);
00503 
00504 
00508     DeclException0 (ExcInvalidBoundaryIndicator);
00512     DeclException0 (ExcInvalidComponentMask);
00516     DeclException0 (ExcInvalidCoefficient);
00520     DeclException0 (ExcInvalidBoundaryFunction);
00524     DeclException2 (ExcIncompatibleNumberOfElements,
00525                     int, int,
00526                     << "The number of elements " << arg1 << " and " << arg2
00527                     << " of the vectors do not match!");
00531     DeclException0 (ExcInvalidSolutionVector);
00535     DeclException0 (ExcNoSolutions);
00536 };
00537 
00538 
00539 
00551 template <int spacedim>
00552 class KellyErrorEstimator<1,spacedim>
00553 {
00554   public:
00602     template <typename InputVector, class DH>
00603     static void estimate (const Mapping<1,spacedim>  &mapping,
00604                           const DH   &dof,
00605                           const Quadrature<0> &quadrature,
00606                           const typename FunctionMap<spacedim>::type &neumann_bc,
00607                           const InputVector       &solution,
00608                           Vector<float>           &error,
00609                           const std::vector<bool> &component_mask = std::vector<bool>(),
00610                           const Function<spacedim>     *coefficients   = 0,
00611                           const unsigned int       n_threads = multithread_info.n_default_threads,
00612                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00613                           const types::material_id_t       material_id = types::invalid_material_id);
00614 
00620     template <typename InputVector, class DH>
00621     static void estimate (const DH   &dof,
00622                           const Quadrature<0> &quadrature,
00623                           const typename FunctionMap<spacedim>::type &neumann_bc,
00624                           const InputVector       &solution,
00625                           Vector<float>           &error,
00626                           const std::vector<bool> &component_mask = std::vector<bool>(),
00627                           const Function<spacedim>     *coefficients   = 0,
00628                           const unsigned int       n_threads = multithread_info.n_default_threads,
00629                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00630                           const types::material_id_t       material_id = types::invalid_material_id);
00631 
00659     template <typename InputVector, class DH>
00660     static void estimate (const Mapping<1,spacedim>          &mapping,
00661                           const DH       &dof,
00662                           const Quadrature<0>     &quadrature,
00663                           const typename FunctionMap<spacedim>::type &neumann_bc,
00664                           const std::vector<const InputVector *> &solutions,
00665                           std::vector<Vector<float>*> &errors,
00666                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00667                           const Function<spacedim>         *coefficients   = 0,
00668                           const unsigned int           n_threads = multithread_info.n_default_threads,
00669                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00670                           const types::material_id_t           material_id = types::invalid_material_id);
00671 
00677     template <typename InputVector, class DH>
00678     static void estimate (const DH       &dof,
00679                           const Quadrature<0>     &quadrature,
00680                           const typename FunctionMap<spacedim>::type &neumann_bc,
00681                           const std::vector<const InputVector *> &solutions,
00682                           std::vector<Vector<float>*> &errors,
00683                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00684                           const Function<spacedim>         *coefficients   = 0,
00685                           const unsigned int           n_threads = multithread_info.n_default_threads,
00686                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00687                           const types::material_id_t           material_id = types::invalid_material_id);
00688 
00689 
00696     template <typename InputVector, class DH>
00697     static void estimate (const Mapping<1,spacedim>      &mapping,
00698                           const DH                &dof,
00699                           const hp::QCollection<0> &quadrature,
00700                           const typename FunctionMap<spacedim>::type &neumann_bc,
00701                           const InputVector       &solution,
00702                           Vector<float>           &error,
00703                           const std::vector<bool> &component_mask = std::vector<bool>(),
00704                           const Function<spacedim>     *coefficients   = 0,
00705                           const unsigned int       n_threads = multithread_info.n_default_threads,
00706                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00707                           const types::material_id_t       material_id = types::invalid_material_id);
00708 
00709 
00716     template <typename InputVector, class DH>
00717     static void estimate (const DH                &dof,
00718                           const hp::QCollection<0> &quadrature,
00719                           const typename FunctionMap<spacedim>::type &neumann_bc,
00720                           const InputVector       &solution,
00721                           Vector<float>           &error,
00722                           const std::vector<bool> &component_mask = std::vector<bool>(),
00723                           const Function<spacedim>     *coefficients   = 0,
00724                           const unsigned int       n_threads = multithread_info.n_default_threads,
00725                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00726                           const types::material_id_t       material_id = types::invalid_material_id);
00727 
00728 
00735     template <typename InputVector, class DH>
00736     static void estimate (const Mapping<1,spacedim>          &mapping,
00737                           const DH                    &dof,
00738                           const hp::QCollection<0> &quadrature,
00739                           const typename FunctionMap<spacedim>::type &neumann_bc,
00740                           const std::vector<const InputVector *> &solutions,
00741                           std::vector<Vector<float>*> &errors,
00742                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00743                           const Function<spacedim>         *coefficients   = 0,
00744                           const unsigned int           n_threads = multithread_info.n_default_threads,
00745                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00746                           const types::material_id_t           material_id = types::invalid_material_id);
00747 
00748 
00755     template <typename InputVector, class DH>
00756     static void estimate (const DH                    &dof,
00757                           const hp::QCollection<0> &quadrature,
00758                           const typename FunctionMap<spacedim>::type &neumann_bc,
00759                           const std::vector<const InputVector *> &solutions,
00760                           std::vector<Vector<float>*> &errors,
00761                           const std::vector<bool>     &component_mask = std::vector<bool>(),
00762                           const Function<spacedim>         *coefficients   = 0,
00763                           const unsigned int           n_threads = multithread_info.n_default_threads,
00764                           const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id,
00765                           const types::material_id_t           material_id = types::invalid_material_id);
00766 
00770     DeclException0 (ExcInvalidBoundaryIndicator);
00774     DeclException0 (ExcInvalidComponentMask);
00778     DeclException0 (ExcInvalidCoefficient);
00782     DeclException0 (ExcInvalidBoundaryFunction);
00786     DeclException2 (ExcIncompatibleNumberOfElements,
00787                     int, int,
00788                     << "The number of elements " << arg1 << " and " << arg2
00789                     << " of the vectors do not match!");
00793     DeclException0 (ExcInvalidSolutionVector);
00797     DeclException0 (ExcNoSolutions);
00798 };
00799 
00800 
00801 
00802 DEAL_II_NAMESPACE_CLOSE
00803 
00804 #endif
00805 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Tue May 22 2012 12:06:18 by doxygen 1.7.3