]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Make the FEValues::reinit functions (much) more efficient in case we repeatedly work...
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 7 Mar 2009 15:16:42 +0000 (15:16 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 7 Mar 2009 15:16:42 +0000 (15:16 +0000)
git-svn-id: https://svn.dealii.org/trunk@18464 0785d39b-7218-0410-832d-ea1e28bc413d

24 files changed:
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_dgp_nonparametric.h
deal.II/deal.II/include/fe/fe_nedelec.h
deal.II/deal.II/include/fe/fe_poly.h
deal.II/deal.II/include/fe/fe_poly.templates.h
deal.II/deal.II/include/fe/fe_poly_tensor.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/include/fe/mapping.h
deal.II/deal.II/include/fe/mapping_cartesian.h
deal.II/deal.II/include/fe/mapping_q.h
deal.II/deal.II/include/fe/mapping_q1.h
deal.II/deal.II/include/fe/mapping_q1_eulerian.h
deal.II/deal.II/include/fe/mapping_q_eulerian.h
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc
deal.II/deal.II/source/fe/fe_nedelec.cc
deal.II/deal.II/source/fe/fe_poly_tensor.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_q.cc
deal.II/deal.II/source/fe/mapping_q1.cc
deal.II/deal.II/source/fe/mapping_q1_eulerian.cc
deal.II/deal.II/source/fe/mapping_q_eulerian.cc

index 8be920321733e460a334a055a7c1f27051871275..7fef35499f37a1f61164da7f876970a97459f61f 100644 (file)
@@ -2588,12 +2588,13 @@ class FiniteElement : public Subscriptor,
                                      * called for the same cell first!
                                      */                                      
     virtual void
-    fill_fe_values (const Mapping<dim,spacedim>                   &mapping,
+    fill_fe_values (const Mapping<dim,spacedim>                               &mapping,
                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                &quadrature,
-                   typename Mapping<dim,spacedim>::InternalDataBase       &mapping_internal,
-                   typename Mapping<dim,spacedim>::InternalDataBase       &fe_internal,
-                   FEValuesData<dim,spacedim>                    &data) const = 0;
+                   const Quadrature<dim>                                     &quadrature,
+                   const enum CellSimilarity::Similarity                      cell_similarity,
+                   typename Mapping<dim,spacedim>::InternalDataBase          &mapping_internal,
+                   typename Mapping<dim,spacedim>::InternalDataBase          &fe_internal,
+                   FEValuesData<dim,spacedim>                                &data) const = 0;
     
                                     /**
                                      * Fill the fields of
index 41e92926eb207194a2455101d27fd28ba00610c9..98b0c4f9e0d38ad86c47474e3ad0359954868d43 100644 (file)
@@ -480,7 +480,8 @@ class FE_DGPNonparametric : public FiniteElement<dim,spacedim>
     virtual void
     fill_fe_values (const Mapping<dim,spacedim> &mapping,
                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                &quadrature,
+                   const Quadrature<dim>                                 &quadrature,
+                   const enum CellSimilarity::Similarity                  cell_similarity,
                    typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
                    FEValuesData<dim,spacedim>& data) const;
index 62cc2404f57772a22a99dabe6d6d5c53df1d2a5c..5270ded777880bf09f3fe92bd199e5ddf7df8667 100644 (file)
@@ -412,11 +412,12 @@ class FE_Nedelec : public FiniteElement<dim,spacedim>
                                      * FiniteElement.
                                      */
     virtual void
-    fill_fe_values (const Mapping<dim,spacedim> &mapping,
+    fill_fe_values (const Mapping<dim,spacedim>                      &mapping,
                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                &quadrature,
-                   typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
-                   typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
+                   const Quadrature<dim>                            &quadrature,
+                   const enum CellSimilarity::Similarity             cell_similarity,
+                   typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
+                   typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
                    FEValuesData<dim,spacedim>& data) const;
     
                                     /**
index ed1b4652ff850cb3e642e4b4554901483d1b4fee..c014e6819a750b605d03c467c2a07dad5de156a1 100644 (file)
@@ -218,12 +218,13 @@ class FE_Poly : public FiniteElement<dim,spacedim>
              const Quadrature<dim>& quadrature) const ;
 
     virtual void
-    fill_fe_values (const Mapping<dim,spacedim> &mapping,
+    fill_fe_values (const Mapping<dim,spacedim>                           &mapping,
                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                &quadrature,
+                   const Quadrature<dim>                                 &quadrature,
+                   const enum CellSimilarity::Similarity                  cell_similarity,
                    typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
-                   FEValuesData<dim,spacedim>& data) const;
+                   FEValuesData<dim,spacedim>&                            data) const;
     
     virtual void
     fill_fe_face_values (const Mapping<dim,spacedim> &mapping,
index 9c54d475260c513317815d06bd04928919a4185e..9f979ba926b2ea1f76e6245dcaf7a72cc5b02410 100644 (file)
@@ -156,8 +156,8 @@ FE_Poly<POLY,dim,spacedim>::update_each (const UpdateFlags flags) const
 template <class POLY, int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
 FE_Poly<POLY,dim,spacedim>::get_data (const UpdateFlags      update_flags,
-                            const Mapping<dim,spacedim>    &mapping,
-                            const Quadrature<dim> &quadrature) const
+                                     const Mapping<dim,spacedim>    &mapping,
+                                     const Quadrature<dim> &quadrature) const
 {
                                   // generate a new data object and
                                   // initialize some fields
@@ -238,12 +238,14 @@ FE_Poly<POLY,dim,spacedim>::get_data (const UpdateFlags      update_flags,
 template <>
 inline
 void
-FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_values (const Mapping<1,2> &mapping,
-                                  const Triangulation<1,2>::cell_iterator   &/*cell*/,
-                                  const Quadrature<1>                       &quadrature,
-                                  Mapping<1,2>::InternalDataBase            &mapping_data,
-                                  Mapping<1,2>::InternalDataBase            &fedata,
-                                  FEValuesData<1,2>                         &data) const
+FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_values 
+  (const Mapping<1,2> &mapping,
+   const Triangulation<1,2>::cell_iterator   &/*cell*/,
+   const Quadrature<1>                       &quadrature,
+   const enum CellSimilarity::Similarity      cell_similarity,
+   Mapping<1,2>::InternalDataBase            &mapping_data,
+   Mapping<1,2>::InternalDataBase            &fedata,
+   FEValuesData<1,2>                         &data) const
 {
                                   // convert data object to internal
                                   // data for this class. fails with
@@ -260,13 +262,13 @@ FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_values (const Mapping<1,2> &ma
        for (unsigned int i=0; i<quadrature.size(); ++i)
          data.shape_values(k,i) = fe_data.shape_values[k][i];
       
-      if (flags & update_gradients)
+      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
        mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k], 
                          mapping_data, mapping_covariant);
     }
   
   //  const typename QProjector<1>::DataSetDescriptor dsd;
