]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove the first_cell/is_first_cell/current_update_flags() mechanism from FiniteElement. 1823/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 1 Nov 2015 15:56:00 +0000 (09:56 -0600)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 2 Nov 2015 08:33:53 +0000 (02:33 -0600)
What this patch does is address the issue that
`FiniteElement::get_data()` (as was previously the case with Mappings
as well) sets a flag in the returned `InternalData` object that
indicates that we are visiting a cell for the first time. `FEValues`
later sets this flag back once we're done with the first cell.

The way this is used is that the various finite element
implementations query `current_update_flags` in
`FE::fill_fe_*_values()`, which equals `update_each|update_once` if we
are on the first cell, and `update_each` on later cells. That means
that on the first cell we update fields in the output object that we
need to initialize only once, e.g., the values of shape functions on
quadrature points, but we won't do this again on later cells.

I don't particularly like this approach. It is opaque and took even me
a long time to reconstruct when I saw it (see #1305 for a
discussion). My preference is to do things we only need to do once in
`FE::get_data()`, and do things we need to do every time in
`FE::fill_fe_*_values()`. This kind of change was already implemented
for mappings in the patches references from #1305, and the current
pull request goes into the same direction.

The current PR only does one half of the necessary changes in order to
keep the patch reasonably self contained and readable. In particular,
what it really does is simply *always* do everything, by storing
`update_once|update_each` in `FE::InternalDataBase::update_each`. In
other words, it is the equivalent to telling finite element
implementations that we are always on the first cell. This requires
more work than in the previous state, although copying shape values
every time is arguably of small expense compared to transforming
gradients.

I will of course fix this up in upcoming patches. It requires
significant shuffling around. I'll open an issue in a minute
discussing the details. I'd simply like to propose this here already
as an incremental step forward, without making patch review too
onerous. The patch is mostly mechanical once you understand what it
does.

15 files changed:
include/deal.II/fe/fe.h
include/deal.II/fe/fe_face.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.cc
source/fe/fe_dgp_nonparametric.cc
source/fe/fe_face.cc
source/fe/fe_poly.cc
source/fe/fe_poly_tensor.cc
source/fe/fe_system.cc
source/fe/fe_values.cc

index c5a990ff67bb69c6e29d8052b4346826e45a7d31..add7dbc80c049771283b02f9e86a6527ecabc99e 100644 (file)
@@ -519,16 +519,6 @@ public:
      */
     virtual ~InternalDataBase ();
 
-    /**
-     * Values updated by the constructor or by reinit.
-     */
-    UpdateFlags          update_flags;
-
-    /**
-     * Values computed by constructor.
-     */
-    UpdateFlags          update_once;
-
     /**
      * A set of update flags specifying the kind of information that
      * an implementation of the FiniteElement interface needs to compute on
@@ -547,30 +537,10 @@ public:
      */
     UpdateFlags          update_each;
 
-    /**
-     * If <tt>first_cell==true</tt> this function returns @p update_flags,
-     * i.e. <tt>update_once|update_each</tt>. If <tt>first_cell==false</tt> it
-     * returns @p update_each.
-     */
-    UpdateFlags  current_update_flags() const;
-
-    /**
-     * Set the @p first_cell flag to @p false. Used by the @p FEValues class
-     * to indicate that we have already done the work on the first cell.
-     */
-    virtual void clear_first_cell ();
-
     /**
      * Return an estimate (in bytes) or the memory consumption of this object.
      */
     virtual std::size_t memory_consumption () const;
-
-  private:
-    /**
-     * Initially set to true, but reset to false when clear_first_cell()
-     * is called.
-     */
-    bool first_cell;
   };
 
 public:
index 352bc36399c8d9cc01b7c0ee1af9347ed45dc1ca..ad9cfa624bf412f13e75b45efee8101671ea2275 100644 (file)
@@ -250,20 +250,14 @@ protected:
     // generate a new data object and initialize some fields
     typename FiniteElement<1,spacedim>::InternalDataBase *data =
       new typename FiniteElement<1,spacedim>::InternalDataBase;
+    data->update_each = update_once(update_flags) | update_each(update_flags);  // FIX: only update_each required
 
-    // check what needs to be initialized only once and what on every
-    // cell/face/subface we visit
-    data->update_once = update_once(update_flags);
-    data->update_each = update_each(update_flags);
-    data->update_flags = data->update_once | data->update_each;
-
-    const UpdateFlags flags(data->update_flags);
     const unsigned int n_q_points = quadrature.size();
     AssertDimension(n_q_points, 1);
     (void)n_q_points;
 
     // No derivatives of this element are implemented.
-    if (flags & update_gradients || flags & update_hessians)
+    if (data->update_each & update_gradients || data->update_each & update_hessians)
       {
         Assert(false, ExcNotImplemented());
       }
index 79331d4e851b89e134d88309f456a75303a84619..c73b58117cd0cc931d880dbfc9eeec0105e379b3 100644 (file)
@@ -222,16 +222,8 @@ protected:
     // generate a new data object and
     // initialize some fields
     InternalData *data = new InternalData;
+    data->update_each = update_each(update_flags) | update_once(update_flags);  // FIX: only update_each required
 
