]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Check in stubs for hp functions. They're presently not implemented, but will be next.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 14 Dec 2006 22:00:23 +0000 (22:00 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 14 Dec 2006 22:00:23 +0000 (22:00 +0000)
git-svn-id: https://svn.dealii.org/trunk@14247 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/error_estimator.h
deal.II/deal.II/source/numerics/error_estimator.cc

index 2d17e5f29aa42dfed0fd744bb5cbb6621363fb43..421583b0030eeb143edb22edb68265f7e2f15791 100644 (file)
@@ -396,6 +396,84 @@ class KellyErrorEstimator
                           const unsigned int           subdomain_id = deal_II_numbers::invalid_unsigned_int,
                           const unsigned int           material_id = deal_II_numbers::invalid_unsigned_int);
 
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    template <typename InputVector, class DH>
+    static void estimate (const Mapping<dim>      &mapping,
+                         const DH                &dof,
+                         const hp::QCollection<dim-1> &quadrature,
+                         const typename FunctionMap<dim>::type &neumann_bc,
+                         const InputVector       &solution,
+                         Vector<float>           &error,
+                         const std::vector<bool> &component_mask = std::vector<bool>(),
+                         const Function<dim>     *coefficients   = 0,
+                         const unsigned int       n_threads = multithread_info.n_default_threads,
+                          const unsigned int       subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int       material_id = deal_II_numbers::invalid_unsigned_int);
+
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    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 InputVector       &solution,
+                         Vector<float>           &error,
+                         const std::vector<bool> &component_mask = std::vector<bool>(),
+                         const Function<dim>     *coefficients   = 0,
+                         const unsigned int       n_threads = multithread_info.n_default_threads,
+                          const unsigned int       subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int       material_id = deal_II_numbers::invalid_unsigned_int);
+
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    template <typename InputVector, class DH>
+    static void estimate (const Mapping<dim>          &mapping,
+                         const DH                    &dof,
+                         const hp::QCollection<dim-1> &quadrature,
+                         const typename FunctionMap<dim>::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 unsigned int           n_threads = multithread_info.n_default_threads,
+                          const unsigned int           subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int           material_id = deal_II_numbers::invalid_unsigned_int);
+
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    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 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 unsigned int           n_threads = multithread_info.n_default_threads,
+                          const unsigned int           subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int           material_id = deal_II_numbers::invalid_unsigned_int);
+    
     
                                     /**
                                      * Exception
@@ -841,6 +919,83 @@ class KellyErrorEstimator<1>
                           const unsigned int           subdomain_id = deal_II_numbers::invalid_unsigned_int,
                           const unsigned int           material_id = deal_II_numbers::invalid_unsigned_int);
 
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    template <typename InputVector, class DH>
+    static void estimate (const Mapping<1>      &mapping,
+                         const DH                &dof,
+                         const hp::QCollection<0> &quadrature,
+                         const typename FunctionMap<1>::type &neumann_bc,
+                         const InputVector       &solution,
+                         Vector<float>           &error,
+                         const std::vector<bool> &component_mask = std::vector<bool>(),
+                         const Function<1>     *coefficients   = 0,
+                         const unsigned int       n_threads = multithread_info.n_default_threads,
+                          const unsigned int       subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int       material_id = deal_II_numbers::invalid_unsigned_int);
+
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    template <typename InputVector, class DH>
+    static void estimate (const DH                &dof,
+                         const hp::QCollection<0> &quadrature,
+                         const typename FunctionMap<1>::type &neumann_bc,
+                         const InputVector       &solution,
+                         Vector<float>           &error,
+                         const std::vector<bool> &component_mask = std::vector<bool>(),
+                         const Function<1>     *coefficients   = 0,
+                         const unsigned int       n_threads = multithread_info.n_default_threads,
+                          const unsigned int       subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int       material_id = deal_II_numbers::invalid_unsigned_int);
+
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    template <typename InputVector, class DH>
+    static void estimate (const Mapping<1>          &mapping,
+                         const DH                    &dof,
+                         const hp::QCollection<0> &quadrature,
+                         const typename FunctionMap<1>::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 unsigned int           n_threads = multithread_info.n_default_threads,
+                          const unsigned int           subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int           material_id = deal_II_numbers::invalid_unsigned_int);
+
+
+                                    /**
+                                     * Equivalent to the set of functions
+                                     * above, except that this one takes a
+                                     * quadrature collection for hp finite
+                                     * element dof handlers.
+                                     */
+    template <typename InputVector, class DH>
+    static void estimate (const DH                    &dof,
+                         const hp::QCollection<0> &quadrature,
+                         const typename FunctionMap<1>::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 unsigned int           n_threads = multithread_info.n_default_threads,
+                          const unsigned int           subdomain_id = deal_II_numbers::invalid_unsigned_int,
+                          const unsigned int           material_id = deal_II_numbers::invalid_unsigned_int);
     
                                     /**
                                      * Exception
index 7c4d86ca39e6dd18358dfa54eae1e7b4695640b8..c224564e410f913132e2b7716d6eb774ab4dacfb 100644 (file)
@@ -134,6 +134,91 @@ estimate (const DH   &dof_handler,
 
 
 
+template <typename InputVector, class DH>
+void
+KellyErrorEstimator<1>::
+estimate (const Mapping<1>      &mapping,
+          const DH   &dof_handler,
+          const hp::QCollection<0> &quadrature,
+          const FunctionMap<1>::type &neumann_bc,
+          const InputVector       &solution,
+          Vector<float>           &error,
+          const std::vector<bool> &component_mask,
+          const Function<1>     *coefficients,
+          const unsigned int       n_threads,
+          const unsigned int       subdomain_id,
+          const unsigned int       material_id)
+{
+                                  // just pass on to the other function
+  const std::vector<const InputVector *> solutions (1, &solution);
+  std::vector<Vector<float>*>              errors (1, &error);
+  estimate (mapping, dof_handler, quadrature, neumann_bc, solutions, errors,
+           component_mask, coefficients, n_threads, subdomain_id, material_id);
+}
+
+
+template <typename InputVector, class DH>
+void
+KellyErrorEstimator<1>::
+estimate (const DH   &dof_handler,
+          const hp::QCollection<0> &quadrature,
+          const FunctionMap<1>::type &neumann_bc,
+          const InputVector       &solution,
+          Vector<float>           &error,
+          const std::vector<bool> &component_mask,
+          const Function<1>     *coefficients,
+          const unsigned int       n_threads,
+          const unsigned int       subdomain_id,
+          const unsigned int       material_id)
+{
+  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
+  estimate(StaticMappingQ1<1>::mapping, dof_handler, quadrature, neumann_bc, solution,
+          error, component_mask, coefficients, n_threads, subdomain_id, material_id);
+}
+
+
+
+template <typename InputVector, class DH>
+void
+KellyErrorEstimator<1>::
+estimate (const DH   &dof_handler,
+          const hp::QCollection<0> &quadrature,
+          const FunctionMap<1>::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 unsigned int       n_threads,
+          const unsigned int       subdomain_id,
+          const unsigned int       material_id)
+{
+  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
+  estimate(StaticMappingQ1<1>::mapping, dof_handler, quadrature, neumann_bc, solutions,
+          errors, component_mask, coefficients, n_threads, subdomain_id, material_id);
+}
+
+
+
+
+template <typename InputVector, class DH>
+void KellyErrorEstimator<1>::
+estimate (const Mapping<1>                    &mapping,
+          const DH                 &dof_handler,
+          const hp::QCollection<0>                 &,
+          const FunctionMap<1>::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 unsigned int,
+          const unsigned int                  subdomain_id,
+          const unsigned int                  material_id)
+{
+  Assert (false, ExcInternalError());
+}
+
+
+
 template <typename InputVector, class DH>
 void KellyErrorEstimator<1>::
 estimate (const Mapping<1>                    &mapping,
@@ -318,7 +403,7 @@ estimate (const Mapping<1>                    &mapping,
                    
                       for (unsigned int s=0; s<n_solution_vectors; ++s)
                         grad_neighbor[s] = v;
-                    };
+                    }
                 }
               else
                                                  // fill with zeroes.
@@ -464,6 +549,51 @@ estimate (const DH                &dof_handler,
 }
 
 
+template <int dim>
+template <typename InputVector, class DH>
+void
+KellyErrorEstimator<dim>::
+estimate (const Mapping<dim>      &mapping,
+          const DH                &dof_handler,
+          const hp::QCollection<dim-1> &quadrature,
+          const typename FunctionMap<dim>::type &neumann_bc,
+          const InputVector       &solution,
+          Vector<float>           &error,
+          const std::vector<bool> &component_mask,
+          const Function<dim>     *coefficients,
+          const unsigned int       n_threads,
+          const unsigned int       subdomain_id,
+          const unsigned int       material_id)
+{
+                                  // just pass on to the other function
+  const std::vector<const InputVector *> solutions (1, &solution);
+  std::vector<Vector<float>*>              errors (1, &error);
+  estimate (mapping, dof_handler, quadrature, neumann_bc, solutions, errors,
+           component_mask, coefficients, n_threads, subdomain_id, material_id);
+}
+
+
+template <int dim>
+template <typename InputVector, class DH>
+void
+KellyErrorEstimator<dim>::
+estimate (const DH                &dof_handler,
+          const hp::QCollection<dim-1> &quadrature,
+          const typename FunctionMap<dim>::type &neumann_bc,
+          const InputVector       &solution,
+          Vector<float>           &error,
+          const std::vector<bool> &component_mask,
+          const Function<dim>     *coefficients,
+          const unsigned int       n_threads,
+          const unsigned int       subdomain_id,
+          const unsigned int       material_id)
+{
+  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
+  estimate(StaticMappingQ1<dim>::mapping, dof_handler, quadrature, neumann_bc, solution,
+          error, component_mask, coefficients, n_threads,
+           subdomain_id, material_id);
+}
+
 
 
 template <int dim>
@@ -720,6 +850,27 @@ estimate_some (const hp::MappingCollection<dim>                  &mapping,
 
 
 
+template <int dim>
+template <typename InputVector, class DH>
+void
+KellyErrorEstimator<dim>::
+estimate (const Mapping<dim>                  &/*mapping*/,
+          const DH                            &/*dof_handler*/,
+          const hp::QCollection<dim-1>        &/*quadrature*/,
+          const typename FunctionMap<dim>::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 unsigned int                   /*n_threads_*/,
+          const unsigned int                   /*subdomain_id*/,
+          const unsigned int                   /*material_id*/)
+{
+  Assert (false, ExcNotImplemented());
+}
+
+
+
 template <int dim>
 template <typename InputVector, class DH>
 void
