From 476b45e5bb6e235c51358980c341a1b89923491b Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Mon, 11 Apr 2011 12:46:49 +0000 Subject: [PATCH] Kelly error estimator and Solution Transfer now work in codim one git-svn-id: https://svn.dealii.org/trunk@23566 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 5 + .../deal.II/numerics/error_estimator.h | 34 ++--- deal.II/source/numerics/error_estimator.cc | 82 +++++------ .../source/numerics/error_estimator.inst.in | 139 ++++++++++++++---- 4 files changed, 170 insertions(+), 90 deletions(-) diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index f62a76a807..dac4e791d5 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -29,6 +29,11 @@ inconvenience this causes.

General

    +
  1. Fixed: Added some instantiations to make KellyErrorEstimator and SolutionTransfer +work in codimension one. Fixed some dim in spacedim. +
    +(Luca Heltai, 2011/04/11) +
  2. Fixed: Added some instantiations to make anisotropic refinement work in codimension one.
    diff --git a/deal.II/include/deal.II/numerics/error_estimator.h b/deal.II/include/deal.II/numerics/error_estimator.h index 9a8bd60395..f3860d25ec 100644 --- a/deal.II/include/deal.II/numerics/error_estimator.h +++ b/deal.II/include/deal.II/numerics/error_estimator.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2009, 2010 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2009, 2010, 2011 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -329,11 +329,11 @@ class KellyErrorEstimator static void estimate (const Mapping &mapping, const DH &dof, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = numbers::invalid_unsigned_int, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -346,11 +346,11 @@ class KellyErrorEstimator template static void estimate (const DH &dof, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -386,11 +386,11 @@ class KellyErrorEstimator static void estimate (const Mapping &mapping, const DH &dof, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -403,11 +403,11 @@ class KellyErrorEstimator template static void estimate (const DH &dof, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -423,11 +423,11 @@ class KellyErrorEstimator static void estimate (const Mapping &mapping, const DH &dof, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -442,11 +442,11 @@ class KellyErrorEstimator template static void estimate (const DH &dof, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -462,11 +462,11 @@ class KellyErrorEstimator static void estimate (const Mapping &mapping, const DH &dof, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function *coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); @@ -481,11 +481,11 @@ class KellyErrorEstimator template static void estimate (const DH &dof, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask = std::vector(), - const Function *coefficients = 0, + const Function*coefficients = 0, const unsigned int n_threads = multithread_info.n_default_threads, const types::subdomain_id_t subdomain_id = types::invalid_subdomain_id, const unsigned int material_id = numbers::invalid_unsigned_int); diff --git a/deal.II/source/numerics/error_estimator.cc b/deal.II/source/numerics/error_estimator.cc index 0c6f4127a3..d97bb1af50 100644 --- a/deal.II/source/numerics/error_estimator.cc +++ b/deal.II/source/numerics/error_estimator.cc @@ -154,9 +154,9 @@ namespace internal * the faces of the current and * potentially of neighbor cells. */ - dealii::hp::FEFaceValues fe_face_values_cell; - dealii::hp::FEFaceValues fe_face_values_neighbor; - dealii::hp::FESubfaceValues fe_subface_values; + dealii::hp::FEFaceValues fe_face_values_cell; + dealii::hp::FEFaceValues fe_face_values_neighbor; + dealii::hp::FESubfaceValues fe_subface_values; /** * A vector to store the jump @@ -192,18 +192,18 @@ namespace internal * index of the solution * vector. */ - std::vector > > > psi; + std::vector > > > psi; /** * The same vector for a neighbor cell */ - std::vector > > > neighbor_psi; + std::vector > > > neighbor_psi; /** * The normal vectors of the finite * element function on one face */ - std::vector > normal_vectors; + std::vector > normal_vectors; /** * Two arrays needed for the @@ -238,9 +238,9 @@ namespace internal * the KellyErrorEstimator::estimate() * function. */ - const typename FunctionMap::type *neumann_bc; + const typename FunctionMap::type *neumann_bc; const std::vector component_mask; - const Function *coefficients; + const Function *coefficients; /** * Constructor. @@ -253,9 +253,9 @@ namespace internal const unsigned int n_solution_vectors, const types::subdomain_id_t subdomain_id, const unsigned int material_id, - const typename FunctionMap::type *neumann_bc, + const typename FunctionMap::type *neumann_bc, const std::vector component_mask, - const Function *coefficients); + const Function *coefficients); /** * Resize the arrays so that they fit the @@ -277,9 +277,9 @@ namespace internal const unsigned int n_solution_vectors, const types::subdomain_id_t subdomain_id, const unsigned int material_id, - const typename FunctionMap::type *neumann_bc, + const typename FunctionMap::type *neumann_bc, const std::vector component_mask, - const Function *coefficients) + const Function *coefficients) : finite_element (fe), face_quadratures (face_quadratures), @@ -305,13 +305,13 @@ namespace internal (face_quadratures.max_n_quadrature_points(), std::vector (fe.n_components()))), psi (n_solution_vectors, - std::vector > > + std::vector > > (face_quadratures.max_n_quadrature_points(), - std::vector > (fe.n_components()))), + std::vector > (fe.n_components()))), neighbor_psi (n_solution_vectors, - std::vector > > + std::vector > > (face_quadratures.max_n_quadrature_points(), - std::vector > (fe.n_components()))), + std::vector > (fe.n_components()))), normal_vectors (face_quadratures.max_n_quadrature_points()), coefficient_values1 (face_quadratures.max_n_quadrature_points()), coefficient_values (face_quadratures.max_n_quadrature_points(), @@ -414,8 +414,8 @@ namespace internal std::map > &local_face_integrals, const typename DH::active_cell_iterator &cell, const unsigned int face_no, - dealii::hp::FEFaceValues &fe_face_values_cell, - dealii::hp::FEFaceValues &fe_face_values_neighbor) + dealii::hp::FEFaceValues &fe_face_values_cell, + dealii::hp::FEFaceValues &fe_face_values_neighbor) { const unsigned int dim = DH::dimension; @@ -625,8 +625,8 @@ namespace internal std::map > &local_face_integrals, const typename DH::active_cell_iterator &cell, const unsigned int face_no, - dealii::hp::FEFaceValues &fe_face_values, - dealii::hp::FESubfaceValues &fe_subface_values) + dealii::hp::FEFaceValues &fe_face_values, + dealii::hp::FESubfaceValues &fe_subface_values) { const unsigned int dim = DH::dimension; @@ -1435,11 +1435,11 @@ KellyErrorEstimator:: estimate (const Mapping &mapping, const DH &dof_handler, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) @@ -1458,17 +1458,17 @@ void KellyErrorEstimator:: estimate (const DH &dof_handler, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); - estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solution, + estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solution, error, component_mask, coefficients, n_threads, subdomain_id, material_id); } @@ -1481,11 +1481,11 @@ KellyErrorEstimator:: estimate (const Mapping &mapping, const DH &dof_handler, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) @@ -1504,17 +1504,17 @@ void KellyErrorEstimator:: estimate (const DH &dof_handler, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const InputVector &solution, Vector &error, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); - estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solution, + estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solution, error, component_mask, coefficients, n_threads, subdomain_id, material_id); } @@ -1529,11 +1529,11 @@ KellyErrorEstimator:: estimate (const Mapping &mapping, const DH &dof_handler, const hp::QCollection &face_quadratures, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask_, - const Function *coefficients, + const Function *coefficients, const unsigned int , const types::subdomain_id_t subdomain_id_, const unsigned int material_id) @@ -1573,7 +1573,7 @@ estimate (const Mapping &mapping, Assert (solutions.size() == errors.size(), ExcIncompatibleNumberOfElements(solutions.size(), errors.size())); - for (typename FunctionMap::type::const_iterator i=neumann_bc.begin(); + for (typename FunctionMap::type::const_iterator i=neumann_bc.begin(); i!=neumann_bc.end(); ++i) Assert (i->second->n_components == n_components, ExcInvalidBoundaryFunction()); @@ -1706,11 +1706,11 @@ KellyErrorEstimator:: estimate (const Mapping &mapping, const DH &dof_handler, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) @@ -1728,17 +1728,17 @@ template template void KellyErrorEstimator::estimate (const DH &dof_handler, const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); - estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solutions, + estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solutions, errors, component_mask, coefficients, n_threads, subdomain_id, material_id); } @@ -1749,17 +1749,17 @@ template template void KellyErrorEstimator::estimate (const DH &dof_handler, const hp::QCollection &quadrature, - const typename FunctionMap::type &neumann_bc, + const typename FunctionMap::type &neumann_bc, const std::vector &solutions, std::vector*> &errors, const std::vector &component_mask, - const Function *coefficients, + const Function *coefficients, const unsigned int n_threads, const types::subdomain_id_t subdomain_id, const unsigned int material_id) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); - estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solutions, + estimate(StaticMappingQ1::mapping, dof_handler, quadrature, neumann_bc, solutions, errors, component_mask, coefficients, n_threads, subdomain_id, material_id); } diff --git a/deal.II/source/numerics/error_estimator.inst.in b/deal.II/source/numerics/error_estimator.inst.in index 2400e5eb5b..387828adad 100644 --- a/deal.II/source/numerics/error_estimator.inst.in +++ b/deal.II/source/numerics/error_estimator.inst.in @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2010 by the deal.II authors +// Copyright (C) 2010, 2011 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -17,6 +17,11 @@ for (deal_II_dimension : DIMENSIONS) template class KellyErrorEstimator; #endif + +#if deal_II_dimension == 2 +template class KellyErrorEstimator; +#endif + // instantiate the externally visible functions. define a list of functions // for vectors, where the vector/matrix can be replaced using a preprocessor // variable VectorType/MatrixType @@ -24,8 +29,8 @@ template class KellyErrorEstimator; template \ void \ KellyErrorEstimator:: \ -estimate (const Mapping &, \ - const DH &, \ +estimate > (const Mapping &, \ + const DH &, \ const Q &, \ const FunctionMap::type &, \ const InputVector &, \ @@ -39,7 +44,7 @@ estimate (const Mapping template \ void \ KellyErrorEstimator:: \ -estimate (const DH &, \ +estimate > (const DH &, \ const Q &, \ const FunctionMap::type &, \ const InputVector &, \ @@ -53,8 +58,8 @@ estimate (const DH &, \ template \ void \ KellyErrorEstimator:: \ -estimate (const Mapping &, \ - const DH &, \ +estimate > (const Mapping &, \ + const DH &, \ const Q &, \ const FunctionMap::type &, \ const std::vector &, \ @@ -68,7 +73,7 @@ estimate (const Mapping template \ void \ KellyErrorEstimator:: \ -estimate (const DH &, \ +estimate > (const DH &, \ const Q &, \ const FunctionMap::type &, \ const std::vector &, \ @@ -77,49 +82,119 @@ estimate (const DH &, \ const Function *, \ const unsigned int , \ const unsigned int , \ - const unsigned int) + const unsigned int);\ + +#if deal_II_dimension == 2 + +#define INSTANTIATE_CODIM(InputVector,DH,Q) \ +template \ +void \ +KellyErrorEstimator:: \ +estimate > (const Mapping &, \ + const DH &, \ + const Q &, \ + const FunctionMap::type &, \ + const InputVector &, \ + Vector &, \ + const std::vector &, \ + const Function *, \ + const unsigned int , \ + const unsigned int , \ + const unsigned int); \ + \ +template \ +void \ +KellyErrorEstimator:: \ +estimate > (const DH &, \ + const Q &, \ + const FunctionMap::type &, \ + const InputVector &, \ + Vector &, \ + const std::vector &, \ + const Function *, \ + const unsigned int , \ + const unsigned int , \ + const unsigned int); \ + \ +template \ +void \ +KellyErrorEstimator:: \ +estimate > (const Mapping &, \ + const DH &, \ + const Q &, \ + const FunctionMap::type &, \ + const std::vector &, \ + std::vector*> &, \ + const std::vector &, \ + const Function *, \ + const unsigned int , \ + const unsigned int , \ + const unsigned int); \ + \ +template \ +void \ +KellyErrorEstimator:: \ +estimate > (const DH &, \ + const Q &, \ + const FunctionMap::type &, \ + const std::vector &, \ + std::vector*> &, \ + const std::vector &, \ + const Function *, \ + const unsigned int , \ + const unsigned int , \ + const unsigned int);\ + +#else +#define INSTANTIATE_CODIM(InputVector,DH,Q) +#endif + #define INSTANTIATE(InputVector,DH) \ - INSTANTIATE_1(InputVector,DH,Quadrature); \ - INSTANTIATE_1(InputVector,DH,hp::QCollection) + INSTANTIATE_1(InputVector,DH,Quadrature) \ + INSTANTIATE_1(InputVector,DH,hp::QCollection) \ + INSTANTIATE_CODIM(InputVector, DH, Quadrature) \ + INSTANTIATE_CODIM(InputVector, DH, hp::QCollection) + -INSTANTIATE(Vector,DoFHandler); -INSTANTIATE(Vector,DoFHandler); -INSTANTIATE(BlockVector,DoFHandler); -INSTANTIATE(BlockVector,DoFHandler); +INSTANTIATE(Vector,DoFHandler) +INSTANTIATE(Vector,DoFHandler) +INSTANTIATE(BlockVector,DoFHandler) +INSTANTIATE(BlockVector,DoFHandler) -INSTANTIATE(Vector,hp::DoFHandler); -INSTANTIATE(Vector,hp::DoFHandler); -INSTANTIATE(BlockVector,hp::DoFHandler); -INSTANTIATE(BlockVector,hp::DoFHandler); +INSTANTIATE(Vector,hp::DoFHandler) +INSTANTIATE(Vector,hp::DoFHandler) +INSTANTIATE(BlockVector,hp::DoFHandler) +INSTANTIATE(BlockVector,hp::DoFHandler) #ifdef DEAL_II_USE_PETSC -INSTANTIATE(PETScWrappers::Vector,DoFHandler); -INSTANTIATE(PETScWrappers::BlockVector,DoFHandler); +INSTANTIATE(PETScWrappers::Vector,DoFHandler) +INSTANTIATE(PETScWrappers::BlockVector,DoFHandler) -INSTANTIATE(PETScWrappers::Vector,hp::DoFHandler); -INSTANTIATE(PETScWrappers::BlockVector,hp::DoFHandler); +INSTANTIATE(PETScWrappers::Vector,hp::DoFHandler) +INSTANTIATE(PETScWrappers::BlockVector,hp::DoFHandler) -INSTANTIATE(PETScWrappers::MPI::Vector,DoFHandler); -INSTANTIATE(PETScWrappers::MPI::BlockVector,DoFHandler); +INSTANTIATE(PETScWrappers::MPI::Vector,DoFHandler) +INSTANTIATE(PETScWrappers::MPI::BlockVector,DoFHandler) #endif #ifdef DEAL_II_USE_TRILINOS -INSTANTIATE(TrilinosWrappers::Vector,DoFHandler); -INSTANTIATE(TrilinosWrappers::BlockVector,DoFHandler); +INSTANTIATE(TrilinosWrappers::Vector,DoFHandler) +INSTANTIATE(TrilinosWrappers::BlockVector,DoFHandler) -INSTANTIATE(TrilinosWrappers::Vector,hp::DoFHandler); -INSTANTIATE(TrilinosWrappers::BlockVector,hp::DoFHandler); +INSTANTIATE(TrilinosWrappers::Vector,hp::DoFHandler) +INSTANTIATE(TrilinosWrappers::BlockVector,hp::DoFHandler) -INSTANTIATE(TrilinosWrappers::MPI::Vector,DoFHandler); -INSTANTIATE(TrilinosWrappers::MPI::BlockVector,DoFHandler); +INSTANTIATE(TrilinosWrappers::MPI::Vector,DoFHandler) +INSTANTIATE(TrilinosWrappers::MPI::BlockVector,DoFHandler) -INSTANTIATE(TrilinosWrappers::MPI::Vector,hp::DoFHandler); -INSTANTIATE(TrilinosWrappers::MPI::BlockVector,hp::DoFHandler); +INSTANTIATE(TrilinosWrappers::MPI::Vector,hp::DoFHandler) +INSTANTIATE(TrilinosWrappers::MPI::BlockVector,hp::DoFHandler) #endif #undef INSTANTIATE #undef INSTANTIATE_1 +#undef INSTANTIATE_CODIM } -- 2.39.5