]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fixed KellyErrorEstimator in codim one
authorLuca Heltai <luca.heltai@sissa.it>
Tue, 12 Apr 2011 17:54:47 +0000 (17:54 +0000)
committerLuca Heltai <luca.heltai@sissa.it>
Tue, 12 Apr 2011 17:54:47 +0000 (17:54 +0000)
git-svn-id: https://svn.dealii.org/trunk@23579 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/fe/fe_poly.templates.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 e397274e52b0d61279fa2a708dde7674bd969cee..ac4f01733d788e82c26107e872935aa1449241d4 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2006, 2007, 2008, 2009, 2010 by the deal.II authors
+//    Copyright (C) 2006, 2007, 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
@@ -508,68 +508,68 @@ fill_fe_face_values (const Mapping<dim,spacedim>                   &mapping,
 
 
 //codimension 1
-template <>
-inline
-void
-FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_subface_values (const Mapping<1,2> &,
-                                                              const Triangulation<1,2>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const unsigned int,
-                                                              const Quadrature<0> &,
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              FEValuesData<1,2> &) const
-{
-  AssertThrow(false, ExcNotImplemented());
-}
-
-
-template <>
-inline
-void
-FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_subface_values (const Mapping<2,3> &,
-                                                              const Triangulation<2,3>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const unsigned int,
-                                                              const Quadrature<1> &,
-                                                              Mapping<2,3>::InternalDataBase &,
-                                                              Mapping<2,3>::InternalDataBase &,
-                                                              FEValuesData<2,3> &) const
-{
-  AssertThrow(false, ExcNotImplemented());
-}
-
-
-template <>
-inline
-void
-FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_subface_values (const Mapping<1,2> &,
-                                                              const Triangulation<1,2>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const unsigned int,
-                                                              const Quadrature<0> &,
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              FEValuesData<1,2> &) const
-{
-  AssertThrow(false, ExcNotImplemented());
-}
-
-
-template <>
-inline
-void
-FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_subface_values (const Mapping<2,3> &,
-                                                              const Triangulation<2,3>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const unsigned int,
-                                                              const Quadrature<1> &,
-                                                              Mapping<2,3>::InternalDataBase &,
-                                                              Mapping<2,3>::InternalDataBase &,
-                                                              FEValuesData<2,3> &) const
-{
-  AssertThrow(false, ExcNotImplemented());
-}
+// template <>
+// inline
+// void
+// FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_subface_values (const Mapping<1,2> &,
+//                                                            const Triangulation<1,2>::cell_iterator &,
+//                                                            const unsigned int,
+//                                                            const unsigned int,
+//                                                            const Quadrature<0> &,
+//                                                            Mapping<1,2>::InternalDataBase &,
+//                                                            Mapping<1,2>::InternalDataBase &,
+//                                                            FEValuesData<1,2> &) const
+// {
+//   AssertThrow(false, ExcNotImplemented());
+// }
+
+
+// template <>
+// inline
+// void
+// FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_subface_values (const Mapping<2,3> &,
+//                                                            const Triangulation<2,3>::cell_iterator &,
+//                                                            const unsigned int,
+//                                                            const unsigned int,
+//                                                            const Quadrature<1> &,
+//                                                            Mapping<2,3>::InternalDataBase &,
+//                                                            Mapping<2,3>::InternalDataBase &,
+//                                                            FEValuesData<2,3> &) const
+// {
+//   AssertThrow(false, ExcNotImplemented());
+// }
+
+
+// template <>
+// inline
+// void
+// FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_subface_values (const Mapping<1,2> &,
+//                                                            const Triangulation<1,2>::cell_iterator &,
+//                                                            const unsigned int,
+//                                                            const unsigned int,
+//                                                            const Quadrature<0> &,
+//                                                            Mapping<1,2>::InternalDataBase &,
+//                                                            Mapping<1,2>::InternalDataBase &,
+//                                                            FEValuesData<1,2> &) const
+// {
+//   AssertThrow(false, ExcNotImplemented());
+// }
+
+
+// template <>
+// inline
+// void
+// FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_subface_values (const Mapping<2,3> &,
+//                                                            const Triangulation<2,3>::cell_iterator &,
+//                                                            const unsigned int,
+//                                                            const unsigned int,
+//                                                            const Quadrature<1> &,
+//                                                            Mapping<2,3>::InternalDataBase &,
+//                                                            Mapping<2,3>::InternalDataBase &,
+//                                                            FEValuesData<2,3> &) const
+// {
+//   AssertThrow(false, ExcNotImplemented());
+// }
 
 
 
index 9381a6bc95d936b41c81d8c596f9bb6fa8407864..af647adcf4f64a5040fa7a93b80a9a222372749d 100644 (file)
@@ -592,11 +592,11 @@ class KellyErrorEstimator<1,spacedim>
     static void estimate (const Mapping<1,spacedim>  &mapping,
                          const DH   &dof,
                          const Quadrature<0> &quadrature,
-                         const FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<1>     *coefficients   = 0,
+                         const Function<spacedim>     *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);
@@ -609,11 +609,11 @@ class KellyErrorEstimator<1,spacedim>
     template <typename InputVector, class DH>
     static void estimate (const DH   &dof,
                          const Quadrature<0> &quadrature,
-                         const FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<1>     *coefficients   = 0,
+                         const Function<spacedim>     *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);
@@ -649,11 +649,11 @@ class KellyErrorEstimator<1,spacedim>
     static void estimate (const Mapping<1,spacedim>          &mapping,
                          const DH       &dof,
                          const Quadrature<0>     &quadrature,
-                         const FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::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<1>         *coefficients   = 0,
+                         const Function<spacedim>         *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);
@@ -666,11 +666,11 @@ class KellyErrorEstimator<1,spacedim>
     template <typename InputVector, class DH>
     static void estimate (const DH       &dof,
                          const Quadrature<0>     &quadrature,
-                         const FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::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<1>         *coefficients   = 0,
+                         const Function<spacedim>         *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);
@@ -686,11 +686,11 @@ class KellyErrorEstimator<1,spacedim>
     static void estimate (const Mapping<1,spacedim>      &mapping,
                          const DH                &dof,
                          const hp::QCollection<0> &quadrature,
-                         const typename FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<1>     *coefficients   = 0,
+                         const Function<spacedim>     *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);
@@ -705,11 +705,11 @@ class KellyErrorEstimator<1,spacedim>
     template <typename InputVector, class DH>
     static void estimate (const DH                &dof,
                          const hp::QCollection<0> &quadrature,
-                         const typename FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::type &neumann_bc,
                          const InputVector       &solution,
                          Vector<float>           &error,
                          const std::vector<bool> &component_mask = std::vector<bool>(),
-                         const Function<1>     *coefficients   = 0,
+                         const Function<spacedim>     *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);
@@ -725,11 +725,11 @@ class KellyErrorEstimator<1,spacedim>
     static void estimate (const Mapping<1,spacedim>          &mapping,
                          const DH                    &dof,
                          const hp::QCollection<0> &quadrature,
-                         const typename FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::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<1>         *coefficients   = 0,
+                         const Function<spacedim>         *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);
@@ -744,11 +744,11 @@ class KellyErrorEstimator<1,spacedim>
     template <typename InputVector, class DH>
     static void estimate (const DH                    &dof,
                          const hp::QCollection<0> &quadrature,
-                         const typename FunctionMap<1>::type &neumann_bc,
+                         const typename FunctionMap<spacedim>::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<1>         *coefficients   = 0,
+                         const Function<spacedim>         *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 534ede17c6d65e7f44c035e0d495c6a36f47732b..1cf9a8304358954c1f42f3ae8c9f8797ec1cfa91 100644 (file)
@@ -1008,11 +1008,11 @@ KellyErrorEstimator<1,spacedim>::
 estimate (const Mapping<1,spacedim>      &mapping,
           const DH   &dof_handler,
           const Quadrature<0> &quadrature,
-          const FunctionMap<1>::type &neumann_bc,
+          const typename FunctionMap<spacedim>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<1>     *coefficients,
+          const Function<spacedim>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1032,11 +1032,11 @@ void
 KellyErrorEstimator<1,spacedim>::
 estimate (const DH   &dof_handler,
           const Quadrature<0> &quadrature,
-          const FunctionMap<1>::type &neumann_bc,
+          const typename FunctionMap<spacedim>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<1>     *coefficients,
+          const Function<spacedim>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1054,11 +1054,11 @@ void
 KellyErrorEstimator<1,spacedim>::
 estimate (const DH   &dof_handler,
           const Quadrature<0> &quadrature,
-          const FunctionMap<1>::type &neumann_bc,
+          const typename FunctionMap<spacedim>::type &neumann_bc,
           const std::vector<const InputVector*> &solutions,
           std::vector<Vector<float>*> &errors,
           const std::vector<bool> &component_mask,
-          const Function<1>     *coefficients,
+          const Function<spacedim>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1077,11 +1077,11 @@ KellyErrorEstimator<1,spacedim>::
 estimate (const Mapping<1,spacedim>      &mapping,
           const DH   &dof_handler,
           const hp::QCollection<0> &quadrature,
-          const FunctionMap<1>::type &neumann_bc,
+          const typename FunctionMap<spacedim>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<1>     *coefficients,
+          const Function<spacedim>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1100,11 +1100,11 @@ void
 KellyErrorEstimator<1,spacedim>::
 estimate (const DH   &dof_handler,
           const hp::QCollection<0> &quadrature,
-          const FunctionMap<1>::type &neumann_bc,
+          const typename FunctionMap<spacedim>::type &neumann_bc,
           const InputVector       &solution,
           Vector<float>           &error,
           const std::vector<bool> &component_mask,
-          const Function<1>     *coefficients,
+          const Function<spacedim>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1122,11 +1122,11 @@ void
 KellyErrorEstimator<1,spacedim>::
 estimate (const DH   &dof_handler,
           const hp::QCollection<0> &quadrature,
-          const FunctionMap<1>::type &neumann_bc,
+          const typename FunctionMap<spacedim>::type &neumann_bc,
           const std::vector<const InputVector*> &solutions,
           std::vector<Vector<float>*> &errors,
           const std::vector<bool> &component_mask,
-          const Function<1>     *coefficients,
+          const Function<spacedim>     *coefficients,
           const unsigned int       n_threads,
           const types::subdomain_id_t subdomain_id,
           const unsigned int       material_id)
@@ -1145,11 +1145,11 @@ void KellyErrorEstimator<1,spacedim>::
 estimate (const Mapping<1,spacedim>                    &/*mapping*/,
           const DH                            &/*dof_handler*/,
           const hp::QCollection<0>            &,
-          const FunctionMap<1>::type          &/*neumann_bc*/,
+          const typename FunctionMap<spacedim>::type          &/*neumann_bc*/,
           const std::vector<const InputVector *> &/*solutions*/,
           std::vector<Vector<float>*>            &/*errors*/,
           const std::vector<bool>                &/*component_mask_*/,
-          const Function<1>                   */*coefficient*/,
+          const Function<spacedim>                   */*coefficient*/,
           const unsigned int,
           const types::subdomain_id_t          /*subdomain_id*/,
           const unsigned int                   /*material_id*/)
@@ -1165,11 +1165,11 @@ void KellyErrorEstimator<1,spacedim>::
 estimate (const Mapping<1,spacedim>                    &mapping,
           const DH                 &dof_handler,
           const Quadrature<0>                 &,
-          const FunctionMap<1>::type          &neumann_bc,
+          const typename FunctionMap<spacedim>::type          &neumann_bc,
           const std::vector<const InputVector *> &solutions,
           std::vector<Vector<float>*>              &errors,
           const std::vector<bool>                  &component_mask_,
-          const Function<1>                   *coefficient,
+          const Function<spacedim>                   *coefficient,
           const unsigned int,
           const types::subdomain_id_t         subdomain_id_,
           const unsigned int                  material_id)
@@ -1208,7 +1208,7 @@ estimate (const Mapping<1,spacedim>                    &mapping,
   Assert (neumann_bc.find(255) == neumann_bc.end(),
          ExcInvalidBoundaryIndicator());
 
-  for (FunctionMap<1>::type::const_iterator i=neumann_bc.begin();
+  for (typename FunctionMap<spacedim>::type::const_iterator i=neumann_bc.begin();
        i!=neumann_bc.end(); ++i)
     Assert (i->second->n_components == n_components, ExcInvalidBoundaryFunction());
 
@@ -1245,7 +1245,7 @@ estimate (const Mapping<1,spacedim>                    &mapping,
          (coefficient->n_components == 1),
          ExcInvalidCoefficient());
 
-  for (FunctionMap<1>::type::const_iterator i=neumann_bc.begin();
+  for (typename FunctionMap<spacedim>::type::const_iterator i=neumann_bc.begin();
        i!=neumann_bc.end(); ++i)
     Assert (i->second->n_components == n_components,
             ExcInvalidBoundaryFunction());
@@ -1262,10 +1262,10 @@ estimate (const Mapping<1,spacedim>                    &mapping,
                                   // need several auxiliary fields,
                                   // depending on the way we get it
                                   // (see below)
-  std::vector<std::vector<std::vector<Tensor<1,1> > > >
+  std::vector<std::vector<std::vector<Tensor<1,spacedim> > > >
     gradients_here (n_solution_vectors,
-                   std::vector<std::vector<Tensor<1,1> > >(2, std::vector<Tensor<1,1> >(n_components)));
-  std::vector<std::vector<std::vector<Tensor<1,1> > > >
+                   std::vector<std::vector<Tensor<1,spacedim> > >(2, std::vector<Tensor<1,spacedim> >(n_components)));
+  std::vector<std::vector<std::vector<Tensor<1,spacedim> > > >
     gradients_neighbor (gradients_here);
   std::vector<Vector<double> >
     grad_neighbor (n_solution_vectors, Vector<double>(n_components));
index 387828adade29bf53c89efcc030d89e58c7a0235..da65ff9c7ad9b6da9d435e64dba33c8c86134c9d 100644 (file)
@@ -84,7 +84,7 @@ estimate<InputVector,DH<deal_II_dimension> > (const DH<deal_II_dimension>
           const unsigned int           , \
          const unsigned int);\
  
-#if deal_II_dimension == 2
+#if deal_II_dimension == 2 || deal_II_dimension == 1
 
 #define INSTANTIATE_CODIM(InputVector,DH,Q)    \
 template    \

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.