-    // check what needs to be
-    // initialized only once and what
-    // on every cell/face/subface we
-    // visit
-    data->update_once = update_once(update_flags);
-    data->update_each = update_each(update_flags);
-    data->update_flags = data->update_once | data->update_each;
-
-    const UpdateFlags flags(data->update_flags);
     const unsigned int n_q_points = quadrature.size();
 
     // some scratch arrays
@@ -244,28 +236,28 @@ protected:
     // initialize fields only if really
     // necessary. otherwise, don't
     // allocate memory
-    if (flags & update_values)
+    if (update_flags & update_values)
       {
         values.resize (this->dofs_per_cell);
         data->shape_values.resize (this->dofs_per_cell,
                                    std::vector<double> (n_q_points));
       }
 
-    if (flags & update_gradients)
+    if (update_flags & update_gradients)
       {
         grads.resize (this->dofs_per_cell);
         data->shape_gradients.resize (this->dofs_per_cell,
                                       std::vector<Tensor<1,dim> > (n_q_points));
       }
 
-    if (flags & update_hessians)
+    if (update_flags & update_hessians)
       {
         grad_grads.resize (this->dofs_per_cell);
         data->shape_hessians.resize (this->dofs_per_cell,
                                      std::vector<Tensor<2,dim> > (n_q_points));
       }
 
