]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Change compute_mapping_support_points to return its data.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 1 Oct 2015 21:45:07 +0000 (16:45 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 1 Oct 2015 21:45:07 +0000 (16:45 -0500)
The previous version returned the list of points via its last argument,
but this is no longer necessary in times of C++11.

include/deal.II/fe/mapping_q1_eulerian.h
include/deal.II/fe/mapping_q_eulerian.h
include/deal.II/fe/mapping_q_generic.h
source/fe/mapping_q1_eulerian.cc
source/fe/mapping_q_eulerian.cc
source/fe/mapping_q_generic.cc

index fc17a1734a8d4825789ccc8cb6a4d177c6a1cbe1..af7dc1771b9a28106904e6cd3c434f5d718a28d0 100644 (file)
@@ -151,11 +151,12 @@ protected:
   /**
    * Compute the support points of the mapping. For the current class, these are
    * the vertices, as obtained by calling Mapping::get_vertices().
+   * See the documentation of MappingQGeneric::compute_mapping_support_points()
+   * for more information.
    */
   virtual
-  void
-  compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                 std::vector<Point<spacedim> > &a) const;
+  std::vector<Point<spacedim> >
+  compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const;
 
   /**
    * Reference to the vector of shifts.
index f4bd91316b074866b2bf18e9401c573a00deb4b0..7be9b15f48bd970b6f22f212a7b0235a83c36a85 100644 (file)
@@ -211,13 +211,13 @@ private:
   mutable Threads::Mutex fe_values_mutex;
 
   /**
-   * Compute the positions of the support points in the current configuration
+   * Compute the positions of the support points in the current configuration.
+   * See the documentation of MappingQGeneric::compute_mapping_support_points()
+   * for more information.
    */
   virtual
-  void
-  compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                 std::vector<Point<spacedim> > &a) const;
-
+  std::vector<Point<spacedim> >
+  compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const;
 };
 
 /*@}*/
index a4049bdc56082de09346972c75a57b1ae028f6fc..a96dbcb9f7c259ed9434befcaefc3bba19071560 100644 (file)
@@ -509,9 +509,12 @@ protected:
   Table<2,double> support_point_weights_on_hex;
 
   /**
-   * An interface that derived classes have to implement and that
-   * computes the locations of support points for the mapping. For
-   * example, for MappingQ1 these are the vertices. However, other
+   * Return the locations of support points for the mapping. For
+   * example, for $Q_1$ mappings these are the vertices, and for higher
+   * order polynomial mappings they are the vertices plus interior
+   * points on edges, faces, and the cell interior that are placed
+   * in consultation with the Manifold description of the domain and
+   * its boundary. However, other
    * classes may override this function differently. In particular,
    * the MappingQ1Eulerian class does exactly this by not computing
    * the support points from the geometry of the current cell but
@@ -537,9 +540,8 @@ protected:
    * the locations of interior points.
    */
   virtual
