]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Kelly error estimator and Solution Transfer now work in codim one
authorLuca Heltai <luca.heltai@sissa.it>
Mon, 11 Apr 2011 12:46:49 +0000 (12:46 +0000)
committerLuca Heltai <luca.heltai@sissa.it>
Mon, 11 Apr 2011 12:46:49 +0000 (12:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@23566 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/numerics/error_estimator.h
deal.II/source/numerics/error_estimator.cc
deal.II/source/numerics/error_estimator.inst.in

index f62a76a807b3bb0262a9ae97df0fe24cf608dd1b..dac4e791d5eaf86f91025d08f2b9707631238dd1 100644 (file)
@@ -29,6 +29,11 @@ inconvenience this causes.
 <h3>General</h3>
 
 <ol>
+<li> Fixed: Added some instantiations to make KellyErrorEstimator and SolutionTransfer
+work in  codimension one. Fixed some dim in spacedim.
+<br>
+(Luca Heltai, 2011/04/11)
+
 <li> Fixed: Added some instantiations to make anisotropic refinement work 
 in codimension one.
 <br>
index 9a8bd603950f698a754653a117946773eb8a8409..f3860d25ec8a2c0ab3d5c037fb944dfab1df5f30 100644 (file)
@@ -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<dim, spacedim>      &mapping,
                          const DH                &dof,
                          const Quadrature<dim-1> &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<dim>     *coefficients   = 0,
+                         const Function<DH::space_dimension>     *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 <typename InputVector, class DH>
     static void estimate (const DH                &dof,
                          const Quadrature<dim-1> &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<dim>     *coefficients   = 0,
+                         const Function<DH::space_dimension>     *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<dim, spacedim>          &mapping,
                          const DH                    &dof,
                          const Quadrature<dim-1>     &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const std::vector<const InputVector *> &solutions,
                          std::vector<Vector<float>*> &errors,
                          const std::vector<bool>     &component_mask = std::vector<bool>(),
-                         const Function<dim>         *coefficients   = 0,
+                         const Function<DH::space_dimension>         *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 <typename InputVector, class DH>
     static void estimate (const DH                    &dof,
                          const Quadrature<dim-1>     &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const std::vector<const InputVector *> &solutions,
                          std::vector<Vector<float>*> &errors,
                          const std::vector<bool>     &component_mask = std::vector<bool>(),
-                         const Function<dim>         *coefficients   = 0,
+                         const Function<DH::space_dimension>         *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<dim, spacedim>      &mapping,
                          const DH                &dof,
                          const hp::QCollection<dim-1> &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<dim>     *coefficients   = 0,
+                         const Function<DH::space_dimension>     *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 <typename InputVector, class DH>
     static void estimate (const DH                &dof,
                          const hp::QCollection<dim-1> &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<dim>     *coefficients   = 0,
+                         const Function<DH::space_dimension>     *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<dim, spacedim>          &mapping,
                          const DH                    &dof,
                          const hp::QCollection<dim-1> &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const std::vector<const InputVector *> &solutions,
                          std::vector<Vector<float>*> &errors,
                          const std::vector<bool>     &component_mask = std::vector<bool>(),
-                         const Function<dim>         *coefficients   = 0,
+                         const Function<DH::space_dimension>         *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 <typename InputVector, class DH>
     static void estimate (const DH                    &dof,
                          const hp::QCollection<dim-1> &quadrature,
-                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                          const std::vector<const InputVector *> &solutions,
                          std::vector<Vector<float>*> &errors,
                          const std::vector<bool>     &component_mask = std::vector<bool>(),
-                         const Function<dim>         *coefficients   = 0,
+                         const Function<DH::space_dimension>*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);
index 0c6f4127a3bd2cc94669a3b0f192c06dfaa8ae3b..d97bb1af5007e3c5c59ffddf84ccc271b495a3cf 100644 (file)
@@ -154,9 +154,9 @@ namespace internal
                                          * the faces of the current and
                                          * potentially of neighbor cells.
                                          */
-       dealii::hp::FEFaceValues<dim>    fe_face_values_cell;
-       dealii::hp::FEFaceValues<dim>    fe_face_values_neighbor;
-       dealii::hp::FESubfaceValues<dim> fe_subface_values;
+       dealii::hp::FEFaceValues<dim,spacedim>    fe_face_values_cell;
+       dealii::hp::FEFaceValues<dim,spacedim>    fe_face_values_neighbor;
+       dealii::hp::FESubfaceValues<dim,spacedim> fe_subface_values;
 
                                         /**
                                          * A vector to store the jump
@@ -192,18 +192,18 @@ namespace internal
                                          * index of the solution
                                          * vector.
                                          */
-       std::vector<std::vector<std::vector<Tensor<1,dim> > > > psi;
+       std::vector<std::vector<std::vector<Tensor<1,spacedim> > > > psi;
 
                                         /**
                                          * The same vector for a neighbor cell
                                          */
-       std::vector<std::vector<std::vector<Tensor<1,dim> > > > neighbor_psi;
+       std::vector<std::vector<std::vector<Tensor<1,spacedim> > > > neighbor_psi;
 
                                         /**
                                          * The normal vectors of the finite
                                          * element function on one face
                                          */
-       std::vector<Point<dim> > normal_vectors;
+       std::vector<Point<spacedim> > normal_vectors;
 
                                         /**
                                          * Two arrays needed for the
@@ -238,9 +238,9 @@ namespace internal
                                          * the KellyErrorEstimator::estimate()
                                          * function.
                                          */
-       const typename FunctionMap<dim>::type *neumann_bc;
+       const typename FunctionMap<DH::space_dimension>::type *neumann_bc;
        const std::vector<bool>                component_mask;
-       const Function<dim>                   *coefficients;
+       const Function<DH::space_dimension>                   *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<DH::dimension>::type *neumann_bc,
+                     const typename FunctionMap<DH::space_dimension>::type *neumann_bc,
                      const std::vector<bool>                component_mask,
-                     const Function<DH::dimension>                   *coefficients);
+                     const Function<DH::space_dimension>                   *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<DH::dimension>::type *neumann_bc,
+                 const typename FunctionMap<DH::space_dimension>::type *neumann_bc,
                  const std::vector<bool>                component_mask,
-                 const Function<DH::dimension>         *coefficients)
+                 const Function<DH::space_dimension>         *coefficients)
                    :
                    finite_element (fe),
                    face_quadratures (face_quadratures),
@@ -305,13 +305,13 @@ namespace internal
                         (face_quadratures.max_n_quadrature_points(),
                          std::vector<double> (fe.n_components()))),
                    psi (n_solution_vectors,
-                        std::vector<std::vector<Tensor<1,dim> > >
+                        std::vector<std::vector<Tensor<1,spacedim> > >
                         (face_quadratures.max_n_quadrature_points(),
-                         std::vector<Tensor<1,dim> > (fe.n_components()))),
+                         std::vector<Tensor<1,spacedim> > (fe.n_components()))),
                    neighbor_psi (n_solution_vectors,
-                                 std::vector<std::vector<Tensor<1,dim> > >
+                                 std::vector<std::vector<Tensor<1,spacedim> > >
                                  (face_quadratures.max_n_quadrature_points(),
-                                  std::vector<Tensor<1,dim> > (fe.n_components()))),
+                                  std::vector<Tensor<1,spacedim> > (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<typename DH::face_iterator,std::vector<double> > &local_face_integrals,
                                 const typename DH::active_cell_iterator &cell,
                                 const unsigned int                       face_no,
-                                dealii::hp::FEFaceValues<DH::dimension> &fe_face_values_cell,
-                                dealii::hp::FEFaceValues<DH::dimension> &fe_face_values_neighbor)
+                                dealii::hp::FEFaceValues<DH::dimension, DH::space_dimension> &fe_face_values_cell,
+                                dealii::hp::FEFaceValues<DH::dimension, DH::space_dimension> &fe_face_values_neighbor)
     {
       const unsigned int dim = DH::dimension;
 
@@ -625,8 +625,8 @@ namespace internal
                                   std::map<typename DH::face_iterator,std::vector<double> > &local_face_integrals,
                                   const typename DH::active_cell_iterator    &cell,
                                   const unsigned int                          face_no,
-                                  dealii::hp::FEFaceValues<DH::dimension>    &fe_face_values,
-                                  dealii::hp::FESubfaceValues<DH::dimension> &fe_subface_values)
+                                  dealii::hp::FEFaceValues<DH::dimension,DH::space_dimension>    &fe_face_values,
+                                  dealii::hp::FESubfaceValues<DH::dimension, DH::space_dimension> &fe_subface_values)
     {
       const unsigned int dim = DH::dimension;
 
@@ -1435,11 +1435,11 @@ KellyErrorEstimator<dim, spacedim>::
 estimate (const Mapping<dim, spacedim>      &mapping,
           const DH                &dof_handler,
           const Quadrature<dim-1> &quadrature,
-          const typename FunctionMap<dim>::type &neumann_bc,
+          const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<dim>     *coefficients,
+          const Function<DH::space_dimension>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1458,17 +1458,17 @@ void
 KellyErrorEstimator<dim,spacedim>::
 estimate (const DH                &dof_handler,
           const Quadrature<dim-1> &quadrature,
-          const typename FunctionMap<dim>::type &neumann_bc,
+          const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<dim>     *coefficients,
+          const Function<DH::space_dimension>     *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<dim>::mapping, dof_handler, quadrature, neumann_bc, solution,
+  estimate(StaticMappingQ1<dim,spacedim>::mapping, dof_handler, quadrature, neumann_bc, solution,
           error, component_mask, coefficients, n_threads,
            subdomain_id, material_id);
 }
@@ -1481,11 +1481,11 @@ KellyErrorEstimator<dim, spacedim>::
 estimate (const Mapping<dim, spacedim>      &mapping,
           const DH                &dof_handler,
           const hp::QCollection<dim-1> &quadrature,
-          const typename FunctionMap<dim>::type &neumann_bc,
+          const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<dim>     *coefficients,
+          const Function<DH::space_dimension>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1504,17 +1504,17 @@ void
 KellyErrorEstimator<dim, spacedim>::
 estimate (const DH                &dof_handler,
           const hp::QCollection<dim-1> &quadrature,
-          const typename FunctionMap<dim>::type &neumann_bc,
+          const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<dim>     *coefficients,
+          const Function<DH::space_dimension>     *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<dim>::mapping, dof_handler, quadrature, neumann_bc, solution,
+  estimate(StaticMappingQ1<dim, spacedim>::mapping, dof_handler, quadrature, neumann_bc, solution,
           error, component_mask, coefficients, n_threads,
            subdomain_id, material_id);
 }
@@ -1529,11 +1529,11 @@ KellyErrorEstimator<dim, spacedim>::
 estimate (const Mapping<dim, spacedim>                  &mapping,
           const DH                            &dof_handler,
           const hp::QCollection<dim-1>        &face_quadratures,
-          const typename FunctionMap<dim>::type &neumann_bc,
+          const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
           const std::vector<const InputVector *> &solutions,
           std::vector<Vector<float>*>              &errors,
           const std::vector<bool>                  &component_mask_,
-          const Function<dim>                 *coefficients,
+          const Function<DH::space_dimension>                 *coefficients,
           const unsigned int                   ,
           const types::subdomain_id_t          subdomain_id_,
           const unsigned int                   material_id)
@@ -1573,7 +1573,7 @@ estimate (const Mapping<dim, spacedim>                  &mapping,
   Assert (solutions.size() == errors.size(),
          ExcIncompatibleNumberOfElements(solutions.size(), errors.size()));
 
-  for (typename FunctionMap<dim>::type::const_iterator i=neumann_bc.begin();
+  for (typename FunctionMap<DH::space_dimension>::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<dim, spacedim>::
 estimate (const Mapping<dim, spacedim>                  &mapping,
           const DH                            &dof_handler,
           const Quadrature<dim-1>             &quadrature,
-          const typename FunctionMap<dim>::type &neumann_bc,
+          const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
           const std::vector<const InputVector *> &solutions,
           std::vector<Vector<float>*>              &errors,
           const std::vector<bool>                  &component_mask,
-          const Function<dim>                 *coefficients,
+          const Function<DH::space_dimension>                 *coefficients,
           const unsigned int                   n_threads,
           const types::subdomain_id_t          subdomain_id,
           const unsigned int                   material_id)
@@ -1728,17 +1728,17 @@ template <int dim, int spacedim>
 template <typename InputVector, class DH>
 void KellyErrorEstimator<dim, spacedim>::estimate (const DH                            &dof_handler,
                                         const Quadrature<dim-1>             &quadrature,
-                                        const typename FunctionMap<dim>::type &neumann_bc,
+                                        const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                                         const std::vector<const InputVector *> &solutions,
                                         std::vector<Vector<float>*>              &errors,
                                         const std::vector<bool>                  &component_mask,
-                                        const Function<dim>                 *coefficients,
+                                        const Function<DH::space_dimension>                 *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<dim>::mapping, dof_handler, quadrature, neumann_bc, solutions,
+  estimate(StaticMappingQ1<dim, spacedim>::mapping, dof_handler, quadrature, neumann_bc, solutions,
           errors, component_mask, coefficients, n_threads,
            subdomain_id, material_id);
 }
@@ -1749,17 +1749,17 @@ template <int dim, int spacedim>
 template <typename InputVector, class DH>
 void KellyErrorEstimator<dim, spacedim>::estimate (const DH                            &dof_handler,
                                         const hp::QCollection<dim-1>             &quadrature,
-                                        const typename FunctionMap<dim>::type &neumann_bc,
+                                        const typename FunctionMap<DH::space_dimension>::type &neumann_bc,
                                         const std::vector<const InputVector *> &solutions,
                                         std::vector<Vector<float>*>              &errors,
                                         const std::vector<bool>                  &component_mask,
-                                        const Function<dim>                 *coefficients,
+                                        const Function<DH::space_dimension>                 *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<dim>::mapping, dof_handler, quadrature, neumann_bc, solutions,
+  estimate(StaticMappingQ1<dim, spacedim>::mapping, dof_handler, quadrature, neumann_bc, solutions,
           errors, component_mask, coefficients, n_threads,
            subdomain_id, material_id);
 }
index 2400e5eb5b5b3bb88fc911c018f8e98f1421bb7e..387828adade29bf53c89efcc030d89e58c7a0235 100644 (file)
@@ -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<deal_II_dimension, deal_II_dimension>;
 #endif
 
+
+#if deal_II_dimension == 2
+template class KellyErrorEstimator<deal_II_dimension, deal_II_dimension+1>;
+#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<deal_II_dimension, deal_II_dimension>;
 template    \
 void    \
 KellyErrorEstimator<deal_II_dimension, deal_II_dimension>::    \
-estimate<InputVector,DH > (const Mapping<deal_II_dimension, deal_II_dimension>      &,    \
-          const DH   &,    \
+estimate<InputVector,DH<deal_II_dimension> > (const Mapping<deal_II_dimension, deal_II_dimension>      &,    \
+          const DH<deal_II_dimension>   &,    \
           const Q<deal_II_dimension-1> &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const InputVector       &,    \
@@ -39,7 +44,7 @@ estimate<InputVector,DH > (const Mapping<deal_II_dimension, deal_II_dimension>
 template    \
 void    \
 KellyErrorEstimator<deal_II_dimension, deal_II_dimension>::    \
-estimate<InputVector,DH > (const DH   &,    \
+estimate<InputVector,DH<deal_II_dimension> > (const DH<deal_II_dimension>   &,    \
           const Q<deal_II_dimension-1> &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const InputVector       &,    \
@@ -53,8 +58,8 @@ estimate<InputVector,DH > (const DH   &,    \
 template    \
 void    \
 KellyErrorEstimator<deal_II_dimension, deal_II_dimension>::    \
-estimate<InputVector,DH > (const Mapping<deal_II_dimension, deal_II_dimension>          &,    \
-          const DH       &,    \
+estimate<InputVector,DH<deal_II_dimension> > (const Mapping<deal_II_dimension, deal_II_dimension>          &,    \
+          const DH<deal_II_dimension>       &,    \
           const Q<deal_II_dimension-1>     &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const std::vector<const InputVector *> &,    \
@@ -68,7 +73,7 @@ estimate<InputVector,DH > (const Mapping<deal_II_dimension, deal_II_dimension>
 template    \
 void    \
 KellyErrorEstimator<deal_II_dimension, deal_II_dimension>::    \
-estimate<InputVector,DH > (const DH       &,    \
+estimate<InputVector,DH<deal_II_dimension> > (const DH<deal_II_dimension>       &,    \
           const Q<deal_II_dimension-1>     &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const std::vector<const InputVector *> &,    \
@@ -77,49 +82,119 @@ estimate<InputVector,DH > (const DH       &,    \
           const Function<deal_II_dimension>         *,    \
           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<deal_II_dimension, deal_II_dimension+1>::    \
+estimate<InputVector,DH<deal_II_dimension,deal_II_dimension+1> > (const Mapping<deal_II_dimension, deal_II_dimension+1>      &,    \
+          const DH<deal_II_dimension,deal_II_dimension+1>   &,    \
+          const Q<deal_II_dimension-1> &,    \
+          const FunctionMap<deal_II_dimension+1>::type &,    \
+          const InputVector       &,    \
+          Vector<float>           &,    \
+          const std::vector<bool> &,    \
+          const Function<deal_II_dimension+1>     *,    \
+          const unsigned int       , \
+          const unsigned int       , \
+          const unsigned int);    \
+    \
+template    \
+void    \
+KellyErrorEstimator<deal_II_dimension, deal_II_dimension+1>::    \
+estimate<InputVector,DH<deal_II_dimension,deal_II_dimension+1> > (const DH<deal_II_dimension,deal_II_dimension+1>   &,    \
+          const Q<deal_II_dimension-1> &,    \
+          const FunctionMap<deal_II_dimension+1>::type &,    \
+          const InputVector       &,    \
+          Vector<float>           &,    \
+          const std::vector<bool> &,    \
+          const Function<deal_II_dimension+1>     *,    \
+          const unsigned int       , \
+          const unsigned int       , \
+          const unsigned int);    \
+        \
+template    \
+void    \
+KellyErrorEstimator<deal_II_dimension, deal_II_dimension+1>::    \
+estimate<InputVector,DH<deal_II_dimension,deal_II_dimension+1> > (const Mapping<deal_II_dimension, deal_II_dimension+1>          &,    \
+          const DH<deal_II_dimension,deal_II_dimension+1>       &,    \
+          const Q<deal_II_dimension-1>     &,    \
+          const FunctionMap<deal_II_dimension+1>::type &,    \
+          const std::vector<const InputVector *> &,    \
+          std::vector<Vector<float>*> &,    \
+          const std::vector<bool>     &,    \
+          const Function<deal_II_dimension+1>         *,    \
+          const unsigned int           , \
+          const unsigned int           , \
+          const unsigned int);    \
+    \
+template    \
+void    \
+KellyErrorEstimator<deal_II_dimension, deal_II_dimension+1>::    \
+estimate<InputVector,DH<deal_II_dimension,deal_II_dimension+1> > (const DH<deal_II_dimension,deal_II_dimension+1>       &,    \
+          const Q<deal_II_dimension-1>     &,    \
+          const FunctionMap<deal_II_dimension+1>::type &,    \
+          const std::vector<const InputVector *> &,    \
+          std::vector<Vector<float>*> &,    \
+          const std::vector<bool>     &,    \
+          const Function<deal_II_dimension+1>         *,    \
+          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<double>,DoFHandler<deal_II_dimension>);
-INSTANTIATE(Vector<float>,DoFHandler<deal_II_dimension>);
-INSTANTIATE(BlockVector<double>,DoFHandler<deal_II_dimension>);
-INSTANTIATE(BlockVector<float>,DoFHandler<deal_II_dimension>);
+INSTANTIATE(Vector<double>,DoFHandler)
+INSTANTIATE(Vector<float>,DoFHandler)
+INSTANTIATE(BlockVector<double>,DoFHandler)
+INSTANTIATE(BlockVector<float>,DoFHandler)
 
-INSTANTIATE(Vector<double>,hp::DoFHandler<deal_II_dimension>);
-INSTANTIATE(Vector<float>,hp::DoFHandler<deal_II_dimension>);
-INSTANTIATE(BlockVector<double>,hp::DoFHandler<deal_II_dimension>);
-INSTANTIATE(BlockVector<float>,hp::DoFHandler<deal_II_dimension>);
+INSTANTIATE(Vector<double>,hp::DoFHandler)
+INSTANTIATE(Vector<float>,hp::DoFHandler)
+INSTANTIATE(BlockVector<double>,hp::DoFHandler)
+INSTANTIATE(BlockVector<float>,hp::DoFHandler)
 
 #ifdef DEAL_II_USE_PETSC
-INSTANTIATE(PETScWrappers::Vector,DoFHandler<deal_II_dimension>);
-INSTANTIATE(PETScWrappers::BlockVector,DoFHandler<deal_II_dimension>);
+INSTANTIATE(PETScWrappers::Vector,DoFHandler)
+INSTANTIATE(PETScWrappers::BlockVector,DoFHandler)
 
-INSTANTIATE(PETScWrappers::Vector,hp::DoFHandler<deal_II_dimension>);
-INSTANTIATE(PETScWrappers::BlockVector,hp::DoFHandler<deal_II_dimension>);
+INSTANTIATE(PETScWrappers::Vector,hp::DoFHandler)
+INSTANTIATE(PETScWrappers::BlockVector,hp::DoFHandler)
 
-INSTANTIATE(PETScWrappers::MPI::Vector,DoFHandler<deal_II_dimension>);
-INSTANTIATE(PETScWrappers::MPI::BlockVector,DoFHandler<deal_II_dimension>);
+INSTANTIATE(PETScWrappers::MPI::Vector,DoFHandler)
+INSTANTIATE(PETScWrappers::MPI::BlockVector,DoFHandler)
 
 #endif
 
 #ifdef DEAL_II_USE_TRILINOS
-INSTANTIATE(TrilinosWrappers::Vector,DoFHandler<deal_II_dimension>);
-INSTANTIATE(TrilinosWrappers::BlockVector,DoFHandler<deal_II_dimension>);
+INSTANTIATE(TrilinosWrappers::Vector,DoFHandler)
+INSTANTIATE(TrilinosWrappers::BlockVector,DoFHandler)
 
-INSTANTIATE(TrilinosWrappers::Vector,hp::DoFHandler<deal_II_dimension>);
-INSTANTIATE(TrilinosWrappers::BlockVector,hp::DoFHandler<deal_II_dimension>);
+INSTANTIATE(TrilinosWrappers::Vector,hp::DoFHandler)
+INSTANTIATE(TrilinosWrappers::BlockVector,hp::DoFHandler)
 
-INSTANTIATE(TrilinosWrappers::MPI::Vector,DoFHandler<deal_II_dimension>);
-INSTANTIATE(TrilinosWrappers::MPI::BlockVector,DoFHandler<deal_II_dimension>);
+INSTANTIATE(TrilinosWrappers::MPI::Vector,DoFHandler)
+INSTANTIATE(TrilinosWrappers::MPI::BlockVector,DoFHandler)
 
-INSTANTIATE(TrilinosWrappers::MPI::Vector,hp::DoFHandler<deal_II_dimension>);
-INSTANTIATE(TrilinosWrappers::MPI::BlockVector,hp::DoFHandler<deal_II_dimension>);
+INSTANTIATE(TrilinosWrappers::MPI::Vector,hp::DoFHandler)
+INSTANTIATE(TrilinosWrappers::MPI::BlockVector,hp::DoFHandler)
 #endif
 
 #undef INSTANTIATE
 #undef INSTANTIATE_1
+#undef INSTANTIATE_CODIM
 }

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.