]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move initialization of MappingQ1::InternalData into this class. 1425/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 24 Aug 2015 03:58:42 +0000 (22:58 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 25 Aug 2015 01:57:18 +0000 (20:57 -0500)
This allows removing the poorly named 'compute_data()' function from
the interface of MappingQ1 (don't all functions somehow compute some
data?).

include/deal.II/fe/mapping_q1.h
source/fe/mapping_q.cc
source/fe/mapping_q1.cc

index c735f065ae8229b4032c9422e7aee30153549e0b..48b9d00ffbce3668c9ef588b914b10138bdce768 100644 (file)
@@ -165,6 +165,24 @@ public:
      */
     InternalData(const unsigned int n_shape_functions);
 
+    /**
+     * Initialize the object's member variables related to cell data
+     * based on the given arguments.
+     */
+    void
+    initialize (const UpdateFlags      update_flags,
+                const Quadrature<dim> &quadrature,
+                const unsigned int     n_original_q_points);
+
+    /**
+     * Initialize the object's member variables related to cell and
+     * face data based on the given arguments.
+     */
+    void
+    initialize_face (const UpdateFlags      update_flags,
+                     const Quadrature<dim> &quadrature,
+                     const unsigned int     n_original_q_points);
+
     /**
      * Shape function at quadrature point. Shape functions are in tensor
      * product order, so vertices must be reordered to obtain transformation.
@@ -414,27 +432,6 @@ protected:
   void compute_shapes (const std::vector<Point<dim> > &unit_points,
                        InternalData &data) const;
 
-  /**
-   * Do the computations for the @p get_data functions. Here, the data vectors
-   * of @p InternalData are reinitialized to proper size and shape values are
-   * computed.
-   */
-  void compute_data (const UpdateFlags flags,
-                     const Quadrature<dim> &quadrature,
-                     const unsigned int n_orig_q_points,
-                     InternalData &data) const;
-
-  /**
-   * Do the computations for the @p get_face_data functions. Here, the data
-   * vectors of @p InternalData are reinitialized to proper size and shape
-   * values and derivatives are computed. Furthermore @p unit_tangential
-   * vectors of the face are computed.
-   */
-  void compute_face_data (const UpdateFlags flags,
-                          const Quadrature<dim> &quadrature,
-                          const unsigned int n_orig_q_points,
-                          InternalData &data) const;
-
   /**
    * Compute shape values and/or derivatives.
    */
index 29fa65d00e1cec13b9fd5f7806308967d973e6e4..ad5784c1de549a787eee665d0cc92afc6bde3d86 100644 (file)
@@ -203,22 +203,26 @@ MappingQ<dim,spacedim>::get_data (const UpdateFlags update_flags,
   // fill the data of both the Q_p and the Q_1 objects in parallel
   Threads::TaskGroup<> tasks;
   tasks += Threads::new_task (std_cxx11::function<void()>
-                              (std_cxx11::bind(&MappingQ<dim,spacedim>::compute_data,
-                                               this,
+                              (std_cxx11::bind(&MappingQ1<dim,spacedim>::InternalData::initialize,
+                                               data,
                                                update_flags,
                                                std_cxx11::cref(quadrature),
-                                               quadrature.size(),
-                                               std_cxx11::ref(*data))));
+                                               quadrature.size())));
   if (!use_mapping_q_on_all_cells)
     tasks += Threads::new_task (std_cxx11::function<void()>
-                                (std_cxx11::bind(&MappingQ<dim,spacedim>::compute_data,
-                                                 this,
+                                (std_cxx11::bind(&MappingQ1<dim,spacedim>::InternalData::initialize,
+                                                 &data->mapping_q1_data,
                                                  update_flags,
                                                  std_cxx11::cref(quadrature),
-                                                 quadrature.size(),
-                                                 std_cxx11::ref(data->mapping_q1_data))));
+                                                 quadrature.size())));
   tasks.join_all ();
 
+  // TODO: parallelize this as well
+  compute_shapes (quadrature.get_points(), *data);
+  if (!use_mapping_q_on_all_cells)
+    compute_shapes (quadrature.get_points(), data->mapping_q1_data);
+
+
   return data;
 }
 