-  void
-  compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                  std::vector<Point<spacedim> > &a) const;
+  std::vector<Point<spacedim> >
+  compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell_iterator &cell) const;
 
   /**
    * Transforms the point @p p on the real cell to the corresponding point on
index 4393d813cb9877dde1f700f3906f97d04c82b789..08dedcf0243752b595c9ce77aee934e706aba593 100644 (file)
@@ -89,17 +89,18 @@ get_vertices
 
 
 template<int dim, class EulerVectorType, int spacedim>
-void
+std::vector<Point<spacedim> >
 MappingQ1Eulerian<dim,EulerVectorType,spacedim>::
-compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                               std::vector<Point<spacedim> > &a) const
+compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
 {
   const std_cxx11::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell>
   vertices = this->get_vertices(cell);
 
-  a.resize(GeometryInfo<dim>::vertices_per_cell);
+  std::vector<Point<spacedim> > a(GeometryInfo<dim>::vertices_per_cell);
   for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
     a[i] = vertices[i];
+
+  return a;
 }
 
 
index 51b86544983a447bb303f6a337d4719460b3721f..56d4cbad0da1822d6e5a86b08c5b5b65c7c042fc 100644 (file)
@@ -133,8 +133,7 @@ get_vertices
 (const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
 {
   // get the vertices as the first 2^dim mapping support points
-  std::vector<Point<spacedim> > a;
-  compute_mapping_support_points(cell, a);
+  const std::vector<Point<spacedim> > a = compute_mapping_support_points(cell);
 
   std_cxx11::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell> vertex_locations;
   std::copy (a.begin(),
@@ -147,10 +146,9 @@ get_vertices
 
 
 template <int dim, class EulerVectorType, int spacedim>
-void
+std::vector<Point<spacedim> >
 MappingQEulerian<dim, EulerVectorType, spacedim>::
-compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                std::vector<Point<spacedim> > &a) const
+compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
 {
   // first, basic assertion with respect to vector size,
 
@@ -197,13 +195,15 @@ compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell
 
   // and finally compute the positions of the support points in the deformed
   // configuration.
-  a.resize(n_support_pts);
+  std::vector<Point<spacedim> > a(n_support_pts);
   for (unsigned int q=0; q<n_support_pts; ++q)
     {
       a[q] = fe_values.quadrature_point(q);
       for (unsigned int d=0; d<spacedim; ++d)
         a[q](d) += shift_vector[q](d);
     }
+
+  return a;
 }
 
 
index fa4e637d0dfbe2e6870e09951915da79135ffaf6..42fea442f7f638f88032d8f9d33ea9ea56216057 100644 (file)
@@ -1040,8 +1040,8 @@ transform_unit_to_real_cell (const typename Triangulation<dim,spacedim>::cell_it
               FiniteElementData<dim> (get_dpo_vector<dim>(polynomial_degree), 1,
                                       polynomial_degree)));
 
-  std::vector<Point<spacedim> > support_points (tensor_pols.n());
-  this->compute_mapping_support_points(cell, support_points);
+  const std::vector<Point<spacedim> > support_points
+    = this->compute_mapping_support_points(cell);
 
   Point<spacedim> mapped_point;
   for (unsigned int i=0; i<tensor_pols.n(); ++i)
@@ -1438,7 +1438,7 @@ transform_real_to_unit_cell_internal
   std_cxx11::unique_ptr<InternalData> mdata (get_data(update_flags,
                                                       point_quadrature));
 
-  this->compute_mapping_support_points (cell, mdata->mapping_support_points);
+  mdata->mapping_support_points = this->compute_mapping_support_points (cell);
 
   // dispatch to the various specializations for spacedim=dim,
   // spacedim=dim+1, etc
@@ -1464,7 +1464,7 @@ transform_real_to_unit_cell_internal
   std_cxx11::unique_ptr<InternalData> mdata (get_data(update_flags,
                                                       point_quadrature));
 
-  this->compute_mapping_support_points (cell, mdata->mapping_support_points);
+  mdata->mapping_support_points = this->compute_mapping_support_points (cell);
 
   // dispatch to the various specializations for spacedim=dim,
   // spacedim=dim+1, etc
@@ -1490,7 +1490,7 @@ transform_real_to_unit_cell_internal
   std_cxx11::unique_ptr<InternalData> mdata (get_data(update_flags,
                                                       point_quadrature));
 
-  this->compute_mapping_support_points (cell, mdata->mapping_support_points);
+  mdata->mapping_support_points = this->compute_mapping_support_points (cell);
 
   // dispatch to the various specializations for spacedim=dim,
   // spacedim=dim+1, etc
@@ -1516,7 +1516,7 @@ transform_real_to_unit_cell_internal
   std_cxx11::unique_ptr<InternalData> mdata (get_data(update_flags,
                                                       point_quadrature));
 
-  this->compute_mapping_support_points (cell, mdata->mapping_support_points);
+  mdata->mapping_support_points = this->compute_mapping_support_points (cell);
 
   // dispatch to the various specializations for spacedim=dim,
   // spacedim=dim+1, etc
@@ -1542,7 +1542,7 @@ transform_real_to_unit_cell_internal
   std_cxx11::unique_ptr<InternalData> mdata (get_data(update_flags,
                                                       point_quadrature));
 
-  this->compute_mapping_support_points (cell, mdata->mapping_support_points);
+  mdata->mapping_support_points = this->compute_mapping_support_points (cell);
 
   // dispatch to the various specializations for spacedim=dim,
   // spacedim=dim+1, etc
@@ -1938,8 +1938,8 @@ transform_real_to_unit_cell (const typename Triangulation<dim,spacedim>::cell_it
       // Find the initial value for the Newton iteration by a normal
       // projection to the least square plane determined by the vertices
       // of the cell
-      std::vector<Point<spacedim> > a;
-      this->compute_mapping_support_points (cell,a);
+      const std::vector<Point<spacedim> > a
+        = this->compute_mapping_support_points (cell);
       Assert(a.size() == GeometryInfo<dim>::vertices_per_cell,
              ExcInternalError());
       initial_p_unit = internal::MappingQ1::transform_real_to_unit_cell_initial_guess<dim,spacedim>(a,p);
@@ -1955,8 +1955,8 @@ transform_real_to_unit_cell (const typename Triangulation<dim,spacedim>::cell_it
           // we do this by first getting all support points, then
           // throwing away all but the vertices, and finally calling
           // the same function as above
-          std::vector<Point<spacedim> > a;
-          this->compute_mapping_support_points (cell,a);
+          std::vector<Point<spacedim> > a
+            = this->compute_mapping_support_points (cell);
           a.resize(GeometryInfo<dim>::vertices_per_cell);
           initial_p_unit = internal::MappingQ1::transform_real_to_unit_cell_initial_guess<dim,spacedim>(a,p);
         }
@@ -2588,7 +2588,7 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
       ||
       (cell != data.cell_of_current_support_points))
     {
-      this->compute_mapping_support_points(cell, data.mapping_support_points);
+      data.mapping_support_points = this->compute_mapping_support_points(cell);
       data.cell_of_current_support_points = cell;
     }
 
@@ -2976,7 +2976,7 @@ fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &
       ||
       (cell != data.cell_of_current_support_points))
     {
-      this->compute_mapping_support_points(cell, data.mapping_support_points);
+      data.mapping_support_points = this->compute_mapping_support_points(cell);
       data.cell_of_current_support_points = cell;
     }
 
@@ -3021,7 +3021,7 @@ fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterato
       ||
       (cell != data.cell_of_current_support_points))
     {
-      this->compute_mapping_support_points(cell, data.mapping_support_points);
+      data.mapping_support_points = this->compute_mapping_support_points(cell);
       data.cell_of_current_support_points = cell;
     }
 
@@ -3862,13 +3862,12 @@ add_quad_support_points(const typename Triangulation<dim,spacedim>::cell_iterato
 
 
 template<int dim, int spacedim>
-void
+std::vector<Point<spacedim> >
 MappingQGeneric<dim,spacedim>::
-compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                               std::vector<Point<spacedim> > &a) const
+compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
 {
   // get the vertices first
-  a.resize(GeometryInfo<dim>::vertices_per_cell);
+  std::vector<Point<spacedim> > a(GeometryInfo<dim>::vertices_per_cell);
   for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
     a[i] = cell->vertex(i);
 
@@ -3906,6 +3905,8 @@ compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_
         Assert(false, ExcNotImplemented());
         break;
       }
+
+  return a;
 }
 
 

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.