]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Instead of setting or clearing the first_cell flag at mostly arbitrary places, do...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 21 Oct 2002 14:19:54 +0000 (14:19 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 21 Oct 2002 14:19:54 +0000 (14:19 +0000)
git-svn-id: https://svn.dealii.org/trunk@6693 0785d39b-7218-0410-832d-ea1e28bc413d

12 files changed:
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/include/fe/mapping.h
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_dgp.cc
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc
deal.II/deal.II/source/fe/fe_dgq.cc
deal.II/deal.II/source/fe/fe_nedelec.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_system.cc
deal.II/deal.II/source/fe/fe_values.cc
deal.II/deal.II/source/fe/mapping_cartesian.cc
deal.II/deal.II/source/fe/mapping_q1.cc

index 2d01d61b33c2ec87f06287779ee69603ed561fba..da53e2df548d96a7b1a6574ff1ac059a79e2e9a5 100644 (file)
@@ -376,6 +376,12 @@ class FESystem : public FiniteElement<dim>
                                      * Implementation of the same
                                      * function in
                                      * @ref{FiniteElement}.
+                                     *
+                                     * Passes on control to
+                                     * @p{compute_fill} that does the
+                                     * work for all three
+                                     * @p{fill_fe*_values}
+                                     * functions.
                                      */
     virtual void
     fill_fe_values (const Mapping<dim>                      &mapping,
@@ -389,6 +395,11 @@ class FESystem : public FiniteElement<dim>
                                      * Implementation of the same
                                      * function in
                                      * @ref{FiniteElement}.
+                                     *
+                                     * Passes on control to
+                                     * @p{compute_fill} that does the
+                                     * work for all three
+                                     * @p{fill_fe*_values} functions.
                                      */    
     virtual void
     fill_fe_face_values (const Mapping<dim>                   &mapping,
@@ -403,6 +414,11 @@ class FESystem : public FiniteElement<dim>
                                      * Implementation of the same
                                      * function in
                                      * @ref{FiniteElement}.
+                                     *
+                                     * Passes on control to
+                                     * @p{compute_fill} that does the
+                                     * work for all three
+                                     * @p{fill_fe*_values} functions.
                                      */
     virtual void
     fill_fe_subface_values (const Mapping<dim>                   &mapping,
@@ -416,18 +432,26 @@ class FESystem : public FiniteElement<dim>
     
 
                                     /**
+                                     * Do the work for the three
+                                     * @p{fill_fe*_values} functions.
+                                     * 
                                      * Calls (among other things)
                                      * @p{fill_fe_([sub]face)_values}
                                      * of the base elements. Calls
                                      * @p{fill_fe_values} if
-                                     * @p{face_no==-1} and
-                                     * @p{sub_no==-1}; calls
-                                     * @p{fill_fe_face_values} if
-                                     * @p{face_no==-1} and
-                                     * @p{sub_no!=-1}; and calls
-                                     * @p{fill_fe_subface_values} if 
-                                     * @p{face_no!=-1} and
-                                     * @p{sub_no!=-1}.
+                                     * @p{face_no==invalid_face_no}
+                                     * and
+                                     * @p{sub_no==invalid_face_no};
+                                     * calls @p{fill_fe_face_values}
+                                     * if
+                                     * @p{face_no==invalid_face_no}
+                                     * and
+                                     * @p{sub_no!=invalid_face_no};
+                                     * and calls
+                                     * @p{fill_fe_subface_values} if
+                                     * @p{face_no!=invalid_face_no}
+                                     * and
+                                     * @p{sub_no!=invalid_face_no}.
                                      */
     template <int dim_1>
     void compute_fill (const Mapping<dim>                   &mapping,
@@ -435,8 +459,8 @@ class FESystem : public FiniteElement<dim>
                       const unsigned int                    face_no,
                       const unsigned int                    sub_no,
                       const Quadrature<dim_1>              &quadrature,
-                      typename Mapping<dim>::InternalDataBase      &mapping_data,
-                      typename Mapping<dim>::InternalDataBase      &fe_data,
+                      typename Mapping<dim>::InternalDataBase &mapping_data,
+                      typename Mapping<dim>::InternalDataBase &fe_data,
                       FEValuesData<dim>                    &data) const ;
 
   private:
@@ -743,6 +767,22 @@ class FESystem : public FiniteElement<dim>
                                          * classes).
                                          */
        void delete_fe_values_data (const unsigned int base_no);
+
+                                         /**
+                                          * 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 52618daa7fc88b0abf2e8f0c4393ab4364493554..b1a02a77506678ec42e0c8d2a100838490c029cb 100644 (file)
@@ -170,20 +170,38 @@ class Mapping : public Subscriptor
                                          */
        UpdateFlags  current_update_flags() const;
 
-                                        /**
-                                         * Determine if this is the
-                                         * first cell visited. Then,
-                                         * we need to use
-                                         * @p{update_once} to fill
-                                         * cell-independent
-                                         * fields. Users of derived
-                                         * classes will have to set
-                                         * this field to @p{false};
-                                         * it is initialized to
-                                         * @p{true}.
-                                         */
-       bool first_cell;
-
+                                         /**
+                                          * Return whether we are
+                                          * presently initializing
+                                          * data for the first
+                                          * cell. The value of the
+                                          * field this function is
+                                          * returning is set to
+                                          * @p{true} in the
+                                          * constructor, and cleared
+                                          * by the @p{FEValues} class
+                                          * after the first cell has
+                                          * been initialized.
+                                          *
+                                          * This function is used to
+                                          * determine whether we need
+                                          * to use the @p{update_once}
+                                          * flags for computing data,
+                                          * or whether we can use the
+                                          * @p{update_each} flags.
+                                          */
+        bool is_first_cell () 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
@@ -191,6 +209,13 @@ class Mapping : public Subscriptor
                                          * object.
                                          */
        virtual unsigned int memory_consumption () const;
+
+      private:
+                                         /**
+                                          * The value returned by
+                                          * @p{is_first_cell}.
+                                          */
+        bool first_cell;
     };
     
                                     /**
@@ -515,19 +540,37 @@ template<> const unsigned int Mapping<3>::normal_directions[6];
 template <int dim>
 inline
 UpdateFlags
-Mapping<dim>::InternalDataBase::current_update_flags() const
+Mapping<dim>::InternalDataBase::current_update_flags () const
 {
   if (first_cell)
     {
-      Assert(update_flags==(update_once|update_each), ExcInternalError());
+      Assert (update_flags==(update_once|update_each),
+              ExcInternalError());
       return update_flags;
     }
   else
     return update_each;
-  
-  return update_default;
-}
+};
+
+
 
+template <int dim>
+inline
+bool
+Mapping<dim>::InternalDataBase::is_first_cell () const
+{
+  return first_cell;
+};
+
+
+
+template <int dim>
+inline
+void
+Mapping<dim>::InternalDataBase::clear_first_cell ()
+{
+  first_cell = false;
+};
 
 
 #endif
index 935e9b692a613edf9e73f15535642bb540f43c6e..7223c5b0b8235fcf73f7f45497e5f5d5a84cef91 100644 (file)
@@ -27,8 +27,7 @@
 #include <numeric>
 
 
-//TODO[GK]: in all FEs, somehow make the handling of first_cell more intuitive. presently, this flag is cleared at the end of the fill_fe_*_values function, but this leads to confusion, if for example in the FESystems class, we ask for the present values of the update flags, using update_current, because the element thinks it is already on the second cell and the FESystem does not want to do some actions it is supposed to do. thus, in the respective FESystem function, we store these update_values once before calling the fill_* functions of the subelements, and reuse the stored values. if we would request their values only after initializing the base element (i.e. requesting them at the point _where we actually need them_!), this would fail, and we would not copy the shape_values, for example.
-//TODO[GK]: the solution would probably to not set and clear the first_cell flag in the FEs themselves, but to have a function which sets or clears this flag, and call it from the FEValues object, which should always know what the present state is. this way, we could be more sure that this flag is actually consistent with the present state, and among the various elements that are involved in a FESystem.
+//TODO[GK]: It is relatively confusing with all these update_each, update_once, and in particular update_current flags. there also does not seem to be much documentation at a higher level, except for the documentation of some member variables on this. It would be nice, if we could at least get rid of update_current, to make things a little bit more clear
 
 // if necessary try to work around a bug in the IBM xlC compiler
 #ifdef XLC_WORK_AROUND_STD_BUG
index b95d42aade8d97e0cabbc79c7042cb168e550f89..71d8111e1696a82e3d1f536dda0ef3c67fd475fd 100644 (file)
@@ -319,8 +319,6 @@ FE_DGP<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
   
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, 0, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -368,8 +366,6 @@ FE_DGP<dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -419,8 +415,6 @@ FE_DGP<dim>::fill_fe_subface_values (const Mapping<dim>                   &mappi
   
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
index 11fd8f9586da714e253d349eecf15e29e87d92cd..068c0b361c95e93663343fa81519b271a69fe2a7 100644 (file)
@@ -286,7 +286,6 @@ FE_DGPNonparametric<dim>::fill_fe_values (
              data.shape_2nd_derivatives[k][i] = fe_data.grad_grads[k];
          }
       }
-  fe_data.first_cell = false;
 }
 
 
@@ -328,7 +327,6 @@ FE_DGPNonparametric<dim>::fill_fe_face_values (
              data.shape_2nd_derivatives[k][i] = fe_data.grad_grads[k];
          }
       }
-  fe_data.first_cell = false;
 }
 
 
@@ -371,7 +369,6 @@ FE_DGPNonparametric<dim>::fill_fe_subface_values (
              data.shape_2nd_derivatives[k][i] = fe_data.grad_grads[k];
          }
       }
-  fe_data.first_cell = false;
 }
 
 
index 3715bf68a368b5372ce165195ae550b30977fed2..08c8ae502768830a550d14f868b4c4bf074cba35 100644 (file)
@@ -534,8 +534,6 @@ FE_DGQ<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
   
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, 0, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -583,8 +581,6 @@ FE_DGQ<dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -634,8 +630,6 @@ FE_DGQ<dim>::fill_fe_subface_values (const Mapping<dim>                   &mappi
   
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
index e4cb3128475566ef9e8c64545d9404af91c50fb4..7372f96276fec7d1d6e7f8fa2dfb5b75f401379b 100644 (file)
@@ -938,8 +938,6 @@ FE_Nedelec<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, 0, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 };
 
 
@@ -1081,8 +1079,6 @@ FE_Nedelec<dim>::fill_fe_face_values (const Mapping<dim>                   &mapp
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -1226,8 +1222,6 @@ FE_Nedelec<dim>::fill_fe_subface_values (const Mapping<dim>                   &m
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 };
 
 
index 9a7e53f6872be12183d133dd2024c6a218894e0b..5d79895182cf5e6bb4bfb2e24e4789c686c51068 100644 (file)
@@ -1199,8 +1199,6 @@ FE_Q<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, 0, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -1248,8 +1246,6 @@ FE_Q<dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
 
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
@@ -1299,8 +1295,6 @@ FE_Q<dim>::fill_fe_subface_values (const Mapping<dim>                   &mapping
   
   if (flags & update_second_derivatives)
     compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-  
-  fe_data.first_cell = false;
 }
 
 
index a4f2e49aaadfd85e99434fdde4150d0aedf15579..cc1fa8f31081ec5a820b78409670393388723670 100644 (file)
@@ -127,6 +127,20 @@ InternalData::delete_fe_values_data (const unsigned int base_no)
 
 
 
+template <int dim>
+void
+FESystem<dim>::InternalData::clear_first_cell ()
+{
+                                   // call respective function of base
+                                   // class
+  FiniteElementBase<dim>::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 ------------------- */
 
 
@@ -529,12 +543,13 @@ FESystem<dim>::get_data (const UpdateFlags      flags_,
 
 template <int dim>
 void
-FESystem<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
-                              const typename DoFHandler<dim>::cell_iterator &cell,
-                              const Quadrature<dim>                &quadrature,
-                              typename Mapping<dim>::InternalDataBase       &mapping_data,
-                              typename Mapping<dim>::InternalDataBase       &fe_data,
-                              FEValuesData<dim>                    &data) const
+FESystem<dim>::
+fill_fe_values (const Mapping<dim>                   &mapping,
+                const typename DoFHandler<dim>::cell_iterator &cell,
+                const Quadrature<dim>                &quadrature,
+                typename Mapping<dim>::InternalDataBase &mapping_data,
+                typename Mapping<dim>::InternalDataBase &fe_data,
+                FEValuesData<dim>                    &data) const
 {
   compute_fill(mapping, cell, invalid_face_number, invalid_face_number,
               quadrature, mapping_data, fe_data, data);
@@ -544,16 +559,17 @@ FESystem<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
 
 template <int dim>
 void
-FESystem<dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
-                                   const typename DoFHandler<dim>::cell_iterator &cell,
-                                   const unsigned int                    face_no,
-                                   const Quadrature<dim-1>              &quadrature,
-                                   typename Mapping<dim>::InternalDataBase       &mapping_data,
-                                   typename Mapping<dim>::InternalDataBase       &fe_data,
-                                   FEValuesData<dim>                    &data) const
+FESystem<dim>::
+fill_fe_face_values (const Mapping<dim>                   &mapping,
+                     const typename DoFHandler<dim>::cell_iterator &cell,
+                     const unsigned int                    face_no,
+                     const Quadrature<dim-1>              &quadrature,
+                     typename Mapping<dim>::InternalDataBase &mapping_data,
+                     typename Mapping<dim>::InternalDataBase &fe_data,
+                     FEValuesData<dim>                    &data) const
 {
-  compute_fill(mapping, cell, face_no, invalid_face_number,
-              quadrature, mapping_data, fe_data, data);
+  compute_fill (mapping, cell, face_no, invalid_face_number,
+                quadrature, mapping_data, fe_data, data);
 };
 
 
@@ -561,17 +577,18 @@ FESystem<dim>::fill_fe_face_values (const Mapping<dim>                   &mappin
 
 template <int dim>
 void
-FESystem<dim>::fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                                      const typename DoFHandler<dim>::cell_iterator &cell,
-                                      const unsigned int                    face_no,
-                                      const unsigned int                    sub_no,
-                                      const Quadrature<dim-1>              &quadrature,
-                                      typename Mapping<dim>::InternalDataBase       &mapping_data,
-                                      typename Mapping<dim>::InternalDataBase       &fe_data,
-                                      FEValuesData<dim>                    &data) const
+FESystem<dim>::
+fill_fe_subface_values (const Mapping<dim>                   &mapping,
+                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const unsigned int                    face_no,
+                        const unsigned int                    sub_no,
+                        const Quadrature<dim-1>              &quadrature,
+                        typename Mapping<dim>::InternalDataBase &mapping_data,
+                        typename Mapping<dim>::InternalDataBase &fe_data,
+                        FEValuesData<dim>                    &data) const
 {
-  compute_fill(mapping, cell, face_no, sub_no,
-              quadrature, mapping_data, fe_data, data);
+  compute_fill (mapping, cell, face_no, sub_no,
+                quadrature, mapping_data, fe_data, data);
 }
 
 
@@ -579,14 +596,15 @@ FESystem<dim>::fill_fe_subface_values (const Mapping<dim>                   &map
 template <int dim>
 template <int dim_1>
 void
-FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
-                            const typename DoFHandler<dim>::cell_iterator &cell,
-                            const unsigned int                    face_no,
-                            const unsigned int                    sub_no,
-                            const Quadrature<dim_1>              &quadrature,
-                            typename Mapping<dim>::InternalDataBase       &mapping_data,
-                            typename Mapping<dim>::InternalDataBase       &fedata,
-                            FEValuesData<dim>                    &data) const
+FESystem<dim>::
+compute_fill (const Mapping<dim>                   &mapping,
+              const typename DoFHandler<dim>::cell_iterator &cell,
+              const unsigned int                    face_no,
+              const unsigned int                    sub_no,
+              const Quadrature<dim_1>              &quadrature,
+              typename Mapping<dim>::InternalDataBase &mapping_data,
+              typename Mapping<dim>::InternalDataBase &fedata,
+              FEValuesData<dim>                    &data) const
 {       
   const unsigned int n_q_points = quadrature.n_quadrature_points;
   
@@ -596,8 +614,9 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
                                   // possible
   InternalData & fe_data = dynamic_cast<InternalData&> (fedata);
   
-                                  // Either dim_1==dim (fill_fe_values)
-                                  // or dim_1==dim-1 (fill_fe_(sub)face_values)
+                                  // Either dim_1==dim
+                                  // (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() :
@@ -606,7 +625,7 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
 
   if (flags & (update_values | update_gradients))
     {
-      if (fe_data.first_cell)
+      if (fe_data.is_first_cell())
        {
                                           // Initialize the
                                           // FEValuesDatas for the
@@ -631,8 +650,8 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
                base_fe_data = fe_data.get_fe_data(base_no);
              
                                               // compute update flags ...
-             const UpdateFlags base_update_flags(mapping_data.update_flags
-                                                 | base_fe_data.update_flags);
+             const UpdateFlags base_update_flags
+                = mapping_data.update_flags | base_fe_data.update_flags;
              
                                               // Initialize the FEValuesDatas
                                               // for the base elements.
@@ -663,7 +682,8 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
          Assert(dim_1==dim, 
                 typename FiniteElementData<dim>::
                 ExcSpaceDimensionMismatch(dim_1,dim));
-         cell_quadrature=dynamic_cast<const Quadrature<dim> *>(quadrature_base_pointer);
+         cell_quadrature
+            = dynamic_cast<const Quadrature<dim> *>(quadrature_base_pointer);
          Assert (cell_quadrature != 0, ExcInternalError());
        }
       else
@@ -671,7 +691,8 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
          Assert(dim_1==dim-1, 
                 typename FiniteElementData<dim>::
                 ExcSpaceDimensionMismatch(dim_1,dim-1));
-         face_quadrature=dynamic_cast<const Quadrature<dim-1> *>(quadrature_base_pointer);
+         face_quadrature
+            = dynamic_cast<const Quadrature<dim-1> *>(quadrature_base_pointer);
          Assert (face_quadrature != 0, ExcInternalError());
        }
 
@@ -701,18 +722,6 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
                                            // data on each face,
                                            // therefore use
                                            // base_fe_data.update_flags.
-                                           //
-                                           // Store these flags into
-                                           // base_flags before
-                                           // calling
-                                           // base_fe.fill_fe_([sub]face_)values
-                                           // as the latter changes
-                                           // the return value of
-                                           // base_fe_data.current_update_flags()
-          const UpdateFlags base_flags(dim_1==dim ?
-                                       base_fe_data.current_update_flags() :
-                                       base_fe_data.update_flags);       
-          
          if (face_no==invalid_face_number)
            base_fe.fill_fe_values(mapping, cell,
                                   *cell_quadrature, mapping_data, base_fe_data, base_data);
@@ -752,6 +761,9 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
                                            // those shape functions
                                            // that belong to a given
                                            // base element
+          const UpdateFlags base_flags(dim_1==dim ?
+                                       base_fe_data.current_update_flags() :
+                                       base_fe_data.update_flags);          
           for (unsigned int system_index=0; system_index<this->dofs_per_cell;
                ++system_index)
             if (this->system_to_base_table[system_index].first.first == base_no)
@@ -816,12 +828,8 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
               };
         };
 
-      if (fe_data.first_cell)
+      if (fe_data.is_first_cell())
        {
-         fe_data.first_cell = false;
-         for (unsigned int base_no=0; base_no<n_base_elements(); ++base_no)
-           Assert(fe_data.get_fe_data(base_no).first_cell==false, ExcInternalError());
-         
                                           // delete FEValuesDatas that
                                           // are not needed any more
          for (unsigned int base_no=0; base_no<n_base_elements(); ++base_no)
index a525ea75c0ce87d7aa96bf36ca7a7f1819f90cd8..3575a51e284a8854b82b6e41417e911c8de60e7b 100644 (file)
@@ -607,6 +607,9 @@ void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell)
                                *this->mapping_data,
                                *this->fe_data,
                                *this);
