]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Merge most of the changes in the fe subdirectory from the hp branch to mainline,...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 27 Jun 2004 21:23:35 +0000 (21:23 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 27 Jun 2004 21:23:35 +0000 (21:23 +0000)
git-svn-id: https://svn.dealii.org/trunk@9470 0785d39b-7218-0410-832d-ea1e28bc413d

29 files changed:
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_base.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_raviart_thomas.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/include/fe/hp_fe_values.h [new file with mode: 0644]
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/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_data.cc
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc
deal.II/deal.II/source/fe/fe_dgq.cc
deal.II/deal.II/source/fe/fe_nedelec.cc
deal.II/deal.II/source/fe/fe_poly.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_raviart_thomas.cc
deal.II/deal.II/source/fe/fe_system.cc
deal.II/deal.II/source/fe/fe_tools.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

index 19779ee73ab37c74774852d74a9835329d8b88c9..e3c9e4f3e95563ad43907163d15351501e8d9ff1 100644 (file)
@@ -23,6 +23,7 @@ template <int dim> class FEValues;
 template <int dim> class FEFaceValues;
 template <int dim> class FESubfaceValues;
 template <int dim> class FESystem;
+template <int dim> class FECollection;
 
 /*!@addtogroup febase */
 /*@{*/
@@ -242,7 +243,8 @@ class FiniteElement : public FiniteElementBase<dim>
                                      * a copy constructor.
                                      *
                                      * This function is needed by the
-                                     * constructors of @p FESystem.
+                                     * constructors of FESystem as well
+                                     * as by the FECollection class.
                                      */
     virtual FiniteElement<dim> *clone() const = 0;
     
@@ -309,7 +311,7 @@ class FiniteElement : public FiniteElementBase<dim>
                                      */                                      
     virtual void
     fill_fe_values (const Mapping<dim>                   &mapping,
-                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                &quadrature,
                    typename Mapping<dim>::InternalDataBase       &mapping_internal,
                    typename Mapping<dim>::InternalDataBase       &fe_internal,
@@ -328,7 +330,7 @@ class FiniteElement : public FiniteElementBase<dim>
                                      */                                      
     virtual void
     fill_fe_face_values (const Mapping<dim>                   &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const Quadrature<dim-1>              &quadrature,
                         typename Mapping<dim>::InternalDataBase       &mapping_internal,
@@ -348,7 +350,7 @@ class FiniteElement : public FiniteElementBase<dim>
                                      */                                      
     virtual void
     fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                    face_no,
                            const unsigned int                    sub_no,
                            const Quadrature<dim-1>              &quadrature,
@@ -365,6 +367,7 @@ class FiniteElement : public FiniteElementBase<dim>
     friend class FEFaceValues<dim>;
     friend class FESubfaceValues<dim>;
     friend class FESystem<dim>;
+    friend class FECollection<dim>;
 };
 
 /*@}*/
index 551f7058df80d380730df50e5c2ba7948931d5ea..df31b9fddf0913708385bc9e3cee53c4b07f1b3f 100644 (file)
@@ -31,7 +31,6 @@ template<int dim> class FESystem;
 /*!@addtogroup febase */
 /*@{*/
 
-
 /**
  * Dimension independent data for finite elements. See the derived
  * class FiniteElementBase class for information on its use. All
@@ -904,7 +903,7 @@ class FiniteElementBase : public Subscriptor,
                                      * @p system_to_base_index
                                      * function.
                                      */
-    std::pair<unsigned int,unsigned int>
+    std::pair<unsigned int, unsigned int>
     system_to_component_index (const unsigned int index) const;    
   
                                     /**
@@ -912,7 +911,7 @@ class FiniteElementBase : public Subscriptor,
                                      * shape functions and their
                                      * indices on a face.
                                      */
-    std::pair<unsigned int,unsigned int>
+    std::pair<unsigned int, unsigned int>
     face_system_to_component_index (const unsigned int index) const;
 
                                      /**
@@ -1661,7 +1660,7 @@ class FiniteElementBase : public Subscriptor,
                                      * gradients.
                                      */
     void compute_2nd (const Mapping<dim>                      &mapping,
-                     const typename DoFHandler<dim>::cell_iterator    &cell,
+                     const typename Triangulation<dim>::cell_iterator    &cell,
                      const unsigned int                       offset,
                      typename Mapping<dim>::InternalDataBase &mapping_internal,
                      InternalDataBase                        &fe_internal,
index 928f1c2426c9c0062dfb83d8a8d5eab0df4f331f..191c368da16a3ce332bbac7e8d57c578dd1229c1 100644 (file)
@@ -310,7 +310,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_values (const Mapping<dim> &mapping,
-                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                &quadrature,
                    typename Mapping<dim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim>::InternalDataBase      &fe_internal,
@@ -323,7 +323,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_face_values (const Mapping<dim> &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const Quadrature<dim-1>                &quadrature,
                         typename Mapping<dim>::InternalDataBase      &mapping_internal,
@@ -337,7 +337,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_subface_values (const Mapping<dim> &mapping,
-                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                    face_no,
                            const unsigned int                    sub_no,
                            const Quadrature<dim-1>                &quadrature,
index 4755163d5a3492da04bb260d19321d75d24d7af2..330da4aae88c11dbeb1f7f8d1a1dca4848f5ce97 100644 (file)
@@ -413,7 +413,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_values (const Mapping<dim> &mapping,
-                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                &quadrature,
                    typename Mapping<dim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim>::InternalDataBase      &fe_internal,
@@ -426,7 +426,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_face_values (const Mapping<dim> &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const Quadrature<dim-1>                &quadrature,
                         typename Mapping<dim>::InternalDataBase      &mapping_internal,
@@ -440,7 +440,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_subface_values (const Mapping<dim> &mapping,
-                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                    face_no,
                            const unsigned int                    sub_no,
                            const Quadrature<dim-1>                &quadrature,
index 39818dcba2603e072a14587c86464ac3f743e265..cc54a3843894185a8c69ae3a96f9e74a07a3477b 100644 (file)
@@ -225,7 +225,7 @@ class FE_Poly : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_values (const Mapping<dim> &mapping,
-                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                &quadrature,
                    typename Mapping<dim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim>::InternalDataBase      &fe_internal,
@@ -237,7 +237,7 @@ class FE_Poly : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_face_values (const Mapping<dim> &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const Quadrature<dim-1>                &quadrature,
                         typename Mapping<dim>::InternalDataBase      &mapping_internal,
@@ -250,7 +250,7 @@ class FE_Poly : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_subface_values (const Mapping<dim> &mapping,
-                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                    face_no,
                            const unsigned int                    sub_no,
                            const Quadrature<dim-1>                &quadrature,
index 9532dc1b4de8cec835ece2771ee9cada681dcb21..e4d3c4a2c2c3cccba25e21582bee5a622cb1cde3 100644 (file)
@@ -318,7 +318,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_values (const Mapping<dim> &mapping,
-                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                &quadrature,
                    typename Mapping<dim>::InternalDataBase      &mapping_internal,
                    typename Mapping<dim>::InternalDataBase      &fe_internal,
@@ -331,7 +331,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_face_values (const Mapping<dim> &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const Quadrature<dim-1>                &quadrature,
                         typename Mapping<dim>::InternalDataBase      &mapping_internal,
@@ -345,7 +345,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_subface_values (const Mapping<dim> &mapping,
-                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                    face_no,
                            const unsigned int                    sub_no,
                            const Quadrature<dim-1>                &quadrature,
index a4e4b870e97c33229505987e67ee3165c7eb5f9a..63d64dc91142c578d28574831fe49ee6151896e8 100644 (file)
@@ -448,7 +448,7 @@ class FESystem : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_values (const Mapping<dim>                      &mapping,
-                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                   &quadrature,
                    typename Mapping<dim>::InternalDataBase &mapping_data,
                    typename Mapping<dim>::InternalDataBase &fe_data,
@@ -466,7 +466,7 @@ class FESystem : public FiniteElement<dim>
                                      */    
     virtual void
     fill_fe_face_values (const Mapping<dim>                   &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const Quadrature<dim-1>              &quadrature,
                         typename Mapping<dim>::InternalDataBase      &mapping_data,
@@ -485,7 +485,7 @@ class FESystem : public FiniteElement<dim>
                                      */
     virtual void
     fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                    face_no,
                            const unsigned int                    sub_no,
                            const Quadrature<dim-1>              &quadrature,
@@ -518,7 +518,7 @@ class FESystem : public FiniteElement<dim>
                                      */
     template <int dim_1>
     void compute_fill (const Mapping<dim>                   &mapping,
-                      const typename DoFHandler<dim>::cell_iterator &cell,
+                      const typename Triangulation<dim>::cell_iterator &cell,
                       const unsigned int                    face_no,
                       const unsigned int                    sub_no,
                       const Quadrature<dim_1>              &quadrature,
index f08a7056840fffa6fd852d4ba91afc101becd7a7..8072d3dbad5063c547ca83c8e85ae8129e32fdbb 100644 (file)
 
 
 
+#include <base/config.h>
+#include <base/exceptions.h>
+
+#include <vector>
+
 template <typename number> class FullMatrix;
+template <typename number> class Vector;
 template <int dim> class FiniteElement;
 template <int dim> class DoFHandler;
-template <typename number> class Vector;
+template <int dim> class hpDoFHandler;
 template <int dim> class FE_Q;
 class ConstraintMatrix;
 
@@ -31,6 +37,7 @@ class ConstraintMatrix;
 /*!@addtogroup febase */
 /*@{*/
 
+
 /**
  * This class performs interpolations and extrapolations of discrete
  * functions of one @p FiniteElement @p fe1 to another @p FiniteElement
@@ -41,7 +48,7 @@ class ConstraintMatrix;
  * $id-I_h$ that is needed for evaluating $(id-I_h)z$ for e.g. the
  * dual solution $z$.
  *
- * @author Ralf Hartmann, Guido Kanschat, 2000, 2004
+ * @author Ralf Hartmann, 2000; Wolfgang Bangerth, 2003, Guido Kanschat, 2000, 2004
  */
 class FETools
 {
@@ -88,9 +95,11 @@ class FETools
                                      * is an embedding matrix.
                                      */
     template <int dim, typename number>
-    static void get_interpolation_matrix(const FiniteElement<dim> &fe1,
-                                        const FiniteElement<dim> &fe2,
-                                        FullMatrix<number> &interpolation_matrix);
+    static
+    void
+    get_interpolation_matrix(const FiniteElement<dim> &fe1,
+                             const FiniteElement<dim> &fe2,
+                             FullMatrix<number> &interpolation_matrix);
     
                                     /**
                                      * Gives the interpolation matrix
@@ -111,9 +120,11 @@ class FETools
                                      * only the unit matrix.
                                      */
     template <int dim, typename number>
-    static void get_back_interpolation_matrix(const FiniteElement<dim> &fe1,
-                                             const FiniteElement<dim> &fe2,
-                                             FullMatrix<number> &interpolation_matrix);
+    static
+    void
+    get_back_interpolation_matrix(const FiniteElement<dim> &fe1,
+                                  const FiniteElement<dim> &fe2,
+                                  FullMatrix<number> &interpolation_matrix);
 
                                     /**
                                      * Gives the unit matrix minus the
@@ -136,8 +147,8 @@ class FETools
     static
     void
     get_interpolation_difference_matrix(const FiniteElement<dim> &fe1,
-                                       const FiniteElement<dim> &fe2,
-                                       FullMatrix<number> &difference_matrix);
+                                        const FiniteElement<dim> &fe2,
+                                        FullMatrix<number> &difference_matrix);
 
                                     /**
                                      * Compute the local
@@ -207,11 +218,16 @@ class FETools
                                      * your hanging node constraints
                                      * object.
                                      */
-    template <int dim, class InVector, class OutVector>
-    static void interpolate (const DoFHandler<dim>& dof1,
-                            const InVector&        u1,
-                            const DoFHandler<dim>& dof2,
-                            OutVector&             u2);
+    template <int dim,
+              template <int> class DH1,
+              template <int> class DH2,
+              class InVector, class OutVector>
+    static
+    void
+    interpolate (const DH1<dim> &dof1,
+                 const InVector &u1,
+                 const DH2<dim> &dof2,
+                 OutVector      &u2);
     
                                     /**
                                      * Gives the interpolation of a
@@ -252,6 +268,23 @@ class FETools
                             const DoFHandler<dim>&  dof2,
                             const ConstraintMatrix& constraints,
                             OutVector&              u2);
+
+
+                                     /**
+                                      * Same as last function, except
+                                      * that one or both of the dof
+                                      * handler objects might be of
+                                      * type @p hpDoFHandler.
+                                      */
+    template <int dim,
+              template <int> class DH1,
+              template <int> class DH2,              
+              class InVector, class OutVector>
+    static void interpolate (const DH1<dim>         &dof1,
+                            const InVector         &u1,
+                            const DH2<dim>         &dof2,
+                            const ConstraintMatrix &constraints,
+                            OutVector&              u2);
     
 
                                     /**
@@ -281,12 +314,26 @@ class FETools
                                      * mapping.
                                      */
     template <int dim, class InVector, class OutVector>
-    static void back_interpolate (const DoFHandler<dim>&    dof1,
-                                 const InVector&           u1,
-                                 const FiniteElement<dim>fe2,
-                                 OutVector&                u1_interpolated);
+    static void back_interpolate (const DoFHandler<dim>    &dof1,
+                                 const InVector           &u1,
+                                 const FiniteElement<dim> &fe2,
+                                 OutVector                &u1_interpolated);
 
-                                    /**
+                                     /**
+                                      * Same as last function, except
+                                      * that the dof handler objects
+                                      * might be of type
+                                      * @p hpDoFHandler.
+                                      */
+    template <int dim,
+              template <int> class DH,
+              class InVector, class OutVector>
+    static void back_interpolate (const DH<dim>            &dof1,
+                                 const InVector           &u1,
+                                 const FiniteElement<dim> &fe2,
+                                 OutVector                &u1_interpolated);
+
+                                     /**
                                      * Gives the interpolation of the
                                      * @p dof1-function @p u1 to a
                                      * @p dof2-function, and
index 0b1968b28c13dc7cfd75bccaa997b7686602c67b..18a620e16142e9391ca079ab106e7bc5a2900219 100644 (file)
 #include <base/tensor.h>
 #include <base/quadrature.h>
 #include <base/table.h>
+#include <lac/vector.h>
+#include <lac/block_vector.h>
+#include <lac/petsc_vector.h>
+#include <lac/petsc_block_vector.h>
 #include <grid/tria.h>
 #include <grid/tria_iterator.h>
 #include <dofs/dof_handler.h>
 #include <fe/fe.h>
 #include <fe/fe_update_flags.h>
 #include <fe/mapping.h>
+#include <multigrid/mg_dof_handler.h>
+#include <multigrid/mg_dof_accessor.h>
 
 #include <algorithm>
+#include <memory>
 
 template <int dim> class Quadrature;
 
+
 /*!@addtogroup febase */
 /*@{*/
 
@@ -307,7 +315,7 @@ class FEValuesData
  * updating. Still, it is recommended to give ALL needed update flags
  * to @p FEValues.
  *
- * @author Wolfgang Bangerth, 1998, Guido Kanschat, 2001
+ * @author Wolfgang Bangerth, 1998, 2003, Guido Kanschat, 2001
  */
 template <int dim>
 class FEValuesBase : protected FEValuesData<dim>
@@ -824,10 +832,10 @@ class FEValuesBase : protected FEValuesData<dim>
     const FiniteElement<dim> & get_fe () const;
 
                                     /**
-                                     * Constant reference to the
-                                     * current cell
+                                     * Return a triangulation
+                                     * iterator to the current cell.
                                      */
-    typename Triangulation<dim>::cell_iterator get_cell () const;
+    const typename Triangulation<dim>::cell_iterator get_cell () const;
 
                                     /**
                                      * Determine an estimate for the
@@ -880,6 +888,469 @@ class FEValuesBase : protected FEValuesData<dim>
     DeclException0 (ExcFENotPrimitive);
     
   protected:
+                                     /**
+                                      * Objects of the @p FEValues
+                                      * class need to store a pointer
+                                      * (i.e. an iterator) to the
+                                      * present cell in order to be
+                                      * able to extract the values of
+                                      * the degrees of freedom on this
+                                      * cell in the
+                                      * @p get_function_values and
+                                      * assorted functions. On the
+                                      * other hand, this class should
+                                      * also work for different
+                                      * iterators, as long as they
+                                      * have the same interface to
+                                      * extract the DoF values (i.e.,
+                                      * for example, they need to have
+                                      * a @p get_interpolated_dof_values
+                                      * function).
+                                      *
+                                      * This calls for a common base
+                                      * class of iterator classes, and
+                                      * making the functions we need
+                                      * here @p virtual. On the other
+                                      * hand, this is the only place
+                                      * in the library where we need
+                                      * this, and introducing a base
+                                      * class of iterators and making
+                                      * a function virtual penalizes
+                                      * @em{all} users of the
+                                      * iterators, which are basically
+                                      * intended as very fast accessor
+                                      * functions. So we do not want
+                                      * to do this. Rather, what we do
+                                      * here is making the functions
+                                      * we need virtual only for use
+                                      * with @em{this class}. The idea
+                                      * is the following: have a
+                                      * common base class which
+                                      * declares some pure virtual
+                                      * functions, and for each
+                                      * possible iterator type, we
+                                      * have a derived class which
+                                      * stores the iterator to the
+                                      * cell and implements these
+                                      * functions. Since the iterator
+                                      * classes have the same
+                                      * interface, we can make the
+                                      * derived classes a template,
+                                      * templatized on the iterator
+                                      * type.
+                                      *
+                                      * This way, the use of virtual
+                                      * functions is restricted to
+                                      * only this class, and other
+                                      * users of iterators do not have
+                                      * to bear the negative effects.
+                                      *
+                                      * @author Wolfgang Bangerth, 2003
+                                      */
+    class CellIteratorBase 
+    {
+      public:
+                                         /**
+                                          * Destructor. Made virtual
+                                          * since we store only
+                                          * pointers to the base
+                                          * class.
+                                          */
+        virtual ~CellIteratorBase ();
+        
+                                         /**
+                                          * Conversion operator to an
+                                          * iterator for
+                                          * triangulations. This
+                                          * conversion is implicit for
+                                          * the original iterators,
+                                          * since they are derived
+                                          * classes. However, since
+                                          * here we have kind of a
+                                          * parallel class hierarchy,
+                                          * we have to have a
+                                          * conversion operator.
+                                          */
+        virtual
+        operator const typename Triangulation<dim>::cell_iterator () const = 0;
+
+                                         /**
+                                          * Return the number of
+                                          * degrees of freedom the DoF
+                                          * handler object has to
+                                          * which the iterator belongs
+                                          * to.
+                                          */
+        virtual
+        unsigned int
+        n_dofs_for_dof_handler () const = 0;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const Vector<double> &in,
+                                     Vector<double>       &out) const = 0;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const Vector<float> &in,
+                                     Vector<float>       &out) const = 0;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const BlockVector<double> &in,
+                                     Vector<double>       &out) const = 0;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const BlockVector<float> &in,
+                                     Vector<float>       &out) const = 0;
+
+#ifdef DEAL_II_USE_PETSC
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const PETScWrappers::Vector &in,
+                                     Vector<PetscScalar>         &out) const = 0;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const PETScWrappers::BlockVector &in,
+                                     Vector<PetscScalar>              &out) const = 0;
+#endif
+    };
+
+
+                                     /**
+                                      * Implementation of derived
+                                      * classes of the
+                                      * CellIteratorBase
+                                      * interface. See there for a
+                                      * description of the use of
+                                      * these classes.
+                                      * 
+                                      * @author Wolfgang Bangerth, 2003
+                                      */
+    template <typename CI>
+    class CellIterator : public CellIteratorBase
+    {
+      public:
+                                         /**
+                                          * Constructor. Take an
+                                          * iterator and store it in
+                                          * this class.
+                                          */
+        CellIterator (const CI &cell);
+        
+                                         /**
+                                          * Conversion operator to an
+                                          * iterator for
+                                          * triangulations. This
+                                          * conversion is implicit for
+                                          * the original iterators,
+                                          * since they are derived
+                                          * classes. However, since
+                                          * here we have kind of a
+                                          * parallel class hierarchy,
+                                          * we have to have a
+                                          * conversion operator.
+                                          */
+        virtual
+        operator const typename Triangulation<dim>::cell_iterator () const;
+        
+                                         /**
+                                          * Return the number of
+                                          * degrees of freedom the DoF
+                                          * handler object has to
+                                          * which the iterator belongs
+                                          * to.
+                                          */
+        virtual
+        unsigned int
+        n_dofs_for_dof_handler () const;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const Vector<double> &in,
+                                     Vector<double>       &out) const;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const Vector<float> &in,
+                                     Vector<float>       &out) const;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const BlockVector<double> &in,
+                                     Vector<double>       &out) const;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const BlockVector<float> &in,
+                                     Vector<float>       &out) const;
+
+#ifdef DEAL_II_USE_PETSC
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const PETScWrappers::Vector &in,
+                                     Vector<PetscScalar>         &out) const;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const PETScWrappers::BlockVector &in,
+                                     Vector<PetscScalar>              &out) const;
+#endif
+      private:
+                                         /**
+                                          * Copy of the iterator which
+                                          * we use in this object.
+                                          */
+        const CI cell;
+    };
+
+
+                                     /**
+                                      * Implementation of a derived
+                                      * class of the
+                                      * CellIteratorBase
+                                      * interface. See there for a
+                                      * description of the use of
+                                      * these classes.
+                                      *
+                                      * This class is basically a
+                                      * specialization of the general
+                                      * template for iterators into
+                                      * @p Triangulation objects (but
+                                      * since C++ does not allow
+                                      * something like this for nested
+                                      * classes, it runs under a
+                                      * separate name). Since these do
+                                      * not implement the interface
+                                      * that we would like to call,
+                                      * the functions of this class
+                                      * cannot be implemented
+                                      * meaningfully. However, most
+                                      * functions of the @p FEValues
+                                      * class do not make any use of
+                                      * degrees of freedom at all, so
+                                      * it should be possible to call
+                                      * @p FEValues::reinit with a tria
+                                      * iterator only; this class
+                                      * makes this possible, but
+                                      * whenever one of the functions
+                                      * of @p FEValues tries to call
+                                      * any of the functions of this
+                                      * class, an exception will be
+                                      * raised reminding the user that
+                                      * if she wants to use these
+                                      * features, then the
+                                      * @p FEValues object has to be
+                                      * reinitialized with a cell
+                                      * iterator that allows to
+                                      * extract degree of freedom
+                                      * information.
+                                      * 
+                                      * @author Wolfgang Bangerth, 2003
+                                      */
+    class TriaCellIterator : public CellIteratorBase
+    {
+      public:
+                                         /**
+                                          * Constructor. Take an
+                                          * iterator and store it in
+                                          * this class.
+                                          */
+        TriaCellIterator (const typename Triangulation<dim>::cell_iterator &cell);
+        
+                                         /**
+                                          * Conversion operator to an
+                                          * iterator for
+                                          * triangulations. This
+                                          * conversion is implicit for
+                                          * the original iterators,
+                                          * since they are derived
+                                          * classes. However, since
+                                          * here we have kind of a
+                                          * parallel class hierarchy,
+                                          * we have to have a
+                                          * conversion operator. Here,
+                                          * the conversion is trivial,
+                                          * from and to the same time.
+                                          */
+        virtual
+        operator const typename Triangulation<dim>::cell_iterator () const;
+
+                                         /**
+                                          * Implement the respective
+                                          * function of the base
+                                          * class. Since this is not
+                                          * possible, we just raise an
+                                          * error.
+                                          */
+        virtual
+        unsigned int
+        n_dofs_for_dof_handler () const;
+
+                                         /**
+                                          * Implement the respective
+                                          * function of the base
+                                          * class. Since this is not
+                                          * possible, we just raise an
+                                          * error.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const Vector<double> &in,
+                                     Vector<double>       &out) const;
+
+                                         /**
+                                          * Implement the respective
+                                          * function of the base
+                                          * class. Since this is not
+                                          * possible, we just raise an
+                                          * error.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const Vector<float> &in,
+                                     Vector<float>       &out) const;
+
+                                         /**
+                                          * Implement the respective
+                                          * function of the base
+                                          * class. Since this is not
+                                          * possible, we just raise an
+                                          * error.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const BlockVector<double> &in,
+                                     Vector<double>       &out) const;
+
+                                         /**
+                                          * Implement the respective
+                                          * function of the base
+                                          * class. Since this is not
+                                          * possible, we just raise an
+                                          * error.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const BlockVector<float> &in,
+                                     Vector<float>       &out) const;
+
+#ifdef DEAL_II_USE_PETSC
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const PETScWrappers::Vector &in,
+                                     Vector<PetscScalar>         &out) const;
+
+                                         /**
+                                          * Call
+                                          * @p get_interpolated_dof_values
+                                          * of the iterator with the
+                                          * given arguments.
+                                          */
+        virtual
+        void
+        get_interpolated_dof_values (const PETScWrappers::BlockVector &in,
+                                     Vector<PetscScalar>              &out) const;
+#endif
+      private:
+                                         /**
+                                          * Copy of the iterator which
+                                          * we use in this object.
+                                          */
+        const typename Triangulation<dim>::cell_iterator cell;
+
+                                         /**
+                                          * String to be displayed
+                                          * whenever one of the
+                                          * functions of this class is
+                                          * called. Make it a static
+                                          * member variable, since we
+                                          * show the same message for
+                                          * all member functions.
+                                          */
+        static const char * const message_string;
+    };
+    
                                     /**
                                      * Store the cell selected last time
                                      * the @p reinit function was called
@@ -887,8 +1358,8 @@ class FEValuesBase : protected FEValuesData<dim>
                                      * to the <tt>get_function_*</tt> functions
                                      * safer.
                                      */
-    typename DoFHandler<dim>::cell_iterator present_cell;
-
+    std::auto_ptr<const CellIteratorBase> present_cell;
+    
                                     /**
                                      * Storage for the mapping object.
                                      */
@@ -937,6 +1408,23 @@ class FEValuesBase : protected FEValuesData<dim>
                                      * @p MappingQ1 implicitly.
                                      */
     static const Mapping<dim> &get_default_mapping();
+
+  private:
+                                     /**
+                                      * Copy constructor. Since
+                                      * objects of this class are not
+                                      * copyable, we make it private,
+                                      * and also do not implement it.
+                                      */
+    FEValuesBase (const FEValuesBase &);
+    
+                                     /**
+                                      * Copy operator. Since
+                                      * objects of this class are not
+                                      * copyable, we make it private,
+                                      * and also do not implement it.
+                                      */
+    FEValuesBase & operator= (const FEValuesBase &);    
 };
 
 
@@ -975,11 +1463,62 @@ class FEValues : public FEValuesBase<dim>
              const UpdateFlags         update_flags);
     
                                     /**
-                                     * Reinitialize the gradients, Jacobi
-                                     * determinants, etc for the given cell
-                                     * and the given finite element.
+                                     * Reinitialize the gradients,
+                                     * Jacobi determinants, etc for
+                                     * the given cell of type
+                                     * "iterator into a DoFHandler
+                                     * object", and the finite
+                                     * element associated with this
+                                     * object. It is assumed that the
+                                     * finite element used by the
+                                     * given cell is also the one
+                                     * used by this @p FEValues
+                                     * object.
                                      */
-    void reinit (const typename DoFHandler<dim>::cell_iterator &);
+    void reinit (const typename DoFHandler<dim>::cell_iterator &cell);
+
+                                    /**
+                                     * Reinitialize the gradients,
+                                     * Jacobi determinants, etc for
+                                     * the given cell of type
+                                     * "iterator into a MGDoFHandler
+                                     * object", and the finite
+                                     * element associated with this
+                                     * object. It is assumed that the
+                                     * finite element used by the
+                                     * given cell is also the one
+                                     * used by this @p FEValues
+                                     * object.
+                                     */
+    void reinit (const typename MGDoFHandler<dim>::cell_iterator &cell);
+
+                                    /**
+                                     * Reinitialize the gradients,
+                                     * Jacobi determinants, etc for
+                                     * the given cell of type
+                                     * "iterator into a Triangulation
+                                     * object", and the given finite
+                                     * element. Since iterators into
+                                     * triangulation alone only
+                                     * convey information about the
+                                     * geometry of a cell, but not
+                                     * about degrees of freedom
+                                     * possibly associated with this
+                                     * cell, you will not be able to
+                                     * call some functions of this
+                                     * class if they need information
+                                     * about degrees of
+                                     * freedom. These functions are,
+                                     * above all, the
+                                     * <tt>get_function_value/grads/2nd_derivatives</tt>
+                                     * functions. If you want to call
+                                     * these functions, you have to
+                                     * call the @p reinit variants
+                                     * that take iterators into
+                                     * @p DoFHandler or other DoF
+                                     * handler type objects.
+                                     */
+    void reinit (const typename Triangulation<dim>::cell_iterator &cell);
 
                                     /**
                                      * Return a reference to the copy
@@ -995,6 +1534,13 @@ class FEValues : public FEValuesBase<dim>
                                      */
     unsigned int memory_consumption () const;
 
+                                     /**
+                                      * Return a reference to this
+                                      * very object.
+                                      */
+//TODO: explain reason    
+    const FEValues<dim> & get_present_fe_values () const;
+    
   private:
                                     /**
                                      * Store a copy of the quadrature
@@ -1007,6 +1553,20 @@ class FEValues : public FEValuesBase<dim>
                                      * constructors.
                                      */
     void initialize (const UpdateFlags update_flags);
+
+                                     /**
+                                      * The reinit() functions do
+                                      * only that part of the work
+                                      * that requires knowledge of the
+                                      * type of iterator. After
+                                      * setting present_cell(),
+                                      * they pass on to this function,
+                                      * which does the real work, and
+                                      * which is independent of the
+                                      * actual type of the cell
+                                      * iterator.
+                                      */
+    void do_reinit ();
 };
 
 
@@ -1139,18 +1699,18 @@ class FEFaceValues : public FEFaceValuesBase<dim>
                                      * matching the quadrature rule
                                      * and update flags.
                                      */
-    FEFaceValues (const Mapping<dim>mapping,
-                 const FiniteElement<dim>fe,
-                 const Quadrature<dim-1>quadrature,
-                 const UpdateFlags update_flags);
+    FEFaceValues (const Mapping<dim>       &mapping,
+                 const FiniteElement<dim> &fe,
+                 const Quadrature<dim-1>  &quadrature,
+                 const UpdateFlags         update_flags);
 
                                      /**
                                      * Constructor. Uses MappingQ1
                                      * implicitly.
                                      */
-    FEFaceValues (const FiniteElement<dim>fe,
-                 const Quadrature<dim-1>quadrature,
-                 const UpdateFlags update_flags);
+    FEFaceValues (const FiniteElement<dim> &fe,
+                 const Quadrature<dim-1>  &quadrature,
+                 const UpdateFlags         update_flags);
 
                                     /**
                                      * Reinitialize the gradients, Jacobi
@@ -1158,9 +1718,44 @@ class FEFaceValues : public FEFaceValuesBase<dim>
                                      * number @p face_no of @p cell
                                      * and the given finite element.
                                      */
-    void reinit (const typename DoFHandler<dim>::cell_iterator& cell,
-                const unsigned int face_no);
-
+    void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
+                const unsigned int                            face_no);
+
+                                    /**
+                                     * Reinitialize the gradients,
+                                     * Jacobi determinants, etc for
+                                     * the given face on given cell
+                                     * of type "iterator into a
+                                     * Triangulation object", and the
+                                     * given finite element. Since
+                                     * iterators into triangulation
+                                     * alone only convey information
+                                     * about the geometry of a cell,
+                                     * but not about degrees of
+                                     * freedom possibly associated
+                                     * with this cell, you will not
+                                     * be able to call some functions
+                                     * of this class if they need
+                                     * information about degrees of
+                                     * freedom. These functions are,
+                                     * above all, the
+                                     * <tt>get_function_value/grads/2nd_derivatives</tt>
+                                     * functions. If you want to call
+                                     * these functions, you have to
+                                     * call the @p reinit variants
+                                     * that take iterators into
+                                     * @p DoFHandler or other DoF
+                                     * handler type objects.
+                                     */
+    void reinit (const typename Triangulation<dim>::cell_iterator &cell,
+                const unsigned int                    face_no);
+    
+                                     /**
+                                      * Return a reference to this
+                                      * very object.
+                                      */
+//TODO: explain reason    
+    const FEFaceValues<dim> & get_present_fe_values () const;
   private:
 
                                     /**
@@ -1168,6 +1763,20 @@ class FEFaceValues : public FEFaceValuesBase<dim>
                                      * constructors.
                                      */
     void initialize (const UpdateFlags update_flags);    
+
+                                     /**
+                                      * The reinit() functions do
+                                      * only that part of the work
+                                      * that requires knowledge of the
+                                      * type of iterator. After
+                                      * setting present_cell(),
+                                      * they pass on to this function,
+                                      * which does the real work, and
+                                      * which is independent of the
+                                      * actual type of the cell
+                                      * iterator.
+                                      */
+    void do_reinit (const unsigned int face_no);
 };
 
 
@@ -1199,30 +1808,74 @@ class FESubfaceValues : public FEFaceValuesBase<dim>
                                      * matching the quadrature rule
                                      * and update flags.
                                      */
-    FESubfaceValues (const Mapping<dim>mapping,
-                    const FiniteElement<dim>fe,
-                    const Quadrature<dim-1>face_quadrature,
-                    const UpdateFlags update_flags);
+    FESubfaceValues (const Mapping<dim>       &mapping,
+                    const FiniteElement<dim> &fe,
+                    const Quadrature<dim-1>  &face_quadrature,
+                    const UpdateFlags         update_flags);
 
                                      /**
                                      * Constructor. Uses MappingQ1
                                      * implicitly.
                                      */
-    FESubfaceValues (const FiniteElement<dim>& fe,
-                    const Quadrature<dim-1>& face_quadrature,
-                    const UpdateFlags update_flags);
-
-                                    /**
-                                     * Reinitialize the gradients, Jacobi
-                                     * determinants, etc for the face with
-                                     * number @p face_no of @p cell
-                                     * and the given finite element.
-                                     */
-    void reinit (const typename DoFHandler<dim>::cell_iterator& cell,
-                const unsigned int face_no,
-                const unsigned int subface_no);
+    FESubfaceValues (const FiniteElement<dim> &fe,
+                    const Quadrature<dim-1>  &face_quadrature,
+                    const UpdateFlags         update_flags);
+
+                                    /**
+                                     * Reinitialize the gradients,
+                                     * Jacobi determinants, etc for
+                                     * the given cell of type
+                                     * "iterator into a DoFHandler
+                                     * object", and the finite
+                                     * element associated with this
+                                     * object. It is assumed that the
+                                     * finite element used by the
+                                     * given cell is also the one
+                                     * used by this
+                                     * @p FESubfaceValues object.
+                                     */
+    void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
+                const unsigned int                    face_no,
+                const unsigned int                    subface_no);
+
+                                    /**
+                                     * Reinitialize the gradients,
+                                     * Jacobi determinants, etc for
+                                     * the given subface on given
+                                     * cell of type "iterator into a
+                                     * Triangulation object", and the
+                                     * given finite element. Since
+                                     * iterators into triangulation
+                                     * alone only convey information
+                                     * about the geometry of a cell,
+                                     * but not about degrees of
+                                     * freedom possibly associated
+                                     * with this cell, you will not
+                                     * be able to call some functions
+                                     * of this class if they need
+                                     * information about degrees of
+                                     * freedom. These functions are,
+                                     * above all, the
+                                     * <tt>get_function_value/grads/2nd_derivatives</tt>
+                                     * functions. If you want to call
+                                     * these functions, you have to
+                                     * call the @p reinit variants
+                                     * that take iterators into
+                                     * @p DoFHandler or other DoF
+                                     * handler type objects.
+                                     */
+    void reinit (const typename Triangulation<dim>::cell_iterator &cell,
+                const unsigned int                    face_no,
+                const unsigned int                    subface_no);
+    
+                                     /**
+                                      * Return a reference to this
+                                      * very object.
+                                      */
+//TODO: explain reason    
+    const FESubfaceValues<dim> & get_present_fe_values () const;
 
-                                    /**
+                                     /**
                                      * Exception
                                      */
     DeclException0 (ExcReinitCalledWithBoundaryFace);
@@ -1239,6 +1892,21 @@ class FESubfaceValues : public FEFaceValuesBase<dim>
                                      * constructors.
                                      */
     void initialize (const UpdateFlags update_flags);
+
+                                     /**
+                                      * The reinit() functions do
+                                      * only that part of the work
+                                      * that requires knowledge of the
+                                      * type of iterator. After
+                                      * setting present_cell(),
+                                      * they pass on to this function,
+                                      * which does the real work, and
+                                      * which is independent of the
+                                      * actual type of the cell
+                                      * iterator.
+                                      */
+    void do_reinit (const unsigned int face_no,
+                    const unsigned int subface_no);
 };
 
 /*@}*/