-    if (flags & update_3rd_derivatives)
+    if (update_flags & update_3rd_derivatives)
       {
         third_derivatives.resize (this->dofs_per_cell);
         data->shape_3rd_derivatives.resize (this->dofs_per_cell,
@@ -279,8 +271,8 @@ protected:
     // unit cell, and need to be
     // transformed when visiting an
     // actual cell
-    if (flags & (update_values | update_gradients
-                 | update_hessians | update_3rd_derivatives) )
+    if (update_flags & (update_values | update_gradients
+                        | update_hessians | update_3rd_derivatives) )
       for (unsigned int i=0; i<n_q_points; ++i)
         {
           poly_space.compute(quadrature.point(i),
@@ -288,19 +280,19 @@ protected:
                              third_derivatives,
                              fourth_derivatives);
 
-          if (flags & update_values)
+          if (update_flags & update_values)
             for (unsigned int k=0; k<this->dofs_per_cell; ++k)
               data->shape_values[k][i] = values[k];
 
-          if (flags & update_gradients)
+          if (update_flags & update_gradients)
             for (unsigned int k=0; k<this->dofs_per_cell; ++k)
               data->shape_gradients[k][i] = grads[k];
 
-          if (flags & update_hessians)
+          if (update_flags & update_hessians)
             for (unsigned int k=0; k<this->dofs_per_cell; ++k)
               data->shape_hessians[k][i] = grad_grads[k];
 
-          if (flags & update_3rd_derivatives)
+          if (update_flags & update_3rd_derivatives)
             for (unsigned int k=0; k<this->dofs_per_cell; ++k)
               data->shape_3rd_derivatives[k][i] = third_derivatives[k];
         }
index bdaa97d2e19109b4f14c5e4df154fc038796d272..74519fda038825cfa346228330ec300df38994f1 100644 (file)
@@ -240,7 +240,7 @@ fill_fe_values (const Mapping<dim,spacedim>                                  &ma
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.current_update_flags());
+  const UpdateFlags flags(fe_data.update_each);
 
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
@@ -312,7 +312,7 @@ fill_fe_face_values (const Mapping<dim,spacedim>
                                                 cell->face_rotation(face),
                                                 quadrature.size());
 
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+  const UpdateFlags flags(fe_data.update_each);
 
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
@@ -388,7 +388,7 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
                                                    quadrature.size(),
                                                    cell->subface_case(face));
 
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+  const UpdateFlags flags(fe_data.update_each);
 
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
index 697d5e801ce4096e5f4dd407ce1edc5f269fb808..31d85ea509289fd98d66ed770f99844f4786d664 100644 (file)
@@ -99,16 +99,8 @@ protected:
     // generate a new data object and
     // initialize some fields
     InternalData *data = new InternalData;
+    data->update_each = update_once(update_flags) | update_each(update_flags);  // FIX: only update_each required
 
-    // check what needs to be
-    // initialized only once and what
-    // on every cell/face/subface we
-    // visit
-    data->update_once = update_once(update_flags);
-    data->update_each = update_each(update_flags);
-    data->update_flags = data->update_once | data->update_each;
-
-    const UpdateFlags flags(data->update_flags);
     const unsigned int n_q_points = quadrature.size();
 
     // some scratch arrays
@@ -121,7 +113,7 @@ protected:
     // initialize fields only if really
     // necessary. otherwise, don't
     // allocate memory
-    if (flags & update_values)
+    if (data->update_each & update_values)
       {
         values.resize (poly_space.n());
         data->shape_values.resize (poly_space.n(),
@@ -139,7 +131,7 @@ protected:
       }
     // No derivatives of this element
     // are implemented.
-    if (flags & update_gradients || flags & update_hessians)
+    if (data->update_each & update_gradients || data->update_each & update_hessians)
       {
         Assert(false, ExcNotImplemented());
       }
index 786f0bd3e935d3c683df135b53d79d4cdb424d65..5497997a91a521326985832dd9f6962446607fd9 100644 (file)
@@ -118,9 +118,7 @@ fill_fe_face_values (const Mapping<dim,spacedim> &,
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
-
-  if (flags & update_values)
+  if (fe_data.update_each & update_values)
     for (unsigned int i=0; i<quadrature.size(); ++i)
       {
         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
@@ -136,7 +134,10 @@ fill_fe_face_values (const Mapping<dim,spacedim> &,
                 for (unsigned int k=0; k<this->dofs_per_quad; ++k)
                   output_data.shape_values(foffset+k,i) = fe_data.shape_values[k+this->first_face_quad_index][i];
               }
+
+            // fall through...
           }
+
           case 2:
           {
             // Fill data for line shape functions
@@ -150,7 +151,10 @@ fill_fe_face_values (const Mapping<dim,spacedim> &,
                         = fe_data.shape_values[k+(line*this->dofs_per_line)+this->first_face_line_index][i];
                   }
               }
+
+            // fall through...
           }
+
           case 1:
           {
             // Fill data for vertex shape functions
@@ -185,12 +189,10 @@ fill_fe_subface_values (const Mapping<dim,spacedim> &,
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
-
   const unsigned int foffset = fe_data.shape_values.size() * face;
   const unsigned int offset = subface*quadrature.size();
 
-  if (flags & update_values)
+  if (fe_data.update_each & update_values)
     for (unsigned int i=0; i<quadrature.size(); ++i)
       {
         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
@@ -199,8 +201,8 @@ fill_fe_subface_values (const Mapping<dim,spacedim> &,
           output_data.shape_values(foffset+k,i) = fe_data.shape_values[k][i+offset];
       }
 
-  Assert (!(flags & update_gradients), ExcNotImplemented());
-  Assert (!(flags & update_hessians), ExcNotImplemented());
+  Assert (!(fe_data.update_each & update_gradients), ExcNotImplemented());
+  Assert (!(fe_data.update_each & update_hessians), ExcNotImplemented());
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 78bb7dbfc12905a55c9585b02418cdb1390b8eb9..fefaba3f587b7e149a8fe24014c2a4c4f922f3ff 100644 (file)
@@ -188,16 +188,8 @@ protected:
     // generate a new data object and
     // initialize some fields
     InternalData *data = new InternalData;
+    data->update_each = update_each(update_flags) | update_once(update_flags);  // FIX: only update_each required
 
-    // check what needs to be
-    // initialized only once and what
-    // on every cell/face/subface we
-    // visit
-    data->update_once = update_once(update_flags);
-    data->update_each = update_each(update_flags);
-    data->update_flags = data->update_once | data->update_each;
-
-    const UpdateFlags flags(data->update_flags);
     const unsigned int n_q_points = quadrature.size();
 
     // some scratch arrays
@@ -207,13 +199,13 @@ protected:
     std::vector<Tensor<4,dim> > third_derivatives(0);
     std::vector<Tensor<5,dim> > fourth_derivatives(0);
 
-    if (flags & (update_values | update_gradients | update_hessians) )
+    if (update_flags & (update_values | update_gradients | update_hessians) )
       data->sign_change.resize (this->dofs_per_cell);
 
     // initialize fields only if really
     // necessary. otherwise, don't
     // allocate memory
-    if (flags & update_values)
+    if (update_flags & update_values)
       {
         values.resize (this->dofs_per_cell);
         data->shape_values.resize (this->dofs_per_cell);
@@ -223,7 +215,7 @@ protected:
           data->transformed_shape_values.resize (n_q_points);
       }
 
-    if (flags & update_gradients)
+    if (update_flags & update_gradients)
       {
         grads.resize (this->dofs_per_cell);
         data->shape_grads.resize (this->dofs_per_cell);
@@ -241,7 +233,7 @@ protected:
           data->untransformed_shape_grads.resize(n_q_points);
       }
 
-    if (flags & update_hessians)
+    if (update_flags & update_hessians)
       {
         grad_grads.resize (this->dofs_per_cell);
         data->shape_grad_grads.resize (this->dofs_per_cell);
@@ -259,7 +251,7 @@ protected:
     // node values N_i holds
     // N_i(v_j)=\delta_ij for all basis
     // functions v_j
-    if (flags & (update_values | update_gradients))
+    if (update_flags & (update_values | update_gradients))
       for (unsigned int k=0; k<n_q_points; ++k)
         {
           poly_space.compute(quadrature.point(k),
@@ -267,7 +259,7 @@ protected:
                              third_derivatives,
                              fourth_derivatives);
 
-          if (flags & update_values)
+          if (update_flags & update_values)
             {
               if (inverse_node_matrix.n_cols() == 0)
                 for (unsigned int i=0; i<this->dofs_per_cell; ++i)
@@ -282,7 +274,7 @@ protected:
                   }
             }
 
-          if (flags & update_gradients)
+          if (update_flags & update_gradients)
             {
               if (inverse_node_matrix.n_cols() == 0)
                 for (unsigned int i=0; i<this->dofs_per_cell; ++i)
@@ -297,7 +289,7 @@ protected:
                   }
             }
 
-          if (flags & update_hessians)
+          if (update_flags & update_hessians)
             {
               if (inverse_node_matrix.n_cols() == 0)
                 for (unsigned int i=0; i<this->dofs_per_cell; ++i)
index 999de4106a61449b185f0584a32743eeafc8a4b6..ff256ef3a6f6651531d9d3795574451567132c75 100644 (file)
@@ -1145,15 +1145,6 @@ private:
     internal::FEValues::FiniteElementRelatedData<dim,spacedim> &
     get_fe_output_object (const unsigned int base_no) const;
 
-    /**
-     * Set the @p first_cell flag to @p false. Used by the @p FEValues class
-     * to indicate that we have already done the work on the first cell.
-     *
-     * In addition to calling the respective function of the base class, this
-     * function also calls the functions of the sub-data objects.
-     */
-    virtual void clear_first_cell ();
-
   private:
 
     /**
index 4bf0d9cc123d00e26e2670adcfe76a31223bb594..8107ac7e5fcd2f7421a6c239f2ce8f55b77b5792 100644 (file)
@@ -36,10 +36,7 @@ DEAL_II_NAMESPACE_OPEN
 
 template <int dim, int spacedim>
 FiniteElement<dim, spacedim>::InternalDataBase::InternalDataBase ():
-  update_flags(update_default),
-  update_once(update_default),
-  update_each(update_default),
-  first_cell(true)
+  update_each(update_default)
 {}
 
 
@@ -59,32 +56,6 @@ FiniteElement<dim, spacedim>::InternalDataBase::memory_consumption () const
 
 
 
-template <int dim, int spacedim>
-UpdateFlags
-FiniteElement<dim,spacedim>::InternalDataBase::current_update_flags () const
-{
-  if (first_cell)
-    {
-      Assert (update_flags==(update_once|update_each),
-              ExcInternalError());
-      return update_flags;
-    }
-  else
-    return update_each;
-}
-
-
-
-template <int dim, int spacedim>
-void
-FiniteElement<dim,spacedim>::InternalDataBase::clear_first_cell ()
-{
-  first_cell = false;
-}
-
-
-
-
 template <int dim, int spacedim>
 FiniteElement<dim,spacedim>::
 FiniteElement (const FiniteElementData<dim> &fe_data,
index 4588b0c4b8eaf6b8425d1a5a7be5c1f450c14577..00045beebd62696f0a86f04ff959fac410fe0f3b 100644 (file)
@@ -255,13 +255,7 @@ FE_DGPNonparametric<dim,spacedim>::get_data (
 {
   // generate a new data object
   typename FiniteElement<dim,spacedim>::InternalDataBase *data = new typename FiniteElement<dim,spacedim>::InternalDataBase;
-  // check what needs to be
-  // initialized only once and what
-  // on every cell/face/subface we
-  // visit
-  data->update_once = update_once(update_flags);
-  data->update_each = update_each(update_flags);
-  data->update_flags = data->update_once | data->update_each;
+  data->update_each = update_once(update_flags) | update_each(update_flags);   // FIX: only update_each required
 
   // other than that, there is nothing we can add here as discussed
   // in the general documentation of this class
@@ -287,18 +281,17 @@ fill_fe_values (const Mapping<dim,spacedim> &,
                 internal::FEValues::FiniteElementRelatedData<dim,spacedim> &output_data,
                 const CellSimilarity::Similarity /*cell_similarity*/) const
 {
-  const UpdateFlags flags(fe_data.current_update_flags());
-  Assert (flags & update_quadrature_points, ExcInternalError());
+  Assert (fe_data.update_each & update_quadrature_points, ExcInternalError());
 
   const unsigned int n_q_points = mapping_data.quadrature_points.size();
 
-  std::vector<double> values(flags & update_values ? this->dofs_per_cell : 0);
-  std::vector<Tensor<1,dim> > grads(flags & update_gradients ? this->dofs_per_cell : 0);
-  std::vector<Tensor<2,dim> > grad_grads(flags & update_hessians ? this->dofs_per_cell : 0);
+  std::vector<double> values(fe_data.update_each & update_values ? this->dofs_per_cell : 0);
+  std::vector<Tensor<1,dim> > grads(fe_data.update_each & update_gradients ? this->dofs_per_cell : 0);
+  std::vector<Tensor<2,dim> > grad_grads(fe_data.update_each & update_hessians ? this->dofs_per_cell : 0);
   std::vector<Tensor<3,dim> > empty_vector_of_3rd_order_tensors;
   std::vector<Tensor<4,dim> > empty_vector_of_4th_order_tensors;
 
-  if (flags & (update_values | update_gradients))
+  if (fe_data.update_each & (update_values | update_gradients))
     for (unsigned int i=0; i<n_q_points; ++i)
       {
         polynomial_space.compute(mapping_data.quadrature_points[i],
@@ -307,11 +300,11 @@ fill_fe_values (const Mapping<dim,spacedim> &,
                                  empty_vector_of_4th_order_tensors);
         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
           {
-            if (flags & update_values)
+            if (fe_data.update_each & update_values)
               output_data.shape_values[k][i] = values[k];
-            if (flags & update_gradients)
+            if (fe_data.update_each & update_gradients)
               output_data.shape_gradients[k][i] = grads[k];
-            if (flags & update_hessians)
+            if (fe_data.update_each & update_hessians)
               output_data.shape_hessians[k][i] = grad_grads[k];
           }
       }
@@ -331,18 +324,17 @@ fill_fe_face_values (const Mapping<dim,spacedim> &,
                      const internal::FEValues::MappingRelatedData<dim,spacedim> &mapping_data,
                      internal::FEValues::FiniteElementRelatedData<dim,spacedim> &output_data) const
 {
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
-  Assert (flags & update_quadrature_points, ExcInternalError());
+  Assert (fe_data.update_each & update_quadrature_points, ExcInternalError());
 
   const unsigned int n_q_points = mapping_data.quadrature_points.size();
 
-  std::vector<double> values(flags & update_values ? this->dofs_per_cell : 0);
-  std::vector<Tensor<1,dim> > grads(flags & update_gradients ? this->dofs_per_cell : 0);
-  std::vector<Tensor<2,dim> > grad_grads(flags & update_hessians ? this->dofs_per_cell : 0);
+  std::vector<double> values(fe_data.update_each & update_values ? this->dofs_per_cell : 0);
+  std::vector<Tensor<1,dim> > grads(fe_data.update_each & update_gradients ? this->dofs_per_cell : 0);
+  std::vector<Tensor<2,dim> > grad_grads(fe_data.update_each & update_hessians ? this->dofs_per_cell : 0);
   std::vector<Tensor<3,dim> > empty_vector_of_3rd_order_tensors;
   std::vector<Tensor<4,dim> > empty_vector_of_4th_order_tensors;
 
-  if (flags & (update_values | update_gradients))
+  if (fe_data.update_each & (update_values | update_gradients))
     for (unsigned int i=0; i<n_q_points; ++i)
       {
         polynomial_space.compute(mapping_data.quadrature_points[i],
@@ -351,11 +343,11 @@ fill_fe_face_values (const Mapping<dim,spacedim> &,
                                  empty_vector_of_4th_order_tensors);
         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
           {
-            if (flags & update_values)
+            if (fe_data.update_each & update_values)
               output_data.shape_values[k][i] = values[k];
-            if (flags & update_gradients)
+            if (fe_data.update_each & update_gradients)
               output_data.shape_gradients[k][i] = grads[k];
-            if (flags & update_hessians)
+            if (fe_data.update_each & update_hessians)
               output_data.shape_hessians[k][i] = grad_grads[k];
           }
       }
@@ -376,18 +368,17 @@ fill_fe_subface_values (const Mapping<dim,spacedim> &,
                         const internal::FEValues::MappingRelatedData<dim,spacedim> &mapping_data,
                         internal::FEValues::FiniteElementRelatedData<dim,spacedim> &output_data) const
 {
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
-  Assert (flags & update_quadrature_points, ExcInternalError());
+  Assert (fe_data.update_each & update_quadrature_points, ExcInternalError());
 
   const unsigned int n_q_points = mapping_data.quadrature_points.size();
 
-  std::vector<double> values(flags & update_values ? this->dofs_per_cell : 0);
-  std::vector<Tensor<1,dim> > grads(flags & update_gradients ? this->dofs_per_cell : 0);
-  std::vector<Tensor<2,dim> > grad_grads(flags & update_hessians ? this->dofs_per_cell : 0);
+  std::vector<double> values(fe_data.update_each & update_values ? this->dofs_per_cell : 0);
+  std::vector<Tensor<1,dim> > grads(fe_data.update_each & update_gradients ? this->dofs_per_cell : 0);
+  std::vector<Tensor<2,dim> > grad_grads(fe_data.update_each & update_hessians ? this->dofs_per_cell : 0);
   std::vector<Tensor<3,dim> > empty_vector_of_3rd_order_tensors;
   std::vector<Tensor<4,dim> > empty_vector_of_4th_order_tensors;
 
-  if (flags & (update_values | update_gradients))
+  if (fe_data.update_each & (update_values | update_gradients))
     for (unsigned int i=0; i<n_q_points; ++i)
       {
         polynomial_space.compute(mapping_data.quadrature_points[i],
@@ -396,11 +387,11 @@ fill_fe_subface_values (const Mapping<dim,spacedim> &,
                                  empty_vector_of_4th_order_tensors);
         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
           {
-            if (flags & update_values)
+            if (fe_data.update_each & update_values)
               output_data.shape_values[k][i] = values[k];
-            if (flags & update_gradients)
+            if (fe_data.update_each & update_gradients)
               output_data.shape_gradients[k][i] = grads[k];
-            if (flags & update_hessians)
+            if (fe_data.update_each & update_hessians)
               output_data.shape_hessians[k][i] = grad_grads[k];
           }
       }
index 9732fab950b1810a41dbbc097d1c2b75fb74d7cc..7f3c8b94584a9b38d6a0b2cbe1d8159e8ac6a3ef 100644 (file)
@@ -494,10 +494,8 @@ fill_fe_face_values (const Mapping<1,spacedim> &,
                      const internal::FEValues::MappingRelatedData<1,spacedim> &,
                      internal::FEValues::FiniteElementRelatedData<1,spacedim> &output_data) const
 {
-  const UpdateFlags flags(fedata.update_once | fedata.update_each);
-
   const unsigned int foffset = face;
-  if (flags & update_values)
+  if (fedata.update_each & update_values)
     {
       for (unsigned int k=0; k<this->dofs_per_cell; ++k)
         output_data.shape_values(k,0) = 0.;
index 7981da494635b73e89d092b52c0b3b0fde572a1d..c4c49b6ddd5a0637d75d444ac00d06ae7490dd25 100644 (file)
@@ -47,21 +47,19 @@ fill_fe_values (const Mapping<1,2>                                &mapping,
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.current_update_flags());
-
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
-      if (flags & update_values)
+      if (fe_data.update_each & update_values)
         for (unsigned int i=0; i<quadrature.size(); ++i)
           output_data.shape_values(k,i) = fe_data.shape_values[k][i];
 
-      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
         mapping.transform (fe_data.shape_gradients[k],
                            mapping_covariant,
                            mapping_internal,
                            output_data.shape_gradients[k]);
 
-      if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_hessians[k],
                              mapping_covariant_gradient,
@@ -75,7 +73,7 @@ fill_fe_values (const Mapping<1,2>                                &mapping,
                 * output_data.shape_gradients[k][i][j];
         }
 
-      if (flags & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_3rd_derivatives[k],
                              mapping_covariant_hessian,
@@ -107,21 +105,19 @@ fill_fe_values (const Mapping<2,3>                                &mapping,
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.current_update_flags());
-
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
-      if (flags & update_values)
+      if (fe_data.update_each & update_values)
         for (unsigned int i=0; i<quadrature.size(); ++i)
           output_data.shape_values(k,i) = fe_data.shape_values[k][i];
 
-      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
         mapping.transform (fe_data.shape_gradients[k],
                            mapping_covariant,
                            mapping_internal,
                            output_data.shape_gradients[k]);
 
-      if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_hessians[k],
                              mapping_covariant_gradient,
@@ -135,7 +131,7 @@ fill_fe_values (const Mapping<2,3>                                &mapping,
                 * output_data.shape_gradients[k][i][j];
         }
 
-      if (flags & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_3rd_derivatives[k],
                              mapping_covariant_hessian,
@@ -168,21 +164,19 @@ fill_fe_values (const Mapping<1,2>                                &mapping,
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.current_update_flags());
-
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
-      if (flags & update_values)
+      if (fe_data.update_each & update_values)
         for (unsigned int i=0; i<quadrature.size(); ++i)
           output_data.shape_values(k,i) = fe_data.shape_values[k][i];
 
-      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
         mapping.transform (fe_data.shape_gradients[k],
                            mapping_covariant,
                            mapping_internal,
                            output_data.shape_gradients[k]);
 
-      if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_hessians[k],
                              mapping_covariant_gradient,
@@ -196,7 +190,7 @@ fill_fe_values (const Mapping<1,2>                                &mapping,
                 * output_data.shape_gradients[k][i][j];
         }
 
-      if (flags & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_3rd_derivatives[k],
                              mapping_covariant_hessian,
@@ -224,22 +218,20 @@ fill_fe_values (const Mapping<2,3>                                &mapping,
   Assert (dynamic_cast<const InternalData *> (&fedata) != 0, ExcInternalError());
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
-  const UpdateFlags flags(fe_data.current_update_flags());
-
   for (unsigned int k=0; k<this->dofs_per_cell; ++k)
     {
-      if (flags & update_values)
+      if (fe_data.update_each & update_values)
         for (unsigned int i=0; i<quadrature.size(); ++i)
           output_data.shape_values(k,i) = fe_data.shape_values[k][i];
 
 
-      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
         mapping.transform (fe_data.shape_gradients[k],
                            mapping_covariant,
                            mapping_internal,
                            output_data.shape_gradients[k]);
 
-      if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_hessians[k],
                              mapping_covariant_gradient,
@@ -253,7 +245,7 @@ fill_fe_values (const Mapping<2,3>                                &mapping,
                 * output_data.shape_gradients[k][i][j];
         }
 
-      if (flags & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
+      if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
         {
           mapping.transform (fe_data.shape_3rd_derivatives[k],
                              mapping_covariant_hessian,
index 8939e9d05d0329d2433fbcd99b18a86667e147f1..cd116087e8922c3dc505a0754ce84471f5f846ad 100644 (file)
@@ -307,11 +307,10 @@ fill_fe_values (const Mapping<dim,spacedim>                                  &ma
   const InternalData &fe_data = static_cast<const InternalData &> (fedata);
 
   const unsigned int n_q_points = quadrature.size();
-  const UpdateFlags flags(fe_data.current_update_flags());
 
-  Assert(!(flags & update_values) || fe_data.shape_values.size() == this->dofs_per_cell,
+  Assert(!(fe_data.update_each & update_values) || fe_data.shape_values.size() == this->dofs_per_cell,
          ExcDimensionMismatch(fe_data.shape_values.size(), this->dofs_per_cell));
-  Assert(!(flags & update_values) || fe_data.shape_values[0].size() == n_q_points,
+  Assert(!(fe_data.update_each & update_values) || fe_data.shape_values[0].size() == n_q_points,
          ExcDimensionMismatch(fe_data.shape_values[0].size(), n_q_points));
 
   // Create table with sign changes, due to the special structure of the RT elements.
@@ -338,7 +337,7 @@ fill_fe_values (const Mapping<dim,spacedim>                                  &ma
       // the previous one; or, even if it is a translation, if we use mappings
       // other than the standard mappings that require us to recompute values
       // and derivatives because of possible sign changes
-      if (flags & update_values &&
+      if (fe_data.update_each & update_values &&
           ((cell_similarity != CellSimilarity::translation)
            ||
            ((mapping_type == mapping_piola) || (mapping_type == mapping_raviart_thomas)
@@ -404,7 +403,7 @@ fill_fe_values (const Mapping<dim,spacedim>                                  &ma
         }
 
       // update gradients. apply the same logic as above
-      if (flags & update_gradients
+      if (fe_data.update_each & update_gradients
           &&
           ((cell_similarity != CellSimilarity::translation)
            ||
@@ -532,7 +531,7 @@ fill_fe_values (const Mapping<dim,spacedim>                                  &ma
         }
 
       // update hessians. apply the same logic as above
-      if (flags & update_hessians
+      if (fe_data.update_each & update_hessians
           &&
           ((cell_similarity != CellSimilarity::translation)
            ||
@@ -742,7 +741,7 @@ fill_fe_values (const Mapping<dim,spacedim>                                  &ma
         }
 
       // third derivatives are not implemented
-      if (flags & update_3rd_derivatives
+      if (fe_data.update_each & update_3rd_derivatives
           &&
           ((cell_similarity != CellSimilarity::translation)
            ||
@@ -788,8 +787,6 @@ fill_fe_face_values (const Mapping<dim,spacedim>
                                                 cell->face_rotation(face),
                                                 n_q_points);
 
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
-
 //TODO: Size assertions
 
 // Create table with sign changes, due to the special structure of the RT elements.
@@ -810,7 +807,7 @@ fill_fe_face_values (const Mapping<dim,spacedim>
       const unsigned int first = output_data.shape_function_to_row_table[i * this->n_components() +
                                  this->get_nonzero_components(i).first_selected_component()];
 
-      if (flags & update_values)
+      if (fe_data.update_each & update_values)
         {
           switch (mapping_type)
             {
@@ -880,7 +877,7 @@ fill_fe_face_values (const Mapping<dim,spacedim>
             }
         }
 
-      if (flags & update_gradients)
+      if (fe_data.update_each & update_gradients)
         {
           VectorSlice< std::vector<Tensor<2,spacedim> > >
           transformed_shape_grads (fe_data.transformed_shape_grads, offset, n_q_points);
@@ -1004,7 +1001,7 @@ fill_fe_face_values (const Mapping<dim,spacedim>
             }
         }
 
-      if (flags & update_hessians)
+      if (fe_data.update_each & update_hessians)
 
         {
 
@@ -1208,7 +1205,7 @@ fill_fe_face_values (const Mapping<dim,spacedim>
         }
 
       // third derivatives are not implemented
-      if (flags & update_3rd_derivatives)
+      if (fe_data.update_each & update_3rd_derivatives)
         {
           Assert(false, ExcNotImplemented())
         }
@@ -1251,8 +1248,6 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
                                                    n_q_points,
                                                    cell->subface_case(face));
 
-  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
-
 //   Assert(mapping_type == independent
 //       || ( mapping_type == independent_on_cartesian
 //            && dynamic_cast<const MappingCartesian<dim>*>(&mapping) != 0),
@@ -1276,7 +1271,7 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
       const unsigned int first = output_data.shape_function_to_row_table[i * this->n_components() +
                                  this->get_nonzero_components(i).first_selected_component()];
 
-      if (flags & update_values)
+      if (fe_data.update_each & update_values)
         {
           switch (mapping_type)
             {
@@ -1343,7 +1338,7 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
             }
         }
 
-      if (flags & update_gradients)
+      if (fe_data.update_each & update_gradients)
         {
           VectorSlice< std::vector<Tensor<2, spacedim > > >
           transformed_shape_grads (fe_data.transformed_shape_grads, offset, n_q_points);
@@ -1469,7 +1464,7 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
             }
         }
 
-      if (flags & update_hessians)
+      if (fe_data.update_each & update_hessians)
 
         {
 
@@ -1672,7 +1667,7 @@ fill_fe_subface_values (const Mapping<dim,spacedim>
         }
 
       // third derivatives are not implemented
-      if (flags & update_3rd_derivatives)
+      if (fe_data.update_each & update_3rd_derivatives)
         {
           Assert(false, ExcNotImplemented())
         }
index afd837eca0bb05816e481768a64c7a4796f55bac..95f097d15620023327414d32a41991f2cade4ebf 100644 (file)
@@ -104,20 +104,6 @@ InternalData::get_fe_output_object (const unsigned int base_no) const
 
 
 
-template <int dim, int spacedim>
-void
-FESystem<dim,spacedim>::InternalData::clear_first_cell ()
-{
-  // call respective function of base
-  // class
-  FiniteElement<dim,spacedim>::InternalDataBase::clear_first_cell ();
-  // then the functions of all the
-  // sub-objects
-  for (unsigned int i=0; i<base_fe_datas.size(); ++i)
-    base_fe_datas[i]->clear_first_cell ();
-}
-
-
 /* ---------------------------------- FESystem ------------------- */
 
 
@@ -893,10 +879,7 @@ FESystem<dim,spacedim>::get_data (const UpdateFlags      flags_,
 {
   UpdateFlags flags = flags_;
   InternalData *data = new InternalData(this->n_base_elements());
-
-  data->update_once = update_once (flags);
-  data->update_each = update_each (flags);
-  flags = data->update_once | data->update_each;
+  data->update_each = update_each(flags) | update_once(flags);
 
   // get data objects from each of the base elements and store
   // them. do the creation of these objects in parallel as their
@@ -922,13 +905,11 @@ FESystem<dim,spacedim>::get_data (const UpdateFlags      flags_,
       internal::FEValues::FiniteElementRelatedData<dim,spacedim> &base_fe_output_object
         = data->get_fe_output_object(base_no);
       base_fe_output_object.initialize (quadrature.size(), base_element(base_no),
-                                        flags | base_fe_data->update_flags);
+                                        flags | base_fe_data->update_each);
 
       // then store the pointer to the base internal object
       data->set_fe_data(base_no, base_fe_data);
     }
-  data->update_flags = data->update_once |
-                       data->update_each;
   return data;
 }
 
@@ -944,10 +925,7 @@ FESystem<dim,spacedim>::get_face_data (
 {
   UpdateFlags flags = flags_;
   InternalData *data = new InternalData(this->n_base_elements());
-
-  data->update_once = update_once (flags);
-  data->update_each = update_each (flags);
-  flags = data->update_once | data->update_each;
+  data->update_each = update_each(flags) | update_once(flags);
 
   // get data objects from each of the base elements and store
   // them. do the creation of these objects in parallel as their
@@ -973,13 +951,11 @@ FESystem<dim,spacedim>::get_face_data (
       internal::FEValues::FiniteElementRelatedData<dim,spacedim> &base_fe_output_object
         = data->get_fe_output_object(base_no);
       base_fe_output_object.initialize (quadrature.size(), base_element(base_no),
-                                        flags | base_fe_data->update_flags);
+                                        flags | base_fe_data->update_each);
 
       // then store the pointer to the base internal object
       data->set_fe_data(base_no, base_fe_data);
     }
-  data->update_flags = data->update_once |
-                       data->update_each;
   return data;
 }
 
@@ -997,10 +973,7 @@ FESystem<dim,spacedim>::get_subface_data (
 {
   UpdateFlags flags = flags_;
   InternalData *data = new InternalData(this->n_base_elements());
-
-  data->update_once = update_once (flags);
-  data->update_each = update_each (flags);
-  flags = data->update_once | data->update_each;
+  data->update_each = update_each(flags) | update_once(flags);
 
   // get data objects from each of the base elements and store
   // them. do the creation of these objects in parallel as their
@@ -1026,13 +999,11 @@ FESystem<dim,spacedim>::get_subface_data (
       internal::FEValues::FiniteElementRelatedData<dim,spacedim> &base_fe_output_object
         = data->get_fe_output_object(base_no);
       base_fe_output_object.initialize (quadrature.size(), base_element(base_no),
-                                        flags | base_fe_data->update_flags);
+                                        flags | base_fe_data->update_each);
 
       // then store the pointer to the base internal object
       data->set_fe_data(base_no, base_fe_data);
     }
-  data->update_flags = data->update_once |
-                       data->update_each;
   return data;
 }
 
@@ -1183,9 +1154,7 @@ compute_fill_one_base (const Mapping<dim,spacedim>                      &mapping
   // all shape functions of the composed element, and here only treat
   // those shape functions that belong to a given base element
   //TODO: Introduce the needed table and loop only over base element shape functions. This here is not efficient at all AND very bad style
-  const UpdateFlags base_flags(dim_1==dim ?
-                               base_fe_data.current_update_flags() :
-                               base_fe_data.update_flags);
+  const UpdateFlags base_flags = base_fe_data.update_each;
 
   // if the current cell is just a translation of the previous one,
   // the underlying data has not changed, and we don't even need to
@@ -1270,9 +1239,7 @@ compute_fill (const Mapping<dim,spacedim>                      &mapping,
   // (fill_fe_values) or dim_1==dim-1
   // (fill_fe_(sub)face_values)
   Assert(dim_1==dim || dim_1==dim-1, ExcInternalError());
-  const UpdateFlags flags(dim_1==dim ?
-                          fe_data.current_update_flags() :
-                          fe_data.update_flags);
+  const UpdateFlags flags = fe_data.update_each;
 
 
   if (flags & (update_values | update_gradients
index f5341ebee6c6b58b1b11c817666f2f30bdf75cc7..2457767b9ba6149fede46df6435b3cf00be323ce 100644 (file)
@@ -3768,8 +3768,6 @@ void FEValues<dim,spacedim>::do_reinit ()
                                 this->mapping_output,
                                 this->finite_element_output,
                                 this->cell_similarity);
-
-  this->fe_data->clear_first_cell ();
 }
 
 
@@ -3974,8 +3972,6 @@ void FEFaceValues<dim,spacedim>::do_reinit (const unsigned int face_no)
                                      *this->fe_data,
                                      this->mapping_output,
                                      this->finite_element_output);
-
-  this->fe_data->clear_first_cell ();
 }
 
 
@@ -4217,8 +4213,6 @@ void FESubfaceValues<dim,spacedim>::do_reinit (const unsigned int face_no,
                                         *this->fe_data,
                                         this->mapping_output,
                                         this->finite_element_output);
-
-  this->fe_data->clear_first_cell ();
 }
 
 

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.