+
+  this->fe_data->clear_first_cell ();
+  this->mapping_data->clear_first_cell ();
 }
 
 
@@ -758,6 +761,9 @@ void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &c
                                     *this->mapping_data,
                                     *this->fe_data,
                                     *this);
+
+  this->fe_data->clear_first_cell ();
+  this->mapping_data->clear_first_cell ();
 };
 
 
@@ -853,6 +859,9 @@ void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator
                                        *this->mapping_data,
                                        *this->fe_data,
                                        *this);
+
+  this->fe_data->clear_first_cell ();
+  this->mapping_data->clear_first_cell ();
 };
 
 
index 7fa2d74579879c65f1109899cc7ef486f5f2367d..af9c98727612335a4f078a26a4dffa60adf7029e 100644 (file)
@@ -268,8 +268,6 @@ MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterat
                                       // vectors on a face are equal
       std::fill (normal_vectors.begin(), normal_vectors.end(), n);
     }
-  
-  data.first_cell = false;
 }
 
 
index 8a31a3129b1dde3c172ac93754de64adcdf5e640..eb5b90d0fe3b8e27454f6eb7b58ce36446a657d7 100644 (file)
@@ -505,8 +505,6 @@ MappingQ1<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterator &cel
   if (update_flags & update_covariant_transformation)
     for (unsigned int point=0; point<npts; ++point)
       data.covariant[point] = invert(data.contravariant[point]);
-
-  data.first_cell = false;
 }
 
 

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.