@@ -1551,10 +2219,10 @@ FEValuesBase<dim>::get_mapping () const
 
 template <int dim>
 inline
-typename Triangulation<dim>::cell_iterator
+const typename Triangulation<dim>::cell_iterator
 FEValuesBase<dim>::get_cell () const
 {
-  return present_cell;
+  return *present_cell;
 }
 
 /*------------------------ Inline functions: FEValues ----------------------------*/
@@ -1569,6 +2237,16 @@ FEValues<dim>::get_quadrature () const
 }
 
 
+
+template <int dim>
+inline
+const FEValues<dim> &
+FEValues<dim>::get_present_fe_values () const 
+{
+  return *this;
+}
+
+
 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
 
 
@@ -1586,6 +2264,39 @@ FEFaceValuesBase<dim>::normal_vector (const unsigned int i) const
 }
 
 
+
+/*------------------------ Inline functions: FE*FaceValues --------------------*/
+
+template <int dim>
+inline
+const Quadrature<dim-1> &
+FEFaceValuesBase<dim>::get_quadrature () const 
+{
+  return quadrature;
+}
+
+
+
+template <int dim>
+inline
+const FEFaceValues<dim> &
+FEFaceValues<dim>::get_present_fe_values () const 
+{
+  return *this;
+}
+
+
+
+template <int dim>
+inline
+const FESubfaceValues<dim> &
+FESubfaceValues<dim>::get_present_fe_values () const 
+{
+  return *this;
+}
+
+
+
 template <int dim>
 inline
 const Tensor<1,dim> &
