]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove Mapping::update_once/each(). 1361/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 14 Aug 2015 01:18:08 +0000 (20:18 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 15 Aug 2015 11:42:04 +0000 (06:42 -0500)
These two functions computed something that, in essence, is only important
for the internal implementation of the mapping classes. It should thus
not be part of the public interface of these classes. Consequently,
remove them and replace them by the only thing that is of interest,
namely to compute transitive closure of the set of flags one needs to
compute -- as now done using the requires_update_flags() function.

doc/news/changes.h
include/deal.II/fe/fe.h
include/deal.II/fe/mapping.h
include/deal.II/fe/mapping_cartesian.h
include/deal.II/fe/mapping_fe_field.h
include/deal.II/fe/mapping_q1.h
source/fe/fe_values.cc
source/fe/mapping_cartesian.cc
source/fe/mapping_fe_field.cc
source/fe/mapping_q1.cc

index ff44c44573687301b594e7fd21687ad0c12a497f..b4f44d535b74b60f55c8c9c70ad93eadb61004da 100644 (file)
@@ -47,7 +47,17 @@ inconvenience this causes.
   class hierarchy. As part of a general overhaul, the FEValuesData class
   has also been removed.
   <br>
-  (Wolfgang Bangerth, 2015/07/20-2015/08/06)
+  (Wolfgang Bangerth, 2015/07/20-2015/08/13)
+  </li>
+
+  <li> Changed: The functions update_once() and update_each() in the
+  Mapping classes computed information that was, in essence, only of use
+  internally. No external piece of code actually needed to know which
+  pieces of information a mapping could compute once and which they needed
+  to compute on every cell. Consequently, these two functions have been
+  removed and have been replaced by Mapping::requires_update_flags().
+  <br>
+  (Wolfgang Bangerth, 2015/07/20-2015/08/13)
   </li>
 
   <li> Changed: The function DoFRenumbering::random() now produces different
index 1268ebfbf4bc3a80cc8c0e0020748de81ccd863b..e7852f39ec14f8b55f4abbc3f77ad32715426bc9 100644 (file)
@@ -419,7 +419,20 @@ public:
     UpdateFlags          update_once;
 
     /**
-     * Values updated on each cell by reinit.
+     * A set of update flags specifying the kind of information that
+     * an implementation of the FiniteElement interface needs to compute on
+     * each cell or face, i.e., in FiniteElement::fill_fe_values() and
+     * friends.
+     *
+     * This set of flags is stored here by implementations of
+     * FiniteElement::get_data(), FiniteElement::get_face_data(), or
+     * FiniteElement::get_subface_data(), and is that subset of the update
+     * flags passed to those functions that require re-computation on
+     * every cell. (The subset of the flags corresponding to
+     * information that can be computed once and for all already at
+     * the time of the call to FiniteElement::get_data() -- or an
+     * implementation of that interface -- need not be stored here
+     * because it has already been taken care of.)
      */
     UpdateFlags          update_each;
 
index c6a5ecfbc58b40a937087fe6752088970444d80e..fcd23e9d69039a922ea1aca5af4ffdfcbcb57882 100644 (file)
@@ -116,15 +116,6 @@ enum MappingType
  * The interface for filling the tables of FEValues is provided. Everything
  * else has to happen in derived classes.
  *
- * The following paragraph applies to the implementation of FEValues. Usage of
- * the class is as follows: first, call the functions @p update_once and @p
- * update_each with the update flags you need. This includes the flags needed
- * by the FiniteElement. Then call <tt>get_*_data</tt> and with the or'd
- * results.  This will initialize and return some internal data structures. On
- * the first cell, call <tt>fill_fe_*_values</tt> with the result of @p
- * update_once. Finally, on each cell, use <tt>fill_fe_*_values</tt> with the
- * result of @p update_each to compute values for a special cell.
- *
  * <h3>Mathematics of the mapping</h3>
  *
  * The mapping is a transformation $\mathbf x = \Phi(\mathbf{\hat x})$ which
@@ -195,13 +186,6 @@ enum MappingType
  * performed through the functions transform(). See the documentation there
  * for possible choices.
  *
- * <h3>Technical notes</h3>
- *
- * A hint to implementors: no function except the two functions @p
- * update_once and @p update_each may add any flags.
- *
- * For more information about the <tt>spacedim</tt> template parameter check
- * the documentation of FiniteElement or the one of Triangulation.
  *
  * <h3>References</h3>
  *
@@ -475,7 +459,20 @@ public:
     virtual ~InternalDataBase ();
 
     /**
-     * Values updated on each cell by Mapping::fill_fe_values() and friends.
+     * A set of update flags specifying the kind of information that
+     * an implementation of the Mapping interface needs to compute on
+     * each cell or face, i.e., in Mapping::fill_fe_values() and
+     * friends.
+     *
+     * This set of flags is stored here by implementations of
+     * Mapping::get_data(), Mapping::get_face_data(), or
+     * Mapping::get_subface_data(), and is that subset of the update
+     * flags passed to those functions that require re-computation on
+     * every cell. (The subset of the flags corresponding to
+     * information that can be computed once and for all already at
+     * the time of the call to Mapping::get_data() -- or an
+     * implementation of that interface -- need not be stored here
+     * because it has already been taken care of.)
      */
     UpdateFlags          update_each;
 
@@ -488,23 +485,28 @@ public:
 
 protected:
   /**
-   * Indicate fields to be updated in the constructor of FEValues. Especially,
-   * fields not asked for by FEValues, but computed for efficiency reasons
-   * will be notified here.
+   * Given a set of update flags, compute which other quantities <i>also</i>
+   * need to be computed in order to satisfy the request by the given flags.
+   * Then return the combination of the original set of flags and those
+   * just computed.
    *
-   * See
-   * @ref UpdateFlagsEssay.
-   */
-  virtual UpdateFlags update_once (const UpdateFlags) const = 0;
-
-  /**
-   * The same as update_once(), but for the flags to be updated for each grid
-   * cell.
+   * As an example, if @p update_flags contains update_JxW_values
+   * (i.e., the product of the determinant of the Jacobian and the
+   * weights provided by the quadrature formula), a mapping may
+   * require the computation of the full Jacobian matrix in order to
+   * compute its determinant. They would then return not just
+   * update_JxW_values, but also update_jacobians. (This is not how it
+   * is actually done internally in the derived classes that compute
+   * the JxW values -- they set update_contravariant_transformation
+   * instead, from which the determinant can also be computed -- but
+   * this does not take away from the instructiveness of the example.)
    *
    * See
    * @ref UpdateFlagsEssay.
    */
-  virtual UpdateFlags update_each (const UpdateFlags) const = 0;
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const = 0;
 
   /**
    * Prepare internal data structures and fill in values independent of the
@@ -513,7 +515,7 @@ protected:
    */
   virtual
   InternalDataBase *
-  get_data (const UpdateFlags,
+  get_data (const UpdateFlags      update_flags,
             const Quadrature<dim> &quadrature) const = 0;
 
   /**
@@ -524,8 +526,8 @@ protected:
    */
   virtual
   InternalDataBase *
-  get_face_data (const UpdateFlags flags,
-                 const Quadrature<dim-1>quadrature) const = 0;
+  get_face_data (const UpdateFlags        update_flags,
+                 const Quadrature<dim-1> &quadrature) const = 0;
 
   /**
    * Prepare internal data structure for transformation of children of faces
@@ -535,8 +537,8 @@ protected:
    */
   virtual
   InternalDataBase *
-  get_subface_data (const UpdateFlags flags,
-                    const Quadrature<dim-1>quadrature) const = 0;
+  get_subface_data (const UpdateFlags        update_flags,
+                    const Quadrature<dim-1> &quadrature) const = 0;
 
   /**
    * Compute information about the mapping from the reference cell
index 69ce9ddddc039325fee2a2eb8e88e91cadd2bc7b..c003cf1193ec0500e6d6ef396d06a5f20d073dad 100644 (file)
@@ -211,8 +211,13 @@ protected:
                      std::vector<Point<dim> > &normal_vectors) const;
 
 private:
-  virtual UpdateFlags update_once (const UpdateFlags) const;
-  virtual UpdateFlags update_each (const UpdateFlags) const;
+  /**
+   * Implementation of the corresponding function in the base class,
+   * Mapping::requires_update_flags(). See there for more information.
+   */
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
 
   /**
    * Value to indicate that a given face or subface number is invalid.
index 2086daf122fb5107870322c6315b62d1414d2e09..341adbc4dc9df003caa53af76a1c03f397816ad6 100644 (file)
@@ -506,18 +506,12 @@ private:
                           typename MappingFEField<dim, spacedim>::InternalData &data) const;
 
   /**
-   * Reimplemented from Mapping. See the documentation of the base class for
-   * detailed information.
-   */
-  virtual UpdateFlags
-  update_once (const UpdateFlags in) const;
-
-  /**
-   * Reimplemented from Mapping. See the documentation of the base class for
-   * detailed information.
+   * Implementation of the corresponding function in the base class,
+   * Mapping::requires_update_flags(). See there for more information.
    */
-  virtual UpdateFlags
-  update_each (const UpdateFlags in) const;
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
 
   /**
    * Reimplemented from Mapping. See the documentation of the base class for
index b3cd90fd66e1fc19d528ec65f2f30502596e63bb..7f050b1ccdc248ab213b78bd1210aa2ce0b48cf5 100644 (file)
@@ -474,49 +474,14 @@ protected:
 
 
 private:
-  /**
-   * Implementation of the interface in Mapping.
-   *
-   * Description of effects:
-   * <ul>
-   * <li> if @p update_quadrature_points is required, the output will contain
-   * @p update_transformation_values. This computes the values of the
-   * transformation basis polynomials at the unit cell quadrature points.
-   * <li> if any of @p update_covariant_transformation, @p
-   * update_contravariant_transformation, @p update_JxW_values, @p
-   * update_boundary_forms, @p update_normal_vectors is required, the output
-   * will contain @p update_transformation_gradients to compute derivatives of
-   * the transformation basis polynomials.
-   * </ul>
-   */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
 
   /**
-   * Implementation of the interface in Mapping.
-   *
-   * Description of effects if @p flags contains:
-   * <ul>
-   * <li> <code>update_quadrature_points</code> is copied to the output to
-   * compute the quadrature points on the real cell.
-   * <li> <code>update_JxW_values</code> is copied and requires @p
-   * update_boundary_forms on faces. The latter, because the surface element
-   * is just the norm of the boundary form.
-   * <li> <code>update_normal_vectors</code> is copied and requires @p
-   * update_boundary_forms. The latter, because the normal vector is the
-   * normalized boundary form.
-   * <li> <code>update_covariant_transformation</code> is copied and requires
-   * @p update_contravariant_transformation, since it is computed as the
-   * inverse of the latter.
-   * <li> <code>update_JxW_values</code> is copied and requires
-   * <code>update_contravariant_transformation</code>, since it is computed as
-   * one over determinant of the latter.
-   * <li> <code>update_boundary_forms</code> is copied and requires
-   * <code>update_contravariant_transformation</code>, since the boundary form
-   * is computed as the contravariant image of the normal vector to the unit
-   * cell.
-   * </ul>
+   * Implementation of the corresponding function in the base class,
+   * Mapping::requires_update_flags(). See there for more information.
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
 
   /**
    * Implementation of the Mapping::get_data() interface.
index 5301c6c116c2dc63ebea3a110e1e378eefc0005b..8fd7bbf61346b7bab17ecea259578f41d671ec80 100644 (file)
@@ -3214,8 +3214,7 @@ FEValuesBase<dim,spacedim>::compute_update_flags (const UpdateFlags update_flags
   UpdateFlags flags = update_flags
                       | fe->update_once (update_flags)
                       | fe->update_each (update_flags);
-  flags |= mapping->update_once (flags)
-           | mapping->update_each (flags);
+  flags |= mapping->requires_update_flags (flags);
 
   return flags;
 }
index 5c04dce90c650d6291b60a76dc6421e1eb0c8d53..8aca482fe4688261ef2ac6ad613ecdc812cd3359 100644 (file)
@@ -59,17 +59,13 @@ MappingCartesian<dim, spacedim>::InternalData::memory_consumption () const
 
 template<int dim, int spacedim>
 UpdateFlags
-MappingCartesian<dim, spacedim>::update_once (const UpdateFlags) const
-{
-  return update_default;
-}
-
-
-
-template<int dim, int spacedim>
-UpdateFlags
-MappingCartesian<dim, spacedim>::update_each (const UpdateFlags in) const
+MappingCartesian<dim, spacedim>::requires_update_flags (const UpdateFlags in) const
 {
+  // this mapping is pretty simple in that it can basically compute
+  // every piece of information wanted by FEValues without requiring
+  // computing any other quantities. boundary forms are one exception
+  // since they can be computed from the normal vectors without much
+  // further ado
   UpdateFlags out = in;
   if (out & update_boundary_forms)
     out |= update_normal_vectors;
@@ -86,7 +82,10 @@ MappingCartesian<dim, spacedim>::get_data (const UpdateFlags      update_flags,
 {
   InternalData *data = new InternalData (q);
 
-  data->update_each = update_each(update_flags);
+  // 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);
 
   return data;
 }
@@ -101,7 +100,14 @@ MappingCartesian<dim, spacedim>::get_face_data (const UpdateFlags update_flags,
   InternalData *data
     = new InternalData (QProjector<dim>::project_to_all_faces(quadrature));
 
-  data->update_each = update_each(update_flags);
+  // verify that we have computed the transitive hull of the required
+  // flags and that FEValues has faithfully passed them on to us
+  Assert (update_flags == requires_update_flags (update_flags),
+          ExcInternalError());
+
+  // store the flags in the internal data object so we can access them
+  // in fill_fe_*_values()
+  data->update_each = update_flags;
 
   return data;
 }
@@ -116,7 +122,14 @@ MappingCartesian<dim, spacedim>::get_subface_data (const UpdateFlags update_flag
   InternalData *data
     = new InternalData (QProjector<dim>::project_to_all_subfaces(quadrature));
 
-  data->update_each = update_each(update_flags);
+  // verify that we have computed the transitive hull of the required
+  // flags and that FEValues has faithfully passed them on to us
+  Assert (update_flags == requires_update_flags (update_flags),
+          ExcInternalError());
+
+  // store the flags in the internal data object so we can access them
+  // in fill_fe_*_values()
+  data->update_each = update_flags;
 
   return data;
 }
index 9d538b8f264e9aba93465240bffeeb0e1336e453..a136a0be412df067ffa2ade4aa43fbcf29ac1463 100644 (file)
@@ -125,43 +125,15 @@ MappingFEField<dim,spacedim,VECTOR,DH>::compute_shapes_virtual (
 
 template<int dim, int spacedim, class VECTOR, class DH>
 UpdateFlags
-MappingFEField<dim,spacedim,VECTOR,DH>::update_once (const UpdateFlags) const
+MappingFEField<dim,spacedim,VECTOR,DH>::requires_update_flags (const UpdateFlags in) const
 {
-  // nothing here
-  return update_default;
-}
-
-
-
-template<int dim, int spacedim, class VECTOR, class DH>
-UpdateFlags
-MappingFEField<dim,spacedim,VECTOR,DH>::update_each (const UpdateFlags in) const
-{
-  // Select flags of concern for the
-  // transformation.
-  UpdateFlags out = UpdateFlags(in & (update_quadrature_points
-                                      | update_covariant_transformation
-                                      | update_contravariant_transformation
-                                      | update_JxW_values
-                                      | update_boundary_forms
-                                      | update_normal_vectors
-                                      | update_volume_elements
-                                      | update_jacobians
-                                      | update_jacobian_grads
-                                      | update_inverse_jacobians));
-
-  // add flags if the respective
-  // quantities are necessary to
-  // compute what we need. note that
-  // some flags appear in both
-  // conditions and in subsequent
-  // set operations. this leads to
-  // some circular logic. the only
-  // way to treat this is to
-  // iterate. since there are 5
-  // if-clauses in the loop, it will
-  // take at most 4 iterations to
+  // add flags if the respective quantities are necessary to compute
+  // what we need. note that some flags appear in both conditions and
+  // in subsequent set operations. this leads to some circular
+  // logic. the only way to treat this is to iterate. since there are
+  // 5 if-clauses in the loop, it will take at most 4 iterations to
   // converge. do them:
+  UpdateFlags out = in;
   for (unsigned int i=0; i<5; ++i)
     {
       // The following is a little incorrect:
@@ -214,36 +186,38 @@ MappingFEField<dim,spacedim,VECTOR,DH>::compute_data (const UpdateFlags      upd
                                                       const unsigned int     n_original_q_points,
                                                       InternalData           &data) const
 {
-  data.update_each = update_each(update_flags);
+  // 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();
-  const UpdateFlags flags = update_once(update_flags) | update_each(update_flags);
 
   // see if we need the (transformation) shape function values
   // and/or gradients and resize the necessary arrays
-  if (flags & update_quadrature_points)
+  if (data.update_each & update_quadrature_points)
     data.shape_values.resize(data.n_shape_functions * n_q_points);
 
-  if (flags & (update_covariant_transformation
-               | update_contravariant_transformation
-               | update_JxW_values
-               | update_boundary_forms
-               | update_normal_vectors
-               | update_jacobians
-               | update_jacobian_grads
-               | update_inverse_jacobians))
+  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 (flags & update_covariant_transformation)
+  if (data.update_each & update_covariant_transformation)
     data.covariant.resize(n_original_q_points);
 
-  if (flags & update_contravariant_transformation)
+  if (data.update_each & update_contravariant_transformation)
     data.contravariant.resize(n_original_q_points);
 
-  if (flags & update_volume_elements)
+  if (data.update_each & update_volume_elements)
     data.volume_elements.resize(n_original_q_points);
 
-  if (flags & update_jacobian_grads)
+  if (data.update_each & update_jacobian_grads)
     data.shape_second_derivatives.resize(data.n_shape_functions * n_q_points);
 
   compute_shapes_virtual (q.get_points(), data);
@@ -261,7 +235,7 @@ MappingFEField<dim,spacedim,VECTOR,DH>::compute_face_data (const UpdateFlags upd
 
   if (dim > 1)
     {
-      if (update_flags & update_boundary_forms)
+      if (data.update_each & update_boundary_forms)
         {
           data.aux.resize (dim-1, std::vector<Tensor<1,spacedim> > (n_original_q_points));
 
index 67357d2e7688adaddc28efcbf968444ff9887611..fb913e9dc9ec34df3683e363aa958bfabdfc974b 100644 (file)
@@ -477,43 +477,15 @@ compute_shapes_virtual (const std::vector<Point<dim> > &unit_points,
 
 template<int dim, int spacedim>
 UpdateFlags
-MappingQ1<dim,spacedim>::update_once (const UpdateFlags) const
+MappingQ1<dim,spacedim>::requires_update_flags (const UpdateFlags in) const
 {
-  // nothing here
-  return update_default;
-}
-
-
-
-template<int dim, int spacedim>
-UpdateFlags
-MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
-{
-  // Select flags of concern for the
-  // transformation.
-  UpdateFlags out = UpdateFlags(in & (update_quadrature_points
-                                      | update_covariant_transformation
-                                      | update_contravariant_transformation
-                                      | update_JxW_values
-                                      | update_boundary_forms
-                                      | update_normal_vectors
-                                      | update_volume_elements
-                                      | update_jacobians
-                                      | update_jacobian_grads
-                                      | update_inverse_jacobians));
-
-  // add flags if the respective
-  // quantities are necessary to
-  // compute what we need. note that
-  // some flags appear in both the
-  // conditions and in subsequent
-  // set operations. this leads to
-  // some circular logic. the only
-  // way to treat this is to
-  // iterate. since there are 5
-  // if-clauses in the loop, it will
-  // take at most 5 iterations to
+  // add flags if the respective quantities are necessary to compute
+  // what we need. note that some flags appear in both the conditions
+  // and in subsequent set operations. this leads to some circular
+  // logic. the only way to treat this is to iterate. since there are
+  // 5 if-clauses in the loop, it will take at most 5 iterations to
   // converge. do them:
+  UpdateFlags out = in;
   for (unsigned int i=0; i<5; ++i)
     {
       // The following is a little incorrect:
@@ -564,36 +536,38 @@ MappingQ1<dim,spacedim>::compute_data (const UpdateFlags      update_flags,
                                        const unsigned int     n_original_q_points,
                                        InternalData          &data) const
 {
-  data.update_each = update_each(update_flags);
+  // 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();
-  const UpdateFlags flags = update_once(update_flags) | update_each(update_flags);
 
   // see if we need the (transformation) shape function values
   // and/or gradients and resize the necessary arrays
-  if (flags & update_quadrature_points)
+  if (data.update_each & update_quadrature_points)
     data.shape_values.resize(data.n_shape_functions * n_q_points);
 
-  if (flags & (update_covariant_transformation
-               | update_contravariant_transformation
-               | update_JxW_values
-               | update_boundary_forms
-               | update_normal_vectors
-               | update_jacobians
-               | update_jacobian_grads
-               | update_inverse_jacobians))
+  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 (flags & update_covariant_transformation)
+  if (data.update_each & update_covariant_transformation)
     data.covariant.resize(n_original_q_points);
 
-  if (flags & update_contravariant_transformation)
+  if (data.update_each & update_contravariant_transformation)
     data.contravariant.resize(n_original_q_points);
 
-  if (flags & update_volume_elements)
+  if (data.update_each & update_volume_elements)
     data.volume_elements.resize(n_original_q_points);
 
-  if (flags & update_jacobian_grads)
+  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);
@@ -624,7 +598,7 @@ MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
 
   if (dim > 1)
     {
-      if (update_flags & update_boundary_forms)
+      if (data.update_each & update_boundary_forms)
         {
           data.aux.resize (dim-1, std::vector<Tensor<1,spacedim> > (n_original_q_points));
 

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.