]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Replace update_once/each by requires_update_flags() in FiniteElement. 1838/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 3 Nov 2015 04:33:52 +0000 (22:33 -0600)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 6 Nov 2015 00:21:06 +0000 (18:21 -0600)
The purpose of this patch is to align the finite element classes with the way the
mapping classes have already been converted. Specifically, there is no need for
any of the users of finite element classes to actually know whether a FE implementation
wants to treat a particular flags as update_once or update_each. This is an internal
decision. Rather, all we need to know is what flags they need overall. This is now
communicated by the new FiniteElement::requires_update_flags() function.

The update_once() and update_each() functions have been retained -- for now -- as
internal functions individual elements can implement, but they are no longer virtual.

16 files changed:
include/deal.II/fe/fe.h
include/deal.II/fe/fe_dgp_nonparametric.h
include/deal.II/fe/fe_face.h
include/deal.II/fe/fe_nothing.h
include/deal.II/fe/fe_poly.h
include/deal.II/fe/fe_poly.templates.h
include/deal.II/fe/fe_poly_face.h
include/deal.II/fe/fe_poly_face.templates.h
include/deal.II/fe/fe_poly_tensor.h
include/deal.II/fe/fe_system.h
source/fe/fe_dgp_nonparametric.cc
source/fe/fe_face.cc
source/fe/fe_nothing.cc
source/fe/fe_poly_tensor.cc
source/fe/fe_system.cc
source/fe/fe_values.cc

index 3f3acf796286127f1d01e8bc8cc67366933ee39c..70f222997d6b20a36ffb9f4b5fd1597dd4178408 100644 (file)
@@ -2169,36 +2169,28 @@ protected:
   compute_n_nonzero_components (const std::vector<ComponentMask> &nonzero_components);
 
   /**
-   * Determine the values a finite element should compute on initialization of
-   * data for FEValues.
+   * 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.
+   *
+   * As an example, if @p update_flags contains update_gradients
+   * a finite element class will typically
+   * require the computation of the inverse of the Jacobian matrix in order to
+   * rotate the gradient of shape functions on the reference cell to the
+   * real cell. It would then return not just
+   * update_gradients, but also update_covariant_transformation, the flag that
+   * makes the mapping class produce the inverse of the Jacobian matrix.
    *
-   * Given a set of flags indicating what quantities are requested from a
-   * FEValues object, update_once() and update_each() compute which values
-   * must really be computed. Then, the <tt>fill_*_values</tt> functions are
-   * called with the result of these.
-   *
-   * Furthermore, values must be computed either on the unit cell or on the
-   * physical cell. For instance, the function values of FE_Q do only depend
-   * on the quadrature points on the unit cell. Therefore, this flags will be
-   * returned by update_once(). The gradients require computation of the
-   * covariant transformation matrix. Therefore, @p
-   * update_covariant_transformation and @p update_gradients will be returned
-   * by update_each().
-   *
-   * For an example see the same function in the derived class FE_Q.
-   */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const = 0;
-
-  /**
-   * Complementary function for update_once().
-   *
-   * While update_once() returns the values to be computed on the unit cell
-   * for yielding the required data, this function determines the values that
-   * must be recomputed on each cell.
+   * An extensive discussion of the interaction between this function and
+   * FEValues can be found in the @ref FE_vs_Mapping_vs_FEValues documentation
+   * module.
    *
-   * Refer to update_once() for more details.
+   * @see UpdateFlags
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const = 0;
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const = 0;
 
   /**
    * Prepare internal data structures and fill in values independent of the
index ce173c0743a1ba1859567d43483c4833d0659557..50556f7097d94f9d64edced6b394b436928bed7d 100644 (file)
@@ -283,6 +283,11 @@ public:
    */
   virtual std::string get_name () const;
 
+  // for documentation, see the FiniteElement base class
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
+
   /**
    * Return the value of the @p ith shape function at the point @p p. See the
    * FiniteElement base class for more information about the semantics of this
@@ -582,7 +587,7 @@ private:
    *
    * For the present kind of finite element, this is exactly the case.
    */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  UpdateFlags update_once (const UpdateFlags flags) const;
 
   /**
    * This is the opposite to the above function: given a set of flags
@@ -593,7 +598,7 @@ private:
    * (for example, we often need the covariant transformation when gradients
    * need to be computed), include this in the result as well.
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  UpdateFlags update_each (const UpdateFlags flags) const;
 
   /**
    * Degree of the polynomials.
index ad9cfa624bf412f13e75b45efee8101671ea2275..acafae2abe753c8d4937e60ee3bd2e184b6279ef 100644 (file)
@@ -168,6 +168,11 @@ public:
    */
   virtual std::string get_name () const;
 