@@ -1601,14 +2312,5 @@ FEFaceValuesBase<dim>::boundary_form (const unsigned int i) const
 
 
 
-template <int dim>
-inline
-const Quadrature<dim-1> &
-FEFaceValuesBase<dim>::get_quadrature () const 
-{
-  return quadrature;
-}
-
-
 
 #endif
diff --git a/deal.II/deal.II/include/fe/hp_fe_values.h b/deal.II/deal.II/include/fe/hp_fe_values.h
new file mode 100644 (file)
index 0000000..09c7cb1
--- /dev/null
@@ -0,0 +1,455 @@
+//----------------------------  hp_fe_values.h  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003, 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  hp_fe_values.h  ---------------------------
+#ifndef __deal2__hp_fe_values_h
+#define __deal2__hp_fe_values_h
+
+#include <base/config.h>
+#include <fe/fe.h>
+#include <fe/fe_collection.h>
+#include <fe/fe_values.h>
+
+#include <map>
+#include <boost/shared_ptr.hpp>
+
+template <int dim> class MappingQ1;
+template <int dim> class FiniteElement;
+
+
+namespace internal 
+{
+/**
+ * Map between finite element objects and @p FEValues (the second
+ * template parameter, which might be <tt>FEValues<dim></tt>,
+ * <tt>FEFaceValues<dim></tt>, ...) objects. The <tt>hpFE*Values</tt> classes
+ * use this to hold an <tt>FE*Values</tt> object for each finite element
+ * that is used in the triangulation that it integrates on.
+ */
+  template <int dim, class FEValues>
+  class FEValuesMap
+  {
+    public:
+                                       /**
+                                        * Make destructor virtual,
+                                        * since we have virtual
+                                        * functions in this class.
+                                        */
+      virtual ~FEValuesMap ();
+      
+                                       /**
+                                        * Return a reference to the
+                                        * @p FEValues object selected
+                                        * by the last call to
+                                        * @p select_fe_values. The
+                                        * returned value is a constant
+                                        * reference since the only
+                                        * non-const function in the
+                                        * underlying class would be
+                                        * @p reinit, which you must
+                                        * not call directly any way;
+                                        * rather, use the @p reinit
+                                        * function of the
+                                        * <tt>hpFE*Values</tt> class.
+                                        */
+      const FEValues & get_present_fe_values () const;
+
+    protected:
+                                       /**
+                                        * Create an object of type
+                                        * @p FEValues for this
+                                        * particular finite
+                                        * element. Since the type of
+                                        * @p FEValues is not known
+                                        * here, we have to leave this
+                                        * to derived classes.
+                                        */
+      virtual
+      FEValues * create_fe_values (const FiniteElement<dim> &fe) const = 0;
+      
+                                       /**
+                                        * Select the @p FEValues
+                                        * object corresponding to the
+                                        * finite element object given
+                                        * as argument. If there is no
+                                        * such @p FEValues object
+                                        * yet, then one is created by
+                                        * calling @p create_fe_values
+                                        * with this finite element
+                                        * object.
+                                        *
+                                        * A non-const reference to
+                                        * this object is returned.
+                                        */
+      FEValues &
+      select_fe_values (const FiniteElement<dim> &fe);
+      
+    private:
+                                       /**
+                                        * Have a map from pointers to
+                                        * finite elements to objects
+                                        * of type @p FEValues.
+                                        *
+                                        * Note that the shared pointer
+                                        * will make sure that the
+                                        * memory previously allocated
+                                        * is released when the
+                                        * destructor of this class is
+                                        * run.
+                                        */
+      std::map<SmartPointer<const FiniteElement<dim> >,
+               boost::shared_ptr<FEValues> > fe_to_fe_values_map;
+
+                                       /**
+                                        * Pointer to the @p FEValues
+                                        * object that is used for the
+                                        * present cell. This always
+                                        * points to one of the objects
+                                        * in the map above, and to
+                                        * which one is determined by
+                                        * the last call to
+                                        * @p select_fe_values.
+                                        */
+      boost::shared_ptr<FEValues> present_fe_values;
+  };
+
+
+/**
+ * Base class for the <tt>hpFE*Values</tt> classes, storing the data that
+ * is common to them. The first template parameter denotes the space
+ * dimension we are in, the second the dimensionality of the object
+ * that we integrate on, i.e. for usual @p hpFEValues it is equal to
+ * the first one, while for face integration it is one less.
+ *
+ * @author Wolfgang Bangerth, 2003
+ */
+  template <int dim, int q_dim>
+  class hpFEValuesBase 
+  {
+    public:
+                                       /**
+                                        * Constructor. Set the fields
+                                        * of this class to the values
+                                        * indicated by the parameters
+                                        * to the constructor.
+                                        */
+      hpFEValuesBase (const Mapping<dim>      &mapping,
+                      const Quadrature<q_dim> &quadrature,
+                      const UpdateFlags        update_flags);
+
+                                       /**
+                                        * Constructor. Set the fields
+                                        * of this class to the values
+                                        * indicated by the parameters
+                                        * to the constructor, and
+                                        * choose a @p MappingQ1
+                                        * object for the mapping
+                                        * object.
+                                        */
+      hpFEValuesBase (const Quadrature<q_dim> &quadrature,
+                      const UpdateFlags        update_flags);
+      
+    protected:
+                                       /**
+                                        * Pointer to a mapping object as
+                                        * provided to the constructor,
+                                        * or a default @p MappingQ1
+                                        * object if none was given.
+                                        */
+      const SmartPointer<const Mapping<dim> > mapping;
+    
+                                       /**
+                                        * Pointer to a quadrature object
+                                        * as provided to the
+                                        * constructor.
+                                        */
+      const SmartPointer<const Quadrature<q_dim> > quadrature;
+
+                                       /**
+                                        * Values of the update flags as
+                                        * given to the constructor.
+                                        */
+      const UpdateFlags update_flags;
+
+    private:
+                                       /**
+                                        * Default mapping, in case
+                                        * none was provided through
+                                        * the call to the constructor.
+                                        */
+      static const MappingQ1<dim> default_mapping;
+  };
+}
+
+
+
+template <int dim>
+class hpFEValues : public internal::FEValuesMap<dim,FEValues<dim> >,
+                   protected internal::hpFEValuesBase<dim,dim>
+{
+  public:
+                                     /**
+                                      * Constructor. Initialize this
+                                      * object with the given
+                                      * parameters.
+                                      *
+                                      * The finite element
+                                      * collection parameter is
+                                      * actually ignored, but is in
+                                      * the signature of this
+                                      * function to make it
+                                      * compatible with the
+                                      * signature of the respective
+                                      * constructor of the usual
+                                      * FEValues object, with
+                                      * the respective parameter in
+                                      * that function also being the
+                                      * return value of the
+                                      * DoFHandler<tt>::get_fe()</tt>
+                                      * function.
+                                      */
+    hpFEValues (const Mapping<dim>      &mapping,
+                const FECollection<dim> &fe_collection,
+                const Quadrature<dim>   &quadrature,
+                const UpdateFlags        update_flags);
+
+                                     /**
+                                      * Constructor. Initialize this
+                                      * object with the given
+                                      * parameters, and choose a
+                                      * @p MappingQ1 object for the
+                                      * mapping object.
+                                      *
+                                      * The finite element
+                                      * collection parameter is
+                                      * actually ignored, but is in
+                                      * the signature of this
+                                      * function to make it
+                                      * compatible with the
+                                      * signature of the respective
+                                      * constructor of the usual
+                                      * FEValues object, with
+                                      * the respective parameter in
+                                      * that function also being the
+                                      * return value of the
+                                      * DoFHandler<tt>::get_fe()</tt>
+                                      * function.
+                                      */
+    hpFEValues (const FECollection<dim> &fe_collection,
+                const Quadrature<dim>   &quadrature,
+                const UpdateFlags        update_flags);
+
+                                     /**
+                                      * Reinitialize the object for
+                                      * the given cell. This selects
+                                      * the right @p FEValues object
+                                      * for the finite element in use
+                                      * by the cell given, and calling
+                                      * the @p reinit function on
+                                      * this object.
+                                      */
+    void
+    reinit (const typename hpDoFHandler<dim>::cell_iterator &cell);
+
+  protected:
+                                     /**
+                                      * Create an object of type
+                                      * @p FEValues for this
+                                      * particular finite element.
+                                      */
+    virtual
+    FEValues<dim> *
+    create_fe_values (const FiniteElement<dim> &fe) const;
+};
+
+
+
+template <int dim>
+class hpFEFaceValues : public internal::FEValuesMap<dim,FEFaceValues<dim> >,
+                       protected internal::hpFEValuesBase<dim,dim-1>
+{
+  public:
+                                     /**
+                                      * Constructor. Initialize this
+                                      * object with the given
+                                      * parameters.
+                                      *
+                                      * The finite element
+                                      * collection parameter is
+                                      * actually ignored, but is in
+                                      * the signature of this
+                                      * function to make it
+                                      * compatible with the
+                                      * signature of the respective
+                                      * constructor of the usual
+                                      * FEValues object, with
+                                      * the respective parameter in
+                                      * that function also being the
+                                      * return value of the
+                                      * DoFHandler<tt>::get_fe()</tt>
+                                      * function.
+                                      */
+    hpFEFaceValues (const Mapping<dim>      &mapping,
+                    const FECollection<dim> &fe_collection,
+                    const Quadrature<dim-1> &quadrature,
+                    const UpdateFlags        update_flags);
+
+                                     /**
+                                      * Constructor. Initialize this
+                                      * object with the given
+                                      * parameters, and choose a
+                                      * @p MappingQ1 object for the
+                                      * mapping object.
+                                      *
+                                      * The finite element
+                                      * collection parameter is
+                                      * actually ignored, but is in
+                                      * the signature of this
+                                      * function to make it
+                                      * compatible with the
+                                      * signature of the respective
+                                      * constructor of the usual
+                                      * FEValues object, with
+                                      * the respective parameter in
+                                      * that function also being the
+                                      * return value of the
+                                      * DoFHandler<tt>::get_fe()</tt>
+                                      * function.
+                                      */
+    hpFEFaceValues (const FECollection<dim> &fe_collection,
+                    const Quadrature<dim-1> &quadrature,
+                    const UpdateFlags        update_flags);
+
+                                     /**
+                                      * Reinitialize the object for
+                                      * the given cell. This selects
+                                      * the right @p FEValues object
+                                      * for the finite element in use
+                                      * by the cell given, and calling
+                                      * the @p reinit function on
+                                      * this object.
+                                      */
+    void
+    reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
+            const unsigned int face_no);
+
+  protected:
+                                     /**
+                                      * Create an object of type
+                                      * @p FEValues for this
+                                      * particular finite element.
+                                      */
+    virtual
+    FEFaceValues<dim> *
+    create_fe_values (const FiniteElement<dim> &fe) const;
+};
+
+
+
+template <int dim>
+class hpFESubfaceValues : public internal::FEValuesMap<dim,FESubfaceValues<dim> >,
+                          protected internal::hpFEValuesBase<dim,dim-1>
+{
+  public:
+                                     /**
+                                      * Constructor. Initialize this
+                                      * object with the given
+                                      * parameters.
+                                      *
+                                      * The finite element
+                                      * collection parameter is
+                                      * actually ignored, but is in
+                                      * the signature of this
+                                      * function to make it
+                                      * compatible with the
+                                      * signature of the respective
+                                      * constructor of the usual
+                                      * FEValues object, with
+                                      * the respective parameter in
+                                      * that function also being the
+                                      * return value of the
+                                      * DoFHandler<tt>::get_fe()</tt>
+                                      * function.
+                                      */
+    hpFESubfaceValues (const Mapping<dim>      &mapping,
+                       const FECollection<dim> &fe_collection,
+                       const Quadrature<dim-1> &quadrature,
+                       const UpdateFlags        update_flags);
+
+                                     /**
+                                      * Constructor. Initialize this
+                                      * object with the given
+                                      * parameters, and choose a
+                                      * @p MappingQ1 object for the
+                                      * mapping object.
+                                      *
+                                      * The finite element
+                                      * collection parameter is
+                                      * actually ignored, but is in
+                                      * the signature of this
+                                      * function to make it
+                                      * compatible with the
+                                      * signature of the respective
+                                      * constructor of the usual
+                                      * FEValues object, with
+                                      * the respective parameter in
+                                      * that function also being the
+                                      * return value of the
+                                      * DoFHandler<tt>::get_fe()</tt>
+                                      * function.
+                                      */
+    hpFESubfaceValues (const FECollection<dim> &fe_collection,
+                       const Quadrature<dim-1> &quadrature,
+                       const UpdateFlags        update_flags);
+
+                                     /**
+                                      * Reinitialize the object for
+                                      * the given cell. This selects
+                                      * the right @p FEValues object
+                                      * for the finite element in use
+                                      * by the cell given, and calling
+                                      * the @p reinit function on
+                                      * this object.
+                                      */
+    void
+    reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
+            const unsigned int face_no,
+            const unsigned int subface_no);
+
+  protected:
+                                     /**
+                                      * Create an object of type
+                                      * @p FEValues for this
+                                      * particular finite element.
+                                      */
+    virtual
+    FESubfaceValues<dim> *
+    create_fe_values (const FiniteElement<dim> &fe) const;
+};
+
+
+// -------------- inline and template functions --------------
+
+namespace internal 
+{
+  template <int dim, class FEValues>
+  inline
+  const FEValues &
+  FEValuesMap<dim,FEValues>::get_present_fe_values () const
+  {
+    return *present_fe_values;
+  }
+}
+
+
+
+
+#endif
index 5613229587e5615c9dde658300592774b7df5bde..dcfbb476f1468adca2e90a19b066ed0b6b944287 100644 (file)
 
 
 #include <base/config.h>