-  if (flags & update_hessians)
+  if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
     {
       AssertThrow(false, ExcNotImplemented());
 /*       this->compute_2nd (mapping, cell, dsd.cell(), */
@@ -279,12 +281,14 @@ FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_values (const Mapping<1,2> &ma
 template <>
 inline
 void
-FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_values (const Mapping<2,3> &mapping,
-                                  const Triangulation<2,3>::cell_iterator &/* cell */,
-                                  const Quadrature<2> &quadrature,
-                                  Mapping<2,3>::InternalDataBase & mapping_data,
-                                  Mapping<2,3>::InternalDataBase &fedata,
-                                  FEValuesData<2,3> &data) const
+FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_values 
+  (const Mapping<2,3> &mapping,
+   const Triangulation<2,3>::cell_iterator &/* cell */,
+   const Quadrature<2> &quadrature,
+   const enum CellSimilarity::Similarity cell_similarity,
+   Mapping<2,3>::InternalDataBase & mapping_data,
+   Mapping<2,3>::InternalDataBase &fedata,
+   FEValuesData<2,3> &data) const
 {
 
                                   // assert that the following dynamics
@@ -300,13 +304,13 @@ FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_values (const Mapping<2,3> &ma
        for (unsigned int i=0; i<quadrature.size(); ++i)
          data.shape_values(k,i) = fe_data.shape_values[k][i];     
 
-       if (flags & update_gradients)
+       if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
         mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
                           mapping_data, mapping_covariant);
     }
   
 /*   const QProjector<2>::DataSetDescriptor dsd; */
-  if (flags & update_hessians)
+  if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
     {
       AssertThrow(false, ExcNotImplemented());
 /*       this->compute_2nd (mapping, cell, dsd.cell(), */
@@ -323,6 +327,7 @@ FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_values (
   const Mapping<1,2>& mapping,
   const Triangulation<1,2>::cell_iterator& /*cell*/,
   const Quadrature<1>& quadrature,
+  const enum CellSimilarity::Similarity cell_similarity,
   Mapping<1,2>::InternalDataBase& mapping_data,
   Mapping<1,2>::InternalDataBase& fedata,
   FEValuesData<1,2>& data) const
@@ -345,13 +350,13 @@ FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_values (
       
 // TODO: I would think this should work. Guido
       
-      if (flags & update_gradients)
+      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
        mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
                          mapping_data, mapping_covariant);
     }
   
                                   //  const typename QProjector<1>::DataSetDescriptor dsd;
-  if (flags & update_hessians)
+  if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
     {
       AssertThrow(false, ExcNotImplemented());
 /*       this->compute_2nd (mapping, cell, dsd.cell(), */
@@ -363,12 +368,14 @@ FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_values (
 template <>
 inline
 void
-FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_values (const Mapping<2,3>& mapping,
-                                  const Triangulation<2,3>::cell_iterator &/* cell */,
-                                  const Quadrature<2>& quadrature,
-                                  Mapping<2,3>::InternalDataBase& mapping_data,
-                                  Mapping<2,3>::InternalDataBase &fedata,
-                                  FEValuesData<2,3> &data) const
+FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_values 
+  (const Mapping<2,3>& mapping,
+   const Triangulation<2,3>::cell_iterator &/* cell */,
+   const Quadrature<2>& quadrature,
+   const enum CellSimilarity::Similarity cell_similarity,
+   Mapping<2,3>::InternalDataBase& mapping_data,
+   Mapping<2,3>::InternalDataBase &fedata,
+   FEValuesData<2,3> &data) const
 {
   Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
   InternalData &fe_data = static_cast<InternalData &> (fedata);
@@ -382,13 +389,13 @@ FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_values (const Mapping<2,3>& mapping,
          data.shape_values(k,i) = fe_data.shape_values[k][i];
       
 
-       if (flags & update_gradients)
+       if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
         mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
                                      mapping_data, mapping_covariant);
     }
   
 /*   const QProjector<2>::DataSetDescriptor dsd; */
-  if (flags & update_hessians)
+  if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
     {
       AssertThrow(false, ExcNotImplemented());
 /*       this->compute_2nd (mapping, cell, dsd.cell(), */
@@ -400,18 +407,19 @@ FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_values (const Mapping<2,3>& mapping,
 
 template <class POLY, int dim, int spacedim>
 void
-FE_Poly<POLY,dim,spacedim>::fill_fe_values (const Mapping<dim,spacedim>                   &mapping,
-                                  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                  const Quadrature<dim>                &quadrature,
-                                  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                                  typename Mapping<dim,spacedim>::InternalDataBase &fedata,
-                                  FEValuesData<dim,spacedim>                    &data) const
+FE_Poly<POLY,dim,spacedim>::fill_fe_values 
+  (const Mapping<dim,spacedim>                      &mapping,
+   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+   const Quadrature<dim>                            &quadrature,
+   const enum CellSimilarity::Similarity             cell_similarity,
+   typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+   typename Mapping<dim,spacedim>::InternalDataBase &fedata,
+   FEValuesData<dim,spacedim>                       &data) const
 {
                                   // convert data object to internal
                                   // data for this class. fails with
                                   // an exception if that is not
-                                  // possible
-        
+                                  // possible   
   Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
   InternalData &fe_data = static_cast<InternalData &> (fedata);
   
@@ -423,13 +431,13 @@ FE_Poly<POLY,dim,spacedim>::fill_fe_values (const Mapping<dim,spacedim>
        for (unsigned int i=0; i<quadrature.size(); ++i)
          data.shape_values(k,i) = fe_data.shape_values[k][i];
       
-      if (flags & update_gradients)
+      if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
        mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
                          mapping_data, mapping_covariant);
     }
   
   const typename QProjector<dim>::DataSetDescriptor dsd;
-  if (flags & update_hessians)
+  if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
     this->compute_2nd (mapping, cell, dsd.cell(),
                       mapping_data, fe_data, data);
 }
@@ -439,13 +447,14 @@ FE_Poly<POLY,dim,spacedim>::fill_fe_values (const Mapping<dim,spacedim>
 template <>
 inline
 void
-FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_face_values (const Mapping<1,2> &,
-                                                              const Triangulation<1,2>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const Quadrature<0> &, 
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              FEValuesData<1,2> &) const
+FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_face_values 
+  (const Mapping<1,2> &,
+   const Triangulation<1,2>::cell_iterator &,
+   const unsigned int,
+   const Quadrature<0> &, 
+   Mapping<1,2>::InternalDataBase &,
+   Mapping<1,2>::InternalDataBase &,
+   FEValuesData<1,2> &) const
 {
   AssertThrow(false, ExcNotImplemented());
 }  
@@ -472,12 +481,12 @@ template <>
 inline
 void
 FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_face_values (const Mapping<1,2> &,
-                                                              const Triangulation<1,2>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const Quadrature<0> &, 
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              Mapping<1,2>::InternalDataBase &,
-                                                              FEValuesData<1,2> &) const
+                                                     const Triangulation<1,2>::cell_iterator &,
+                                                     const unsigned int,
+                                                     const Quadrature<0> &, 
+                                                     Mapping<1,2>::InternalDataBase &,
+                                                     Mapping<1,2>::InternalDataBase &,
+                                                     FEValuesData<1,2> &) const
 {
   AssertThrow(false, ExcNotImplemented());
 }  
@@ -488,12 +497,12 @@ template <>
 inline
 void
 FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_face_values (const Mapping<2,3> &,
-                                                              const Triangulation<2,3>::cell_iterator &,
-                                                              const unsigned int,
-                                                              const Quadrature<1> &, 
-                                                              Mapping<2,3>::InternalDataBase &,
-                                                              Mapping<2,3>::InternalDataBase &,
-                                                              FEValuesData<2,3> &) const
+                                                     const Triangulation<2,3>::cell_iterator &,
+                                                     const unsigned int,
+                                                     const Quadrature<1> &, 
+                                                     Mapping<2,3>::InternalDataBase &,
+                                                     Mapping<2,3>::InternalDataBase &,
+                                                     FEValuesData<2,3> &) const
 {
   AssertThrow(false, ExcNotImplemented());
 }  
@@ -501,13 +510,14 @@ FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_face_values (const Mapping<2,3> &,
 
 template <class POLY, int dim, int spacedim>
 void
-FE_Poly<POLY,dim,spacedim>::fill_fe_face_values (const Mapping<dim,spacedim>                   &mapping,
-                               const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                               const unsigned int                    face,
-                               const Quadrature<dim-1>              &quadrature,
-                               typename Mapping<dim,spacedim>::InternalDataBase       &mapping_data,
-                               typename Mapping<dim,spacedim>::InternalDataBase       &fedata,
-                               FEValuesData<dim,spacedim>                    &data) const
+FE_Poly<POLY,dim,spacedim>::
+fill_fe_face_values (const Mapping<dim,spacedim>                   &mapping,
+                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                    const unsigned int                    face,
+                    const Quadrature<dim-1>              &quadrature,
+                    typename Mapping<dim,spacedim>::InternalDataBase       &mapping_data,
+                    typename Mapping<dim,spacedim>::InternalDataBase       &fedata,
+                    FEValuesData<dim,spacedim>                    &data) const
 {
                                   // convert data object to internal
                                   // data for this class. fails with
index e4dfabb33abfc53bd6006f4465634a4d20ab661d..7d641c4450612a112e0318e2ec9a9dbc39766aa3 100644 (file)
@@ -234,9 +234,10 @@ class FE_PolyTensor : public FiniteElement<dim,spacedim>
              const Quadrature<dim>& quadrature) const ;
 
     virtual void
-    fill_fe_values (const Mapping<dim,spacedim> &mapping,
+    fill_fe_values (const Mapping<dim,spacedim>                           &mapping,
                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                &quadrature,
+                   const Quadrature<dim>                                 &quadrature,
+                   const enum CellSimilarity::Similarity                  cell_similarity,
                    typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
                    FEValuesData<dim,spacedim>& data) const;
index 28a1f90548b7c4d33ab499f149c6ced15341c7a8..336f74e7ccaef87b920b1b686bec4fad7652824e 100644 (file)
@@ -611,7 +611,8 @@ class FESystem : public FiniteElement<dim,spacedim>
     virtual void
     fill_fe_values (const Mapping<dim,spacedim>                      &mapping,
                    const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                   &quadrature,
+                   const Quadrature<dim>                            &quadrature,
+                   const enum CellSimilarity::Similarity             cell_similarity,
                    typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
                    typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
                    FEValuesData<dim,spacedim>                       &data) const;
@@ -679,14 +680,15 @@ class FESystem : public FiniteElement<dim,spacedim>
                                      * <tt>sub_no!=invalid_face_no</tt>.
                                      */
     template <int dim_1>
-    void compute_fill (const Mapping<dim,spacedim>                   &mapping,
+    void compute_fill (const Mapping<dim,spacedim>                      &mapping,
                       const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                      const unsigned int                    face_no,
-                      const unsigned int                    sub_no,
-                      const Quadrature<dim_1>              &quadrature,
+                      const unsigned int                                face_no,
+                      const unsigned int                                sub_no,
+                      const Quadrature<dim_1>                          &quadrature,
+                      const enum CellSimilarity::Similarity             cell_similarity,
                       typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
                       typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
-                      FEValuesData<dim,spacedim>                    &data) const ;
+                      FEValuesData<dim,spacedim>                       &data) const ;
 
   private:
 
index 3e8b5e70c11e9318c555da1a4bad149e17aba4aa..4d77e951775ce5b7500851abf3663317cc214279 100644 (file)
@@ -2224,6 +2224,12 @@ class FEValuesBase : protected FEValuesData<dim,spacedim>,
                                      */
     const FiniteElement<dim,spacedim> & get_fe () const;
 
+                                    /**
+                                     * Return the update flags set
+                                     * for this object.
+                                     */
+    UpdateFlags get_update_flags () const;
+    
                                     /**
                                      * Return a triangulation
                                      * iterator to the current cell.
@@ -2809,6 +2815,28 @@ class FEValues : public FEValuesBase<dim,spacedim>
                                      */
     void initialize (const UpdateFlags update_flags);
 
+                                    /**
+                                     * An enum variable that can store
+                                     * different states of the current cell
+                                     * in comparison to the previously
+                                     * visited cell. If wanted, additional
+                                     * states can be checked here and used
+                                     * in one of the methods used during
+                                     * reinit.
+                                     */
+    enum CellSimilarity::Similarity cell_similarity;
+
+                                    /**
+                                     * A function that checks whether the
+                                     * new cell is similar to the one
+                                     * previously used. Then, a significant
+                                     * amount of the data can be reused,
+                                     * e.g. the derivatives of the basis
+                                     * functions in real space, shape_grad.
+                                     */
+    void
+    check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
+
                                      /**
                                       * The reinit() functions do
                                       * only that part of the work
@@ -4083,6 +4111,15 @@ FEValuesBase<dim,spacedim>::get_mapping () const
 
 template <int dim, int spacedim>  
 inline
+UpdateFlags
+FEValuesBase<dim,spacedim>::get_update_flags () const
+{
+  return this->update_flags;
+}
+
+
+
+template <int dim, int spacedim>  
 const typename Triangulation<dim,spacedim>::cell_iterator
 FEValuesBase<dim,spacedim>::get_cell () const
 {
index 442bd5441ba8aa00779839ba9c97b7e6128b881e..cfbaebbc29ac4e4c61649916044960b45c449de9 100644 (file)
@@ -32,6 +32,24 @@ template <int dim, int spacedim> class FEValues;
 template <int dim, int spacedim> class FEFaceValues;
 template <int dim, int spacedim> class FESubfaceValues;
 
+                                    /**
+                                     * An enum variable that can store
+                                     * different states of the current cell
+                                     * in comparison to the previously
+                                     * visited cell. If wanted, additional
+                                     * states can be checked here and used
+                                     * in one of the methods used during
+                                     * reinit.
+                                     */
+namespace CellSimilarity
+{
+  enum Similarity 
+    {
+      no_similarity, 
+      translation
+    };
+}
+
                                      /**
                                       * The transformation type used
                                       * for the Mapping::transform() functions.
@@ -545,14 +563,15 @@ class Mapping : public Subscriptor
                                       */
     virtual void
     fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                         &quadrature,
-                   InternalDataBase                              &internal,
-                   std::vector<Point<spacedim> >                 &quadrature_points,
-                   std::vector<double>                           &JxW_values,
-                   std::vector<Tensor<2,spacedim> >              &jacobians,
-                   std::vector<Tensor<3,spacedim> >              &jacobian_grads,
-                   std::vector<Tensor<2,spacedim> >              &inverse_jacobians,
-                   std::vector<Point<spacedim> >                 &cell_normal_vectors
+                   const Quadrature<dim>                                     &quadrature,
+                   const enum CellSimilarity::Similarity                                  cell_similarity,
+                   InternalDataBase                                          &internal,
+                   std::vector<Point<spacedim> >                             &quadrature_points,
+                   std::vector<double>                                       &JxW_values,
+                   std::vector<Tensor<2,spacedim> >                          &jacobians,
+                   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+                   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
+                   std::vector<Point<spacedim> >                             &cell_normal_vectors
                   ) const=0;
 
 
index dff1212763e8bc221e6dc3bc90cf272406d01423..9b094b6ef862bbb5b7eb67b7db4fe6c446a6ecb1 100644 (file)
@@ -60,13 +60,14 @@ class MappingCartesian : public Mapping<dim,spacedim>
 
     virtual void
     fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>& quadrature,
-                   typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-                   std::vector<Point<spacedim> >                     &quadrature_points,
-                   std::vector<double>                               &JxW_values,
-                   std::vector<Tensor<2,spacedim> >                  &jacobians,
-                   std::vector<Tensor<3,spacedim> >                  &jacobian_grads,
-                   std::vector<Tensor<2,spacedim> >                  &inverse_jacobians,
+                   const Quadrature<dim>                                     &quadrature,
+                   const enum CellSimilarity::Similarity      cell_similarity,
+                   typename Mapping<dim, spacedim>::InternalDataBase         &mapping_data,
+                   std::vector<Point<spacedim> >                             &quadrature_points,
+                   std::vector<double>                                       &JxW_values,
+                   std::vector<Tensor<2,spacedim> >                          &jacobians,
+                   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+                   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
                    std::vector<Point<spacedim> >   &) const ;
 
 
@@ -180,6 +181,7 @@ class MappingCartesian : public Mapping<dim,spacedim>
     void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                       const unsigned int face_no,
                       const unsigned int sub_no,
+                      const enum CellSimilarity::Similarity cell_similarity,
                       InternalData& data,
                       std::vector<Point<dim> > &quadrature_points,
                       std::vector<Point<dim> >& normal_vectors) const;
index 24c8d11e6019fe89844970ec9112191bd6984869..038144041b26479d39ff5a7709b1fe638fef95ab 100644 (file)
@@ -205,6 +205,7 @@ class MappingQ : public MappingQ1<dim,spacedim>
     virtual void
     fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                    const Quadrature<dim>                                     &quadrature,
+                   enum CellSimilarity::Similarity                            cell_similarity,
                    typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
                    typename std::vector<Point<spacedim> >                    &quadrature_points,
                    std::vector<double>                                       &JxW_values,
index cfc67ab395da72a10c34e3660c85b47e109e9e9f..f21b8eee9a8c56e051fa8373ff6287b0a834da1e 100644 (file)
@@ -294,7 +294,8 @@ class MappingQ1 : public Mapping<dim,spacedim>
                                      */
     virtual void
     fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                   const Quadrature<dim>                                     & quadrature,
+                   const Quadrature<dim>                                     &quadrature,
+                   const enum CellSimilarity::Similarity           cell_similarity,
                    typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
                    typename std::vector<Point<spacedim> >                    &quadrature_points,
                    std::vector<double>                                       &JxW_values,
@@ -384,6 +385,7 @@ class MappingQ1 : public Mapping<dim,spacedim>
     void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                       const unsigned int      npts,
                       const DataSetDescriptor data_set,
+                      const enum CellSimilarity::Similarity cell_similarity,
                       InternalData           &data,
                       std::vector<Point<spacedim> > &quadrature_points) const;
     
index 949f37c599a272449bb2ef656bd2d3801029937f..540c377bfcaca782edcef3954589c31baab84e7d 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2009 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -143,6 +143,23 @@ class MappingQ1Eulerian : public MappingQ1<dim,spacedim>
 
 
   protected:
+                                    /**
+                                     * Implementation of the interface in
+                                     * MappingQ1. Overrides the function in
+                                     * the base class, since we cannot use
+                                     * any cell similarity for this class.
+                                     */
+    virtual void
+    fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                   const Quadrature<dim>                                     &quadrature,
+                   enum CellSimilarity::Similarity                            cell_similarity,
+                   typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
+                   typename std::vector<Point<spacedim> >                    &quadrature_points,
+                   std::vector<double>                                       &JxW_values,
+                   std::vector<Tensor<2,spacedim> >                          &jacobians,
+                   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+                   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
+                   std::vector<Point<spacedim> >                             &cell_normal_vectors) const;
 
                                     /**
                                      * Reference to the vector of
index c5517c31f6d29ebfa5ac65529c67718fef9a8846..809bab67f8d0ea74909c1394919cd9d0247bce82 100644 (file)
@@ -146,6 +146,24 @@ class MappingQEulerian : public MappingQ<dim, spacedim>
                    << "-- expected size " << arg2);
 
   protected:
+                                    /**
+                                     * Implementation of the interface in
+                                     * MappingQ. Overrides the function in
+                                     * the base class, since we cannot use
+                                     * any cell similarity for this class.
+                                     */
+    virtual void
+    fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                   const Quadrature<dim>                                     &quadrature,
+                   enum CellSimilarity::Similarity                            cell_similarity,
+                   typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
+                   typename std::vector<Point<spacedim> >                    &quadrature_points,
+                   std::vector<double>                                       &JxW_values,
+                   std::vector<Tensor<2,spacedim> >                          &jacobians,
+                   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+                   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
+                   std::vector<Point<spacedim> >                             &cell_normal_vectors) const;
+
                                      /**
                                       * Reference to the vector of
                                       * shifts.
index 5cf951458a17b1cc4380223a6d0f46a9abbc36ba..e80074c73a96778de9d86d13edd3c517a62c9de9 100644 (file)
@@ -289,6 +289,7 @@ FE_DGPNonparametric<dim,spacedim>::fill_fe_values (
   const Mapping<dim,spacedim>&,
   const typename Triangulation<dim,spacedim>::cell_iterator&,
   const Quadrature<dim>&,
+  const enum CellSimilarity::Similarity cell_similarity,
   typename Mapping<dim,spacedim>::InternalDataBase&,
   typename Mapping<dim,spacedim>::InternalDataBase& fedata,
   FEValuesData<dim,spacedim>& data) const
@@ -297,13 +298,15 @@ FE_DGPNonparametric<dim,spacedim>::fill_fe_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
-  
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
+
   const UpdateFlags flags(fe_data.current_update_flags());
   Assert (flags & update_quadrature_points, ExcInternalError());
-  
+
   const unsigned int n_q_points = data.quadrature_points.size();
-  
+
   if (flags & (update_values | update_gradients))
     for (unsigned int i=0; i<n_q_points; ++i)
       {
@@ -338,7 +341,9 @@ FE_DGPNonparametric<dim,spacedim>::fill_fe_face_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
 
   const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
   Assert (flags & update_quadrature_points, ExcInternalError());
@@ -380,7 +385,9 @@ FE_DGPNonparametric<dim,spacedim>::fill_fe_subface_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
   
   const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
   Assert (flags & update_quadrature_points, ExcInternalError());
index 71901c7b08ecad7b452f57ee5e3525a9abbb4274..75ad0b8d5a3590698e2343f1e6b77a5ce366fe62 100644 (file)
@@ -990,18 +990,22 @@ FE_Nedelec<dim,spacedim>::get_data (const UpdateFlags      update_flags,
 
 template <int dim, int spacedim>
 void
-FE_Nedelec<dim,spacedim>::fill_fe_values (const Mapping<dim,spacedim>                   &mapping,
-                                const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                const Quadrature<dim>                &quadrature,
-                                typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                                typename Mapping<dim,spacedim>::InternalDataBase &fedata,
-                                FEValuesData<dim,spacedim>                    &data) const
+FE_Nedelec<dim,spacedim>::fill_fe_values 
+  (const Mapping<dim,spacedim>                      &mapping,
+   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+   const Quadrature<dim>                            &quadrature,
+   const enum CellSimilarity::Similarity             cell_similarity,
+   typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+   typename Mapping<dim,spacedim>::InternalDataBase &fedata,
+   FEValuesData<dim,spacedim>                       &data) const
 {
                                   // convert data object to internal
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
 
                                   // get the flags indicating the
                                   // fields that have to be filled
index 04aa19b5016c3cdb9e0c3cecf35da8e663b24b1a..ebc257d682f79ede7c9c1ebd84129d7163b420e8 100644 (file)
@@ -359,18 +359,21 @@ FE_PolyTensor<POLY,dim,spacedim>::get_data (
 template <class POLY, int dim, int spacedim>
 void
 FE_PolyTensor<POLY,dim,spacedim>::fill_fe_values (
-  const Mapping<dim,spacedim>                   &mapping,
+  const Mapping<dim,spacedim>                      &mapping,
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const Quadrature<dim>                &quadrature,
+  const Quadrature<dim>                            &quadrature,
+  const enum CellSimilarity::Similarity             cell_similarity,
   typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
   typename Mapping<dim,spacedim>::InternalDataBase &fedata,
-  FEValuesData<dim,spacedim>                    &data) const
+  FEValuesData<dim,spacedim>                       &data) const
 {
                                   // convert data object to internal
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
 
   const unsigned int n_q_points = quadrature.size();
   const UpdateFlags flags(fe_data.current_update_flags());
@@ -513,7 +516,9 @@ FE_PolyTensor<POLY,dim,spacedim>::fill_fe_face_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
 
   const unsigned int n_q_points = quadrature.size();
                                   // offset determines which data set
@@ -666,7 +671,9 @@ FE_PolyTensor<POLY,dim,spacedim>::fill_fe_subface_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+         ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
 
   const unsigned int n_q_points = quadrature.size();
                                   // offset determines which data set
index 0290b28e697036c574746ae68a61d29eab79705c..9989fb6ea830a2df759fd0d8c1dfb97331b60c30 100644 (file)
@@ -146,7 +146,7 @@ const unsigned int FESystem<dim,spacedim>::invalid_face_number;
 
 template <int dim, int spacedim>
 FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe,
-                        const unsigned int n_elements) :
+                                 const unsigned int n_elements) :
                FiniteElement<dim,spacedim> (multiply_dof_numbers(fe, n_elements),
                                    compute_restriction_is_additive_flags (fe, n_elements),
                                    compute_nonzero_components(fe, n_elements)),
@@ -161,9 +161,9 @@ FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe,
 
 template <int dim, int spacedim>
 FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
-                        const unsigned int        n1,
-                        const FiniteElement<dim,spacedim> &fe2,
-                        const unsigned int        n2) :
+                                 const unsigned int        n1,
+                                 const FiniteElement<dim,spacedim> &fe2,
+                                 const unsigned int        n2) :
                FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1, fe2, n2),
                                    compute_restriction_is_additive_flags (fe1, n1,
                                                                           fe2, n2),
@@ -183,12 +183,12 @@ FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
 
 template <int dim, int spacedim>
 FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
-                        const unsigned int        n1,
-                        const FiniteElement<dim,spacedim> &fe2,
-                        const unsigned int        n2,
-                        const FiniteElement<dim,spacedim> &fe3,
-                        const unsigned int        n3) :
-               FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1,
+                                 const unsigned int        n1,
+                                 const FiniteElement<dim,spacedim> &fe2,
+                                 const unsigned int        n2,
+                                 const FiniteElement<dim,spacedim> &fe3,
+                                 const unsigned int        n3) :
+       FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1,
                                                         fe2, n2,
                                                         fe3, n3),
                                    compute_restriction_is_additive_flags (fe1, n1,
@@ -288,7 +288,7 @@ FESystem<dim,spacedim>::clone() const
 template <int dim, int spacedim>
 double
 FESystem<dim,spacedim>::shape_value (const unsigned int i,
-                           const Point<dim> &p) const
+                                    const Point<dim> &p) const
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (this->is_primitive(i), 
@@ -303,8 +303,8 @@ FESystem<dim,spacedim>::shape_value (const unsigned int i,
 template <int dim, int spacedim>
 double
 FESystem<dim,spacedim>::shape_value_component (const unsigned int i,
-                                     const Point<dim>  &p,
-                                     const unsigned int component) const
+                                              const Point<dim>  &p,
+                                              const unsigned int component) const
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (component < this->n_components(),
@@ -340,7 +340,7 @@ FESystem<dim,spacedim>::shape_value_component (const unsigned int i,
 template <int dim, int spacedim>
 Tensor<1,dim>
 FESystem<dim,spacedim>::shape_grad (const unsigned int i,
-                          const Point<dim> &p) const
+                                   const Point<dim> &p) const
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (this->is_primitive(i),
@@ -392,7 +392,7 @@ FESystem<dim,spacedim>::shape_grad_component (const unsigned int i,
 template <int dim, int spacedim>
 Tensor<2,dim>
 FESystem<dim,spacedim>::shape_grad_grad (const unsigned int i,
-                               const Point<dim> &p) const
+                                        const Point<dim> &p) const
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (this->is_primitive(i), 
@@ -577,8 +577,8 @@ FESystem<dim,spacedim>::update_each (const UpdateFlags flags) const
 template <int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
 FESystem<dim,spacedim>::get_data (const UpdateFlags      flags_,
-                        const Mapping<dim,spacedim>    &mapping,
-                        const Quadrature<dim> &quadrature) const
+                                 const Mapping<dim,spacedim>    &mapping,
+                                 const Quadrature<dim> &quadrature) const
 {
   UpdateFlags flags = flags_;
   
@@ -664,15 +664,16 @@ FESystem<dim,spacedim>::get_data (const UpdateFlags      flags_,
 template <int dim, int spacedim>
 void
 FESystem<dim,spacedim>::
-fill_fe_values (const Mapping<dim,spacedim>                   &mapping,
+fill_fe_values (const Mapping<dim,spacedim>                      &mapping,
                 const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                const Quadrature<dim>                &quadrature,
+                const Quadrature<dim>                            &quadrature,
+               const enum CellSimilarity::Similarity             cell_similarity,
                 typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
                 typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
-                FEValuesData<dim,spacedim>                    &data) const
+                FEValuesData<dim,spacedim>                       &data) const
 {
   compute_fill(mapping, cell, invalid_face_number, invalid_face_number,
-              quadrature, mapping_data, fe_data, data);
+              quadrature, cell_similarity, mapping_data, fe_data, data);
 }
 
 
@@ -688,8 +689,8 @@ fill_fe_face_values (const Mapping<dim,spacedim>                   &mapping,
                      typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
                      FEValuesData<dim,spacedim>                    &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,
+                CellSimilarity::no_similarity, mapping_data, fe_data, data);
 }
 
 
@@ -698,17 +699,17 @@ fill_fe_face_values (const Mapping<dim,spacedim>                   &mapping,
 template <int dim, int spacedim>
 void
 FESystem<dim,spacedim>::
-fill_fe_subface_values (const Mapping<dim,spacedim>                   &mapping,
+fill_fe_subface_values (const Mapping<dim,spacedim>                      &mapping,
                         const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                        const unsigned int                    face_no,
-                        const unsigned int                    sub_no,
-                        const Quadrature<dim-1>              &quadrature,
+                        const unsigned int                                face_no,
+                        const unsigned int                                sub_no,
+                        const Quadrature<dim-1>                          &quadrature,
                         typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
                         typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
-                        FEValuesData<dim,spacedim>                    &data) const
+                        FEValuesData<dim,spacedim>                       &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,
+                CellSimilarity::no_similarity, mapping_data, fe_data, data);
 }
 
 
@@ -717,14 +718,15 @@ template <int dim, int spacedim>
 template <int dim_1>
 void
 FESystem<dim,spacedim>::
-compute_fill (const Mapping<dim,spacedim>                   &mapping,
+compute_fill (const Mapping<dim,spacedim>                      &mapping,
               const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-              const unsigned int                    face_no,
-              const unsigned int                    sub_no,
-              const Quadrature<dim_1>              &quadrature,
+              const unsigned int                                face_no,
+              const unsigned int                                sub_no,
+              const Quadrature<dim_1>                          &quadrature,
+             const enum CellSimilarity::Similarity             cell_similarity,
               typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
               typename Mapping<dim,spacedim>::InternalDataBase &fedata,
-              FEValuesData<dim,spacedim>                    &data) const
+              FEValuesData<dim,spacedim>                       &data) const
 {       
   const unsigned int n_q_points = quadrature.size();
   
@@ -732,7 +734,8 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData & fe_data = dynamic_cast<InternalData&> (fedata);
+  Assert (dynamic_cast<InternalData*> (&fedata) != 0, ExcInternalError());
+  InternalData & fe_data = static_cast<InternalData&> (fedata);
   
                                   // Either dim_1==dim
                                   // (fill_fe_values) or dim_1==dim-1
@@ -800,16 +803,18 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
       if (face_no==invalid_face_number)
        {
          Assert(dim_1==dim, ExcDimensionMismatch(dim_1,dim));
+         Assert (dynamic_cast<const Quadrature<dim> *>(quadrature_base_pointer) != 0, 
+                 ExcInternalError());
          cell_quadrature
-            = dynamic_cast<const Quadrature<dim> *>(quadrature_base_pointer);
-         Assert (cell_quadrature != 0, ExcInternalError());
+            = static_cast<const Quadrature<dim> *>(quadrature_base_pointer);
        }
       else
        {
          Assert(dim_1==dim-1, ExcDimensionMismatch(dim_1,dim-1));
+         Assert (dynamic_cast<const Quadrature<dim-1> *>(quadrature_base_pointer) != 0, 
+                 ExcInternalError());
          face_quadrature
-            = dynamic_cast<const Quadrature<dim-1> *>(quadrature_base_pointer);
-         Assert (face_quadrature != 0, ExcInternalError());
+            = static_cast<const Quadrature<dim-1> *>(quadrature_base_pointer);
        }
 
                                        // let base elements update the
@@ -852,8 +857,8 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
                                            // therefore use
                                            // 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);
+           base_fe.fill_fe_values(mapping, cell, *cell_quadrature, cell_similarity, 
+                                  mapping_data, base_fe_data, base_data);
          else if (sub_no==invalid_face_number)
            base_fe.fill_fe_face_values(mapping, cell, face_no,
                                        *face_quadrature, mapping_data, base_fe_data, base_data);
@@ -892,14 +897,21 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
                                            // 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)
-              {
-                const unsigned int
-                  base_index = this->system_to_base_table[system_index].second;
-                Assert (base_index<base_fe.dofs_per_cell, ExcInternalError());
+                                       base_fe_data.update_flags);
+
+                                          // 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 enter this section
+         if (cell_similarity != CellSimilarity::translation)
+           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)
+               {
+                 const unsigned int
+                   base_index = this->system_to_base_table[system_index].second;
+                 Assert (base_index<base_fe.dofs_per_cell, ExcInternalError());
 
                                                  // now copy. if the
                                                  // shape function is
@@ -917,30 +929,30 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
                                                  // this one value, and
                                                  // to which index to
                                                  // put
-                unsigned int out_index = 0;
-                for (unsigned int i=0; i<system_index; ++i)
-                  out_index += this->n_nonzero_components(i);
-                unsigned int in_index = 0;
-                for (unsigned int i=0; i<base_index; ++i)
-                  in_index += base_fe.n_nonzero_components(i);
+                 unsigned int out_index = 0;
+                 for (unsigned int i=0; i<system_index; ++i)
+                   out_index += this->n_nonzero_components(i);
+                 unsigned int in_index = 0;
+                 for (unsigned int i=0; i<base_index; ++i)
+                   in_index += base_fe.n_nonzero_components(i);
                 
                                                  // then loop over the
                                                  // number of components
                                                  // to be copied
-                Assert (this->n_nonzero_components(system_index) ==
-                        base_fe.n_nonzero_components(base_index),
-                        ExcInternalError());
-                for (unsigned int s=0; s<this->n_nonzero_components(system_index); ++s)
-                  {
-                    if (base_flags & update_values)
-                      for (unsigned int q=0; q<n_q_points; ++q)
-                        data.shape_values[out_index+s][q] =
-                          base_data.shape_values(in_index+s,q);
+                 Assert (this->n_nonzero_components(system_index) ==
+                         base_fe.n_nonzero_components(base_index),
+                         ExcInternalError());
+                 for (unsigned int s=0; s<this->n_nonzero_components(system_index); ++s)
+                   {
+                     if (base_flags & update_values)
+                       for (unsigned int q=0; q<n_q_points; ++q)
+                         data.shape_values[out_index+s][q] =
+                           base_data.shape_values(in_index+s,q);
                     
-                    if (base_flags & update_gradients)
-                      for (unsigned int q=0; q<n_q_points; ++q)
-                        data.shape_gradients[out_index+s][q]=
-                          base_data.shape_gradients[in_index+s][q];
+                     if (base_flags & update_gradients)
+                       for (unsigned int q=0; q<n_q_points; ++q)
+                         data.shape_gradients[out_index+s][q]=
+                           base_data.shape_gradients[in_index+s][q];
 
                                                      // _we_ handle
                                                      // computation of
@@ -951,10 +963,10 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
                                                      // should not
                                                      // have computed
                                                      // them!
-                    Assert (!(base_flags & update_hessians),
-                            ExcInternalError());
-                  };
-              };
+                     Assert (!(base_flags & update_hessians),
+                             ExcInternalError());
+                   };
+               };
         };
 
       if (fe_data.is_first_cell())
@@ -981,7 +993,7 @@ compute_fill (const Mapping<dim,spacedim>                   &mapping,
        }
     }
   
-  if (fe_data.compute_hessians)
+  if (fe_data.compute_hessians && cell_similarity != CellSimilarity::translation)
     {
       unsigned int offset = 0;
       if (face_no != invalid_face_number)
index f70ba8aa1c4208c51b711d1daeb62690c0c7679d..0cb2bd6f11bc6e312e01d2f355b08709d6317c58 100644 (file)
@@ -1690,6 +1690,7 @@ FEValuesData<dim,spacedim>::initialize (const unsigned int        n_quadrature_p
 }
 
 
+
 /*------------------------------- FEValuesBase ---------------------------*/
 
 
