]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Extend all functions of the KellyErrorEstimator to use an arbitrary mapping. For...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 5 Apr 2001 16:53:25 +0000 (16:53 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 5 Apr 2001 16:53:25 +0000 (16:53 +0000)
git-svn-id: https://svn.dealii.org/trunk@4383 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 96f53add2fc76574fc02b5575a229dbe1430622b..f58b695698f673336aa7f55bb4ada700cc463282 100644 (file)
@@ -269,7 +269,8 @@ class KellyErrorEstimator
                                      * implemented in two and three
                                      * dimensions.
                                      */
-    static void estimate (const DoFHandler<dim>   &dof,
+    static void estimate (const Mapping<dim>      &mapping,
+                         const DoFHandler<dim>   &dof,
                          const Quadrature<dim-1> &quadrature,
                          const FunctionMap       &neumann_bc,
                          const Vector<double>    &solution,
@@ -278,6 +279,19 @@ class KellyErrorEstimator
                          const Function<dim>     *coefficients   = 0,
                          unsigned int             n_threads = multithread_info.n_default_threads);
 
+                                    /**
+                                     * Calls the @p{estimate}
+                                     * function, see above, with
+                                     * @p{mapping=MappingQ1<dim>()}.
+                                     */    
+    static void estimate (const DoFHandler<dim>   &dof,
+                         const Quadrature<dim-1> &quadrature,
+                         const FunctionMap       &neumann_bc,
+                         const Vector<double>    &solution,
+                         Vector<float>           &error,
+                         const std::vector<bool> &component_mask = std::vector<bool>(),
+                         const Function<dim>     *coefficients   = 0,
+                         unsigned int             n_threads = multithread_info.n_default_threads);
                                     /**
                                      * Same function as above, but
                                      * accepts more than one solution
@@ -305,14 +319,29 @@ class KellyErrorEstimator
                                      * references, so we had to use a
                                      * vector of pointers.)
                                      */
-    static void estimate (const DoFHandler<dim>               &dof,
-                         const Quadrature<dim-1>             &quadrature,
-                         const FunctionMap                   &neumann_bc,
+    static void estimate (const Mapping<dim>          &mapping,
+                         const DoFHandler<dim>       &dof,
+                         const Quadrature<dim-1>     &quadrature,
+                         const FunctionMap           &neumann_bc,
+                         const std::vector<const Vector<double>*> &solutions,
+                         std::vector<Vector<float>*> &errors,
+                         const std::vector<bool>     &component_mask = std::vector<bool>(),
+                         const Function<dim>         *coefficients   = 0,
+                         unsigned int                 n_threads = multithread_info.n_default_threads);
+
+                                    /**
+                                     * Calls the @p{estimate}
+                                     * function, see above, with
+                                     * @p{mapping=MappingQ1<dim>()}.
+                                     */    
+    static void estimate (const DoFHandler<dim>       &dof,
+                         const Quadrature<dim-1>     &quadrature,
+                         const FunctionMap           &neumann_bc,
                          const std::vector<const Vector<double>*> &solutions,
-                         std::vector<Vector<float>*>              &errors,
-                         const std::vector<bool>                  &component_mask = std::vector<bool>(),
-                         const Function<dim>                 *coefficients   = 0,
-                         unsigned int                         n_threads = multithread_info.n_default_threads);
+                         std::vector<Vector<float>*> &errors,
+                         const std::vector<bool>     &component_mask = std::vector<bool>(),
+                         const Function<dim>         *coefficients   = 0,
+                         unsigned int                 n_threads = multithread_info.n_default_threads);
 
     
                                     /**
@@ -374,7 +403,7 @@ class KellyErrorEstimator
                                      * functions of the error
                                      * estimator is gathered in this
                                      * struct. It is passed as a
-                                     * reference to the seperate
+                                     * reference to the separate
                                      * functions in this class.
                                      *
                                      * The reason for invention of
@@ -415,6 +444,7 @@ class KellyErrorEstimator
                                      */
     struct Data
     {
+       const Mapping<dim>                  &mapping;
        const DoFHandler<dim>               &dof_handler;
        const Quadrature<dim-1>             &quadrature;
        const FunctionMap                   &neumann_bc;
@@ -500,7 +530,8 @@ class KellyErrorEstimator
                                          * class Data. All variables are
                                          * passed as references.
                                          */
-       Data(const DoFHandler<dim>               &dof,
+       Data(const Mapping<dim>                  &mapping,
+            const DoFHandler<dim>               &dof,
             const Quadrature<dim-1>             &quadrature,
             const FunctionMap                   &neumann_bc,
             const std::vector<const Vector<double>*> &solutions,
index c4c7cef0046f4a46903b16e076ce321a14ce9bd1..13c1f0b3d9cdcc108fea2cac41d65e8f3b66544c 100644 (file)
@@ -39,10 +39,6 @@ using namespace std;
 
 
 
-//TODO:[RH,GK] Replace global by local object; better: have two functions, or by default arg
-static const MappingQ1<deal_II_dimension> mapping;
-
-
 
 static
 inline
@@ -55,7 +51,8 @@ double sqr (const double x)
 #if deal_II_dimension == 1
 
 template <>
-KellyErrorEstimator<1>::Data::Data(const DoFHandler<1>                 &,
+KellyErrorEstimator<1>::Data::Data(const Mapping<1>                    &,
+                                  const DoFHandler<1>                 &,
                                   const Quadrature<0>                 &,
                                   const FunctionMap                   &,
                                   const std::vector<const Vector<double>*> &,
@@ -63,6 +60,7 @@ KellyErrorEstimator<1>::Data::Data(const DoFHandler<1>                 &,
                                   const Function<1>                   *,
                                   const unsigned int                   ,
                                   FaceIntegrals                       &):
+               mapping(* static_cast <const Mapping<1> *> (0)),
                dof_handler(* static_cast <const DoFHandler<1> *> (0)),
                quadrature(* static_cast <const Quadrature<0> *> (0)),
                neumann_bc(* static_cast <const FunctionMap *> (0)),
@@ -75,7 +73,8 @@ KellyErrorEstimator<1>::Data::Data(const DoFHandler<1>                 &,
 #else
 
 template <int dim>
-KellyErrorEstimator<dim>::Data::Data(const DoFHandler<dim>               &dof_handler,
+KellyErrorEstimator<dim>::Data::Data(const Mapping<dim>                  &mapping,
+                                    const DoFHandler<dim>               &dof_handler,
                                     const Quadrature<dim-1>             &quadrature,
                                     const FunctionMap                   &neumann_bc,
                                     const std::vector<const Vector<double>*> &solutions,
@@ -83,6 +82,7 @@ KellyErrorEstimator<dim>::Data::Data(const DoFHandler<dim>               &dof_ha
                                     const Function<dim>                 *coefficients,
                                     const unsigned int                   n_threads,
                                     FaceIntegrals                       &face_integrals):
+               mapping (mapping),
                dof_handler (dof_handler),
                quadrature (quadrature),
                neumann_bc (neumann_bc),
@@ -149,7 +149,8 @@ KellyErrorEstimator<dim>::Data::Data(const DoFHandler<dim>               &dof_ha
 // the following function is still independent of dimension, but it
 // calls dimension dependent functions
 template <int dim>
-void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof_handler,
+void KellyErrorEstimator<dim>::estimate (const Mapping<dim>      &mapping,
+                                        const DoFHandler<dim>   &dof_handler,
                                         const Quadrature<dim-1> &quadrature,
                                         const FunctionMap       &neumann_bc,
                                         const Vector<double>    &solution,
@@ -161,11 +162,26 @@ void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof_handler,
                                   // just pass on to the other function
   const std::vector<const Vector<double>*> solutions (1, &solution);
   std::vector<Vector<float>*>              errors (1, &error);
-  estimate (dof_handler, quadrature, neumann_bc, solutions, errors,
+  estimate (mapping, dof_handler, quadrature, neumann_bc, solutions, errors,
            component_mask, coefficients, n_threads);
 };
 
 
+template <int dim>
+void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof_handler,
+                                        const Quadrature<dim-1> &quadrature,
+                                        const FunctionMap       &neumann_bc,
+                                        const Vector<double>    &solution,
+                                        Vector<float>           &error,
+                                        const std::vector<bool> &component_mask,
+                                        const Function<dim>     *coefficients,
+                                        unsigned int             n_threads)
+{
+  static const MappingQ1<dim> mapping;
+  estimate(mapping, dof_handler, quadrature, neumann_bc, solution,
+          error, component_mask, coefficients, n_threads);
+};
+
 
 
 #if deal_II_dimension == 1
@@ -181,7 +197,8 @@ void KellyErrorEstimator<1>::estimate_some (Data &, const unsigned int)
 
 
 template <>
-void KellyErrorEstimator<1>::estimate (const DoFHandler<1>                 &dof_handler,
+void KellyErrorEstimator<1>::estimate (const Mapping<1>                    &mapping,
+                                      const DoFHandler<1>                 &dof_handler,
                                       const Quadrature<0>                 &,
                                       const FunctionMap                   &neumann_bc,
                                       const std::vector<const Vector<double>*> &solutions,
@@ -376,7 +393,7 @@ void KellyErrorEstimator<dim>::estimate_some (Data               &data,
                                   // need not compute all values on the
                                   // neighbor cells, so using two objects
                                   // gives us a performance gain).
-  FEFaceValues<dim> fe_face_values_cell (mapping,
+  FEFaceValues<dim> fe_face_values_cell (data.mapping,
                                         data.dof_handler.get_fe(),
                                         data.quadrature,
                                         UpdateFlags(update_gradients      |
@@ -385,11 +402,11 @@ void KellyErrorEstimator<dim>::estimate_some (Data               &data,
                                                       (data.coefficients != 0))  ?
                                                      update_q_points : 0) |
                                                     update_normal_vectors)); 
-  FEFaceValues<dim> fe_face_values_neighbor (mapping,
+  FEFaceValues<dim> fe_face_values_neighbor (data.mapping,
                                             data.dof_handler.get_fe(),
                                             data.quadrature,
                                             update_gradients);
-  FESubfaceValues<dim> fe_subface_values (mapping,
+  FESubfaceValues<dim> fe_subface_values (data.mapping,
                                          data.dof_handler.get_fe(),
                                          data.quadrature,
                                          update_gradients);
@@ -503,7 +520,8 @@ void KellyErrorEstimator<dim>::estimate_some (Data               &data,
 
 
 template <int dim>
-void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>               &dof_handler,
+void KellyErrorEstimator<dim>::estimate (const Mapping<dim>                  &mapping,
+                                        const DoFHandler<dim>               &dof_handler,
                                         const Quadrature<dim-1>             &quadrature,
                                         const FunctionMap                   &neumann_bc,
                                         const std::vector<const Vector<double>*> &solutions,
@@ -564,7 +582,8 @@ void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>               &do
                                   // components
   std::vector<Data*> data_structures (n_threads);
   for (unsigned int i=0; i<n_threads; ++i)
-    data_structures[i] = new Data (dof_handler,
+    data_structures[i] = new Data (mapping,
+                                  dof_handler,
                                   quadrature,
                                   neumann_bc,
                                   solutions,
@@ -635,6 +654,22 @@ void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>               &do
 
 #endif
 
+template <int dim>
+void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>               &dof_handler,
+                                        const Quadrature<dim-1>             &quadrature,
+                                        const FunctionMap                   &neumann_bc,
+                                        const std::vector<const Vector<double>*> &solutions,
+                                        std::vector<Vector<float>*>              &errors,
+                                        const std::vector<bool>                  &component_mask,
+                                        const Function<dim>                 *coefficients,
+                                        unsigned int                         n_threads)
+{
+  static const MappingQ1<dim> mapping;
+  estimate(mapping, dof_handler, quadrature, neumann_bc, solutions,
+          errors, component_mask, coefficients, n_threads);
+}
+
+
 
 #if deal_II_dimension == 1
 

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.