+  // for documentation, see the FiniteElement base class
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
+
   /**
    * Return the matrix interpolating from a face of of one element to the face
    * of the neighboring element.  The size of the matrix is then
@@ -322,7 +327,7 @@ protected:
    * All other flags of the result are cleared, since everything else must be
    * computed for each cell.
    */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  UpdateFlags update_once (const UpdateFlags flags) const;
 
   /**
    * Determine the values that need to be computed on every cell to be able to
@@ -352,7 +357,7 @@ protected:
    *
    * </ul>
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  UpdateFlags update_each (const UpdateFlags flags) const;
 
 private:
   /**
index 65bf2297aa0deda5e0c2bdb72a644790b26463a5..3985b5dac6cd703f8c24364426efbbe7c70ab3a5 100644 (file)
@@ -113,36 +113,10 @@ public:
   std::string
   get_name() const;
 
-  /**
-   * Determine the values a finite element should compute on initialization of
-   * data for FEValues.
-   *
-   * Given a set of flags indicating what quantities are requested from a
-   * FEValues object, update_once() and update_each() compute which values
-   * must really be computed. Then, the <tt>fill_*_values</tt> functions are
-   * called with the result of these.
-   *
-   * In this case, since the element has zero degrees of freedom and no
-   * information can be computed on it, this function simply returns the
-   * default (empty) set of update flags.
-   */
-
-  virtual
-  UpdateFlags
-  update_once (const UpdateFlags flags) const;
-
-  /**
-   * Complementary function for update_once().
-   *
-   * While update_once() returns the values to be computed on the unit cell
-   * for yielding the required data, this function determines the values that
-   * must be recomputed on each cell.
-   *
-   * Refer to update_once() for more details.
-   */
+  // for documentation, see the FiniteElement base class
   virtual
   UpdateFlags
-  update_each (const UpdateFlags flags) const;
+  requires_update_flags (const UpdateFlags update_flags) const;
 
   /**
    * Return the value of the @p ith shape function at the point @p p. @p p is
@@ -276,11 +250,11 @@ public:
   bool is_dominating() const;
 
 private:
+
   /**
    * If true, this element will dominate any other apart from itself in compare_for_face_domination();
    */
   const bool dominate;
-
 };
 
 
index c73b58117cd0cc931d880dbfc9eeec0105e379b3..a2f3bc3489e9021a8b1637a67a00697581bf1c86 100644 (file)
@@ -83,6 +83,11 @@ public:
    */
   unsigned int get_degree () const;
 
+  // for documentation, see the FiniteElement base class
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
+
   /**
    * Return the numbering of the underlying polynomial space compared to
    * lexicographic ordering of the basis functions. Returns
@@ -347,7 +352,7 @@ protected:
    * All other flags of the result are cleared, since everything else must be
    * computed for each cell.
    */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  UpdateFlags update_once (const UpdateFlags flags) const;
 
   /**
    * Determine the values that need to be computed on every cell to be able to
@@ -377,7 +382,7 @@ protected:
    *
    * </ul>
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  UpdateFlags update_each (const UpdateFlags flags) const;
 
   /**
    * Fields of cell-independent data.
index 74519fda038825cfa346228330ec300df38994f1..da43087d1cf38cc5cc744b108f98877c5273b4c8 100644 (file)
@@ -176,6 +176,13 @@ FE_Poly<POLY,dim,spacedim>::shape_4th_derivative_component (const unsigned int i
 //---------------------------------------------------------------------------
 
 
+template <class POLY, int dim, int spacedim>
+UpdateFlags
+FE_Poly<POLY,dim,spacedim>::requires_update_flags (const UpdateFlags flags) const
+{
+  return update_once(flags) | update_each(flags);
+}
+
 
 
 template <class POLY, int dim, int spacedim>
index 31d85ea509289fd98d66ed770f99844f4786d664..294fcee3679b99e9b86e6b7720ebc159a5a2f593 100644 (file)
@@ -74,6 +74,11 @@ public:
    */
   unsigned int get_degree () const;
 