-#include <cmath>
-#include <base/point.h>
-#include <base/subscriptor.h>
 #include <grid/tria.h>
 #include <dofs/dof_handler.h>
 #include <fe/fe_update_flags.h>
 
+#include <cmath>
+
 template <int dim> class Quadrature;
 template <int dim> class FEValuesData;
 template <int dim> class FEValues;
@@ -445,7 +444,7 @@ class Mapping : public Subscriptor
                                      *
                                      */
     virtual void
-    fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_values (const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                         &quadrature,
                    InternalDataBase                              &internal,
                    std::vector<Point<dim> >                      &quadrature_points,
@@ -473,7 +472,7 @@ class Mapping : public Subscriptor
                                      * transformed quadrature weight.
                                      */
     virtual void
-    fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                        face_no,
                         const Quadrature<dim-1>                  &quadrature,
                         InternalDataBase                         &internal,
@@ -486,7 +485,7 @@ class Mapping : public Subscriptor
                                      * See above.
                                      */
     virtual void
-    fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int                        face_no,
                            const unsigned int                        sub_no,
                            const Quadrature<dim-1>                  &quadrature,
index b174e0969f583b6add4d649de1b668a1e6448d42..0ad66dce357b713ed865afd6f8c840b1512bdb34 100644 (file)
@@ -54,14 +54,14 @@ class MappingCartesian : public Mapping<dim>
                      const Quadrature<dim-1>& quadrature) const;
 
     virtual void
-    fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_values (const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>& quadrature,
                    typename Mapping<dim>::InternalDataBase &mapping_data,
                    std::vector<Point<dim> >        &quadrature_points,
                    std::vector<double>             &JxW_values) const ;
 
     virtual void
-    fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int face_no,
                         const Quadrature<dim-1>& quadrature,
                         typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -70,7 +70,7 @@ class MappingCartesian : public Mapping<dim>
                         std::vector<Tensor<1,dim> >        &boundary_form,
                         std::vector<Point<dim> >        &normal_vectors) const ;
     virtual void
-    fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int face_no,
                            const unsigned int sub_no,
                            const Quadrature<dim-1>& quadrature,
@@ -193,7 +193,7 @@ class MappingCartesian : public Mapping<dim>
                                      * Do the computation for the
                                      * <tt>fill_*</tt> functions.
                                      */