@@ -882,7 +1033,7 @@ estimate (const Mapping<dim>                  &mapping,
       (*errors[n]).reinit (dof_handler.get_tria().n_active_cells());
       for (unsigned int i=0;i<dof_handler.get_tria().n_active_cells();++i)
        (*errors[n])(i)=0;
-    };
+    }
 
   unsigned int present_cell=0;
 
@@ -950,6 +1101,27 @@ void KellyErrorEstimator<dim>::estimate (const DH                            &do
 
 
 
+template <int dim>
+template <typename InputVector, class DH>
+void KellyErrorEstimator<dim>::estimate (const DH                            &dof_handler,
+                                        const hp::QCollection<dim-1>             &quadrature,
+                                        const typename FunctionMap<dim>::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 unsigned int                   n_threads,
+                                         const unsigned int       subdomain_id,
+                                         const unsigned int       material_id)
+{
+  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));  
+  estimate(StaticMappingQ1<dim>::mapping, dof_handler, quadrature, neumann_bc, solutions,
+          errors, component_mask, coefficients, n_threads,
+           subdomain_id, material_id);
+}
+
+
+
 template <int dim>
 template <typename InputVector, class DH>
 void KellyErrorEstimator<dim>::
@@ -1076,8 +1248,8 @@ integrate_over_regular_face (const DH                                 &dof_handl
              for (unsigned int point=0; point<n_q_points; ++point)
                per_thread_data.phi[n][point][component] *=
                  per_thread_data.coefficient_values[point](component);
-       };
-    };
+       }
+    }
 
 
   if (face->at_boundary() == true)