+  // for documentation, see the FiniteElement base class
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
+
 protected:
   /*
    * NOTE: The following functions have their definitions inlined into the class declaration
@@ -196,7 +201,7 @@ protected:
    * All other flags of the result are cleared, since everything else must be
    * computed for each cell.
    */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  UpdateFlags update_once (const UpdateFlags flags) const;
 
   /**
    * Determine the values that need to be computed on every cell to be able to
@@ -226,7 +231,7 @@ protected:
    *
    * </ul>
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  UpdateFlags update_each (const UpdateFlags flags) const;
 
 
   /**
index 5497997a91a521326985832dd9f6962446607fd9..972118fcd84115902cd5859d8b9c1dda47f0b72d 100644 (file)
@@ -49,14 +49,18 @@ FE_PolyFace<POLY,dim,spacedim>::get_degree () const
 //---------------------------------------------------------------------------
 
 
+template <class POLY, int dim, int spacedim>
+UpdateFlags
+FE_PolyFace<POLY,dim,spacedim>::requires_update_flags (const UpdateFlags flags) const
+{
+  return update_once(flags) | update_each(flags);
+}
 
 
 template <class POLY, int dim, int spacedim>
 UpdateFlags
 FE_PolyFace<POLY,dim,spacedim>::update_once (const UpdateFlags) const
 {
-  // for this kind of elements, only the values can be precomputed once and
-  // for all. set this flag if the values are requested at all
   return update_default;
 }
 
index fefaba3f587b7e149a8fe24014c2a4c4f922f3ff..f1a0a303c2c4ef0c9e580f31b1d856c6174d7399 100644 (file)
@@ -126,6 +126,11 @@ public:
                  const std::vector<bool> &restriction_is_additive_flags,
                  const std::vector<ComponentMask> &nonzero_components);
 
+  // for documentation, see the FiniteElement base class
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
+
   /**
    * Since these elements are vector valued, an exception is thrown.
    */
@@ -156,29 +161,30 @@ public:
                                                    const Point<dim> &p,
                                                    const unsigned int component) const;
 
+protected:
+  /**
+   * The mapping type to be used to map shape functions from the reference
+   * cell to the mesh cell.
+   */
+  MappingType mapping_type;
+
+
   /**
    * Given <tt>flags</tt>, determines the values which must be computed only
    * for the reference cell. Make sure, that #mapping_type is set by the
    * derived class, such that this function can operate correctly.
    */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  UpdateFlags update_once (const UpdateFlags flags) const;