-    void compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
+    void compute_fill (const typename Triangulation<dim>::cell_iterator &cell,
                       const unsigned int face_no,
                       const unsigned int sub_no,
                       InternalData& data,
@@ -216,7 +216,7 @@ class MappingCartesian : public Mapping<dim>
 /// @if NoDoc
 
 template <> void MappingCartesian<1>::fill_fe_face_values (
-  const DoFHandler<1>::cell_iterator &,
+  const Triangulation<1>::cell_iterator &,
   const unsigned,
   const Quadrature<0>&,
   Mapping<1>::InternalDataBase&,
@@ -226,7 +226,7 @@ template <> void MappingCartesian<1>::fill_fe_face_values (
   std::vector<Point<1> >&) const;
 
 template <> void MappingCartesian<1>::fill_fe_subface_values (
-  const DoFHandler<1>::cell_iterator &,
+  const Triangulation<1>::cell_iterator &,
   const unsigned,
   const unsigned,
   const Quadrature<0>&,
index 100b4fc95285a56ba44881e723d531fa62391c9f..16b7d79d9f2780670884bebdf5938a52d5859389 100644 (file)
@@ -184,7 +184,7 @@ class MappingQ : public MappingQ1<dim>
                                      * Mapping.
                                      */
     virtual void
-    fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_values (const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>                &quadrature,
                    typename Mapping<dim>::InternalDataBase &mapping_data,
                    typename std::vector<Point<dim> >             &quadrature_points,
@@ -195,7 +195,7 @@ class MappingQ : public MappingQ1<dim>
                                      * Mapping.
                                      */
     virtual void
-    fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int face_no,
                         const Quadrature<dim-1>& quadrature,
                         typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -209,7 +209,7 @@ class MappingQ : public MappingQ1<dim>
                                      * Mapping.
                                      */
     virtual void
-    fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int face_no,
                            const unsigned int sub_no,
                            const Quadrature<dim-1>& quadrature,
index c67cd38f121a8ab6403b9a19fd0cb053aef35925..3dfe781fbc4d98e1ad61f1af1275d299c4e4df09 100644 (file)
@@ -324,7 +324,7 @@ class MappingQ1 : public Mapping<dim>
                                          * @p mapping_support_points are
                                          * stored.
                                          */
-       typename DoFHandler<dim>::cell_iterator cell_of_current_support_points;
+       typename Triangulation<dim>::cell_iterator cell_of_current_support_points;
        
                                         /**
                                          * Default value of this flag
@@ -374,7 +374,7 @@ class MappingQ1 : public Mapping<dim>
                                      * Mapping.
                                      */
     virtual void
-    fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_values (const typename Triangulation<dim>::cell_iterator &cell,
                    const Quadrature<dim>& quadrature,
                    typename Mapping<dim>::InternalDataBase &mapping_data,
                    typename std::vector<Point<dim> >       &quadrature_points,
@@ -385,7 +385,7 @@ class MappingQ1 : public Mapping<dim>
                                      * Mapping.
                                      */
     virtual void
-    fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int face_no,
                         const Quadrature<dim-1>& quadrature,
                         typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -399,7 +399,7 @@ class MappingQ1 : public Mapping<dim>
                                      * Mapping.
                                      */
     virtual void
-    fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+    fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int face_no,
                            const unsigned int sub_no,
                            const Quadrature<dim-1>& quadrature,
@@ -458,7 +458,7 @@ class MappingQ1 : public Mapping<dim>
                                      * Do the computation for the
                                      * <tt>fill_*</tt> functions.
                                      */
-    void compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
+    void compute_fill (const typename Triangulation<dim>::cell_iterator &cell,
                       const unsigned int      npts,
                       const DataSetDescriptor data_set,
                       InternalData           &data,
@@ -468,7 +468,7 @@ class MappingQ1 : public Mapping<dim>
                                      * Do the computation for the
                                      * <tt>fill_*</tt> functions.
                                      */
-    void compute_fill_face (const typename DoFHandler<dim>::cell_iterator &cell,
+    void compute_fill_face (const typename Triangulation<dim>::cell_iterator &cell,
                            const unsigned int      face_no,
                            const bool              is_subface,
                            const unsigned int      npts,
@@ -697,7 +697,7 @@ MappingQ1<1>::compute_face_data (const UpdateFlags,
                                  InternalData &) const;
 
 template <> void MappingQ1<1>::compute_fill_face (
-  const DoFHandler<1>::cell_iterator &,
+  const Triangulation<1>::cell_iterator &,
   const unsigned int,
   const bool,
   const unsigned int,
@@ -710,7 +710,7 @@ template <> void MappingQ1<1>::compute_fill_face (
   std::vector<Point<1> > &) const;
 
 template <> void MappingQ1<1>::fill_fe_face_values (
-  const DoFHandler<1>::cell_iterator &,
+  const Triangulation<1>::cell_iterator &,
   const unsigned,
   const Quadrature<0>&,
   Mapping<1>::InternalDataBase&,
@@ -720,7 +720,7 @@ template <> void MappingQ1<1>::fill_fe_face_values (
   std::vector<Point<1> >&) const;
 
 template <> void MappingQ1<1>::fill_fe_subface_values (
-  const DoFHandler<1>::cell_iterator &,
+  const Triangulation<1>::cell_iterator &,
   const unsigned,
   const unsigned,
   const Quadrature<0>&,
index 93b87e4aab760e79d8fd1ef8daeaf20bffb54f77..0a839120561b98434d48b5e89e8ae75309f83e58 100644 (file)
@@ -12,9 +12,9 @@
 //----------------------------  fe.cc  ---------------------------
 
 
-#include <fe/fe.h>
 #include <base/memory_consumption.h>
 #include <fe/mapping.h>
+#include <fe/fe.h>
 #include <fe/fe_values.h>
 #include <base/quadrature.h>
 #include <grid/tria.h>
@@ -462,7 +462,7 @@ template <int dim>
 void
 FiniteElementBase<dim>::
 compute_2nd (const Mapping<dim>                   &mapping,
-            const typename DoFHandler<dim>::cell_iterator &cell,
+            const typename Triangulation<dim>::cell_iterator &cell,
             const unsigned int,
             typename Mapping<dim>::InternalDataBase &mapping_internal,
             InternalDataBase                     &fe_internal,
index 5769609076acb98a0bddc64bcd8df29191ef2ede..5d3f987a5b55fa01994b103678e1359973039454 100644 (file)
@@ -12,9 +12,8 @@
 //----------------------------------------------------------------
 
 
-#include <fe/fe.h>
 #include <grid/geometry_info.h>
-
+#include <fe/fe.h>
 
 
 template <int dim>
index 32ba0c7bbb2ab8641368444aa5ee36bd4681eb5e..e9111b188caeacc10008ea7c28d8f855b3c8b160 100644 (file)
@@ -286,7 +286,7 @@ template <int dim>
 void
 FE_DGPNonparametric<dim>::fill_fe_values (
   const Mapping<dim>&,
-  const typename DoFHandler<dim>::cell_iterator&,
+  const typename Triangulation<dim>::cell_iterator&,
   const Quadrature<dim>&,
   typename Mapping<dim>::InternalDataBase&,
   typename Mapping<dim>::InternalDataBase& fedata,
@@ -326,7 +326,7 @@ template <int dim>
 void
 FE_DGPNonparametric<dim>::fill_fe_face_values (
   const Mapping<dim>&,
-  const typename DoFHandler<dim>::cell_iterator&,
+  const typename Triangulation<dim>::cell_iterator&,
   const unsigned int,
   const Quadrature<dim-1>&,
   typename Mapping<dim>::InternalDataBase&,
@@ -367,7 +367,7 @@ template <int dim>
 void
 FE_DGPNonparametric<dim>::fill_fe_subface_values (
   const Mapping<dim>&,
-  const typename DoFHandler<dim>::cell_iterator&,
+  const typename Triangulation<dim>::cell_iterator&,
   const unsigned int,
   const unsigned int,
   const Quadrature<dim-1>&,
index c679cc49bf4bf0d981d59420830d53c444c49a78..bd7dc5fc4c0979833f48d809f75c49adc9e790b4 100644 (file)
@@ -537,7 +537,6 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
 // Fill data of FEValues
 //----------------------------------------------------------------------
 
-
 template <int dim>
 bool
 FE_DGQ<dim>::has_support_on_face (const unsigned int shape_index,
index 4ae45426a6af70fbb872475795c89bc5092e00f3..92c7788f45a9dee56577cd66375f7a833d9eec84 100644 (file)
@@ -979,7 +979,7 @@ FE_Nedelec<dim>::get_data (const UpdateFlags      update_flags,
 template <int dim>
 void
 FE_Nedelec<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
-                                const typename DoFHandler<dim>::cell_iterator &cell,
+                                const typename Triangulation<dim>::cell_iterator &cell,
                                 const Quadrature<dim>                &quadrature,
                                 typename Mapping<dim>::InternalDataBase &mapping_data,
                                 typename Mapping<dim>::InternalDataBase &fedata,
@@ -1104,7 +1104,7 @@ FE_Nedelec<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
 template <int dim>
 void
 FE_Nedelec<dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
-                                     const typename DoFHandler<dim>::cell_iterator &cell,
+                                     const typename Triangulation<dim>::cell_iterator &cell,
                                      const unsigned int                    face,
                                      const Quadrature<dim-1>              &quadrature,
                                      typename Mapping<dim>::InternalDataBase       &mapping_data,
@@ -1249,7 +1249,7 @@ FE_Nedelec<dim>::fill_fe_face_values (const Mapping<dim>                   &mapp
 template <int dim>
 void
 FE_Nedelec<dim>::fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                                        const typename DoFHandler<dim>::cell_iterator &cell,
+                                        const typename Triangulation<dim>::cell_iterator &cell,
                                         const unsigned int                    face,
                                         const unsigned int                    subface,
                                         const Quadrature<dim-1>              &quadrature,
index 83e55ce5faa92a9f3434c1aad734fe17a601f013..38e8f75c1a9bdffc9169dda7c34002e46dd5b6e3 100644 (file)
@@ -235,11 +235,11 @@ FE_Poly<POLY,dim>::get_data (const UpdateFlags      update_flags,
 template <class POLY, int dim>
 void
 FE_Poly<POLY,dim>::fill_fe_values (const Mapping<dim>                   &mapping,
-                                   const typename DoFHandler<dim>::cell_iterator &cell,
-                                   const Quadrature<dim>                &quadrature,
-                                   typename Mapping<dim>::InternalDataBase &mapping_data,
-                                   typename Mapping<dim>::InternalDataBase &fedata,
-                                   FEValuesData<dim>                    &data) const
+                                  const typename Triangulation<dim>::cell_iterator &cell,
+                                  const Quadrature<dim>                &quadrature,
+                                  typename Mapping<dim>::InternalDataBase &mapping_data,
+                                  typename Mapping<dim>::InternalDataBase &fedata,
+                                  FEValuesData<dim>                    &data) const
 {
                                   // convert data object to internal
                                   // data for this class. fails with
@@ -277,7 +277,7 @@ FE_Poly<POLY,dim>::fill_fe_values (const Mapping<dim>                   &mapping
 template <class POLY, int dim>
 void
 FE_Poly<POLY,dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
-                               const typename DoFHandler<dim>::cell_iterator &cell,
+                               const typename Triangulation<dim>::cell_iterator &cell,
                                const unsigned int                    face,
                                const Quadrature<dim-1>              &quadrature,
                                typename Mapping<dim>::InternalDataBase       &mapping_data,
@@ -327,7 +327,7 @@ FE_Poly<POLY,dim>::fill_fe_face_values (const Mapping<dim>                   &ma
 template <class POLY, int dim>
 void
 FE_Poly<POLY,dim>::fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                                  const typename DoFHandler<dim>::cell_iterator &cell,
+                                  const typename Triangulation<dim>::cell_iterator &cell,
                                   const unsigned int                    face,
                                   const unsigned int                    subface,
                                   const Quadrature<dim-1>              &quadrature,
index d925e7cde30f9ec8080215526545da77ba862765..8c820606b6fc2bb9b84512b73bf59327a2c362cf 100644 (file)
@@ -930,7 +930,6 @@ FE_Q<1>::initialize_constraints ()
 
 #endif
 
-
 #if deal_II_dimension == 2
 
 template <>
index e045278fe1e9b5cce599f5b9574531e99320cf42..1bc800305a89d8638cace8a23f2ed51991150e4c 100644 (file)
@@ -1493,7 +1493,7 @@ FE_RaviartThomas<dim>::get_data (const UpdateFlags      update_flags,
 template <int dim>
 void
 FE_RaviartThomas<dim>::fill_fe_values (const Mapping<dim>                   &mapping,
-                                       const typename DoFHandler<dim>::cell_iterator &cell,
+                                       const typename Triangulation<dim>::cell_iterator &cell,
                                        const Quadrature<dim>                &quadrature,
                                        typename Mapping<dim>::InternalDataBase &mapping_data,
                                        typename Mapping<dim>::InternalDataBase &fedata,
@@ -1618,7 +1618,7 @@ FE_RaviartThomas<dim>::fill_fe_values (const Mapping<dim>                   &map
 template <int dim>
 void
 FE_RaviartThomas<dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
-                                            const typename DoFHandler<dim>::cell_iterator &cell,
+                                            const typename Triangulation<dim>::cell_iterator &cell,
                                             const unsigned int                    face,
                                             const Quadrature<dim-1>              &quadrature,
                                             typename Mapping<dim>::InternalDataBase       &mapping_data,
@@ -1755,7 +1755,7 @@ FE_RaviartThomas<dim>::fill_fe_face_values (const Mapping<dim>
 template <int dim>
 void
 FE_RaviartThomas<dim>::fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                                               const typename DoFHandler<dim>::cell_iterator &cell,
+                                               const typename Triangulation<dim>::cell_iterator &cell,
                                                const unsigned int                    face,
                                                const unsigned int                    subface,
                                                const Quadrature<dim-1>              &quadrature,
index ad1acc86d257ebbc860d2719a611f33af7bcc67f..4a21e06a252bce71057753c2de4a0fd3b305e584 100644 (file)
@@ -670,7 +670,7 @@ template <int dim>
 void
 FESystem<dim>::
 fill_fe_values (const Mapping<dim>                   &mapping,
-                const typename DoFHandler<dim>::cell_iterator &cell,
+                const typename Triangulation<dim>::cell_iterator &cell,
                 const Quadrature<dim>                &quadrature,
                 typename Mapping<dim>::InternalDataBase &mapping_data,
                 typename Mapping<dim>::InternalDataBase &fe_data,
@@ -686,7 +686,7 @@ template <int dim>
 void
 FESystem<dim>::
 fill_fe_face_values (const Mapping<dim>                   &mapping,
-                     const typename DoFHandler<dim>::cell_iterator &cell,
+                     const typename Triangulation<dim>::cell_iterator &cell,
                      const unsigned int                    face_no,
                      const Quadrature<dim-1>              &quadrature,
                      typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -704,7 +704,7 @@ template <int dim>
 void
 FESystem<dim>::
 fill_fe_subface_values (const Mapping<dim>                   &mapping,
-                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const typename Triangulation<dim>::cell_iterator &cell,
                         const unsigned int                    face_no,
                         const unsigned int                    sub_no,
                         const Quadrature<dim-1>              &quadrature,
@@ -723,7 +723,7 @@ template <int dim_1>
 void
 FESystem<dim>::
 compute_fill (const Mapping<dim>                   &mapping,
-              const typename DoFHandler<dim>::cell_iterator &cell,
+              const typename Triangulation<dim>::cell_iterator &cell,
               const unsigned int                    face_no,
               const unsigned int                    sub_no,
               const Quadrature<dim_1>              &quadrature,
index 28eeb32df805b3cf88a67cfc49332efb632f045d..c6117959e6db2c19642e776f19289cc0d197e035 100644 (file)
@@ -34,6 +34,8 @@
 #include <dofs/dof_handler.h>
 #include <dofs/dof_accessor.h>
 
+#include <boost/shared_ptr.hpp>
+
 #ifdef HAVE_STD_STRINGSTREAM
 #  include <sstream>
 #else
 
 namespace 
 {
+  template <int dim>
+  unsigned int
+  max_dofs_per_cell (const DoFHandler<dim> &dh) 
+  {
+    return dh.get_fe().dofs_per_cell;
+  }
+
+
+  template <int dim>
+  unsigned int
+  max_dofs_per_face (const DoFHandler<dim> &dh) 
+  {
+    return dh.get_fe().dofs_per_face;
+  }
+
+
+  template <int dim>
+  unsigned int
+  n_components (const DoFHandler<dim> &dh) 
+  {
+    return dh.get_fe().n_components();
+  }
+
+
                                    // forwarder function for
                                    // FE::get_interpolation_matrix. we
                                    // will want to call that function
@@ -76,8 +102,6 @@ namespace
     interpolation_matrix = tmp;
   }
 
-
-
                                   // return true if the given pattern
                                   // string matches the given name at
                                   // the first position of the string
@@ -352,6 +376,7 @@ void FETools::get_interpolation_difference_matrix (const FiniteElement<dim> &fe1
 }
 
 
+
 template <int dim, typename number>
 void FETools::get_projection_matrix (const FiniteElement<dim> &fe1,
                                     const FiniteElement<dim> &fe2,
@@ -444,11 +469,16 @@ void FETools::get_projection_matrix (const FiniteElement<dim> &fe1,
 }
 
 
-template <int dim, class InVector, class OutVector>
-void FETools::interpolate(const DoFHandler<dim> &dof1,
-                         const InVector &u1,
-                         const DoFHandler<dim> &dof2,
-                         OutVector &u2)
+
+template <int dim,
+          template <int> class DH1,
+          template <int> class DH2,
+          class InVector, class OutVector>
+void
+FETools::interpolate(const DH1<dim> &dof1,
+                     const InVector &u1,
+                     const DH2<dim> &dof2,
+                     OutVector      &u2)
 {
   ConstraintMatrix dummy;
   dummy.close();
@@ -457,17 +487,21 @@ void FETools::interpolate(const DoFHandler<dim> &dof1,
 
   
 template <int dim, class InVector, class OutVector>
-void FETools::interpolate(const DoFHandler<dim> &dof1,
-                         const InVector &u1,
-                         const DoFHandler<dim> &dof2,
-                         const ConstraintMatrix &constraints,
-                         OutVector &u2)
+void
+FETools::interpolate(const DoFHandler<dim>  &dof1,
+                     const InVector         &u1,
+                     const DoFHandler<dim>  &dof2,
+                     const ConstraintMatrix &constraints,
+                     OutVector              &u2)
 {
-  Assert(&dof1.get_tria()==&dof2.get_tria(), ExcTriangulationMismatch());
+  Assert(&dof1.get_tria() == &dof2.get_tria(), ExcTriangulationMismatch());
   Assert(dof1.get_fe().n_components() == dof2.get_fe().n_components(),
-        ExcDimensionMismatch(dof1.get_fe().n_components(), dof2.get_fe().n_components()));
-  Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
-  Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs()));
+        ExcDimensionMismatch(dof1.get_fe().n_components(),
+                              dof2.get_fe().n_components()));
+  Assert(u1.size()==dof1.n_dofs(),
+         ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
+  Assert(u2.size()==dof2.n_dofs(),
+         ExcDimensionMismatch(u2.size(), dof2.n_dofs()));
 
                                   // for continuous elements on grids
                                   // with hanging nodes we need
@@ -475,9 +509,9 @@ void FETools::interpolate(const DoFHandler<dim> &dof1,
                                   // constraints. Consequentely, if
                                   // there are no constraints then
                                   // hanging nodes are not allowed.
-  const bool hanging_nodes_not_allowed=
-    (dof2.get_fe().dofs_per_vertex != 0)
-    && (constraints.n_constraints()==0);
+  const bool hanging_nodes_not_allowed
+    = ((dof2.get_fe().dofs_per_vertex != 0) &&
+       (constraints.n_constraints() == 0));
   
   const unsigned int dofs_per_cell1=dof1.get_fe().dofs_per_cell;
   const unsigned int dofs_per_cell2=dof2.get_fe().dofs_per_cell;
@@ -495,11 +529,11 @@ void FETools::interpolate(const DoFHandler<dim> &dof1,
                                                 cell2 = dof2.begin_active(),
                                                 endc2 = dof2.end();
 
-  std::vector<unsigned int> index_multiplicity(dof2.n_dofs(),0);
+  std::vector<unsigned int> touch_count(dof2.n_dofs(),0);
   std::vector<unsigned int> dofs (dofs_per_cell2);
   u2 = 0;
 
-  for (; cell1!=endc1, cell2!=endc2; ++cell1, ++cell2) 
+  for (; cell1!=endc1; ++cell1, ++cell2) 
     {
       if (hanging_nodes_not_allowed)
        for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
@@ -513,17 +547,142 @@ void FETools::interpolate(const DoFHandler<dim> &dof1,
       for (unsigned int i=0; i<dofs_per_cell2; ++i)
        {
          u2(dofs[i])+=u2_local(i);
-         ++index_multiplicity[dofs[i]];
+         ++touch_count[dofs[i]];
+       }
+    }
+                                   // cell1 is at the end, so should
+                                   // be cell2
+  Assert (cell2 == endc2, ExcInternalError());
+  
+                                  // when a discontinuous element is
+                                  // interpolated to a continuous
+                                  // one, we take the mean values.
+  for (unsigned int i=0; i<dof2.n_dofs(); ++i)
+    {
+      Assert(touch_count[i]!=0, ExcInternalError());
+      u2(i) /= touch_count[i];
+    }
+
+                                  // Apply hanging node constraints.
+  constraints.distribute(u2);
+}
+
+
+
+template <int dim,
+          template <int> class DH1,
+          template <int> class DH2,              
+          class InVector, class OutVector>
+void
+FETools::interpolate (const DH1<dim>         &dof1,
+                      const InVector         &u1,
+                      const DH2<dim>         &dof2,
+                      const ConstraintMatrix &constraints,
+                      OutVector              &u2)
+{
+  Assert(&dof1.get_tria()==&dof2.get_tria(), ExcTriangulationMismatch());
+
+  Assert(u1.size()==dof1.n_dofs(),
+         ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
+  Assert(u2.size()==dof2.n_dofs(),
+         ExcDimensionMismatch(u2.size(), dof2.n_dofs()));
+
+                                   // allocate vectors at maximal
+                                   // size. will be reinited in inner
+                                   // cell, but Vector makes sure that
+                                   // this does not lead to
+                                   // reallocation of memory
+  Vector<typename OutVector::value_type> u1_local(max_dofs_per_cell(dof1));
+  Vector<typename OutVector::value_type> u2_local(max_dofs_per_cell(dof2));
+
+                                   // have a map for interpolation
+                                   // matrices. shared_ptr make sure
+                                   // that memory is released again
+  std::map<const FiniteElement<dim> *,
+           std::map<const FiniteElement<dim> *,
+                    boost::shared_ptr<FullMatrix<double> > > >
+    interpolation_matrices;
+  
+  typename DH1<dim>::active_cell_iterator cell1 = dof1.begin_active(),
+                                          endc1 = dof1.end();
+  typename DH2<dim>::active_cell_iterator cell2 = dof2.begin_active(),
+                                          endc2 = dof2.end();
+
+  std::vector<unsigned int> touch_count(dof2.n_dofs(),0);
+  std::vector<unsigned int> dofs;
+  dofs.reserve (max_dofs_per_cell (dof2));
+  u2.clear ();
+
+  for (; cell1!=endc1; ++cell1, ++cell2) 
+    {
+      Assert(cell1->get_fe().n_components() == cell2->get_fe().n_components(),
+             ExcDimensionMismatch (cell1->get_fe().n_components(),
+                                   cell2->get_fe().n_components()));
+
+                                       // for continuous elements on
+                                       // grids with hanging nodes we
+                                       // need hanging node
+                                       // constraints. Consequentely,
+                                       // if there are no constraints
+                                       // then hanging nodes are not
+                                       // allowed.
+      const bool hanging_nodes_not_allowed
+        = ((cell2->get_fe().dofs_per_vertex != 0) &&
+           (constraints.n_constraints() == 0));
+  
+      if (hanging_nodes_not_allowed)
+       for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+         Assert (cell1->at_boundary(face) ||
+                 cell1->neighbor(face)->level() == cell1->level(),
+                 ExcHangingNodesNotAllowed(0));
+      
+      
+      const unsigned int dofs_per_cell1 = cell1->get_fe().dofs_per_cell;
+      const unsigned int dofs_per_cell2 = cell2->get_fe().dofs_per_cell;
+      u1_local.reinit (dofs_per_cell1);
+      u2_local.reinit (dofs_per_cell2);
+
+                                       // check if interpolation
+                                       // matrix for this particular
+                                       // pair of elements is already
+                                       // there
+      if (interpolation_matrices[&cell1->get_fe()][&cell2->get_fe()] == 0)
+        {
+          boost::shared_ptr<FullMatrix<double> >
+            interpolation_matrix (new FullMatrix<double> (dofs_per_cell2,
+                                                          dofs_per_cell1));
+          interpolation_matrices[&cell1->get_fe()][&cell2->get_fe()]
+            = interpolation_matrix;
+            
+          FETools::get_interpolation_matrix(cell1->get_fe(),
+                                            cell2->get_fe(),
+                                            *interpolation_matrix);
+        }
+      
+      cell1->get_dof_values(u1, u1_local);
+      interpolation_matrices[&cell1->get_fe()][&cell2->get_fe()]
+        ->vmult(u2_local, u1_local);
+
+      dofs.resize (dofs_per_cell2);
+      cell2->get_dof_indices(dofs);
+
+      for (unsigned int i=0; i<dofs_per_cell2; ++i)
+       {
+         u2(dofs[i])+=u2_local(i);
+         ++touch_count[dofs[i]];
        }
     }
+                                   // cell1 is at the end, so should
+                                   // be cell2
+  Assert (cell2 == endc2, ExcInternalError());
 
                                   // when a discontinuous element is
                                   // interpolated to a continuous
                                   // one, we take the mean values.
   for (unsigned int i=0; i<dof2.n_dofs(); ++i)
     {
-      Assert(index_multiplicity[i]!=0, ExcInternalError());
-      u2(i) /= index_multiplicity[i];
+      Assert(touch_count[i]!=0, ExcInternalError());
+      u2(i) /= touch_count[i];
     }
 
                                   // Apply hanging node constraints.
@@ -531,11 +690,13 @@ void FETools::interpolate(const DoFHandler<dim> &dof1,
 }
 
 
+
 template <int dim, class InVector, class OutVector>
-void FETools::back_interpolate(const DoFHandler<dim> &dof1,
-                              const InVector &u1,
-                              const FiniteElement<dim> &fe2,
-                              OutVector &u1_interpolated)
+void
+FETools::back_interpolate(const DoFHandler<dim>    &dof1,
+                          const InVector           &u1,
+                          const FiniteElement<dim> &fe2,
+                          OutVector                &u1_interpolated)
 {
   Assert(dof1.get_fe().n_components() == fe2.n_components(),
         ExcDimensionMismatch(dof1.get_fe().n_components(), fe2.n_components()));
@@ -545,7 +706,7 @@ void FETools::back_interpolate(const DoFHandler<dim> &dof1,
   
                                   // For continuous elements on grids
                                   // with hanging nodes we need
-                                  // hnaging node
+                                  // hanging node
                                   // constraints. Consequently, when
                                   // the elements are continuous no
                                   // hanging node constraints are
@@ -580,6 +741,78 @@ void FETools::back_interpolate(const DoFHandler<dim> &dof1,
 
 
   
+template <int dim,
+          template <int> class DH,
+          class InVector, class OutVector>
+void
+FETools::back_interpolate(const DH<dim>            &dof1,
+                          const InVector           &u1,
+                          const FiniteElement<dim> &fe2,
+                          OutVector                &u1_interpolated)
+{
+  Assert(u1.size() == dof1.n_dofs(),
+         ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
+  Assert(u1_interpolated.size() == dof1.n_dofs(),
+        ExcDimensionMismatch(u1_interpolated.size(), dof1.n_dofs()));
+  
+  Vector<typename OutVector::value_type> u1_local(max_dofs_per_cell(dof1));
+  Vector<typename OutVector::value_type> u1_int_local(max_dofs_per_cell(dof1));
+  
+  typename DH<dim>::active_cell_iterator cell = dof1.begin_active(),
+                                         endc = dof1.end();
+
+                                   // map from possible fe objects in
+                                   // dof1 to the back_interpolation
+                                   // matrices
+  std::map<const FiniteElement<dim> *,
+           boost::shared_ptr<FullMatrix<double> > > interpolation_matrices;
+  
+  for (; cell!=endc; ++cell) 
+    {
+      Assert(cell->get_fe().n_components() == fe2.n_components(),
+             ExcDimensionMismatch(cell->get_fe().n_components(),
+                                  fe2.n_components()));
+      
+                                       // For continuous elements on
+                                       // grids with hanging nodes we
+                                       // need hanging node
+                                       // constraints. Consequently,
+                                       // when the elements are
+                                       // continuous no hanging node
+                                       // constraints are allowed.
+      const bool hanging_nodes_not_allowed=
+        (cell->get_fe().dofs_per_vertex != 0) || (fe2.dofs_per_vertex != 0);
+      
+      if (hanging_nodes_not_allowed)
+       for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+         Assert (cell->at_boundary(face) ||
+                 cell->neighbor(face)->level() == cell->level(),
+                 ExcHangingNodesNotAllowed(0));
+
+      const unsigned int dofs_per_cell1 = cell->get_fe().dofs_per_cell;
+      
+                                       // make sure back_interpolation
+                                       // matrix is available
+      if (interpolation_matrices[&cell->get_fe()] != 0)
+        {
+          interpolation_matrices[&cell->get_fe()] =
+            boost::shared_ptr<FullMatrix<double> >
+            (new FullMatrix<double>(dofs_per_cell1, dofs_per_cell1));
+          get_back_interpolation_matrix(dof1.get_fe(), fe2,
+                                        *interpolation_matrices[&cell->get_fe()]);
+        }      
+      
+      u1_local.reinit (dofs_per_cell1);
+      u1_int_local.reinit (dofs_per_cell1);
+      
+      cell->get_dof_values(u1, u1_local);
+      interpolation_matrices[&cell->get_fe()]->vmult(u1_int_local, u1_local);
+      cell->set_dof_values(u1_int_local, u1_interpolated);
+    }
+}
+
+
+  
 template <int dim, class InVector, class OutVector>
 void FETools::back_interpolate(const DoFHandler<dim> &dof1,
                               const ConstraintMatrix &constraints1,
index bc65ab2486c5d21536ca39aa71ba393e25d817ff..6ae075c98e3fe39d1ac7ac01df2af07aa4159f60 100644 (file)
@@ -12,9 +12,6 @@
 //----------------------------  fe_values.cc  ---------------------------
 
 
-#include <fe/fe.h>
-#include <fe/mapping_q1.h>
-#include <fe/fe_values.h>
 #include <base/memory_consumption.h>
 #include <base/quadrature.h>
 #include <lac/vector.h>
 #include <grid/tria_accessor.h>
 #include <grid/tria_boundary.h>
 #include <dofs/dof_accessor.h>
+#include <fe/mapping_q1.h>
+#include <fe/fe_values.h>
+#include <fe/fe.h>
 
 #include <iomanip>
 
 
+/* ---------------- FEValuesBase<dim>::CellIteratorBase --------- */
+
+template <int dim>
+inline
+FEValuesBase<dim>::CellIteratorBase::~CellIteratorBase ()
+{}
+
+
+
+/* ---------------- FEValuesBase<dim>::CellIterator<CI> --------- */
+
+
+template <int dim>
+template <typename CI>
+inline
+FEValuesBase<dim>::CellIterator<CI>::CellIterator (const CI &cell)
+                :
+                cell(cell)
+{}
+
+
+
+template <int dim>
+template <typename CI>
+FEValuesBase<dim>::CellIterator<CI>::
+operator const typename Triangulation<dim>::cell_iterator () const
+{
+  return cell;
+}
+
+
+
+template <int dim>
+template <typename CI>
+unsigned int
+FEValuesBase<dim>::CellIterator<CI>::n_dofs_for_dof_handler () const
+{
+  return cell->get_dof_handler().n_dofs();
+}
+
+
+
+template <int dim>
+template <typename CI>
+void
+FEValuesBase<dim>::CellIterator<CI>::
+get_interpolated_dof_values (const Vector<double> &in,
+                             Vector<double>       &out) const
+{
+  cell->get_interpolated_dof_values (in, out);
+}
+
+
+template <int dim>
+template <typename CI>
+void
+FEValuesBase<dim>::CellIterator<CI>::
+get_interpolated_dof_values (const Vector<float> &in,
+                             Vector<float>       &out) const
+{
+  cell->get_interpolated_dof_values (in, out);
+}
+
+
+template <int dim>
+template <typename CI>
+void
+FEValuesBase<dim>::CellIterator<CI>::
+get_interpolated_dof_values (const BlockVector<double> &in,
+                             Vector<double>            &out) const
+{
+  cell->get_interpolated_dof_values (in, out);
+}
+
+
+template <int dim>
+template <typename CI>
+void
+FEValuesBase<dim>::CellIterator<CI>::
+get_interpolated_dof_values (const BlockVector<float> &in,
+                             Vector<float>            &out) const
+{
+  cell->get_interpolated_dof_values (in, out);
+}
+
+
+
+#ifdef DEAL_II_USE_PETSC
+
+template <int dim>
+template <typename CI>
+void
+FEValuesBase<dim>::CellIterator<CI>::
+get_interpolated_dof_values (const PETScWrappers::Vector &in,
+                             Vector<PetscScalar>         &out) const
+{
+  cell->get_interpolated_dof_values (in, out);
+}
+
+
+
+template <int dim>
+template <typename CI>
+void
+FEValuesBase<dim>::CellIterator<CI>::
+get_interpolated_dof_values (const PETScWrappers::BlockVector &in,
+                             Vector<PetscScalar>              &out) const
+{
+  cell->get_interpolated_dof_values (in, out);
+}
+
+#endif
+
+
+/* ---------------- FEValuesBase<dim>::TriaCellIterator --------- */
+
+template <int dim>
+const char * const
+FEValuesBase<dim>::TriaCellIterator::message_string
+= ("You have previously called the FEValues::reinit function with a\n"
+   "cell iterator of type Triangulation<dim>::cell_iterator. However,\n"
+   "when you do this, you cannot call some functions in the FEValues\n"
+   "class, such as the get_function_values/grads/2nd_derivatives\n"
+   "functions. If you need these functions, then you need to call\n"
+   "FEValues::reinit with an iterator type that allows to extract\n"
+   "degrees of freedom, such as DoFHandler<dim>::cell_iterator.");
+
+
+template <int dim>
+inline
+FEValuesBase<dim>::TriaCellIterator::
+TriaCellIterator (const typename Triangulation<dim>::cell_iterator &cell)
+                :
+                cell(cell)
+{}
+
+
+
+template <int dim>
+FEValuesBase<dim>::TriaCellIterator::
+operator const typename Triangulation<dim>::cell_iterator () const
+{
+  return cell;
+}
+
+
+
+template <int dim>
+unsigned int
+FEValuesBase<dim>::TriaCellIterator::n_dofs_for_dof_handler () const
+{
+  Assert (false, ExcMessage (message_string));
+  return 0;
+}
+
+
+template <int dim>
+void
+FEValuesBase<dim>::TriaCellIterator::
+get_interpolated_dof_values (const Vector<double> &,
+                             Vector<double>       &) const
+{
+  Assert (false, ExcMessage (message_string));
+}
+
+
+template <int dim>
+void
+FEValuesBase<dim>::TriaCellIterator::
+get_interpolated_dof_values (const Vector<float> &,
+                             Vector<float>       &) const
+{
+  Assert (false, ExcMessage (message_string));
+}
+
+
+template <int dim>
+void
+FEValuesBase<dim>::TriaCellIterator::
+get_interpolated_dof_values (const BlockVector<double> &,
+                             Vector<double>            &) const
+{
+  Assert (false, ExcMessage (message_string));
+}
+
+
+
+template <int dim>
+void
+FEValuesBase<dim>::TriaCellIterator::
+get_interpolated_dof_values (const BlockVector<float> &,
+                             Vector<float>            &) const
+{
+  Assert (false, ExcMessage (message_string));
+}
+
+
+
+#ifdef DEAL_II_USE_PETSC
+
+template <int dim>
+void
+FEValuesBase<dim>::TriaCellIterator::
+get_interpolated_dof_values (const PETScWrappers::Vector &,
+                             Vector<PetscScalar>         &) const
+{
+  Assert (false, ExcMessage (message_string));
+}
+
+
+
+template <int dim>
+void
+FEValuesBase<dim>::TriaCellIterator::
+get_interpolated_dof_values (const PETScWrappers::BlockVector &,
+                             Vector<PetscScalar>              &) const
+{
+  Assert (false, ExcMessage (message_string));
+}
+
+#endif
+
+
+
+
+/* --------------------- FEValuesData ----------------- */
+
+
 template <int dim>
 void
 FEValuesData<dim>::initialize (const unsigned int        n_quadrature_points,
@@ -144,16 +372,13 @@ void FEValuesBase<dim>::get_function_values (const InputVector            &fe_fu
          ExcWrongNoOfComponents());
   Assert (values.size() == n_quadrature_points,
          ExcWrongVectorSize(values.size(), n_quadrature_points));
-  Assert (fe_function.size() == present_cell->get_dof_handler().n_dofs(),
-         ExcWrongVectorSize(fe_function.size(), present_cell->get_dof_handler().n_dofs()));
+  Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(),
+         ExcWrongVectorSize(fe_function.size(), present_cell->n_dofs_for_dof_handler()));
 
                                   // get function values of dofs
                                   // on this cell
   Vector<typename InputVector::value_type> dof_values (dofs_per_cell);
-  if (present_cell->active())
-    present_cell->get_dof_values (fe_function, dof_values);
-  else
-    present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell->get_interpolated_dof_values(fe_function, dof_values);
 
                                   // initialize with zero
   std::fill_n (values.begin(), n_quadrature_points, 0);
@@ -185,16 +410,13 @@ void FEValuesBase<dim>::get_function_values (const InputVector
     Assert (values[i].size() == n_components, ExcWrongNoOfComponents());
 
   Assert (this->update_flags & update_values, ExcAccessToUninitializedField());
-  Assert (fe_function.size() == present_cell->get_dof_handler().n_dofs(),
-         ExcWrongVectorSize(fe_function.size(), present_cell->get_dof_handler().n_dofs()));
+  Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(),
+         ExcWrongVectorSize(fe_function.size(), present_cell->n_dofs_for_dof_handler()));
     
                                   // get function values of dofs
                                   // on this cell
   Vector<typename InputVector::value_type> dof_values (dofs_per_cell);
-  if (present_cell->active())
-    present_cell->get_dof_values (fe_function, dof_values);
-  else
-    present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell->get_interpolated_dof_values(fe_function, dof_values);
   
                                   // initialize with zero
   for (unsigned i=0;i<values.size();++i)
@@ -252,16 +474,13 @@ get_function_grads (const InputVector           &fe_function,
          ExcWrongNoOfComponents());
   Assert (gradients.size() == n_quadrature_points,
          ExcWrongVectorSize(gradients.size(), n_quadrature_points));
-  Assert (fe_function.size() == present_cell->get_dof_handler().n_dofs(),
-         ExcWrongVectorSize(fe_function.size(), present_cell->get_dof_handler().n_dofs()));
+  Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(),
+         ExcWrongVectorSize(fe_function.size(), present_cell->n_dofs_for_dof_handler()));
 
                                   // get function values of dofs
                                   // on this cell
   Vector<typename InputVector::value_type> dof_values (dofs_per_cell);
-  if (present_cell->active())
-    present_cell->get_dof_values (fe_function, dof_values);
-  else
-    present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell->get_interpolated_dof_values(fe_function, dof_values);
 
                                   // initialize with zero
   std::fill_n (gradients.begin(), n_quadrature_points, Tensor<1,dim>());
@@ -298,16 +517,13 @@ get_function_grads (const InputVector                         &fe_function,
     Assert (gradients[i].size() == n_components, ExcWrongNoOfComponents());
 
   Assert (this->update_flags & update_gradients, ExcAccessToUninitializedField());
-  Assert (fe_function.size() == present_cell->get_dof_handler().n_dofs(),
-         ExcWrongVectorSize(fe_function.size(), present_cell->get_dof_handler().n_dofs()));
+  Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(),
+         ExcWrongVectorSize(fe_function.size(), present_cell->n_dofs_for_dof_handler()));
 
                                   // get function values of dofs
                                   // on this cell
   Vector<typename InputVector::value_type> dof_values (dofs_per_cell);
-  if (present_cell->active())
-    present_cell->get_dof_values (fe_function, dof_values);
-  else
-    present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell->get_interpolated_dof_values(fe_function, dof_values);
 
                                   // initialize with zero
   for (unsigned i=0;i<gradients.size();++i)
@@ -351,16 +567,13 @@ get_function_2nd_derivatives (const InputVector           &fe_function,
   Assert (second_derivatives.size() == n_quadrature_points,
          ExcWrongVectorSize(second_derivatives.size(), n_quadrature_points));
   Assert (this->update_flags & update_second_derivatives, ExcAccessToUninitializedField());
-  Assert (fe_function.size() == present_cell->get_dof_handler().n_dofs(),
-         ExcWrongVectorSize(fe_function.size(), present_cell->get_dof_handler().n_dofs()));
+  Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(),
+         ExcWrongVectorSize(fe_function.size(), present_cell->n_dofs_for_dof_handler()));
 
                                   // get function values of dofs
                                   // on this cell
   Vector<typename InputVector::value_type> dof_values (dofs_per_cell);
-  if (present_cell->active())
-    present_cell->get_dof_values (fe_function, dof_values);
-  else
-    present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell->get_interpolated_dof_values(fe_function, dof_values);
 
                                   // initialize with zero
   std::fill_n (second_derivatives.begin(), n_quadrature_points, Tensor<2,dim>());
@@ -397,16 +610,13 @@ get_function_2nd_derivatives (const InputVector                         &fe_func
     Assert (second_derivs[i].size() == n_components, ExcWrongNoOfComponents());
 
   Assert (this->update_flags & update_second_derivatives, ExcAccessToUninitializedField());
-  Assert (fe_function.size() == present_cell->get_dof_handler().n_dofs(),
-         ExcWrongVectorSize(fe_function.size(), present_cell->get_dof_handler().n_dofs()));
+  Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(),
+         ExcWrongVectorSize(fe_function.size(), present_cell->n_dofs_for_dof_handler()));
 
                                   // get function values of dofs
                                   // on this cell
   Vector<typename InputVector::value_type> dof_values (dofs_per_cell);
-  if (present_cell->active())
-    present_cell->get_dof_values (fe_function, dof_values);
-  else
-    present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell->get_interpolated_dof_values(fe_function, dof_values);
 
                                   // initialize with zero
   for (unsigned i=0;i<second_derivs.size();++i)
@@ -470,7 +680,6 @@ FEValuesBase<dim>::memory_consumption () const
          MemoryConsumption::memory_consumption (this->normal_vectors) +
          MemoryConsumption::memory_consumption (this->boundary_forms) +
          sizeof(this->update_flags) +
-         MemoryConsumption::memory_consumption (present_cell) +
          MemoryConsumption::memory_consumption (n_quadrature_points) +
          MemoryConsumption::memory_consumption (dofs_per_cell) +
          MemoryConsumption::memory_consumption (mapping) +
@@ -580,26 +789,101 @@ FEValues<dim>::initialize (const UpdateFlags update_flags)
 
 
 template <int dim>
-void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell)
+void
+FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell)
 {
                                   // assert that the finite elements
                                   // passed to the constructor and
                                   // used by the DoFHandler used by
                                   // this cell, are the same
   Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
-         static_cast<const FiniteElementData<dim>&>(cell->get_dof_handler().get_fe()),
+         static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
+         typename FEValuesBase<dim>::ExcFEDontMatch());
+
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::template
+     CellIterator<typename DoFHandler<dim>::cell_iterator> (cell));
+
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit ();
+}
+
+
+
+template <int dim>
+void
+FEValues<dim>::reinit (const typename MGDoFHandler<dim>::cell_iterator &cell)
+{
+                                  // assert that the finite elements
+                                  // passed to the constructor and
+                                  // used by the DoFHandler used by
+                                  // this cell, are the same
+  Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
+         static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
          typename FEValuesBase<dim>::ExcFEDontMatch());
 
-  this->present_cell = cell;
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::template
+     CellIterator<typename MGDoFHandler<dim>::cell_iterator> (cell));
+
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit ();
+}
+
+
+
+template <int dim>
+void FEValues<dim>::reinit (const typename Triangulation<dim>::cell_iterator &cell)
+{
+                                   // no FE in this cell, so no check
+                                   // necessary here
+
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::TriaCellIterator (cell));
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit ();
+}
+
+
 
-  this->get_mapping().fill_fe_values(cell,
+template <int dim>
+void FEValues<dim>::do_reinit ()
+{
+  this->get_mapping().fill_fe_values(*this->present_cell,
                                     quadrature,
                                     *this->mapping_data,
                                     this->quadrature_points,
                                     this->JxW_values);
   
   this->get_fe().fill_fe_values(this->get_mapping(),
-                               cell,
+                               *this->present_cell,
                                quadrature,
                                *this->mapping_data,
                                *this->fe_data,
@@ -740,9 +1024,57 @@ void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &c
          static_cast<const FiniteElementData<dim>&>(cell->get_dof_handler().get_fe()),
          typename FEValuesBase<dim>::ExcFEDontMatch());
 
-  this->present_cell = cell;
+  Assert (face_no < GeometryInfo<dim>::faces_per_cell,
+         ExcIndexRange (face_no, 0, GeometryInfo<dim>::faces_per_cell));
   
-  this->get_mapping().fill_fe_face_values(cell, face_no,
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::template
+     CellIterator<typename DoFHandler<dim>::cell_iterator> (cell));
+
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit (face_no);
+}
+
+
+
+template <int dim>
+void FEFaceValues<dim>::reinit (const typename Triangulation<dim>::cell_iterator &cell,
+                               const unsigned int              face_no)
+{
+  Assert (face_no < GeometryInfo<dim>::faces_per_cell,
+         ExcIndexRange (face_no, 0, GeometryInfo<dim>::faces_per_cell));
+
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::TriaCellIterator (cell));
+
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit (face_no);
+}
+
+
+
+template <int dim>
+void FEFaceValues<dim>::do_reinit (const unsigned int face_no)
+{
+  this->get_mapping().fill_fe_face_values(*this->present_cell, face_no,
                                          this->quadrature,
                                          *this->mapping_data,
                                          this->quadrature_points,
@@ -751,7 +1083,7 @@ void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &c
                                          this->normal_vectors);
   
   this->get_fe().fill_fe_face_values(this->get_mapping(),
-                                    cell, face_no,
+                                    *this->present_cell, face_no,
                                     this->quadrature,
                                     *this->mapping_data,
                                     *this->fe_data,
@@ -841,9 +1173,59 @@ void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator
                       "already refined. Iterate over their children "
                       "instead in these cases."));
   
-  this->present_cell  = cell;
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::template
+     CellIterator<typename DoFHandler<dim>::cell_iterator> (cell));
+
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit (face_no, subface_no);
+}
+
+
 
-  this->get_mapping().fill_fe_subface_values(cell, face_no, subface_no,
+template <int dim>
+void FESubfaceValues<dim>::reinit (const typename Triangulation<dim>::cell_iterator &cell,
+                                  const unsigned int         face_no,
+                                  const unsigned int         subface_no)
+{
+  Assert (face_no < GeometryInfo<dim>::faces_per_cell,
+         ExcIndexRange (face_no, 0, GeometryInfo<dim>::faces_per_cell));
+  Assert (subface_no < GeometryInfo<dim>::subfaces_per_face,
+         ExcIndexRange (subface_no, 0, GeometryInfo<dim>::subfaces_per_face));
+  
+                                   // set new cell. auto_ptr will take
+                                   // care that old object gets
+                                   // destroyed and also that this
+                                   // object gets destroyed in the
+                                   // destruction of this class
+  this->present_cell.reset 
+    (new typename FEValuesBase<dim>::TriaCellIterator (cell));
+
+                                   // this was the part of the work
+                                   // that is dependent on the actual
+                                   // data type of the iterator. now
+                                   // pass on to the function doing
+                                   // the real work.
+  do_reinit (face_no, subface_no);
+}
+
+
+
+template <int dim>
+void FESubfaceValues<dim>::do_reinit (const unsigned int face_no,
+                                      const unsigned int subface_no)
+{
+  this->get_mapping().fill_fe_subface_values(*this->present_cell,
+                                             face_no, subface_no,
                                             this->quadrature,
                                             *this->mapping_data,
                                             this->quadrature_points,
@@ -852,7 +1234,8 @@ void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator
                                             this->normal_vectors);
   
   this->get_fe().fill_fe_subface_values(this->get_mapping(),
-                                       cell, face_no, subface_no,
+                                       *this->present_cell,
+                                        face_no, subface_no,
                                        this->quadrature,
                                        *this->mapping_data,
                                        *this->fe_data,
index 587cfa98ad74a5434a859437d04153c535cdd59e..3c86eba3443ffc86a56c7b862a198af71fce621b 100644 (file)
@@ -127,7 +127,7 @@ MappingCartesian<dim>::get_subface_data (const UpdateFlags update_flags,
 
 template <int dim>
 void
-MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingCartesian<dim>::compute_fill (const typename Triangulation<dim>::cell_iterator &cell,
                                     const unsigned int        face_no,
                                     const unsigned int        sub_no,
                                     InternalData             &data,
@@ -276,7 +276,7 @@ MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterat
 template <int dim>
 void
 MappingCartesian<dim>::
-fill_fe_values (const typename DoFHandler<dim>::cell_iterator& cell,
+fill_fe_values (const typename Triangulation<dim>::cell_iterator& cell,
                 const Quadrature<dim>& q, 
                 typename Mapping<dim>::InternalDataBase& mapping_data,
                 std::vector<Point<dim> >& quadrature_points,
@@ -314,7 +314,7 @@ fill_fe_values (const typename DoFHandler<dim>::cell_iterator& cell,
 
 template <int dim>
 void
-MappingCartesian<dim>::fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingCartesian<dim>::fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                                            const unsigned int            face_no,
                                            const Quadrature<dim-1>      &q,
                                            typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -356,7 +356,7 @@ MappingCartesian<dim>::fill_fe_face_values (const typename DoFHandler<dim>::cell
 
 template <int dim>
 void
-MappingCartesian<dim>::fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingCartesian<dim>::fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                                               const unsigned int       face_no,
                                               const unsigned int       sub_no,
                                               const Quadrature<dim-1> &q,
@@ -400,7 +400,7 @@ MappingCartesian<dim>::fill_fe_subface_values (const typename DoFHandler<dim>::c
 
 template <>
 void
-MappingCartesian<1>::fill_fe_face_values (const DoFHandler<1>::cell_iterator &,
+MappingCartesian<1>::fill_fe_face_values (const Triangulation<1>::cell_iterator &,
                                          const unsigned,
                                          const Quadrature<0>&,
                                          Mapping<1>::InternalDataBase&,
@@ -415,7 +415,7 @@ MappingCartesian<1>::fill_fe_face_values (const DoFHandler<1>::cell_iterator &,
 
 template <>
 void
-MappingCartesian<1>::fill_fe_subface_values (const DoFHandler<1>::cell_iterator &,
+MappingCartesian<1>::fill_fe_subface_values (const Triangulation<1>::cell_iterator &,
                                             const unsigned,
                                             const unsigned,
                                             const Quadrature<0>&,
index 43456f694b78497d1ec699b7f5c1607951a6cd2b..112196bd3c8e73f51e97f77bd6823e576c6fa6fe 100644 (file)
@@ -258,7 +258,7 @@ MappingQ<dim>::get_subface_data (const UpdateFlags update_flags,
 
 template <int dim>
 void
-MappingQ<dim>::fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ<dim>::fill_fe_values (const typename Triangulation<dim>::cell_iterator &cell,
                               const Quadrature<dim>                &q,
                               typename Mapping<dim>::InternalDataBase       &mapping_data,
                               std::vector<Point<dim> >             &quadrature_points,
@@ -295,7 +295,7 @@ MappingQ<dim>::fill_fe_values (const typename DoFHandler<dim>::cell_iterator &ce
 
 template <int dim>
 void
-MappingQ<dim>::fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ<dim>::fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                                    const unsigned int       face_no,
                                    const Quadrature<dim-1> &q,
                                    typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -351,7 +351,7 @@ MappingQ<dim>::fill_fe_face_values (const typename DoFHandler<dim>::cell_iterato
 
 template <int dim>
 void
-MappingQ<dim>::fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ<dim>::fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                                       const unsigned int       face_no,
                                       const unsigned int       sub_no,
                                       const Quadrature<dim-1> &q,
index 22abf10b9f789058451f61871727f4d6d5a82ced..f9ceb4bb6d6491bbeb4e182d3fbd9e29454dfb7d 100644 (file)
@@ -439,7 +439,7 @@ MappingQ1<dim>::get_subface_data (const UpdateFlags update_flags,
 
 template <int dim>
 void
-MappingQ1<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ1<dim>::compute_fill (const typename Triangulation<dim>::cell_iterator &cell,
                              const unsigned int       n_q_points,
                              const DataSetDescriptor  data_set,
                              InternalData            &data,
@@ -542,7 +542,7 @@ MappingQ1<dim>::compute_mapping_support_points(
 
 template <int dim>
 void
-MappingQ1<dim>::fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ1<dim>::fill_fe_values (const typename Triangulation<dim>::cell_iterator &cell,
                                const Quadrature<dim>                &q,
                                typename Mapping<dim>::InternalDataBase      &mapping_data,
                                std::vector<Point<dim> >                  &quadrature_points,
@@ -577,7 +577,7 @@ MappingQ1<dim>::fill_fe_values (const typename DoFHandler<dim>::cell_iterator &c
 
 template <int dim>
 void
-MappingQ1<dim>::compute_fill_face (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ1<dim>::compute_fill_face (const typename Triangulation<dim>::cell_iterator &cell,
                                   const unsigned int      face_no,
                                   const bool              is_subface,
                                   const unsigned int      n_q_points,
@@ -668,7 +668,7 @@ MappingQ1<dim>::compute_fill_face (const typename DoFHandler<dim>::cell_iterator
 
 template <int dim>
 void
-MappingQ1<dim>::fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ1<dim>::fill_fe_face_values (const typename Triangulation<dim>::cell_iterator &cell,
                                     const unsigned int       face_no,
                                     const Quadrature<dim-1> &q,
                                     typename Mapping<dim>::InternalDataBase &mapping_data,
@@ -699,7 +699,7 @@ MappingQ1<dim>::fill_fe_face_values (const typename DoFHandler<dim>::cell_iterat
 
 template <int dim>
 void
-MappingQ1<dim>::fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
+MappingQ1<dim>::fill_fe_subface_values (const typename Triangulation<dim>::cell_iterator &cell,
                                        const unsigned int       face_no,
                                        const unsigned int       sub_no,
                                        const Quadrature<dim-1> &q,
@@ -733,7 +733,7 @@ MappingQ1<dim>::fill_fe_subface_values (const typename DoFHandler<dim>::cell_ite
 
 template <>
 void
-MappingQ1<1>::compute_fill_face (const DoFHandler<1>::cell_iterator &,
+MappingQ1<1>::compute_fill_face (const Triangulation<1>::cell_iterator &,
                                 const unsigned int,
                                 const bool,
                                 const unsigned int,
@@ -751,7 +751,7 @@ MappingQ1<1>::compute_fill_face (const DoFHandler<1>::cell_iterator &,
 
 template <>
 void
-MappingQ1<1>::fill_fe_face_values (const DoFHandler<1>::cell_iterator &,
+MappingQ1<1>::fill_fe_face_values (const Triangulation<1>::cell_iterator &,
                                   const unsigned,
                                   const Quadrature<0>&,
                                   Mapping<1>::InternalDataBase&,
@@ -766,7 +766,7 @@ MappingQ1<1>::fill_fe_face_values (const DoFHandler<1>::cell_iterator &,
 
 template <>
 void
-MappingQ1<1>::fill_fe_subface_values (const DoFHandler<1>::cell_iterator &,
+MappingQ1<1>::fill_fe_subface_values (const Triangulation<1>::cell_iterator &,
                                      const unsigned,
                                      const unsigned,
                                      const Quadrature<0>&,
index 40c83268a68e0457e2e1ed60a35499cd60a2d334..67c5a66373d21a730117efeda6c92a17ebc49ef4 100644 (file)
@@ -23,7 +23,7 @@
 
 template <int dim>
 MappingQ1Eulerian<dim>::MappingQ1Eulerian (const Vector<double>  &euler_transform_vectors,
-                                           const DoFHandler<dim> &shiftmap_dof_handler)
+                                          const DoFHandler<dim> &shiftmap_dof_handler)
                      :
                     euler_transform_vectors(euler_transform_vectors),
                     shiftmap_dof_handler(&shiftmap_dof_handler)

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.