@@ -235,22 +239,25 @@ MappingQ<dim,spacedim>::get_face_data (const UpdateFlags update_flags,
   // fill the data of both the Q_p and the Q_1 objects in parallel
   Threads::TaskGroup<> tasks;
   tasks += Threads::new_task (std_cxx11::function<void()>
-                              (std_cxx11::bind(&MappingQ<dim,spacedim>::compute_face_data,
-                                               this,
+                              (std_cxx11::bind(&MappingQ1<dim,spacedim>::InternalData::initialize_face,
+                                               data,
                                                update_flags,
                                                std_cxx11::cref(q),
-                                               quadrature.size(),
-                                               std_cxx11::ref(*data))));
+                                               quadrature.size())));
   if (!use_mapping_q_on_all_cells)
     tasks += Threads::new_task (std_cxx11::function<void()>
-                                (std_cxx11::bind(&MappingQ<dim,spacedim>::compute_face_data,
-                                                 this,
+                                (std_cxx11::bind(&MappingQ1<dim,spacedim>::InternalData::initialize_face,
+                                                 &data->mapping_q1_data,
                                                  update_flags,
                                                  std_cxx11::cref(q),
-                                                 quadrature.size(),
-                                                 std_cxx11::ref(data->mapping_q1_data))));
+                                                 quadrature.size())));
   tasks.join_all ();
 
+  // TODO: parallelize this as well
+  compute_shapes (q.get_points(), *data);
+  if (!use_mapping_q_on_all_cells)
+    compute_shapes (q.get_points(), data->mapping_q1_data);
+
   return data;
 }
 
@@ -267,22 +274,26 @@ MappingQ<dim,spacedim>::get_subface_data (const UpdateFlags update_flags,
   // fill the data of both the Q_p and the Q_1 objects in parallel
   Threads::TaskGroup<> tasks;
   tasks += Threads::new_task (std_cxx11::function<void()>
-                              (std_cxx11::bind(&MappingQ<dim,spacedim>::compute_face_data,
-                                               this,
+                              (std_cxx11::bind(&MappingQ1<dim,spacedim>::InternalData::initialize_face,
+                                               data,
                                                update_flags,
                                                std_cxx11::cref(q),
-                                               quadrature.size(),
-                                               std_cxx11::ref(*data))));
+                                               quadrature.size())));
   if (!use_mapping_q_on_all_cells)
     tasks += Threads::new_task (std_cxx11::function<void()>
-                                (std_cxx11::bind(&MappingQ<dim,spacedim>::compute_face_data,
-                                                 this,
+                                (std_cxx11::bind(&MappingQ1<dim,spacedim>::InternalData::initialize_face,
+                                                 &data->mapping_q1_data,
                                                  update_flags,
                                                  std_cxx11::cref(q),
-                                                 quadrature.size(),
-                                                 std_cxx11::ref(data->mapping_q1_data))));
+                                                 quadrature.size())));
   tasks.join_all ();
 
+  // TODO: parallelize this as well
+  compute_shapes (q.get_points(), *data);
+  if (!use_mapping_q_on_all_cells)
+    compute_shapes (q.get_points(), data->mapping_q1_data);
+
+
   return data;
 }
 
index b1e4ec2a6e1c14e8701016b74e94b07a95356e4f..4d0d1a1eb75812062e15305c564c4ddc410b5639 100644 (file)
@@ -68,6 +68,116 @@ MappingQ1<dim,spacedim>::InternalData::memory_consumption () const
 }
 
 