+
   /**
    * Given <tt>flags</tt>, determines the values which must be computed in
    * each cell cell. Make sure, that #mapping_type is set by the derived
    * class, such that this function can operate correctly.
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  UpdateFlags update_each (const UpdateFlags flags) const;
 
-protected:
-  /**
-   * The mapping type to be used to map shape functions from the reference
-   * cell to the mesh cell.
-   */
-  MappingType mapping_type;
-  /**
-  NOTE: The following function has its definition inlined into the class declaration
-    * because we otherwise run into a compiler error with MS Visual Studio.
-    */
+  /* NOTE: The following function has its definition inlined into the class declaration
+     because we otherwise run into a compiler error with MS Visual Studio. */
   virtual
   typename FiniteElement<dim,spacedim>::InternalDataBase *
   get_data(const UpdateFlags update_flags,
index ff256ef3a6f6651531d9d3795574451567132c75..5547211bd2afb0045a9749a1562ef931f8b036c4 100644 (file)
@@ -430,6 +430,11 @@ public:
    */
   virtual std::string get_name () const;
 
+  // for documentation, see the FiniteElement base class
+  virtual
+  UpdateFlags
+  requires_update_flags (const UpdateFlags update_flags) const;
+
   /**
    * Return the value of the @p ith shape function at the point @p p.  @p p is
    * a point on the reference element. Since this finite element is always
@@ -842,12 +847,12 @@ protected:
   /**
    * Compute flags for initial update only.
    */
-  virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  UpdateFlags update_once (const UpdateFlags flags) const;
 
   /**
    * Compute flags for update on each cell.
    */
-  virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  UpdateFlags update_each (const UpdateFlags flags) const;
 
   /**
    * @p clone function instead of a copy constructor.
index 00045beebd62696f0a86f04ff959fac410fe0f3b..8bc0ff473962bd23b81bc942f0bc51ba34f5a8e1 100644 (file)
@@ -216,6 +216,16 @@ FE_DGPNonparametric<dim,spacedim>::get_dpo_vector (const unsigned int deg)
 }
 
 
+
+template <int dim, int spacedim>
+UpdateFlags
+FE_DGPNonparametric<dim,spacedim>::requires_update_flags (const UpdateFlags flags) const
+{
+  return update_once(flags) | update_each(flags);
+}
+
+
+
 template <int dim, int spacedim>
 UpdateFlags
 FE_DGPNonparametric<dim,spacedim>::update_once (const UpdateFlags) const
index 7f3c8b94584a9b38d6a0b2cbe1d8159e8ac6a3ef..5297500ea5ec529075b0da8ecdf68b83c8d1ae8f 100644 (file)
@@ -449,6 +449,15 @@ FE_FaceQ<1,spacedim>::update_once (const UpdateFlags) const
 
 
 
+template <int spacedim>
+UpdateFlags
+FE_FaceQ<1,spacedim>::requires_update_flags (const UpdateFlags flags) const
+{
+  return update_once(flags) | update_each(flags);
+}
+
+
+
 template <int spacedim>
 UpdateFlags
 FE_FaceQ<1,spacedim>::update_each (const UpdateFlags flags) const
index d122101c4d97c174249c6afcca3cd574fc878d4f..3d5d9733c2451974f8ca5e9ccf3d314c0e0ddfc1 100644 (file)
@@ -72,18 +72,9 @@ FE_Nothing<dim,spacedim>::get_name () const
 
 template <int dim, int spacedim>
 UpdateFlags
-FE_Nothing<dim,spacedim>::update_once (const UpdateFlags /*flags*/) const
+FE_Nothing<dim,spacedim>::requires_update_flags (const UpdateFlags flags) const
 {
-  return update_default;
-}
-
-
-
-template <int dim, int spacedim>
-UpdateFlags
-FE_Nothing<dim,spacedim>::update_each (const UpdateFlags /*flags*/) const
-{
-  return update_default;
+  return flags;
 }
 
 
index cd116087e8922c3dc505a0754ce84471f5f846ad..dfc24c3bc70c6557c408f2143bc351730524779d 100644 (file)
@@ -1675,6 +1675,15 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
 }
 
 
+
+template <class POLY, int dim, int spacedim>
+UpdateFlags
+FE_PolyTensor<POLY,dim,spacedim>::requires_update_flags(const UpdateFlags flags) const
+{
+  return update_once(flags) | update_each(flags);
+}
+
+
 template <class POLY, int dim, int spacedim>
 UpdateFlags
 FE_PolyTensor<POLY,dim,spacedim>::update_once (const UpdateFlags flags) const
