]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Further clean-ups.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Mar 2001 15:08:38 +0000 (15:08 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Mar 2001 15:08:38 +0000 (15:08 +0000)
git-svn-id: https://svn.dealii.org/trunk@4148 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_base.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_dgq.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_system.cc

index cf32082d61f73954296750203b883b5c1c4e7197..28e40381478503eabfd983f556adf364d48444dc 100644 (file)
@@ -48,6 +48,7 @@ template <int dim> class MatrixCreator;
 template <int dim>
 class FiniteElement : public FiniteElementBase<dim>
 {
+  private:
                                     /**
                                      * Copy constructor prohibited.
                                      */
@@ -67,21 +68,6 @@ class FiniteElement : public FiniteElementBase<dim>
                                      */
     virtual ~FiniteElement ();
 
-                                    /**
-                                     * Compute flags for initial
-                                     * update only.
-                                     * @see{FEValuesBase}
-                                     */
-    virtual UpdateFlags update_once (UpdateFlags flags) const = 0;
-  
-                                    /**
-                                     * Compute flags for update on
-                                     * each cell.
-                                     * @see{FEValuesBase}
-                                     */
-    virtual UpdateFlags update_each (UpdateFlags flags) const = 0;
-  
-
                                     /**
                                      * Return the support points of the
                                      * trial functions on the unit cell.
@@ -116,16 +102,73 @@ class FiniteElement : public FiniteElementBase<dim>
                                      * will be zero. This is the standard behavior,
                                      * if the function is not overloaded.
                                      */
-    virtual void get_unit_face_support_points (std::vector<Point<dim-1> > &) const;    
+    virtual void get_unit_face_support_points (std::vector<Point<dim-1> > &points) const;    
     
-    friend class FEValues<dim>;
-    friend class FEFaceValues<dim>;
-    friend class FESubfaceValues<dim>;
-    friend class FESystem<dim>;
-    friend class MatrixCreator<dim>;
-    friend class VectorTools;
+                                    /**
+                                     * Number of base elements in a mixed
+                                     * discretization. This function returns
+                                     * 1 for primitive elements.
+                                     */
+    virtual unsigned int n_base_elements () const;
+    
+                                    /**
+                                     * Access to base element
+                                     * objects.  By default,
+                                     * #base_element(0)# is #this#.
+                                     * This function is overloaded by
+                                     * system elements to allow
+                                     * access to the different
+                                     * components of mixed
+                                     * discretizations.
+                                     */
+    virtual const FiniteElement<dim> & base_element (const unsigned int index) const;
     
+                                    /**
+                                     * Determine an estimate for the
+                                     * memory consumption (in bytes)
+                                     * of this object.
+                                     *
+                                     * This function is made virtual,
+                                     * since finite element objects
+                                     * are usually accessed through
+                                     * pointers to their base class,
+                                     * rather than the class itself.
+                                     */
+    virtual unsigned int memory_consumption () const;
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcBoundaryFaceUsed);
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcJacobiDeterminantHasWrongSign);
+                                    /**
+                                     * Exception
+                                     */
+    DeclException1 (ExcComputationNotUseful,
+                   int,
+                   << "The computation you required from this function is not "
+                   << "feasible or not probable in the present dimension ("
+                   << arg1 << ") because it would be prohibitively expensive.");
+
   protected:
+
+                                    /**
+                                     * Compute flags for initial
+                                     * update only.
+                                     * @see{FEValuesBase}
+                                     */
+    virtual UpdateFlags update_once (UpdateFlags flags) const = 0;
+  
+                                    /**
+                                     * Compute flags for update on
+                                     * each cell.
+                                     * @see{FEValuesBase}
+                                     */
+    virtual UpdateFlags update_each (UpdateFlags flags) const = 0;
+  
                                     /**
                                      * @p{clone} function instead of
                                      * a copy constructor.
@@ -228,58 +271,14 @@ class FiniteElement : public FiniteElementBase<dim>
                            typename Mapping<dim>::InternalDataBase      &fe_internal,
                            FEValuesData<dim> &data) const = 0;
 
-
-  public:
-
-                                    /**
-                                     * Number of base elements in a mixed
-                                     * discretization. This function returns
-                                     * 1 for simple elements.
-                                     */
-    virtual unsigned int n_base_elements () const;
-    
-                                    /**
-                                     * Access to base element
-                                     * objects.  By default,
-                                     * #base_element(0)# is #this#.
-                                     * This function is overloaded by
-                                     * system elements to allow
-                                     * access to the different
-                                     * components of mixed
-                                     * discretizations.
-                                     */
-    virtual const FiniteElement<dim>& base_element (const unsigned int index) const;
-    
-                                    /**
-                                     * Determine an estimate for the
-                                     * memory consumption (in bytes)
-                                     * of this object.
-                                     *
-                                     * This function is made virtual,
-                                     * since finite element objects
-                                     * are usually accessed through
-                                     * pointers to their base class,
-                                     * rather than the class itself.
-                                     */
-    virtual unsigned int memory_consumption () const;
-
                                     /**
-                                     * Exception
+                                     * Declare some other classes as
+                                     * friends of this class.
                                      */