@@ -1115,8 +1287,8 @@ integrate_over_regular_face (const DH                                 &dof_handl
            for (unsigned int component=0; component<n_components; ++component)
              for (unsigned int point=0; point<n_q_points; ++point)
                per_thread_data.phi[n][point][component] -= g[point](component);
-       };
-    };
+       }
+    }
 
 
                                   // now phi contains the following:
@@ -1292,8 +1464,8 @@ integrate_over_irregular_face (const DH                                 &dof_han
                  for (unsigned int point=0; point<n_q_points; ++point)
                    per_thread_data.phi[n][point][component] *=
                      per_thread_data.coefficient_values[point](component);
-           };
-       };
+           }
+       }
 
                                       // get the weights for the
                                       // integration. note that it
@@ -1319,7 +1491,7 @@ integrate_over_irregular_face (const DH                                 &dof_han
                                  per_thread_data.JxW_values[p];
 
       face_integrals[neighbor_child->face(neighbor_neighbor)] = face_integral;
-    };
+    }
 
 
                                   // finally loop over all subfaces to
@@ -1337,7 +1509,7 @@ integrate_over_irregular_face (const DH                                 &dof_han
       
       for (unsigned int n=0; n<n_solution_vectors; ++n)
        sum[n] += face_integrals[face->child(subface_no)][n];
-    };
+    }
 
   face_integrals[face] = sum;
 }