index 95f097d15620023327414d32a41991f2cade4ebf..bbcf05c3ca089caa8b936f70d60229b5631e3471 100644 (file)
@@ -844,28 +844,13 @@ face_to_cell_index (const unsigned int face_dof_index,
 
 template <int dim, int spacedim>
 UpdateFlags
-FESystem<dim,spacedim>::update_once (const UpdateFlags flags) const
+FESystem<dim,spacedim>::requires_update_flags (const UpdateFlags flags) const
 {
   UpdateFlags out = update_default;
   // generate maximal set of flags
   // that are necessary
   for (unsigned int base_no=0; base_no<this->n_base_elements(); ++base_no)
-    out |= base_element(base_no).update_once(flags);
-  return out;
-}
-
-
-
-template <int dim, int spacedim>
-UpdateFlags
-FESystem<dim,spacedim>::update_each (const UpdateFlags flags) const
-{
-  UpdateFlags out = update_default;
-  // generate maximal set of flags
-  // that are necessary
-  for (unsigned int base_no=0; base_no<this->n_base_elements(); ++base_no)
-    out |= base_element(base_no).update_each(flags);
-
+    out |= base_element(base_no).requires_update_flags (flags);
   return out;
 }
 
@@ -873,13 +858,19 @@ FESystem<dim,spacedim>::update_each (const UpdateFlags flags) const
 
 template <int dim, int spacedim>
 typename FiniteElement<dim,spacedim>::InternalDataBase *
-FESystem<dim,spacedim>::get_data (const UpdateFlags      flags_,
+FESystem<dim,spacedim>::get_data (const UpdateFlags      flags,
                                   const Mapping<dim,spacedim>    &mapping,
                                   const Quadrature<dim> &quadrature) const
 {
-  UpdateFlags flags = flags_;
+  // create an internal data object and set the update flags we will need
+  // to deal with. the current object does not make use of these flags,
+  // but we need to nevertheless set them correctly since we look
+  // into the update_each flag of base elements in fill_fe_values,
+  // and so the current object's update_each flag needs to be
+  // correct in case the current FESystem is a base element for another,
+  // higher-level FESystem itself.
   InternalData *data = new InternalData(this->n_base_elements());
-  data->update_each = update_each(flags) | update_once(flags);
+  data->update_each = requires_update_flags(flags);
 
   // get data objects from each of the base elements and store
   // them. do the creation of these objects in parallel as their
@@ -918,14 +909,19 @@ FESystem<dim,spacedim>::get_data (const UpdateFlags      flags_,
 
 template <int dim, int spacedim>
 typename FiniteElement<dim,spacedim>::InternalDataBase *
-FESystem<dim,spacedim>::get_face_data (
-  const UpdateFlags      flags_,
-  const Mapping<dim,spacedim>    &mapping,
-  const Quadrature<dim-1> &quadrature) const
-{
-  UpdateFlags flags = flags_;
+FESystem<dim,spacedim>::get_face_data (const UpdateFlags      flags,
+                                       const Mapping<dim,spacedim>    &mapping,
+                                       const Quadrature<dim-1> &quadrature) const
+{
+  // create an internal data object and set the update flags we will need
+  // to deal with. the current object does not make use of these flags,
+  // but we need to nevertheless set them correctly since we look
+  // into the update_each flag of base elements in fill_fe_values,
+  // and so the current object's update_each flag needs to be
+  // correct in case the current FESystem is a base element for another,
+  // higher-level FESystem itself.
   InternalData *data = new InternalData(this->n_base_elements());
-  data->update_each = update_each(flags) | update_once(flags);
+  data->update_each = requires_update_flags(flags);
 
   // get data objects from each of the base elements and store
   // them. do the creation of these objects in parallel as their
@@ -966,14 +962,19 @@ FESystem<dim,spacedim>::get_face_data (
 
 template <int dim, int spacedim>
 typename FiniteElement<dim,spacedim>::InternalDataBase *
-FESystem<dim,spacedim>::get_subface_data (
-  const UpdateFlags      flags_,
-  const Mapping<dim,spacedim>    &mapping,
-  const Quadrature<dim-1> &quadrature) const
-{
-  UpdateFlags flags = flags_;
+FESystem<dim,spacedim>::get_subface_data (const UpdateFlags      flags,
+                                          const Mapping<dim,spacedim>    &mapping,
+                                          const Quadrature<dim-1> &quadrature) const
+{
+  // create an internal data object and set the update flags we will need
+  // to deal with. the current object does not make use of these flags,
+  // but we need to nevertheless set them correctly since we look
+  // into the update_each flag of base elements in fill_fe_values,
+  // and so the current object's update_each flag needs to be
+  // correct in case the current FESystem is a base element for another,
+  // higher-level FESystem itself.
   InternalData *data = new InternalData(this->n_base_elements());
-  data->update_each = update_each(flags) | update_once(flags);
+  data->update_each = requires_update_flags(flags);
 
   // get data objects from each of the base elements and store
   // them. do the creation of these objects in parallel as their
index 2457767b9ba6149fede46df6435b3cf00be323ce..3bfb04426597c2e064bb032a8df39d2e6ca52a78 100644 (file)
@@ -3436,16 +3436,14 @@ template <int dim, int spacedim>
 UpdateFlags
 FEValuesBase<dim,spacedim>::compute_update_flags (const UpdateFlags update_flags) const
 {
-
-  // first find out which objects
-  // need to be recomputed on each
-  // cell we visit. this we have to
-  // ask the finite element and mapping.
-  // elements are first since they
-  // might require update in mapping
+  // first find out which objects need to be recomputed on each
+  // cell we visit. this we have to ask the finite element and mapping.
+  // elements are first since they might require update in mapping
+  //
+  // there is no need to iterate since mappings will never require
+  // the finite element to compute something for them
   UpdateFlags flags = update_flags
-                      | fe->update_once (update_flags)
-                      | fe->update_each (update_flags);
+                      | fe->requires_update_flags (update_flags);
   flags |= mapping->requires_update_flags (flags);
 
   return flags;

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.