@@ -3238,7 +3239,6 @@ FEValuesBase<dim,spacedim>::compute_update_flags (const UpdateFlags update_flags
 }
 
 
-
 /*------------------------------- FEValues -------------------------------*/
 
 
@@ -3252,9 +3252,9 @@ const unsigned int FEValues<dim,spacedim>::integral_dimension;
 
 template <int dim, int spacedim>
 FEValues<dim,spacedim>::FEValues (const Mapping<dim,spacedim>       &mapping,
-                        const FiniteElement<dim,spacedim> &fe,
-                        const Quadrature<dim>    &q,
-                        const UpdateFlags         update_flags)
+                                 const FiniteElement<dim,spacedim> &fe,
+                                 const Quadrature<dim>             &q,
+                                 const UpdateFlags                  update_flags)
                :
                FEValuesBase<dim,spacedim> (q.size(),
                                            fe.dofs_per_cell,
@@ -3270,15 +3270,15 @@ FEValues<dim,spacedim>::FEValues (const Mapping<dim,spacedim>       &mapping,
 
 template <int dim, int spacedim>
 FEValues<dim,spacedim>::FEValues (const FiniteElement<dim,spacedim> &fe,
-                        const Quadrature<dim>    &q,
-                        const UpdateFlags         update_flags)
+                                 const Quadrature<dim>             &q,
+                                 const UpdateFlags                  update_flags)
                :
                FEValuesBase<dim,spacedim> (q.size(),
-                                  fe.dofs_per_cell,
-                                  update_default,
+                                           fe.dofs_per_cell,
+                                           update_default,
                                            StaticMappingQ1<dim,spacedim>::mapping,
-                                  fe),
-  quadrature (q)
+                                           fe),
+               quadrature (q)
 {
   Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
   initialize (update_flags);
@@ -3330,9 +3330,11 @@ FEValues<dim,spacedim>::reinit (const typename DoFHandler<dim,spacedim>::cell_it
 
   typedef FEValuesBase<dim,spacedim> FEVB;
   Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
-         static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
+         static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
          typename FEVB::ExcFEDontMatch());
 
+  check_cell_similarity(cell);
+
                                    // set new cell. auto_ptr will take
                                    // care that old object gets
                                    // destroyed and also that this
@@ -3360,11 +3362,13 @@ FEValues<dim,spacedim>::reinit (const typename hp::DoFHandler<dim,spacedim>::cel
                                   // passed to the constructor and
                                   // used by the DoFHandler used by
                                   // this cell, are the same
-    typedef FEValuesBase<dim,spacedim> FEVB;
+  typedef FEValuesBase<dim,spacedim> FEVB;
   Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
          static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
          typename FEVB::ExcFEDontMatch());
 
+  check_cell_similarity(cell);
+
                                    // set new cell. auto_ptr will take
                                    // care that old object gets
                                    // destroyed and also that this
@@ -3420,6 +3424,8 @@ FEValues<dim,spacedim>::reinit (const typename MGDoFHandler<dim,spacedim>::cell_
 //       static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
 //       typename FEValuesBase<dim,spacedim>::ExcFEDontMatch());
 
+  check_cell_similarity(cell);
+
                                    // set new cell. auto_ptr will take
                                    // care that old object gets
                                    // destroyed and also that this
@@ -3442,8 +3448,9 @@ FEValues<dim,spacedim>::reinit (const typename MGDoFHandler<dim,spacedim>::cell_
 template <int dim, int spacedim>
 void FEValues<dim,spacedim>::reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell)
 {
-                                   // no FE in this cell, so no check
+                                   // no FE in this cell, so no assertion
                                    // necessary here
+  check_cell_similarity(cell);
 
                                    // set new cell. auto_ptr will take
                                    // care that old object gets
@@ -3462,11 +3469,59 @@ void FEValues<dim,spacedim>::reinit (const typename Triangulation<dim,spacedim>:
 
 
 
+template <int dim, int spacedim>
+inline
+void
+FEValues<dim,spacedim>::check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell)
+{
+                                  // case that there has not been any cell
+                                  // before
+  if (&*this->present_cell == 0)
+    {
+      cell_similarity = CellSimilarity::no_similarity;
+      return;
+    }
+
+  const typename Triangulation<dim,spacedim>::cell_iterator & present_cell = 
+    *this->present_cell;
+
+                                  // test for translation
+  {
+                                  // otherwise, go through the vertices and
+                                  // check...
+    bool is_translation = true;
+    const Point<spacedim> dist = cell->vertex(0) - present_cell->vertex(0);
+    for (unsigned int i=1; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+      {
+       Point<spacedim> dist_new = cell->vertex(i) - present_cell->vertex(i) - dist;
+       if (dist_new.norm_square() > 1e-28)
+         {
+           is_translation = false;
+           break;
+         }
+      }
+
+    cell_similarity = (is_translation == true 
+                      ? 
+                      CellSimilarity::translation 
+                      : 
+                      CellSimilarity::no_similarity);
+    return;
+  }
+
+                                  // TODO: here, one could implement other
+                                  // checks for similarity, e.g. for
+                                  // children of a parallelogram.
+}
+
+
+
 template <int dim, int spacedim>
 void FEValues<dim,spacedim>::do_reinit ()
 {
   this->get_mapping().fill_fe_values(*this->present_cell,
                                     quadrature,
+                                    this->cell_similarity,
                                     *this->mapping_data,
                                     this->quadrature_points,
                                     this->JxW_values,
@@ -3478,6 +3533,7 @@ void FEValues<dim,spacedim>::do_reinit ()
   this->get_fe().fill_fe_values(this->get_mapping(),
                                *this->present_cell,
                                quadrature,
+                               this->cell_similarity,
                                *this->mapping_data,
                                *this->fe_data,
                                *this);
@@ -3502,11 +3558,11 @@ FEValues<dim,spacedim>::memory_consumption () const
 
 template <int dim, int spacedim>
 FEFaceValuesBase<dim,spacedim>::FEFaceValuesBase (const unsigned int n_q_points,
-                                        const unsigned int dofs_per_cell,
-                                        const UpdateFlags,
-                                        const Mapping<dim,spacedim> &mapping,      
-                                        const FiniteElement<dim,spacedim> &fe,
-                                        const Quadrature<dim-1>& quadrature)
+                                                 const unsigned int dofs_per_cell,
+                                                 const UpdateFlags,
+                                                 const Mapping<dim,spacedim> &mapping,      
+                                                 const FiniteElement<dim,spacedim> &fe,
+                                                 const Quadrature<dim-1>& quadrature)
                :
                FEValuesBase<dim,spacedim> (n_q_points,
                                   dofs_per_cell,
index f908b2bb8bdf422087ed2389e7e28c209f4268a5..ec021ded6c2d5fb7fac6a32b5e1d3137b1f4b074 100644 (file)
@@ -79,7 +79,7 @@ MappingCartesian<dim, spacedim>::update_each (const UpdateFlags in) const
 template<int dim, int spacedim>
 typename Mapping<dim, spacedim>::InternalDataBase *
 MappingCartesian<dim, spacedim>::get_data (const UpdateFlags      update_flags,
-                                const Quadrature<dim> &q) const
+                                          const Quadrature<dim> &q) const
 {
   InternalData* data = new InternalData (q);
 
@@ -95,7 +95,7 @@ MappingCartesian<dim, spacedim>::get_data (const UpdateFlags      update_flags,
 template<int dim, int spacedim>
 typename Mapping<dim, spacedim>::InternalDataBase *
 MappingCartesian<dim, spacedim>::get_face_data (const UpdateFlags update_flags,
-                                     const Quadrature<dim-1>& quadrature) const
+                                               const Quadrature<dim-1>& quadrature) const
 {
   InternalData* data
     = new InternalData (QProjector<dim>::project_to_all_faces(quadrature));
@@ -112,7 +112,7 @@ MappingCartesian<dim, spacedim>::get_face_data (const UpdateFlags update_flags,
 template<int dim, int spacedim>
 typename Mapping<dim, spacedim>::InternalDataBase *
 MappingCartesian<dim, spacedim>::get_subface_data (const UpdateFlags update_flags,
-                                        const Quadrature<dim-1> &quadrature) const
+                                                  const Quadrature<dim-1> &quadrature) const
 {
   InternalData* data
     = new InternalData (QProjector<dim>::project_to_all_subfaces(quadrature));
@@ -130,11 +130,12 @@ MappingCartesian<dim, spacedim>::get_subface_data (const UpdateFlags update_flag
 template<int dim, int spacedim>
 void
 MappingCartesian<dim, spacedim>::compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                                    const unsigned int        face_no,
-                                    const unsigned int        sub_no,
-                                    InternalData             &data,
-                                    std::vector<Point<dim> > &quadrature_points,
-                                    std::vector<Point<dim> > &normal_vectors) const
+                                              const unsigned int        face_no,
+                                              const unsigned int        sub_no,
+                                              const enum CellSimilarity::Similarity cell_similarity,
+                                              InternalData             &data,
+                                              std::vector<Point<dim> > &quadrature_points,
+                                              std::vector<Point<dim> > &normal_vectors) const
 {
   const UpdateFlags update_flags(data.current_update_flags());
 
@@ -299,6 +300,7 @@ void
 MappingCartesian<dim, spacedim>::
 fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator& cell,
                 const Quadrature<dim>& q, 
+               const enum CellSimilarity::Similarity cell_similarity,
                typename Mapping<dim,spacedim>::InternalDataBase& mapping_data,
                 std::vector<Point<spacedim> >& quadrature_points,
                 std::vector<double>& JxW_values,
@@ -311,11 +313,12 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator& cell,
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+  Assert (dynamic_cast<InternalData*> (&mapping_data) != 0, ExcInternalError());
+  InternalData &data = static_cast<InternalData&> (mapping_data);
 
   std::vector<Point<dim> > dummy;
   
-  compute_fill (cell, invalid_face_number, invalid_face_number,
+  compute_fill (cell, invalid_face_number, invalid_face_number, cell_similarity,
                data,
                quadrature_points,
                dummy);
@@ -326,40 +329,44 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator& cell,
                                   // local lengths in each coordinate
                                   // direction
   if (data.current_update_flags() & (update_JxW_values | update_volume_elements))
-    {
-      double J = data.length[0];
-      for (unsigned int d=1;d<dim;++d)
-       J *= data.length[d];
-      data.volume_element = J;
-      if (data.current_update_flags() & update_JxW_values)
-       for (unsigned int i=0; i<JxW_values.size();++i)
-         JxW_values[i] = J * q.weight(i);
-    }
+    if (cell_similarity != CellSimilarity::translation)
+      {
+       double J = data.length[0];
+       for (unsigned int d=1;d<dim;++d)
+         J *= data.length[d];
+       data.volume_element = J;
+       if (data.current_update_flags() & update_JxW_values)
+         for (unsigned int i=0; i<JxW_values.size();++i)
+           JxW_values[i] = J * q.weight(i);
+      }
                                   // "compute" Jacobian at the quadrature
                                   // points, which are all the same
   if (data.current_update_flags() & update_jacobians)
-    for (unsigned int i=0; i<jacobians.size();++i)
-      {
-       jacobians[i]=Tensor<2,dim>();
-       for (unsigned int j=0; j<dim; ++j)
-         jacobians[j][j]=data.length[j];
-      }
+    if (cell_similarity != CellSimilarity::translation)
+      for (unsigned int i=0; i<jacobians.size();++i)
+       {
+         jacobians[i]=Tensor<2,dim>();
+         for (unsigned int j=0; j<dim; ++j)
+           jacobians[j][j]=data.length[j];
+       }
                                   // "compute" the derivative of the Jacobian
                                   // at the quadrature points, which are all
                                   // zero of course
   if (data.current_update_flags() & update_jacobian_grads)
-    for (unsigned int i=0; i<jacobian_grads.size();++i)
-      jacobian_grads[i]=Tensor<3,dim>();
+    if (cell_similarity != CellSimilarity::translation)
+      for (unsigned int i=0; i<jacobian_grads.size();++i)
+       jacobian_grads[i]=Tensor<3,dim>();
                                   // "compute" inverse Jacobian at the 
                                   // quadrature points, which are all 
                                   // the same
   if (data.current_update_flags() & update_inverse_jacobians)
-    for (unsigned int i=0; i<inverse_jacobians.size();++i)
-      {
-       inverse_jacobians[i]=Tensor<2,dim>();
-       for (unsigned int j=0; j<dim; ++j)
-         inverse_jacobians[j][j]=1./data.length[j];
-      }
+    if (cell_similarity != CellSimilarity::translation)
+      for (unsigned int i=0; i<inverse_jacobians.size();++i)
+       {
+         inverse_jacobians[i]=Tensor<2,dim>();
+         for (unsigned int j=0; j<dim; ++j)
+           inverse_jacobians[j][j]=1./data.length[j];
+       }
 }
 
 
@@ -419,9 +426,12 @@ MappingCartesian<dim, spacedim>::fill_fe_face_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+  Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
+         ExcInternalError());
+  InternalData &data = static_cast<InternalData&> (mapping_data);
 
-  compute_fill (cell, face_no, invalid_face_number,
+  compute_fill (cell, face_no, invalid_face_number, 
+               CellSimilarity::no_similarity,
                data,
                quadrature_points,
                normal_vectors);
@@ -471,9 +481,10 @@ MappingCartesian<dim, spacedim>::fill_fe_subface_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+  Assert (dynamic_cast<InternalData*> (&mapping_data) != 0, ExcInternalError());
+  InternalData &data = static_cast<InternalData&> (mapping_data);
 
-  compute_fill (cell, face_no, sub_no,
+  compute_fill (cell, face_no, sub_no, CellSimilarity::no_similarity,
                data,
                quadrature_points,
                normal_vectors);
@@ -535,7 +546,7 @@ MappingCartesian<dim,spacedim>::transform (
   AssertDimension (input.size(), output.size());
   Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0, 
          ExcInternalError());  
-  const InternalData &data = dynamic_cast<const InternalData&> (mapping_data);
+  const InternalData &data = static_cast<const InternalData&> (mapping_data);
   
   switch (mapping_type)
     {
index 658eb90f329b9bbeb4c9ffaa9c7c6884605e4497..fc6c6336b55ce60378aa13d8f19cd4a5d070ea81 100644 (file)
@@ -198,7 +198,7 @@ MappingQ<1>::compute_shapes_virtual (const std::vector<Point<1> > &unit_points,
 template<int dim, int spacedim>
 void
 MappingQ<dim,spacedim>::compute_shapes_virtual (const std::vector<Point<dim> > &unit_points,
-                                      typename MappingQ1<dim,spacedim>::InternalData &data) const
+                                               typename MappingQ1<dim,spacedim>::InternalData &data) const
 {
   const unsigned int n_points=unit_points.size();
   std::vector<double> values;
@@ -250,7 +250,7 @@ MappingQ<dim,spacedim>::compute_shapes_virtual (const std::vector<Point<dim> > &
 template<int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
 MappingQ<dim,spacedim>::get_data (const UpdateFlags update_flags,
-                        const Quadrature<dim> &quadrature) const
+                                 const Quadrature<dim> &quadrature) const
 {
   InternalData *data = new InternalData(n_shape_functions);
   this->compute_data (update_flags, quadrature,
@@ -266,7 +266,7 @@ MappingQ<dim,spacedim>::get_data (const UpdateFlags update_flags,
 template<int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
 MappingQ<dim,spacedim>::get_face_data (const UpdateFlags update_flags,
-                             const Quadrature<dim-1>& quadrature) const
+                                      const Quadrature<dim-1>& quadrature) const
 {
   InternalData *data = new InternalData(n_shape_functions);
   const Quadrature<dim> q (QProjector<dim>::project_to_all_faces(quadrature));
@@ -284,7 +284,7 @@ MappingQ<dim,spacedim>::get_face_data (const UpdateFlags update_flags,
 template<int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
 MappingQ<dim,spacedim>::get_subface_data (const UpdateFlags update_flags,
-                                const Quadrature<dim-1>& quadrature) const
+                                         const Quadrature<dim-1>& quadrature) const
 {
   InternalData *data = new InternalData(n_shape_functions);
   const Quadrature<dim> q (QProjector<dim>::project_to_all_subfaces(quadrature));
@@ -298,12 +298,16 @@ MappingQ<dim,spacedim>::get_subface_data (const UpdateFlags update_flags,
 }
 
 
-
+                                   // Note that the CellSimilarity flag is
+                                   // modifyable, since MappingQ can need to
+                                   // recalculate data even when cells are
+                                   // similar.
 template<int dim, int spacedim>
 void
 MappingQ<dim,spacedim>::fill_fe_values (
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
   const Quadrature<dim>                                     &q,
+  enum CellSimilarity::Similarity                            cell_similarity,
   typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
   std::vector<Point<spacedim> >                             &quadrature_points,
   std::vector<double>                                       &JxW_values,
@@ -316,7 +320,8 @@ MappingQ<dim,spacedim>::fill_fe_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+  Assert (dynamic_cast<InternalData*> (&mapping_data) != 0, ExcInternalError());
+  InternalData &data = static_cast<InternalData&> (mapping_data);
 
                                   // check whether this cell needs
                                   // the full mapping or can be
@@ -326,16 +331,30 @@ MappingQ<dim,spacedim>::fill_fe_values (
   data.use_mapping_q1_on_current_cell = !(use_mapping_q_on_all_cells
                                          || cell->has_boundary_lines());
 
-                                  // depending on this result, use
-                                  // this or the other data object
-                                  // for the mapping
+                                  // depending on this result, use this or
+                                  // the other data object for the
+                                  // mapping. furthermore, we need to
+                                  // ensure that the flag indicating
+                                  // whether we can use some similarity has
+                                  // to be modified - for a general
+                                  // MappingQ, the data needs to be
+                                  // recomputed anyway since then the
+                                  // mapping changes the data. this needs
+                                  // to be known also for later operations,
+                                  // so modify the variable here. TODO:
+                                  // still need to check how this will work
+                                  // when the previous cell disabled this
+                                  // flag.
   typename MappingQ1<dim,spacedim>::InternalData *p_data=0;
   if (data.use_mapping_q1_on_current_cell)
     p_data=&data.mapping_q1_data;
   else
-    p_data=&data;
+    {
+      p_data=&data;
+      cell_similarity = CellSimilarity::no_similarity;
+    }
   
-  MappingQ1<dim,spacedim>::fill_fe_values(cell, q, *p_data,
+  MappingQ1<dim,spacedim>::fill_fe_values(cell, q, cell_similarity, *p_data,
                                          quadrature_points, JxW_values,
                                          jacobians, jacobian_grads, inverse_jacobians,
                                          cell_normal_vectors);
@@ -359,7 +378,9 @@ MappingQ<dim,spacedim>::fill_fe_face_values (
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+  Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
+         ExcInternalError());
+  InternalData &data = static_cast<InternalData&> (mapping_data);
   
                                   // check whether this cell needs
                                   // the full mapping or can be
@@ -419,7 +440,9 @@ MappingQ<dim,spacedim>::fill_fe_subface_values (const typename Triangulation<dim
                                   // data for this class. fails with
                                   // an exception if that is not
                                   // possible
-  InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+  Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
+         ExcInternalError());
+  InternalData &data = static_cast<InternalData&> (mapping_data);
 
                                   // check whether this cell needs
                                   // the full mapping or can be
index 6ac2276306614960997417d34e05d081c4ca8720..d7b3aa9f0000c98a758049153bfa54a133867931 100644 (file)
@@ -73,7 +73,7 @@ MappingQ1<dim,spacedim>::MappingQ1 ()
 template<int dim, int spacedim>
 void
 MappingQ1<dim,spacedim>::compute_shapes (const std::vector<Point<dim> > &unit_points,
-                               InternalData &data) const
+                                        InternalData &data) const
 {
                                   // choose either the function implemented
                                   // in this class, or whatever a virtual
@@ -454,9 +454,9 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
 template<int dim, int spacedim>
 void
 MappingQ1<dim,spacedim>::compute_data (const UpdateFlags      update_flags,
-                             const Quadrature<dim> &q,
-                             const unsigned int     n_original_q_points,
-                             InternalData          &data) const
+                                      const Quadrature<dim> &q,
+                                      const unsigned int     n_original_q_points,
+                                      InternalData          &data) const
 {
   const unsigned int n_q_points = q.size();
 
@@ -504,9 +504,9 @@ MappingQ1<dim,spacedim>::get_data (const UpdateFlags update_flags,
 template<int dim, int spacedim>
 void
 MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
-                                  const Quadrature<dim>& q,
-                                  const unsigned int n_original_q_points,
-                                  InternalData& data) const
+                                           const Quadrature<dim>& q,
+                                           const unsigned int n_original_q_points,
+                                           InternalData& data) const
 {
   compute_data (update_flags, q, n_original_q_points, data);
 
@@ -572,7 +572,7 @@ MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
 template<int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
 MappingQ1<dim,spacedim>::get_face_data (const UpdateFlags        update_flags,
-                              const Quadrature<dim-1> &quadrature) const
+                                       const Quadrature<dim-1> &quadrature) const
 {
   InternalData* data = new InternalData(n_shape_functions);
   compute_face_data (update_flags,
@@ -606,33 +606,13 @@ void
 MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                                       const unsigned int  n_q_points,
                                               const DataSetDescriptor  data_set,
+                                      const enum CellSimilarity::Similarity cell_similarity,
                                       InternalData  &data,
                                       std::vector<Point<spacedim> > &quadrature_points) const
 {
 
   const UpdateFlags update_flags(data.current_update_flags());
 
-  if (update_flags & update_quadrature_points)
-    {
-      Assert (quadrature_points.size() == n_q_points,
-             ExcDimensionMismatch(quadrature_points.size(), n_q_points));
-      std::fill(quadrature_points.begin(), quadrature_points.end(),
-               Point<spacedim>());
-    }
-
-  if (update_flags & update_contravariant_transformation)
-    {
-      Assert (data.contravariant.size() == n_q_points,
-             ExcDimensionMismatch(data.contravariant.size(), n_q_points));
-      std::fill(data.contravariant.begin(), data.contravariant.end(),
-               Tensor<2,spacedim>());
-    }
-
-  if (update_flags & update_covariant_transformation)
-    {
-      Assert (data.covariant.size() == n_q_points,
-             ExcDimensionMismatch(data.covariant.size(), n_q_points));
-    }
                                   // if necessary, recompute the
                                   // support points of the
                                   // transformation of this cell
@@ -655,61 +635,89 @@ MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim
 
                                    // first compute quadrature points
   if (update_flags & update_quadrature_points)
-    for (unsigned int point=0; point<n_q_points; ++point)
-      for (unsigned int k=0; k<data.n_shape_functions; ++k)
-       quadrature_points[point] 
-         += data.shape(point+data_set,k) * data.mapping_support_points[k];
-  
+    {
+      Assert (quadrature_points.size() == n_q_points,
+             ExcDimensionMismatch(quadrature_points.size(), n_q_points));
+      std::fill(quadrature_points.begin(), quadrature_points.end(),
+               Point<spacedim>());
+
+      for (unsigned int point=0; point<n_q_points; ++point)
+       for (unsigned int k=0; k<data.n_shape_functions; ++k)
+         quadrature_points[point] 
+           += data.shape(point+data_set,k) * data.mapping_support_points[k];
+    }
+
                                    // then Jacobians
   if (update_flags & update_contravariant_transformation)
-    for (unsigned int point=0; point<n_q_points; ++point)
-      for (unsigned int k=0; k<data.n_shape_functions; ++k)
+    {
+      Assert (data.contravariant.size() == n_q_points,
+             ExcDimensionMismatch(data.contravariant.size(), n_q_points));
+
+                                  // if the current cell is just a
+                                  // translation of the previous one, no
+                                  // need to recompute jacobians...
+      if (cell_similarity != CellSimilarity::translation)
        {
+         std::fill(data.contravariant.begin(), data.contravariant.end(),
+                   Tensor<2,spacedim>());
+         for (unsigned int point=0; point<n_q_points; ++point)
+           for (unsigned int k=0; k<data.n_shape_functions; ++k)
+             {
                                   // some compilers seem to have problems
                                   // to detected that the two innermost
                                   // loops just use the data of the same
                                   // tensors, so get a reference to them by
                                   // hand
-         const Tensor<1,dim> &data_derv = data.derivative(point+data_set, k);
-         const Tensor<1,spacedim> &supp_pts = data.mapping_support_points[k];
+               const Tensor<1,dim> &data_derv = data.derivative(point+data_set, k);
+               const Tensor<1,spacedim> &supp_pts = data.mapping_support_points[k];
 
-         for (unsigned int i=0; i<spacedim; ++i)
-           for (unsigned int j=0; j<dim; ++j)
-             data.contravariant[point][i][j] += data_derv[j] * supp_pts[i];
+               for (unsigned int i=0; i<spacedim; ++i)
+                 for (unsigned int j=0; j<dim; ++j)
+                   data.contravariant[point][i][j] += data_derv[j] * supp_pts[i];
+             }
        }
+    }
 
-  if (update_flags & update_volume_elements)
-    for (unsigned int point=0; point<n_q_points; ++point)
-      data.volume_elements[point] = determinant(data.contravariant[point]);
-
+  if (update_flags & update_covariant_transformation)
+    {
+      Assert (data.covariant.size() == n_q_points,
+             ExcDimensionMismatch(data.covariant.size(), n_q_points));
+      if (cell_similarity != CellSimilarity::translation)
+       {
+         if (dim == spacedim)
                                    // invert contravariant for
                                    // covariant transformation
                                    // matrices
-  if ((update_flags & update_covariant_transformation)&&(dim==spacedim))
-    for (unsigned int point=0; point<n_q_points; ++point)
-      data.covariant[point] = invert(data.contravariant[point]);
+           for (unsigned int point=0; point<n_q_points; ++point)
+             data.covariant[point] = invert(data.contravariant[point]);
 
+         else if (dim == spacedim - 1)
+           {
                                    // CODIMENSION 1 
                                    // calculate left-inversion of the
                                    // contravariant matrices to obtain
                                    // covariant ones (auxiliary
                                    // rectangular fullmatrices are used)
-  else if ((update_flags & update_covariant_transformation)&&(dim==spacedim-1))
-    {
-      FullMatrix<double> contravariant_matrix(spacedim,dim);
-      FullMatrix<double> covariant_matrix(dim,spacedim);
+             FullMatrix<double> contravariant_matrix(spacedim,dim);
+             FullMatrix<double> covariant_matrix(dim,spacedim);
 
-      for (unsigned int point=0; point<n_q_points; ++point)
-       {
-         contravariant_matrix.
-           copy_from(data.contravariant[point],0,spacedim-1,0,dim-1);  
-         covariant_matrix.
-           left_invert(contravariant_matrix);
-         covariant_matrix.
-           copy_to(data.covariant[point],0,dim-1,0,spacedim-1);
+             for (unsigned int point=0; point<n_q_points; ++point)
+               {
+                 contravariant_matrix.
+                   copy_from(data.contravariant[point],0,spacedim-1,0,dim-1);  
+                 covariant_matrix.
+                   left_invert(contravariant_matrix);
+                 covariant_matrix.
+                   copy_to(data.covariant[point],0,dim-1,0,spacedim-1);
+               }
+           }
        }
     }
 
+  if (update_flags & update_volume_elements)
+    if (cell_similarity != CellSimilarity::translation)
+      for (unsigned int point=0; point<n_q_points; ++point)
+       data.volume_elements[point] = determinant(data.contravariant[point]);
   
 }
 
@@ -730,14 +738,15 @@ template<int dim, int spacedim>
 void
 MappingQ1<dim,spacedim>::fill_fe_values (
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const Quadrature<dim>                     &q,
-  typename Mapping<dim,spacedim>::InternalDataBase   &mapping_data,
-  std::vector<Point<spacedim> >             &quadrature_points,
-  std::vector<double>                       &JxW_values,
-  std::vector<Tensor<2,spacedim> >          &jacobians,
-  std::vector<Tensor<3,spacedim> >          &jacobian_grads,
-  std::vector<Tensor<2,spacedim> >          &inverse_jacobians,
-  std::vector<Point<spacedim> >             &cell_normal_vectors) const
+  const Quadrature<dim>                                     &q,
+  const enum CellSimilarity::Similarity                      cell_similarity,
+  typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
+  std::vector<Point<spacedim> >                             &quadrature_points,
+  std::vector<double>                                       &JxW_values,
+  std::vector<Tensor<2,spacedim> >                          &jacobians,
+  std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+  std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
+  std::vector<Point<spacedim> >                             &cell_normal_vectors) const
 {
   // ensure that the following cast is really correct:
   Assert (dynamic_cast<InternalData *>(&mapping_data) != 0, 
@@ -746,14 +755,13 @@ MappingQ1<dim,spacedim>::fill_fe_values (
 
   const unsigned int n_q_points=q.size();
 
-  compute_fill (cell, n_q_points, DataSetDescriptor::cell (),
+  compute_fill (cell, n_q_points, DataSetDescriptor::cell (), cell_similarity,
                 data, quadrature_points);
 
   
   const UpdateFlags update_flags(data.current_update_flags());
   const std::vector<double> &weights=q.get_weights();
 
-
                                   // Multiply quadrature weights by
                                   // Jacobian determinants or by norm
                                   // of the cross product of the
@@ -769,6 +777,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
       Assert( !(update_flags & update_cell_normal_vectors ) ||
              (cell_normal_vectors.size() == n_q_points),
             ExcDimensionMismatch(cell_normal_vectors.size(), n_q_points));
+      if (cell_similarity != CellSimilarity::translation)
        for (unsigned int point=0; point<n_q_points; ++point)
          { 
 
@@ -816,6 +825,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
     {      
       Assert (jacobians.size() == n_q_points,
              ExcDimensionMismatch(jacobians.size(), n_q_points));
+      if (cell_similarity != CellSimilarity::translation)
       for (unsigned int point=0; point<n_q_points; ++point)
        jacobians[point]
          = data.contravariant[point];
@@ -827,19 +837,23 @@ MappingQ1<dim,spacedim>::fill_fe_values (
     {
       Assert (jacobian_grads.size() == n_q_points,
              ExcDimensionMismatch(jacobian_grads.size(), n_q_points));
-      std::fill(jacobian_grads.begin(),
-               jacobian_grads.end(),
-               Tensor<3,spacedim>());
 
-      for (unsigned int point=0; point<n_q_points; ++point)
-       for (unsigned int k=0; k<data.n_shape_functions; ++k)
-         for (unsigned int i=0; i<dim; ++i)
-           for (unsigned int j=0; j<dim; ++j)
-             for (unsigned int l=0; l<dim; ++l)
-               jacobian_grads[point][i][j][l]
-                 += (data.second_derivative(point+DataSetDescriptor::cell (), k)[j][l]
-                     *
-                     data.mapping_support_points[k][i]);
+      if (cell_similarity != CellSimilarity::translation)
+       {
+         std::fill(jacobian_grads.begin(),
+                   jacobian_grads.end(),
+                   Tensor<3,spacedim>());
+
+         for (unsigned int point=0; point<n_q_points; ++point)
+           for (unsigned int k=0; k<data.n_shape_functions; ++k)
+             for (unsigned int i=0; i<dim; ++i)
+               for (unsigned int j=0; j<dim; ++j)
+                 for (unsigned int l=0; l<dim; ++l)
+                   jacobian_grads[point][i][j][l]
+                     += (data.second_derivative(point+DataSetDescriptor::cell (), k)[j][l]
+                         *
+                         data.mapping_support_points[k][i]);
+       }
     }
                                   // copy values from InternalData to vector
                                   // given by reference
@@ -847,9 +861,10 @@ MappingQ1<dim,spacedim>::fill_fe_values (
     {      
       Assert (inverse_jacobians.size() == n_q_points,
              ExcDimensionMismatch(inverse_jacobians.size(), n_q_points));
-      for (unsigned int point=0; point<n_q_points; ++point)
-       inverse_jacobians[point]
-         = transpose(data.covariant[point]);
+      if (cell_similarity != CellSimilarity::translation)
+       for (unsigned int point=0; point<n_q_points; ++point)
+         inverse_jacobians[point]
+           = transpose(data.covariant[point]);
     }
 }
 
@@ -947,7 +962,8 @@ MappingQ1<dim,spacedim>::compute_fill_face (
   std::vector<Tensor<1,dim> > &boundary_forms,
   std::vector<Point<spacedim> >    &normal_vectors) const
 {
-  compute_fill (cell, n_q_points, data_set, data, quadrature_points);
+  compute_fill (cell, n_q_points, data_set, CellSimilarity::no_similarity,
+               data, quadrature_points);
 
   const UpdateFlags update_flags(data.current_update_flags());
   
index 7f4f1d5ac139f7d96d32bd0980393a5d10e2a76e..287a809e10695a4acb5a5030e9ba6249b15d4864 100644 (file)
@@ -27,9 +27,9 @@ template <int dim, class EulerVectorType, int spacedim>
 MappingQ1Eulerian<dim, EulerVectorType, spacedim>::
 MappingQ1Eulerian (const EulerVectorType  &euler_transform_vectors,
                   const DoFHandler<dim> &shiftmap_dof_handler)
-                     :
-                    euler_transform_vectors(euler_transform_vectors),
-                    shiftmap_dof_handler(&shiftmap_dof_handler)
+                   :
+                  euler_transform_vectors(euler_transform_vectors),
+                  shiftmap_dof_handler(&shiftmap_dof_handler)
 {}
 
 
@@ -114,6 +114,33 @@ MappingQ1Eulerian<dim, EulerVectorType, spacedim>::clone () const
 }
 
 
+
+template<int dim, class EulerVectorType, int spacedim>
+void
+MappingQ1Eulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
+  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+  const Quadrature<dim>                                     &q,
+  enum CellSimilarity::Similarity                            cell_similarity,
+  typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
+  std::vector<Point<spacedim> >                             &quadrature_points,
+  std::vector<double>                                       &JxW_values,
+  std::vector<Tensor<2,spacedim> >                          &jacobians,
+  std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+  std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
+  std::vector<Point<spacedim> >                             &cell_normal_vectors) const
+{
+                                  // disable any previously detected
+                                  // similarity and hand on to the
+                                  // respective function of the base class.
+  cell_similarity = CellSimilarity::no_similarity;
+  MappingQ1<dim,spacedim>::fill_fe_values (cell, q, cell_similarity, mapping_data,
+                                          quadrature_points, JxW_values, jacobians,
+                                          jacobian_grads, inverse_jacobians,
+                                          cell_normal_vectors);
+}
+
+
+
 // explicit instantiation
 template class MappingQ1Eulerian<deal_II_dimension, Vector<double> >;
 #ifdef DEAL_II_USE_PETSC
index f12842a2cfdab404a4ab5484ebc1655f78bcf095..cfa9ce19301d7e8f64f12766b2c59b62e2cdb476 100644 (file)
@@ -182,6 +182,32 @@ compute_mapping_support_points
 
 
 
+template<int dim, class EulerVectorType, int spacedim>
+void
+MappingQEulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
+  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+  const Quadrature<dim>                                     &q,
+  enum CellSimilarity::Similarity                            cell_similarity,
+  typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
+  std::vector<Point<spacedim> >                             &quadrature_points,
+  std::vector<double>                                       &JxW_values,
+  std::vector<Tensor<2,spacedim> >                          &jacobians,
+  std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
+  std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
+  std::vector<Point<spacedim> >                             &cell_normal_vectors) const
+{
+                                  // disable any previously detected
+                                  // similarity and hand on to the respective
+                                  // function of the base class.
+  cell_similarity = CellSimilarity::no_similarity;
+  MappingQ<dim,spacedim>::fill_fe_values (cell, q, cell_similarity, mapping_data,
+                                         quadrature_points, JxW_values, jacobians,
+                                         jacobian_grads, inverse_jacobians,
+                                         cell_normal_vectors);
+}
+
+
+
 // explicit instantiation
 template class MappingQEulerian<deal_II_dimension, Vector<double> >;
 #ifdef DEAL_II_USE_PETSC

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.