+template <int dim, int spacedim>
+void
+MappingQ1<dim,spacedim>::InternalData::
+initialize (const UpdateFlags      update_flags,
+            const Quadrature<dim> &q,
+            const unsigned int     n_original_q_points)
+{
+  // store the flags in the internal data object so we can access them
+  // in fill_fe_*_values()
+  this->update_each = update_flags;
+
+  const unsigned int n_q_points = q.size();
+
+  // see if we need the (transformation) shape function values
+  // and/or gradients and resize the necessary arrays
+  if (this->update_each & update_quadrature_points)
+    shape_values.resize(n_shape_functions * n_q_points);
+
+  if (this->update_each & (update_covariant_transformation
+                           | update_contravariant_transformation
+                           | update_JxW_values
+                           | update_boundary_forms
+                           | update_normal_vectors
+                           | update_jacobians
+                           | update_jacobian_grads
+                           | update_inverse_jacobians))
+    shape_derivatives.resize(n_shape_functions * n_q_points);
+
+  if (this->update_each & update_covariant_transformation)
+    covariant.resize(n_original_q_points);
+
+  if (this->update_each & update_contravariant_transformation)
+    contravariant.resize(n_original_q_points);
+
+  if (this->update_each & update_volume_elements)
+    volume_elements.resize(n_original_q_points);
+
+  if (this->update_each & update_jacobian_grads)
+    shape_second_derivatives.resize(n_shape_functions * n_q_points);
+}
+
+
+template <int dim, int spacedim>
+void
+MappingQ1<dim,spacedim>::InternalData::
+initialize_face (const UpdateFlags      update_flags,
+                 const Quadrature<dim> &q,
+                 const unsigned int     n_original_q_points)
+{
+  initialize (update_flags, q, n_original_q_points);
+
+  if (dim > 1)
+    {
+      if (this->update_each & update_boundary_forms)
+        {
+          aux.resize (dim-1, std::vector<Tensor<1,spacedim> > (n_original_q_points));
+
+          // Compute tangentials to the
+          // unit cell.
+          const unsigned int nfaces = GeometryInfo<dim>::faces_per_cell;
+          unit_tangentials.resize (nfaces*(dim-1),
+                                   std::vector<Tensor<1,dim> > (n_original_q_points));
+          if (dim==2)
+            {
+              // ensure a counterclockwise
+              // orientation of tangentials
+              static const int tangential_orientation[4]= {-1,1,1,-1};
+              for (unsigned int i=0; i<nfaces; ++i)
+                {
+                  Tensor<1,dim> tang;
+                  tang[1-i/2]=tangential_orientation[i];
+                  std::fill (unit_tangentials[i].begin(),
+                             unit_tangentials[i].end(), tang);
+                }
+            }
+          else if (dim==3)
+            {
+              for (unsigned int i=0; i<nfaces; ++i)
+                {
+                  Tensor<1,dim> tang1, tang2;
+
+                  const unsigned int nd=
+                    GeometryInfo<dim>::unit_normal_direction[i];
+
+                  // first tangential
+                  // vector in direction
+                  // of the (nd+1)%3 axis
+                  // and inverted in case
+                  // of unit inward normal
+                  tang1[(nd+1)%dim]=GeometryInfo<dim>::unit_normal_orientation[i];
+                  // second tangential
+                  // vector in direction
+                  // of the (nd+2)%3 axis
+                  tang2[(nd+2)%dim]=1.;
+
+                  // same unit tangents
+                  // for all quadrature
+                  // points on this face
+                  std::fill (unit_tangentials[i].begin(),
+                             unit_tangentials[i].end(), tang1);
+                  std::fill (unit_tangentials[nfaces+i].begin(),
+                             unit_tangentials[nfaces+i].end(), tang2);
+                }
+            }
+        }
+    }
+}
+
+
+
 
 template<int dim, int spacedim>
 MappingQ1<dim,spacedim>::MappingQ1 ()
@@ -531,51 +641,6 @@ MappingQ1<dim,spacedim>::requires_update_flags (const UpdateFlags in) const
 }
 
 
-template<int dim, int spacedim>
-void
-MappingQ1<dim,spacedim>::compute_data (const UpdateFlags      update_flags,
-                                       const Quadrature<dim> &q,
-                                       const unsigned int     n_original_q_points,
-                                       InternalData          &data) const
-{
-  // store the flags in the internal data object so we can access them
-  // in fill_fe_*_values(). use the transitive hull of the required
-  // flags
-  data.update_each = requires_update_flags(update_flags);
-
-  const unsigned int n_q_points = q.size();
-
-  // see if we need the (transformation) shape function values
-  // and/or gradients and resize the necessary arrays
-  if (data.update_each & update_quadrature_points)
-    data.shape_values.resize(data.n_shape_functions * n_q_points);
-
-  if (data.update_each & (update_covariant_transformation
-                          | update_contravariant_transformation
-                          | update_JxW_values
-                          | update_boundary_forms
-                          | update_normal_vectors
-                          | update_jacobians
-                          | update_jacobian_grads
-                          | update_inverse_jacobians))
-    data.shape_derivatives.resize(data.n_shape_functions * n_q_points);
-
-  if (data.update_each & update_covariant_transformation)
-    data.covariant.resize(n_original_q_points);
-
-  if (data.update_each & update_contravariant_transformation)
-    data.contravariant.resize(n_original_q_points);
-
-  if (data.update_each & update_volume_elements)
-    data.volume_elements.resize(n_original_q_points);
-
-  if (data.update_each & update_jacobian_grads)
-    data.shape_second_derivatives.resize(data.n_shape_functions * n_q_points);
-
-  compute_shapes (q.get_points(), data);
-}
-
-
 
 template<int dim, int spacedim>
 typename MappingQ1<dim,spacedim>::InternalData *