@@ -1353,16 +1525,14 @@ template class KellyErrorEstimator<deal_II_dimension>;
 
 // 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. note that we need a space between
-// "VectorType" and ">" to disambiguate ">>" when VectorType trails in an
-// angle bracket
-#define INSTANTIATE(InputVector,DH) \
+// variable VectorType/MatrixType
+#define INSTANTIATE_1(InputVector,DH,Q) \
 template    \
 void    \
 KellyErrorEstimator<deal_II_dimension>::    \
 estimate<InputVector,DH > (const Mapping<deal_II_dimension>      &,    \
           const DH   &,    \
-          const Quadrature<deal_II_dimension-1> &,    \
+          const Q<deal_II_dimension-1> &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const InputVector       &,    \
           Vector<float>           &,    \
@@ -1376,7 +1546,7 @@ template    \
 void    \
 KellyErrorEstimator<deal_II_dimension>::    \
 estimate<InputVector,DH > (const DH   &,    \
-          const Quadrature<deal_II_dimension-1> &,    \
+          const Q<deal_II_dimension-1> &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const InputVector       &,    \
           Vector<float>           &,    \
@@ -1391,7 +1561,7 @@ void    \
 KellyErrorEstimator<deal_II_dimension>::    \
 estimate<InputVector,DH > (const Mapping<deal_II_dimension>          &,    \
           const DH       &,    \
-          const Quadrature<deal_II_dimension-1>     &,    \
+          const Q<deal_II_dimension-1>     &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const std::vector<const InputVector *> &,    \
           std::vector<Vector<float>*> &,    \
@@ -1405,7 +1575,7 @@ template    \
 void    \
 KellyErrorEstimator<deal_II_dimension>::    \
 estimate<InputVector,DH > (const DH       &,    \
-          const Quadrature<deal_II_dimension-1>     &,    \
+          const Q<deal_II_dimension-1>     &,    \
           const FunctionMap<deal_II_dimension>::type &,    \
           const std::vector<const InputVector *> &,    \
           std::vector<Vector<float>*> &,    \
@@ -1415,6 +1585,11 @@ estimate<InputVector,DH > (const DH       &,    \
           const unsigned int           , \
           const unsigned int)
 
+#define INSTANTIATE(InputVector,DH) \
+     INSTANTIATE_1(InputVector,DH,Quadrature); \
+     INSTANTIATE_1(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>);

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.