-    DeclException0 (ExcBoundaryFaceUsed);
-                                    /**
-                                     * Exception
-                                     */
-    DeclException0 (ExcJacobiDeterminantHasWrongSign);
-                                    /**
-                                     * Exception
-                                     */
-    DeclException1 (ExcComputationNotUseful,
-                   int,
-                   << "The computation you required from this function is not "
-                   << "feasible or not probable in the present dimension ("
-                   << arg1 << ") because it would be prohibitively expensive.");
-
+    friend class FEValues<dim>;
+    friend class FEFaceValues<dim>;
+    friend class FESubfaceValues<dim>;
+    friend class FESystem<dim>;
 };
 
 
index 1b0ce3ea8f12aa9fc2d8dc6de206a2fa38ae2a22..1eb6f2eddae0c652b4bbf43bac97225a4db124ea 100644 (file)
@@ -335,12 +335,14 @@ class FiniteElementBase : public Subscriptor,
     {
       public:
                                       /**
-                                       * Initialize @ref{FEValues}
-                                       * pointers for given element.
+                                       * Initialize some pointers
+                                       * used in the computation of
+                                       * second derivatives by finite
+                                       * differencing of gradients.
                                        */
-      void initialize (const FiniteElement<dim>* element,
-                      const Mapping<dim>& mapping,
-                      const Quadrature<dim>& quadrature);
+      void initialize_2nd (const FiniteElement<dim> *element,
+                          const Mapping<dim>       &mapping,
+                          const Quadrature<dim>    &quadrature);
       
                                       /**
                                        * Destructor. Needed to avoid
@@ -377,14 +379,16 @@ class FiniteElementBase : public Subscriptor,
                       const std::vector<bool> &restriction_is_additive_flags);
     
                                      /**
-                                     * Compute second differences.
+                                     * Compute second derivatives by
+                                     * finite differences of
+                                     * gradients.
                                      */
-    void compute_2nd (const Mapping<dim> &mapping,
-                     const DoFHandler<dim>::cell_iterator &cell,
-                     const unsigned int offset,
+    void compute_2nd (const Mapping<dim>                      &mapping,
+                     const DoFHandler<dim>::cell_iterator    &cell,
+                     const unsigned int                       offset,
                      typename Mapping<dim>::InternalDataBase &mapping_internal,
-                     InternalDataBasefe_internal,
-                     FEValuesData<dim>data) const;
+                     InternalDataBase                        &fe_internal,
+                     FEValuesData<dim>                       &data) const;
 
                                     /**
                                      * Projection from a fine grid
index 087ea67dff2eafeeb227ca9358e1e7b7600cf25b..295ebc49622151ae04c82336cea7c778d08dd8eb 100644 (file)
@@ -466,10 +466,11 @@ class FESystem : public FiniteElement<dim>
        ~InternalData();
 
                                         /**
-                                         * Flag for computation of
-                                         * second derivatives.
+                                         * Flag indicating whether
+                                         * second derivatives shall
+                                         * be computed.
                                          */
-       bool second_flag;
+       bool compute_second_derivatives;
        
                                         /**
                                          * Gives write-access to the
index 7718537308a57c62643983ef2a3a86ac1b6a26f6..01eda10290284aa4777871616527c175bf883a52 100644 (file)
@@ -34,35 +34,51 @@ using namespace std;
 template <int dim>
 void
 FiniteElementBase<dim>::
-InternalDataBase::initialize (const FiniteElement<dim>* element,
-                             const Mapping<dim>& mapping,
-                             const Quadrature<dim>& quadrature)
+InternalDataBase::initialize_2nd (const FiniteElement<dim> *element,
+                                 const Mapping<dim>       &mapping,
+                                 const Quadrature<dim>    &quadrature)
 {
-                                  // We compute difference
-                                  // quotients of gradients
-  UpdateFlags diff_flags = update_gradients;
-  
-                                  // We will need shifted
-                                  // quadrature formulae
+                                  // if we shall compute second
+                                  // derivatives, then we do so by
+                                  // finite differencing the
+                                  // gradients. that we do by
+                                  // evaluating the gradients of
+                                  // shape values at points shifted
+                                  // star-like a little in each
+                                  // coordinate direction around each
+                                  // quadrature point.
+                                  //
+                                  // therefore generate 2*dim (the
+                                  // number of evaluation points)
+                                  // FEValues objects with slightly
+                                  // shifted positions
   std::vector<Point<dim> > diff_points (quadrature.n_quadrature_points);
-  std::vector<double> diff_weights (quadrature.n_quadrature_points, 0);
+  std::vector<double>      diff_weights (quadrature.n_quadrature_points, 0);
   
-                                  // The star has 2dim points
   differences.resize(2*dim);
-  for (unsigned int d=0;d<dim;++d)
+  for (unsigned int d=0; d<dim; ++d)
     {
       Point<dim> shift;
+//TODO: unify the places where the finite differencing step length is used      
       shift (d) = 1.e-6;
-      for (unsigned int i=0;i<diff_points.size();++i)
+
+                                      // generate points and FEValues
+                                      // objects shifted in
+                                      // plus-direction. note that
+                                      // they only need to compute
+                                      // gradients, not more
+      for (unsigned int i=0; i<diff_points.size(); ++i)
        diff_points[i] = quadrature.point(i) + shift;
-      Quadrature<dim> plus_quad (diff_points, diff_weights);
-      differences[d] =
-       new FEValues<dim> (mapping, *element, plus_quad, diff_flags);
-      for (unsigned int i=0;i<diff_points.size();++i)
+      const Quadrature<dim> plus_quad (diff_points, diff_weights);
+      differences[d] = new FEValues<dim> (mapping, *element,
+                                         plus_quad, update_gradients);
+
+                                      // now same in minus-direction
+      for (unsigned int i=0; i<diff_points.size(); ++i)
        diff_points[i] = quadrature.point(i) - shift;
-      Quadrature<dim> minus_quad (diff_points, diff_weights);
-      differences[d+dim] =
-       new FEValues<dim> (mapping, *element, minus_quad, diff_flags);    
+      const Quadrature<dim> minus_quad (diff_points, diff_weights);
+      differences[d+dim] = new FEValues<dim> (mapping, *element,
+                                             minus_quad, update_gradients); 
     }
 }
 
@@ -72,9 +88,15 @@ InternalDataBase::initialize (const FiniteElement<dim>* element,
 template <int dim>
 FiniteElementBase<dim>::InternalDataBase::~InternalDataBase ()
 {
-  for (unsigned int i=0;i<differences.size ();++i)
+  for (unsigned int i=0; i<differences.size (); ++i)
     if (differences[i] != 0)
-      delete differences[i];
+      {
+                                        // delete pointer and set it
+                                        // to zero to avoid
+                                        // inadvertant use
+       delete differences[i];
+       differences[i] = 0;
+      };
 }
 
 
@@ -82,14 +104,15 @@ FiniteElementBase<dim>::InternalDataBase::~InternalDataBase ()
 
 template <int dim>
 FiniteElementBase<dim>::FiniteElementBase (const FiniteElementData<dim> &fe_data,
-                                          const std::vector<bool> &restriction_is_additive_flags) :
-  FiniteElementData<dim> (fe_data),
-  system_to_component_table(dofs_per_cell),
-  face_system_to_component_table(dofs_per_face),
-  component_to_system_table(components, std::vector<unsigned>(dofs_per_cell)),
-                             face_component_to_system_table(components, std::vector<unsigned>(dofs_per_face)),
-                                                              component_to_base_table(components),
-                                                              restriction_is_additive_flags(restriction_is_additive_flags)
+                                          const std::vector<bool> &restriction_is_additive_flags)
+               :
+               FiniteElementData<dim> (fe_data),
+                system_to_component_table(dofs_per_cell),
+                face_system_to_component_table(dofs_per_face),
+                component_to_system_table(components, std::vector<unsigned>(dofs_per_cell)),
+               face_component_to_system_table(components, std::vector<unsigned>(dofs_per_face)),
+               component_to_base_table(components),
+               restriction_is_additive_flags(restriction_is_additive_flags)
 {
   Assert(restriction_is_additive_flags.size()==fe_data.components,
         ExcDimensionMismatch(restriction_is_additive_flags.size(),fe_data.components));
@@ -100,13 +123,16 @@ FiniteElementBase<dim>::FiniteElementBase (const FiniteElementData<dim> &fe_data
       prolongation[i].reinit (dofs_per_cell, dofs_per_cell);
     };
 
+                                  // first set sizes of some
+                                  // matrices. they will be filled by
+                                  // derived classes, or re-set to
+                                  // zero size of not defined
   switch (dim)
     {
       case 1:
            Assert ((interface_constraints.m() == 0) &&
                    (interface_constraints.n() == 0),
                    ExcInternalError());
-           
            break;
            
       case 2:
@@ -142,6 +168,7 @@ FiniteElementBase<dim>::FiniteElementBase (const FiniteElementData<dim> &fe_data
 };
 
 
+
 template <int dim>
 const FullMatrix<double> &
 FiniteElementBase<dim>::restrict (const unsigned int child) const
@@ -153,6 +180,7 @@ FiniteElementBase<dim>::restrict (const unsigned int child) const
 };
 
 
+
 template <int dim>
 const FullMatrix<double> &
 FiniteElementBase<dim>::prolongate (const unsigned int child) const
@@ -164,6 +192,7 @@ FiniteElementBase<dim>::prolongate (const unsigned int child) const
 };
 
 
+
 template <int dim>
 const FullMatrix<double> &
 FiniteElementBase<dim>::constraints () const
@@ -180,6 +209,7 @@ FiniteElementBase<dim>::constraints () const
 };
 
 
+
 template <int dim>
 bool FiniteElementBase<dim>::operator == (const FiniteElementBase<dim> &f) const
 {
@@ -211,52 +241,76 @@ FiniteElementBase<dim>::memory_consumption () const
 };
 
 
+
 template <int dim>
 void
 FiniteElementBase<dim>::
-compute_2nd (const Mapping<dim> &mapping,
+compute_2nd (const Mapping<dim>                   &mapping,
             const DoFHandler<dim>::cell_iterator &cell,
-            const unsigned int offset,
-            Mapping<dim>::InternalDataBase &mapping_internal,
-            InternalDataBasefe_internal,
-            FEValuesData<dim>data) const
+            const unsigned int                    offset,
+            Mapping<dim>::InternalDataBase       &mapping_internal,
+            InternalDataBase                     &fe_internal,
+            FEValuesData<dim>                    &data) const
 {
                                   // Number of quadrature points
-  const unsigned int n = data.shape_2nd_derivatives[0].size();
-  
-  for (unsigned int d=0;d<dim;++d)
+  const unsigned int n_q_points = data.shape_2nd_derivatives[0].size();
+
+                                  // first reinit the fe_values
+                                  // objects used for the finite
+                                  // differencing stuff
+  for (unsigned int d=0; d<dim; ++d)
     {
       fe_internal.differences[d]->reinit(cell);
       fe_internal.differences[d+dim]->reinit(cell);
     }
 
-  std::vector<std::vector<Tensor<1,dim> > > diff_quot (dim, std::vector<Tensor<1,dim> >(n));
-  std::vector<Tensor<1,dim> > diff_quot2 (n);
-                                  // Loop over shape functions
+                                  // collection of difference
+                                  // quotients of gradients in each
+                                  // direction (first index) and at
+                                  // all q-points (second index)
+  std::vector<std::vector<Tensor<1,dim> > > diff_quot (dim, std::vector<Tensor<1,dim> >(n_q_points));
+  std::vector<Tensor<1,dim> >               diff_quot2 (n_q_points);
+
+                                  // for all shape functions at all
+                                  // quadrature points and difference
+                                  // quotients in all directions:
   for (unsigned int shape=0; shape<dofs_per_cell; ++shape)
     {
-      
-                                      // Fill difference quotients
-      for (unsigned int d1=0;d1<dim;++d1)
-                                        // Loop over quadrature points
-       for (unsigned int k=0;k<n;++k)
+      for (unsigned int d1=0; d1<dim; ++d1)
+       for (unsigned int q=0; q<n_q_points; ++q)
          {
+                                            // get gradient at points
+                                            // shifted slightly to
+                                            // the right and to the
+                                            // left in the present
+                                            // coordinate direction
            const Tensor<1,dim>& right
-             = fe_internal.differences[d1]->shape_grad(shape, k);
+             = fe_internal.differences[d1]->shape_grad(shape, q);
            const Tensor<1,dim>& left
-             = fe_internal.differences[d1+dim]->shape_grad(shape, k);
-           for (unsigned int d=0;d<dim;++d)
-             diff_quot[d][k][d1] = (.5/1.e-6) * (right[d]-left[d]);
+             = fe_internal.differences[d1+dim]->shape_grad(shape, q);
+
+//TODO: unify the places where the finite differencing step length is used      
+                                            // compute the second
+                                            // derivative from a
+                                            // symmetric difference
+                                            // approximation
+           for (unsigned int d=0; d<dim; ++d)
+             diff_quot[d][q][d1] = 1./(2*1.e-6) * (right[d]-left[d]);
          }
-      
-      for (unsigned int d=0;d<dim;++d)
+
+                                      // up to now we still have
+                                      // difference quotients on the
+                                      // unit cell, so transform it
+                                      // to something on the real
+                                      // cell
+      for (unsigned int d=0; d<dim; ++d)
        {
          mapping.transform_covariant (diff_quot2, diff_quot[d],
                                       mapping_internal, offset);
 
-         for (unsigned int k=0;k<n;++k)
-           for (unsigned int d1=0;d1<dim;++d1)
-             data.shape_2nd_derivatives[shape][k][d][d1] = diff_quot2[k][d1];
+         for (unsigned int q=0; q<n_q_points; ++q)
+           for (unsigned int d1=0; d1<dim; ++d1)
+             data.shape_2nd_derivatives[shape][q][d][d1] = diff_quot2[q][d1];
        }
     }
 }
@@ -283,14 +337,17 @@ template <int dim>
 void
 FiniteElement<dim>::get_unit_support_points (std::vector<Point<dim> > &points) const
 {
+                                  // default implementation
   points.resize(0);
 }
 
+
     
 template <int dim>
 void
 FiniteElement<dim>::get_unit_face_support_points (std::vector<Point<dim-1> > &points) const
 {
+                                  // default implementation
   points.resize(0);
 }
 
@@ -298,8 +355,8 @@ FiniteElement<dim>::get_unit_face_support_points (std::vector<Point<dim-1> > &po
 
 template <int dim>
 Mapping<dim>::InternalDataBase*
-FiniteElement<dim>::get_face_data (const UpdateFlags flags,
-                                  const Mapping<dim>mapping,
+FiniteElement<dim>::get_face_data (const UpdateFlags       flags,
+                                  const Mapping<dim>      &mapping,
                                   const Quadrature<dim-1> &quadrature) const
 {
   QProjector<dim> q(quadrature, false);
@@ -307,25 +364,30 @@ FiniteElement<dim>::get_face_data (const UpdateFlags flags,
 }
 
 
+
 template <int dim>
 Mapping<dim>::InternalDataBase*
-FiniteElement<dim>::get_subface_data (const UpdateFlags flags,
-                                     const Mapping<dim>mapping,
+FiniteElement<dim>::get_subface_data (const UpdateFlags        flags,
+                                     const Mapping<dim>      &mapping,
                                      const Quadrature<dim-1> &quadrature) const
 {
   QProjector<dim> q(quadrature, true);
   return get_data (flags, mapping, q);
-  
 }
 
 
+
 template <int dim>
 unsigned int
 FiniteElement<dim>::n_base_elements() const
 {
+                                  // default implementation
   return 1;
 }
 
+
+
+
 template <int dim>
 unsigned int
 FiniteElement<dim>::memory_consumption () const
index 9e9c75693b495d62c9f09e00527cde528f8357d8..6e8055717ef9bca1a73c8769d65b181b0269ca9b 100644 (file)
@@ -394,10 +394,12 @@ FE_DGQ<dim>::get_data (const UpdateFlags update_flags,
                                   std::vector<Tensor<1,dim> >(quadrature.n_quadrature_points));
     }
 
+                                  // if second derivatives through
+                                  // finite differencing is required,
+                                  // then initialize some objects for
+                                  // that
   if (flags & update_second_derivatives)
-    {
-      data->initialize (this, mapping, quadrature);
-    }
+    data->initialize_2nd (this, mapping, quadrature);
   
   
   if (flags & (update_values | update_gradients))
index 29cb80fd3ab7de8a21d470b22e5ba9a86afc8786..0bc632c5ba1c49a5da3533ca86aee863fa2af611 100644 (file)
@@ -989,8 +989,12 @@ FE_Q<dim>::get_data (const UpdateFlags      update_flags,
                                   std::vector<Tensor<1,dim> >(n_q_points));
     }
 
+                                  // if second derivatives through
+                                  // finite differencing is required,
+                                  // then initialize some objects for
+                                  // that
   if (flags & update_second_derivatives)
-    data->initialize (this, mapping, quadrature);
+    data->initialize_2nd (this, mapping, quadrature);
   
   if (flags & (update_values | update_gradients))
     for (unsigned int i=0; i<n_q_points; ++i)
index b5ba860806b068b00f5920c78377f56bae881fb5..18e0c872960e9a8e5242b6bba3ed3b3888c0d72a 100644 (file)
@@ -156,20 +156,24 @@ FESystem<dim>::update_each (UpdateFlags flags) const
 
 template <int dim>
 Mapping<dim>::InternalDataBase*
-FESystem<dim>::get_data (UpdateFlags flags,
-                        const Mapping<dim>mapping,
-                        const Quadrature<dim>quadrature) const
+FESystem<dim>::get_data ( UpdateFlags           flags,
+                        const Mapping<dim>    &mapping,
+                        const Quadrature<dim> &quadrature) const
 {
   InternalData* data = new InternalData(n_base_elements());
 
-  data->second_flag = flags & update_second_derivatives;
-
-                                  // Make sure that this object
-                                  // computes 2nd derivatives itself
-  if (data->second_flag)
+                                  // if second derivatives through
+                                  // finite differencing is required,
+                                  // then initialize some objects for
+                                  // that
+  data->compute_second_derivatives = flags & update_second_derivatives;
+  if (data->compute_second_derivatives)
     {
+                                      // delete
+                                      // update_second_derivatives
+                                      // from flags list
       flags = UpdateFlags (flags ^ update_second_derivatives);
-      data->initialize (this, mapping, quadrature);
+      data->initialize_2nd (this, mapping, quadrature);
     }
   
   
@@ -177,13 +181,17 @@ FESystem<dim>::get_data (UpdateFlags flags,
     {
       typename Mapping<dim>::InternalDataBase *base_fe_data_base =
        base_element(base_no).get_data(flags, mapping, quadrature);
+
       FiniteElementBase<dim>::InternalDataBase *base_fe_data =
        dynamic_cast<typename FiniteElementBase<dim>::InternalDataBase *>
        (base_fe_data_base);
       
       data->set_fe_data(base_no, base_fe_data);
-      data->update_once|=base_fe_data->update_once;
-      data->update_each|=base_fe_data->update_each;
+
+                                      // collect requirements of base
+                                      // elements
+      data->update_once |= base_fe_data->update_once;
+      data->update_each |= base_fe_data->update_each;
       
                                       // The FEValuesData @p{data}
                                       // given to the
@@ -210,7 +218,7 @@ FESystem<dim>::get_data (UpdateFlags flags,
                                       // @p{data}, similar to the
                                       // storing of the
                                       // @p{base_fe_data}s.
-      FEValuesData<dim> *base_data=new FEValuesData<dim>();
+      FEValuesData<dim> *base_data = new FEValuesData<dim>();
       data->set_fe_values_data(base_no, base_data);
     }
   data->update_flags=data->update_once | data->update_each;
@@ -451,7 +459,7 @@ FESystem<dim>::compute_fill (const Mapping<dim>                   &mapping,
            }
        }
     }
-  if (fe_data.second_flag)
+  if (fe_data.compute_second_derivatives)
     {
       unsigned int offset = 0;
       if (face_no != static_cast<unsigned int> (-1))
@@ -1116,14 +1124,22 @@ FESystem<dim>::InternalData::InternalData(const unsigned int n_base_elements):
 template <int dim>
 FESystem<dim>::InternalData::~InternalData()
 {
+                                  // delete pointers and set them to
+                                  // zero to avoid inadvertant use
   for (unsigned int i=0; i<base_fe_datas.size(); ++i)
     if (base_fe_datas[i])
-      delete base_fe_datas[i];
+      {
+       delete base_fe_datas[i];
+       base_fe_datas[i] = 0;
+      };
   
   for (unsigned int i=0; i<base_fe_values_datas.size(); ++i)
     if (base_fe_values_datas[i])
-      delete base_fe_values_datas[i];
-}
+      {
+       delete base_fe_values_datas[i];
+       base_fe_values_datas[i] = 0;
+      };
+};
 
 
 

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.