@@ -583,76 +648,10 @@ MappingQ1<dim,spacedim>::get_data (const UpdateFlags update_flags,
                                    const Quadrature<dim> &q) const
 {
   InternalData *data = new InternalData(n_shape_functions);
-  compute_data (update_flags, q, q.size(), *data);
-  return data;
-}
-
-
+  data->initialize (requires_update_flags(update_flags), q, q.size());
+  compute_shapes (q.get_points(), *data);
 
-template<int dim, int spacedim>
-void
-MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
-                                            const Quadrature<dim> &q,
-                                            const unsigned int n_original_q_points,
-                                            InternalData &data) const
-{
-  compute_data (update_flags, q, n_original_q_points, data);
-
-  if (dim > 1)
-    {
-      if (data.update_each & update_boundary_forms)
-        {
-          data.aux.resize (dim-1, std::vector<Tensor<1,spacedim> > (n_original_q_points));
-
-          // Compute tangentials to the
-          // unit cell.
-          const unsigned int nfaces = GeometryInfo<dim>::faces_per_cell;
-          data.unit_tangentials.resize (nfaces*(dim-1),
-                                        std::vector<Tensor<1,dim> > (n_original_q_points));
-          if (dim==2)
-            {
-              // ensure a counterclockwise
-              // orientation of tangentials
-              static const int tangential_orientation[4]= {-1,1,1,-1};
-              for (unsigned int i=0; i<nfaces; ++i)
-                {
-                  Tensor<1,dim> tang;
-                  tang[1-i/2]=tangential_orientation[i];
-                  std::fill (data.unit_tangentials[i].begin(),
-                             data.unit_tangentials[i].end(), tang);
-                }
-            }
-          else if (dim==3)
-            {
-              for (unsigned int i=0; i<nfaces; ++i)
-                {
-                  Tensor<1,dim> tang1, tang2;
-
-                  const unsigned int nd=
-                    GeometryInfo<dim>::unit_normal_direction[i];
-
-                  // first tangential
-                  // vector in direction
-                  // of the (nd+1)%3 axis
-                  // and inverted in case
-                  // of unit inward normal
-                  tang1[(nd+1)%dim]=GeometryInfo<dim>::unit_normal_orientation[i];
-                  // second tangential
-                  // vector in direction
-                  // of the (nd+2)%3 axis
-                  tang2[(nd+2)%dim]=1.;
-
-                  // same unit tangents
-                  // for all quadrature
-                  // points on this face
-                  std::fill (data.unit_tangentials[i].begin(),
-                             data.unit_tangentials[i].end(), tang1);
-                  std::fill (data.unit_tangentials[nfaces+i].begin(),
-                             data.unit_tangentials[nfaces+i].end(), tang2);
-                }
-            }
-        }
-    }
+  return data;
 }
 
 
@@ -663,10 +662,11 @@ MappingQ1<dim,spacedim>::get_face_data (const UpdateFlags        update_flags,
                                         const Quadrature<dim-1> &quadrature) const
 {
   InternalData *data = new InternalData(n_shape_functions);
-  compute_face_data (update_flags,
-                     QProjector<dim>::project_to_all_faces(quadrature),
-                     quadrature.size(),
-                     *data);
+  data->initialize_face (requires_update_flags(update_flags),
+                         QProjector<dim>::project_to_all_faces(quadrature),
+                         quadrature.size());
+  compute_shapes (QProjector<dim>::project_to_all_faces(quadrature).get_points(),
+                  *data);
 
   return data;
 }
@@ -679,10 +679,12 @@ MappingQ1<dim,spacedim>::get_subface_data (const UpdateFlags update_flags,
                                            const Quadrature<dim-1>& quadrature) const
 {
   InternalData *data = new InternalData(n_shape_functions);
-  compute_face_data (update_flags,
-                     QProjector<dim>::project_to_all_subfaces(quadrature),
-                     quadrature.size(),
-                     *data);
+  data->initialize_face (requires_update_flags(update_flags),
+                         QProjector<dim>::project_to_all_subfaces(quadrature),
+                         quadrature.size());
+  compute_shapes (QProjector<dim>::project_to_all_subfaces(quadrature).get_points(),
+                  *data);
+
 
   return data;
 }

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.