]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Re-apply previous. It now compiles, links and should even work with some luck.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 3 Nov 2010 03:34:27 +0000 (03:34 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 3 Nov 2010 03:34:27 +0000 (03:34 +0000)
git-svn-id: https://svn.dealii.org/trunk@22592 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/fe/mapping.h
deal.II/include/deal.II/fe/mapping_q.h
deal.II/include/deal.II/fe/mapping_q1.h
deal.II/include/deal.II/hp/fe_values.h
deal.II/include/deal.II/numerics/vectors.templates.h
deal.II/source/fe/fe_values.cc
deal.II/source/fe/mapping_q.cc
deal.II/source/fe/mapping_q1.cc
deal.II/source/hp/fe_values.cc
deal.II/source/numerics/vectors.cc

index ea5e92646e51d32d578890d52877c20a0dd082f0..c6e2146b5811ad594f4add4da9252971c7e23cf0 100644 (file)
@@ -332,6 +332,15 @@ through DoFHandler::get_tria() and DoFHandler::get_fe(), respectively.
 <h3>deal.II</h3>
 
 <ol>
+  <li><p>New: The class hp::FEFaceValues and hp::FESubfaceValues were not
+  previously available for meshes that were embedded in a higher
+  dimensional space (i.e. if the codimension was greater than 1). This is
+  now fixed. As a consequence, the VectorTools::interpolate_boundary_values
+  function is now also available for such meshes.
+  <br>
+  (WB, 2010/11/02)
+  </p></li>
+
   <li><p>Fixed: The FEValuesExtractors::Vector class did not work when the dimension
   of the domain was not equal to the dimension of the space in which it is
   embedded. This is now fixed.
index 4cc31d3ecc815132c18234dbfba1f082009b6555..d3dc65cb7473576b292df01ce5afc79ee8080968 100644 (file)
@@ -126,7 +126,7 @@ enum MappingType
  *
  * A hint to implementators: no function except the two functions
  * @p update_once and @p update_each may add any flags.
- * 
+ *
  * For more information about the <tt>spacedim</tt> template parameter
  * check the documentation of FiniteElement or the one of
  * Triangulation.
@@ -138,12 +138,12 @@ template <int dim, int spacedim=dim>
 class Mapping : public Subscriptor
 {
   public:
-    
+
                                     /**
                                      * Virtual destructor.
                                      */
     virtual ~Mapping ();
-    
+
                                     /**
                                      * Transforms the point @p p on
                                      * the unit cell to the point
@@ -154,7 +154,7 @@ class Mapping : public Subscriptor
     transform_unit_to_real_cell (
       const typename Triangulation<dim,spacedim>::cell_iterator &cell,
       const Point<dim>                                 &p) const = 0;
-    
+
                                     /**
                                      * Transforms the point @p p on
                                      * the real cell to the point
@@ -165,7 +165,7 @@ class Mapping : public Subscriptor
     transform_real_to_unit_cell (
       const typename Triangulation<dim,spacedim>::cell_iterator &cell,
       const Point<spacedim>                            &p) const = 0;
-    
+
                                     /**
                                      * Base class for internal data
                                      * of finite element and mapping
@@ -219,7 +219,7 @@ class Mapping : public Subscriptor
                                          * @p first_cell to @p true.
                                          */
         InternalDataBase ();
-       
+
                                         /**
                                          * Virtual destructor for
                                          * derived classes
@@ -231,7 +231,7 @@ class Mapping : public Subscriptor
                                          * by reinit.
                                          */
        UpdateFlags          update_flags;
-       
+
                                         /**
                                          * Values computed by
                                          * constructor.
@@ -286,7 +286,7 @@ class Mapping : public Subscriptor
                                           * the first cell.
                                           */
         virtual void clear_first_cell ();
-        
+
                                         /**
                                          * Return an estimate (in
                                          * bytes) or the memory
@@ -302,14 +302,14 @@ class Mapping : public Subscriptor
                                          * if #update_volume_elements.
                                          */
        std::vector<double> volume_elements;
-       
+
                                         /**
                                          * The positions of the
                                          * mapped (generalized)
                                          * support points.
                                          */
        std::vector<Point<spacedim> > support_point_values;
-       
+
                                         /*
                                          * The Jacobian of the
                                          * transformation in the
@@ -317,7 +317,7 @@ class Mapping : public Subscriptor
                                          * points.
                                          */
        std::vector<Tensor<2,spacedim> > support_point_gradients;
-       
+
                                         /*
                                          * The inverse of the
                                          * Jacobian of the
@@ -326,8 +326,8 @@ class Mapping : public Subscriptor
                                          * points.
                                          */
        std::vector<Tensor<2,spacedim> > support_point_inverse_gradients;
-       
-       
+
+
       private:
                                          /**
                                           * The value returned by
@@ -335,7 +335,7 @@ class Mapping : public Subscriptor
                                           */
         bool first_cell;
     };
-    
+
                                      /**
                                       * Transform a field of
                                       * vectors accorsing to
@@ -360,7 +360,7 @@ class Mapping : public Subscriptor
                VectorSlice<std::vector<Tensor<2,spacedim> > >             output,
                const InternalDataBase &internal,
                const MappingType type) const = 0;
-    
+
                                     /**
                                      * @deprecated Use transform() instead.
                                      */
@@ -391,7 +391,7 @@ class Mapping : public Subscriptor
                                      /**
                                      * @deprecated Use transform() instead.
                                       */
-    
+
     void
     transform_contravariant (const VectorSlice<const std::vector<Tensor<2,dim> > > intput,
                              const unsigned int                 offset,
@@ -405,7 +405,7 @@ class Mapping : public Subscriptor
     const Point<spacedim>& support_point_value(
       const unsigned int index,
       const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
-    
+
                                     /**
                                      * The Jacobian
                                      * matrix of the transformation
@@ -415,7 +415,7 @@ class Mapping : public Subscriptor
     const Tensor<2,spacedim>& support_point_gradient(
       const unsigned int index,
       const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
-    
+
                                     /**
                                      * The inverse Jacobian
                                      * matrix of the transformation
@@ -425,7 +425,7 @@ class Mapping : public Subscriptor
     const Tensor<2,spacedim>& support_point_inverse_gradient(
       const unsigned int index,
       const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
-    
+
                                      /**
                                       * Return a pointer to a copy of the
                                       * present object. The caller of this
@@ -441,7 +441,7 @@ class Mapping : public Subscriptor
                                       */
     virtual
     Mapping<dim,spacedim> * clone () const = 0;
-    
+
                                     /**
                                      * Returns whether the mapping preserves
                                      * vertex locations, i.e. whether the
@@ -466,7 +466,7 @@ class Mapping : public Subscriptor
     DeclException0 (ExcInvalidData);
 
   private:
-    
+
                                     /**
                                      * Indicate fields to be updated
                                      * in the constructor of
@@ -479,7 +479,7 @@ class Mapping : public Subscriptor
                                      * See @ref UpdateFlagsEssay.
                                      */
     virtual UpdateFlags update_once (const UpdateFlags) const = 0;
-    
+
                                     /**
                                      * The same as update_once(),
                                      * but for the flags to be updated for
@@ -507,7 +507,7 @@ class Mapping : public Subscriptor
     virtual InternalDataBase*
     get_face_data (const UpdateFlags flags,
                   const Quadrature<dim-1>& quadrature) const = 0;
-    
+
                                     /**
                                      * Prepare internal data
                                      * structure for transformation
@@ -601,9 +601,9 @@ class Mapping : public Subscriptor
                         const unsigned int                                        face_no,
                         const Quadrature<dim-1>                                   &quadrature,
                         InternalDataBase                                          &internal,
-                        std::vector<Point<dim> >                                  &quadrature_points,
+                        std::vector<Point<spacedim> >                             &quadrature_points,
                         std::vector<double>                                       &JxW_values,
-                        std::vector<Tensor<1,dim> >                               &boundary_form,
+                        std::vector<Tensor<1,spacedim> >                          &boundary_form,
                         std::vector<Point<spacedim> >                             &normal_vectors) const = 0;
 
                                     /**
@@ -615,9 +615,9 @@ class Mapping : public Subscriptor
                            const unsigned int                        sub_no,
                            const Quadrature<dim-1>                  &quadrature,
                            InternalDataBase                         &internal,
-                           std::vector<Point<dim> >        &quadrature_points,
+                           std::vector<Point<spacedim> >        &quadrature_points,
                            std::vector<double>                      &JxW_values,
-                           std::vector<Tensor<1,dim> >     &boundary_form,
+                           std::vector<Tensor<1,spacedim> >     &boundary_form,
                            std::vector<Point<spacedim> >        &normal_vectors) const = 0;
 
                                     /**
@@ -682,7 +682,7 @@ Mapping<dim,spacedim>::support_point_value(
   const typename Mapping<dim,spacedim>::InternalDataBase& internal) const
 {
   AssertIndexRange(index, internal.support_point_values.size());
-  return internal.support_point_values[index];  
+  return internal.support_point_values[index];
 }
 
 
index 4b063526b3cda20402c45f2ed0dd77c62d83eef1..75538439c1f9d976296a7ffc5c6eecb76026df71 100644 (file)
@@ -222,9 +222,9 @@ class MappingQ : public MappingQ1<dim,spacedim>
                         const unsigned int face_no,
                         const Quadrature<dim-1>& quadrature,
                         typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                        typename std::vector<Point<dim> >        &quadrature_points,
+                        typename std::vector<Point<spacedim> >        &quadrature_points,
                         std::vector<double>             &JxW_values,
-                        typename std::vector<Tensor<1,dim> >        &exterior_form,
+                        typename std::vector<Tensor<1,spacedim> >        &exterior_form,
                         typename std::vector<Point<spacedim> >        &normal_vectors) const ;
 
                                     /**
@@ -237,9 +237,9 @@ class MappingQ : public MappingQ1<dim,spacedim>
                            const unsigned int sub_no,
                            const Quadrature<dim-1>& quadrature,
                            typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                           typename std::vector<Point<dim> >        &quadrature_points,
+                           typename std::vector<Point<spacedim> >        &quadrature_points,
                            std::vector<double>             &JxW_values,
-                           typename std::vector<Tensor<1,dim> >        &exterior_form,
+                           typename std::vector<Tensor<1,spacedim> >        &exterior_form,
                            typename std::vector<Point<spacedim> >        &normal_vectors) const ;
 
                                     /**
index 616d535adec69210ecbeecca3cf710d211478faf..33ff0e7e7e6317a51dfa64e1e649b17cd5693e6c 100644 (file)
@@ -312,9 +312,9 @@ class MappingQ1 : public Mapping<dim,spacedim>
                         const unsigned int                               face_no,
                         const Quadrature<dim-1>                          &quadrature,
                         typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                        typename std::vector<Point<dim> >                &quadrature_points,
+                        typename std::vector<Point<spacedim> >                &quadrature_points,
                         std::vector<double>                              &JxW_values,
-                        typename std::vector<Tensor<1,dim> >             &boundary_form,
+                        typename std::vector<Tensor<1,spacedim> >             &boundary_form,
                         typename std::vector<Point<spacedim> >           &normal_vectors) const ;
 
                                     /**
@@ -327,9 +327,9 @@ class MappingQ1 : public Mapping<dim,spacedim>
                            const unsigned int sub_no,
                            const Quadrature<dim-1>& quadrature,
                            typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                           typename std::vector<Point<dim> >        &quadrature_points,
+                           typename std::vector<Point<spacedim> >        &quadrature_points,
                            std::vector<double>             &JxW_values,
-                           typename std::vector<Tensor<1,dim> >        &boundary_form,
+                           typename std::vector<Tensor<1,spacedim> >        &boundary_form,
                            typename std::vector<Point<spacedim> >        &normal_vectors) const ;
 
                                     /**
@@ -399,9 +399,9 @@ class MappingQ1 : public Mapping<dim,spacedim>
                            const DataSetDescriptor data_set,
                            const std::vector<double>   &weights,
                            InternalData           &mapping_data,
-                           std::vector<Point<dim> >    &quadrature_points,
+                           std::vector<Point<spacedim> >    &quadrature_points,
                            std::vector<double>         &JxW_values,
-                           std::vector<Tensor<1,dim> > &boundary_form,
+                           std::vector<Tensor<1,spacedim> > &boundary_form,
                            std::vector<Point<spacedim> > &normal_vectors) const;
 
                                     /**
index f39cbc0874ef0595fb1ff37b5c21f94a28d3548d..308f9946b0d0c090b1ff58fdded7554c730e4263 100644 (file)
@@ -546,8 +546,8 @@ namespace hp
  * @ingroup hp hpcollection
  * @author Wolfgang Bangerth, 2003
  */
-  template <int dim>
-  class FEFaceValues : public internal::hp::FEValuesBase<dim,dim-1,dealii::FEFaceValues<dim> >
+  template <int dim, int spacedim=dim>
+  class FEFaceValues : public internal::hp::FEValuesBase<dim,dim-1,dealii::FEFaceValues<dim,spacedim> >
   {
     public:
                                        /**
@@ -570,8 +570,8 @@ namespace hp
                                         * <tt>DoFHandler::get_fe()</tt>
                                         * function.
                                         */
-      FEFaceValues (const hp::MappingCollection<dim> &mapping_collection,
-                    const hp::FECollection<dim>  &fe_collection,
+      FEFaceValues (const hp::MappingCollection<dim,spacedim> &mapping_collection,
+                    const hp::FECollection<dim,spacedim>  &fe_collection,
                     const hp::QCollection<dim-1>     &q_collection,
                     const UpdateFlags             update_flags);
 
@@ -598,7 +598,7 @@ namespace hp
                                         * <tt>DoFHandler::get_fe()</tt>
                                         * function.
                                         */
-      FEFaceValues (const hp::FECollection<dim>  &fe_collection,
+      FEFaceValues (const hp::FECollection<dim,spacedim>  &fe_collection,
                     const hp::QCollection<dim-1> &q_collection,
                     const UpdateFlags             update_flags);
 
@@ -702,7 +702,7 @@ namespace hp
                                         * this argument is specified.
                                         */
       void
-      reinit (const typename hp::DoFHandler<dim>::cell_iterator &cell,
+      reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
               const unsigned int mapping_index = numbers::invalid_unsigned_int,
@@ -738,7 +738,7 @@ namespace hp
                                         * these last three arguments.
                                         */
       void
-      reinit (const typename dealii::DoFHandler<dim>::cell_iterator &cell,
+      reinit (const typename dealii::DoFHandler<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
               const unsigned int mapping_index = numbers::invalid_unsigned_int,
@@ -774,7 +774,7 @@ namespace hp
                                         * these last three arguments.
                                         */
       void
-      reinit (const typename MGDoFHandler<dim>::cell_iterator &cell,
+      reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
               const unsigned int mapping_index = numbers::invalid_unsigned_int,
@@ -811,7 +811,7 @@ namespace hp
                                         * three arguments.
                                         */
       void
-      reinit (const typename Triangulation<dim>::cell_iterator &cell,
+      reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
               const unsigned int mapping_index = numbers::invalid_unsigned_int,
@@ -827,8 +827,8 @@ namespace hp
  * @ingroup hp hpcollection
  * @author Wolfgang Bangerth, 2003
  */
-  template <int dim>
-  class FESubfaceValues : public internal::hp::FEValuesBase<dim,dim-1,dealii::FESubfaceValues<dim> >
+  template <int dim, int spacedim=dim>
+  class FESubfaceValues : public internal::hp::FEValuesBase<dim,dim-1,dealii::FESubfaceValues<dim,spacedim> >
   {
     public:
                                        /**
@@ -851,8 +851,8 @@ namespace hp
                                         * <tt>DoFHandler::get_fe()</tt>
                                         * function.
                                         */
-      FESubfaceValues (const hp::MappingCollection<dim> &mapping_collection,
-                       const hp::FECollection<dim>  &fe_collection,
+      FESubfaceValues (const hp::MappingCollection<dim,spacedim> &mapping_collection,
+                       const hp::FECollection<dim,spacedim>  &fe_collection,
                        const hp::QCollection<dim-1>     &q_collection,
                        const UpdateFlags             update_flags);
 
@@ -879,7 +879,7 @@ namespace hp
                                         * <tt>DoFHandler::get_fe()</tt>
                                         * function.
                                         */
-      FESubfaceValues (const hp::FECollection<dim> &fe_collection,
+      FESubfaceValues (const hp::FECollection<dim,spacedim> &fe_collection,
                        const hp::QCollection<dim-1>    &q_collection,
                        const UpdateFlags            update_flags);
 
@@ -962,7 +962,7 @@ namespace hp
                                         * this argument is specified.
                                         */
       void
-      reinit (const typename hp::DoFHandler<dim>::cell_iterator &cell,
+      reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int subface_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
@@ -999,7 +999,7 @@ namespace hp
                                         * these last three arguments.
                                         */
       void
-      reinit (const typename dealii::DoFHandler<dim>::cell_iterator &cell,
+      reinit (const typename dealii::DoFHandler<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int subface_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
@@ -1036,7 +1036,7 @@ namespace hp
                                         * these last three arguments.
                                         */
       void
-      reinit (const typename MGDoFHandler<dim>::cell_iterator &cell,
+      reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int subface_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
@@ -1074,7 +1074,7 @@ namespace hp
                                         * three arguments.
                                         */
       void
-      reinit (const typename Triangulation<dim>::cell_iterator &cell,
+      reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
               const unsigned int face_no,
               const unsigned int subface_no,
               const unsigned int q_index = numbers::invalid_unsigned_int,
index e0494237be9a9b360329dc8ebcee6061d51344c6..b27f75adf4f685d47da4523caeb7c70c630eaf01 100644 (file)
@@ -1353,6 +1353,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
                              const std::vector<bool>       &component_mask_)
 {
   const unsigned int dim=DH::dimension;
+  const unsigned int spacedim=DH::space_dimension;
 
   Assert ((component_mask_.size() == 0) ||
          (component_mask_.size() == dof.get_fe().n_components()),
@@ -1373,7 +1374,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
   const unsigned int        n_components = DoFTools::n_components(dof);
   const bool                fe_is_system = (n_components != 1);
 
-  for (typename FunctionMap<DH::space_dimension>::type::const_iterator i=function_map.begin();
+  for (typename FunctionMap<spacedim>::type::const_iterator i=function_map.begin();
        i!=function_map.end(); ++i)
     Assert (n_components == i->second->n_components,
            ExcDimensionMismatch(n_components, i->second->n_components));
@@ -1391,7 +1392,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
   std::vector<unsigned int> face_dofs;
   face_dofs.reserve (DoFTools::max_dofs_per_face(dof));
 
-  std::vector<Point<DH::space_dimension> >  dof_locations;
+  std::vector<Point<spacedim> >  dof_locations;
   dof_locations.reserve (DoFTools::max_dofs_per_face(dof));
 
                                   // array to store the values of
@@ -1411,11 +1412,11 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
                                   // the interpolation points of all
                                   // finite elements that may ever be
                                   // in use
-  hp::FECollection<dim> finite_elements (dof.get_fe());
+  hp::FECollection<dim,spacedim> finite_elements (dof.get_fe());
   hp::QCollection<dim-1>  q_collection;
   for (unsigned int f=0; f<finite_elements.size(); ++f)
     {
-      const FiniteElement<dim> &fe = finite_elements[f];
+      const FiniteElement<dim,spacedim> &fe = finite_elements[f];
 
                                       // generate a quadrature rule
                                       // on the face from the unit
@@ -1483,9 +1484,9 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
                                   // hp::FEFaceValues object that we
                                   // can use to evaluate the boundary
                                   // values at
-  hp::MappingCollection<dim> mapping_collection (mapping);
-  hp::FEFaceValues<dim> x_fe_values (mapping_collection, finite_elements, q_collection,
-                                    update_quadrature_points);
+  hp::MappingCollection<dim,spacedim> mapping_collection (mapping);
+  hp::FEFaceValues<dim,spacedim> x_fe_values (mapping_collection, finite_elements, q_collection,
+                                             update_quadrature_points);
 
   typename DH::active_cell_iterator cell = dof.begin_active(),
                                    endc = dof.end();
@@ -1494,7 +1495,7 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
       for (unsigned int face_no = 0; face_no < GeometryInfo<dim>::faces_per_cell;
           ++face_no)
        {
-         const FiniteElement<dim,DH::space_dimension> &fe = cell->get_fe();
+         const FiniteElement<dim,spacedim> &fe = cell->get_fe();
 
                                           // we can presently deal only with
                                           // primitive elements for boundary
@@ -1526,14 +1527,14 @@ interpolate_boundary_values (const Mapping<DH::dimension, DH::space_dimension>
            {
                                               // face is of the right component
              x_fe_values.reinit(cell, face_no);
-             const FEFaceValues<dim> &fe_values = x_fe_values.get_present_fe_values();
+             const FEFaceValues<dim,spacedim> &fe_values = x_fe_values.get_present_fe_values();
 
                                               // get indices, physical location and
                                               // boundary values of dofs on this
                                               // face
              face_dofs.resize (fe.dofs_per_face);
              face->get_dof_indices (face_dofs, cell->active_fe_index());
-             const std::vector<Point<DH::space_dimension> > &dof_locations
+             const std::vector<Point<spacedim> > &dof_locations
                = fe_values.get_quadrature_points ();
 
              if (fe_is_system)
@@ -1997,7 +1998,7 @@ VectorTools::project_boundary_values (const Mapping<dim, spacedim>       &mappin
     && ! excluded_dofs[dof_to_boundary_mapping[i]])
       {
        Assert(numbers::is_finite(boundary_projection(dof_to_boundary_mapping[i])), ExcNumberNotFinite());
-       
+
                                       // this dof is on one of the
                                       // interesting boundary parts
                                       //
index e89b8a4911bddadf4130feeebbaf8b8b144d0b46..8481733c2f609785e87758162f8ce6eb8dcd6b8e 100644 (file)
@@ -3607,7 +3607,7 @@ FEFaceValues<dim,spacedim>::FEFaceValues (const FiniteElement<dim,spacedim> &fe,
                FEFaceValuesBase<dim,spacedim> (quadrature.size(),
                                                fe.dofs_per_cell,
                                                update_flags,
-                                               StaticMappingQ1<dim>::mapping,
+                                               StaticMappingQ1<dim,spacedim>::mapping,
                                                fe, quadrature)
 {
   Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
@@ -3850,7 +3850,7 @@ FESubfaceValues<dim,spacedim>::FESubfaceValues (const FiniteElement<dim,spacedim
                FEFaceValuesBase<dim,spacedim> (quadrature.size(),
                                                fe.dofs_per_cell,
                                                update_flags,
-                                               StaticMappingQ1<dim>::mapping,
+                                               StaticMappingQ1<dim,spacedim>::mapping,
                                                fe, quadrature)
 {
   Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
@@ -4158,14 +4158,14 @@ template class FEValuesBase<deal_II_dimension,deal_II_dimension+1>;
 template class FEValues<deal_II_dimension,deal_II_dimension+1>;
 template class FEValuesBase<deal_II_dimension,deal_II_dimension+1>::
 CellIterator<DoFHandler<deal_II_dimension,deal_II_dimension+1>::cell_iterator>;
-//template class FEValuesBase<deal_II_dimension,deal_II_dimension+1>::
-//  CellIterator<MGDoFHandler<deal_II_dimension,deal_II_dimension+1>::cell_iterator>;
-
-// #if deal_II_dimension == 2
-// template class FEFaceValuesBase<deal_II_dimension,deal_II_dimension+1>;
-// template class FEFaceValues<deal_II_dimension,deal_II_dimension+1>;
-// template class FESubfaceValues<deal_II_dimension,deal_II_dimension+1>;
-// #endif
+template class FEValuesBase<deal_II_dimension,deal_II_dimension+1>::
+ CellIterator<MGDoFHandler<deal_II_dimension,deal_II_dimension+1>::cell_iterator>;
+
+#if deal_II_dimension == 2
+template class FEFaceValuesBase<deal_II_dimension,deal_II_dimension+1>;
+template class FEFaceValues<deal_II_dimension,deal_II_dimension+1>;
+template class FESubfaceValues<deal_II_dimension,deal_II_dimension+1>;
+#endif
 #endif
 
 
index 82bc61a809945b18b1b8973772103d26f3e535ff..33253296db93c9d1cd07aaeaed012ab299d1a444 100644 (file)
@@ -46,7 +46,7 @@ MappingQ<dim,spacedim>::InternalData::InternalData (const unsigned int n_shape_f
 
 template<int dim, int spacedim>
 unsigned int
-MappingQ<dim,spacedim>::InternalData::memory_consumption () const 
+MappingQ<dim,spacedim>::InternalData::memory_consumption () const
 {
   return (MappingQ1<dim,spacedim>::InternalData::memory_consumption () +
          MemoryConsumption::memory_consumption (unit_normals) +
@@ -143,7 +143,7 @@ MappingQ<dim,spacedim>::MappingQ (const unsigned int p,
   Assert (n_shape_functions==tensor_pols->n(),
          ExcInternalError());
   Assert(n_inner+n_outer==n_shape_functions, ExcInternalError());
-  
+
                                   // build laplace_on_quad_vector
   if (degree>1)
     {
@@ -216,7 +216,7 @@ MappingQ<dim,spacedim>::compute_shapes_virtual (const std::vector<Point<dim> > &
             ExcInternalError());
       grads.resize(n_shape_functions);
     }
-  
+
 //                                // dummy variable of size 0
   std::vector<Tensor<2,dim> > grad2;
   if (data.shape_second_derivatives.size()!=0)
@@ -226,16 +226,16 @@ MappingQ<dim,spacedim>::compute_shapes_virtual (const std::vector<Point<dim> > &
       grad2.resize(n_shape_functions);
     }
 
-  
+
   if (data.shape_values.size()!=0 || data.shape_derivatives.size()!=0)
     for (unsigned int point=0; point<n_points; ++point)
       {
        tensor_pols->compute(unit_points[point], values, grads, grad2);
-       
+
        if (data.shape_values.size()!=0)
          for (unsigned int i=0; i<n_shape_functions; ++i)
            data.shape(point,renumber[i]) = values[i];
-       
+
        if (data.shape_derivatives.size()!=0)
          for (unsigned int i=0; i<n_shape_functions; ++i)
            data.derivative(point,renumber[i]) = grads[i];
@@ -355,7 +355,7 @@ MappingQ<dim,spacedim>::fill_fe_values (
       if (get_degree() > 1)
        cell_similarity = CellSimilarity::invalid_next_cell;
     }
-  
+
   MappingQ1<dim,spacedim>::fill_fe_values(cell, q, *p_data,
                                          quadrature_points, JxW_values,
                                          jacobians, jacobian_grads, inverse_jacobians,
@@ -371,9 +371,9 @@ MappingQ<dim,spacedim>::fill_fe_face_values (
   const unsigned int       face_no,
   const Quadrature<dim-1> &q,
   typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-  std::vector<Point<dim> >     &quadrature_points,
+  std::vector<Point<spacedim> >     &quadrature_points,
   std::vector<double>          &JxW_values,
-  std::vector<Tensor<1,dim> >  &exterior_forms,
+  std::vector<Tensor<1,spacedim> >  &exterior_forms,
   std::vector<Point<spacedim> >     &normal_vectors) const
 {
                                   // convert data object to internal
@@ -383,7 +383,7 @@ MappingQ<dim,spacedim>::fill_fe_face_values (
   Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
          ExcInternalError());
   InternalData &data = static_cast<InternalData&> (mapping_data);
-  
+
                                   // check whether this cell needs
                                   // the full mapping or can be
                                   // treated by a reduced Q1 mapping,
@@ -433,9 +433,9 @@ MappingQ<dim,spacedim>::fill_fe_subface_values (const typename Triangulation<dim
                                       const unsigned int       sub_no,
                                       const Quadrature<dim-1> &q,
                                       typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                                      std::vector<Point<dim> >     &quadrature_points,
+                                      std::vector<Point<spacedim> >     &quadrature_points,
                                       std::vector<double>          &JxW_values,
-                                      std::vector<Tensor<1,dim> >  &exterior_forms,
+                                      std::vector<Tensor<1,spacedim> >  &exterior_forms,
                                       std::vector<Point<spacedim> >     &normal_vectors) const
 {
                                   // convert data object to internal
@@ -518,7 +518,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
                                        // for degree==1, we shouldn't have to
                                        // compute any support points, since
                                        // all of them are on the vertices
-      
+
       case 2:
       {
                                          // (checked these values against the
@@ -533,7 +533,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
         Assert (sizeof(loqv2)/sizeof(loqv2[0]) ==
                 n_inner_2d * n_outer_2d,
                 ExcInternalError());
-       
+
        break;
       }
 
@@ -542,7 +542,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
                                          // (same as above)
        static const double loqv3[4*12]
          ={80/1053., 1/81., 1/81., 11/1053., 25/117., 44/351.,
-           7/117., 16/351., 25/117., 44/351., 7/117., 16/351., 
+           7/117., 16/351., 25/117., 44/351., 7/117., 16/351.,
            1/81., 80/1053., 11/1053., 1/81., 7/117., 16/351.,
            25/117., 44/351., 44/351., 25/117., 16/351., 7/117.,
            1/81., 11/1053., 80/1053., 1/81., 44/351., 25/117.,
@@ -552,9 +552,9 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
         Assert (sizeof(loqv3)/sizeof(loqv3[0]) ==
                 n_inner_2d * n_outer_2d,
                 ExcInternalError());
-        
+
        loqv_ptr=&loqv3[0];
-       
+
        break;
       }
 
@@ -570,7 +570,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.2231273865431891, 0.1346851306015187,
            0.03812914216116723, 0.02913160002633253,
            0.02200737428129391, 0.01600835564431222,
-           
+
            0.00664803151334206, 0.006648031513342719,
            0.002873452861657458, 0.002873452861657626,
            0.07903572682584378, 0.05969238281250031,
@@ -588,7 +588,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.03812914216116729, 0.1346851306015185,
            0.2231273865431898, 0.01600835564431217,
            0.02200737428129394, 0.02913160002633262,
-           
+
            0.006648031513342073, 0.002873452861657473,
            0.006648031513342726, 0.002873452861657636,
            0.1527716818820238, 0.2348152760709273,
@@ -597,7 +597,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.07903572682584376, 0.05969238281250026,
            0.03619864817415824, 0.07903572682584187,
            0.0596923828124998, 0.0361986481741581,
-           
+
            0.01106770833333302, 0.01106770833333336,
            0.01106770833333337, 0.01106770833333374,
            0.06770833333333424, 0.1035156250000011,
@@ -606,7 +606,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.06770833333333422, 0.1035156250000009,
            0.06770833333333436, 0.0677083333333337,
            0.1035156249999988, 0.0677083333333339,
-           
+
            0.002873452861657185, 0.006648031513342362,
            0.002873452861657334, 0.006648031513343038,
            0.02496269311797779, 0.04081948955407401,
@@ -615,7 +615,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.03619864817415819, 0.05969238281250028,
            0.07903572682584407, 0.03619864817415804,
            0.05969238281249986, 0.0790357268258422,
-           
+
            -0.001075744628906913, 0.00191429223907134,
            0.07405921850311592, -0.001075744628905865,
            0.03812914216116729, 0.1346851306015185,
@@ -633,7 +633,7 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.02496269311797776, 0.04081948955407392,
            0.02496269311797785, 0.1527716818820233,
            0.2348152760709258, 0.1527716818820236,
-           
+
            0.001914292239071237, -0.001075744628906803,
            -0.001075744628906778, 0.07405921850311617,
            0.01600835564431228, 0.02200737428129401,
@@ -641,24 +641,24 @@ MappingQ<dim,spacedim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
            0.1346851306015182, 0.2231273865431886,
            0.01600835564431228, 0.02200737428129397,
            0.02913160002633523, 0.03812914216116726,
-           0.1346851306015181, 0.2231273865431886,    
+           0.1346851306015181, 0.2231273865431886,
           };
-        
+
         Assert (sizeof(loqv4)/sizeof(loqv4[0]) ==
                 n_inner_2d * n_outer_2d,
                 ExcInternalError());
-        
+
        loqv_ptr=&loqv4[0];
-       
+
        break;
       }
-      
+
                                       // no other cases implemented,
                                       // so simply fall through
       default:
             break;
     }
-  
+
   if (loqv_ptr!=0)
     {
                                       // precomputed. copy values to
@@ -726,10 +726,10 @@ MappingQ<3>::set_laplace_on_hex_vector(Table<2,double> &lohvs) const
          7/192., 7/192., 7/192., 7/192., 7/192., 7/192., 7/192., 7/192.,
          7/192., 7/192., 7/192., 7/192.,
          1/12., 1/12., 1/12., 1/12., 1/12., 1/12.};
-      
+
       lohv_ptr=&loqv2[0];
     }
-  
+
   if (lohv_ptr!=0)
     {
                                       // precomputed. copy values to
@@ -742,7 +742,7 @@ MappingQ<3>::set_laplace_on_hex_vector(Table<2,double> &lohvs) const
   else
                                     // not precomputed, then do so now
     compute_laplace_vector(lohvs);
-    
+
                                   // the sum of weights of the points
                                   // at the outer rim should be
                                   // one. check this
@@ -794,11 +794,11 @@ MappingQ<dim,spacedim>::compute_laplace_vector(Table<2,double> &lvs) const
                                   // points on the unit cell
   const QGauss<dim> quadrature(degree+1);
   const unsigned int n_q_points=quadrature.size();
-  
+
   InternalData quadrature_data(n_shape_functions);
   quadrature_data.shape_derivatives.resize(n_shape_functions * n_q_points);
   this->compute_shapes(quadrature.get_points(), quadrature_data);
-  
+
                                   // Compute the stiffness matrix of
                                   // the inner dofs
   FullMatrix<double> S(n_inner);
@@ -808,7 +808,7 @@ MappingQ<dim,spacedim>::compute_laplace_vector(Table<2,double> &lvs) const
        S(i,j) += contract(quadrature_data.derivative(point, n_outer+i),
                           quadrature_data.derivative(point, n_outer+j))
                  * quadrature.weight(point);
-  
+
                                   // Compute the components of T to be the
                                   // product of gradients of inner and
                                   // outer shape functions.
@@ -819,15 +819,15 @@ MappingQ<dim,spacedim>::compute_laplace_vector(Table<2,double> &lvs) const
        T(i,k) += contract(quadrature_data.derivative(point, n_outer+i),
                           quadrature_data.derivative(point, k))
                  *quadrature.weight(point);
-  
+
   FullMatrix<double> S_1(n_inner);
   S_1.invert(S);
-  
+
   FullMatrix<double> S_1_T(n_inner, n_outer);
-  
+
                                   // S:=S_1*T
   S_1.mmult(S_1_T,T);
-  
+
                                   // Resize and initialize the
                                   // lvs
   lvs.reinit (n_inner, n_outer);
@@ -857,7 +857,7 @@ MappingQ<dim,spacedim>::apply_laplace_vector(const Table<2,double> &lvs,
                                    // reason for not aborting there
                                    // any more [WB]
   Assert(lvs.n_rows()!=0, ExcLaplaceVectorNotSet(degree));
-  
+
   const unsigned int n_inner_apply=lvs.n_rows();
   Assert(n_inner_apply==n_inner || n_inner_apply==(degree-1)*(degree-1),
         ExcInternalError());
@@ -903,13 +903,13 @@ MappingQ<dim,spacedim>::compute_mapping_support_points(
                                     // cell
     {
       a.resize(GeometryInfo<dim>::vertices_per_cell);
-      
+
       for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
        a[i] = cell->vertex(i);
     }
 }
 
-  
+
 template<int dim, int spacedim>
 void
 MappingQ<dim,spacedim>::compute_support_points_laplace(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
@@ -920,7 +920,7 @@ MappingQ<dim,spacedim>::compute_support_points_laplace(const typename Triangulat
   a.resize(GeometryInfo<dim>::vertices_per_cell);
   for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
     a[i] = cell->vertex(i);
-  
+
   if (degree>1)
     switch (dim)
       {
@@ -996,10 +996,10 @@ MappingQ<1,2>::add_line_support_points (const Triangulation<1,2>::cell_iterator
                                     // boundary description
     {
       std::vector<Point<spacedim> > line_points (degree-1);
-      
+
       const Boundary<dim,spacedim> * const boundary
        = &(cell->get_triangulation().get_boundary(cell->material_id()));
-      
+
       boundary->get_intermediate_points_on_line (cell, line_points);
                                       // Append all points
       a.insert (a.end(), line_points.begin(), line_points.end());
@@ -1028,7 +1028,7 @@ MappingQ<dim,spacedim>::add_line_support_points (const typename Triangulation<di
               (dim != spacedim) ?
               &line->get_triangulation().get_boundary(cell->material_id()):
               &straight_boundary);
-         
+
          a.push_back(boundary->get_new_point_on_line(line));
        };
     }
@@ -1039,7 +1039,7 @@ MappingQ<dim,spacedim>::add_line_support_points (const typename Triangulation<di
                                     // boundary description
     {
       std::vector<Point<spacedim> > line_points (degree-1);
-      
+
                                       // loop over each of the lines,
                                       // and if it is at the
                                       // boundary, then first get the
@@ -1049,14 +1049,14 @@ MappingQ<dim,spacedim>::add_line_support_points (const typename Triangulation<di
       for (unsigned int line_no=0; line_no<GeometryInfo<dim>::lines_per_cell; ++line_no)
        {
          const typename Triangulation<dim,spacedim>::line_iterator line = cell->line(line_no);
-         
+
          const Boundary<dim,spacedim> * const boundary
            = (line->at_boundary()?
               &line->get_triangulation().get_boundary(line->boundary_indicator()) :
               (dim != spacedim) ?
               &line->get_triangulation().get_boundary(cell->material_id()) :
               &straight_boundary);
-         
+
          boundary->get_intermediate_points_on_line (line, line_points);
          if (dim==3)
            {
@@ -1073,7 +1073,7 @@ MappingQ<dim,spacedim>::add_line_support_points (const typename Triangulation<di
                                             // correct orientation. simply
                                             // append all points
            a.insert (a.end(), line_points.begin(), line_points.end());
-         
+
        }
     }
 }
@@ -1103,8 +1103,8 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
                                   // used if only one line of face
                                   // quad is at boundary
   std::vector<Point<3> > b(4*degree);
-  
-  
+
+
                                   // loop over all faces and collect
                                   // points on them
   for (unsigned int face_no=0; face_no<faces_per_cell; ++face_no)
@@ -1117,7 +1117,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
                 face_flip        = cell->face_flip       (face_no),
                 face_rotation    = cell->face_rotation   (face_no);
 
-#ifdef DEBUG      
+#ifdef DEBUG
                                        // some sanity checks up front
       for (unsigned int i=0; i<vertices_per_face; ++i)
         Assert(face->vertex_index(i)==cell->vertex_index(
@@ -1136,7 +1136,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
          face_no, i, face_orientation, face_flip, face_rotation)),
               ExcInternalError());
 #endif
-      
+
                                       // if face at boundary, then
                                       // ask boundary object to
                                       // return intermediate points
@@ -1168,7 +1168,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
          for (unsigned int i=0; i<lines_per_face; ++i)
            if (face->line(i)->at_boundary())
              ++lines_at_boundary;
-         
+
          Assert(lines_at_boundary<=lines_per_face, ExcInternalError());
 
                                           // if at least one of the
@@ -1187,7 +1187,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
                                               // function was already
                                               // called.
              b.resize(4*degree);
-             
+
                                               // b is of size
                                               // 4*degree, make sure
                                               // that this is the
@@ -1195,7 +1195,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
              Assert(b.size()==vertices_per_face+lines_per_face*(degree-1),
                     ExcDimensionMismatch(b.size(),
                                           vertices_per_face+lines_per_face*(degree-1)));
-             
+
                                               // sort the points into b. We
                                               // used access from the cell (not
                                               // from the face) to fill b, so
@@ -1205,7 +1205,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
                                               // standard orientation as well.
               for (unsigned int i=0; i<vertices_per_face; ++i)
                 b[i]=a[GeometryInfo<3>::face_to_cell_vertices(face_no, i)];
-                     
+
               for (unsigned int i=0; i<lines_per_face; ++i)
                 for (unsigned int j=0; j<degree-1; ++j)
                   b[vertices_per_face+i*(degree-1)+j]=
@@ -1218,7 +1218,7 @@ add_quad_support_points(const Triangulation<3>::cell_iterator &cell,
              apply_laplace_vector(laplace_on_quad_vector, b);
              Assert(b.size()==4*degree+(degree-1)*(degree-1),
                     ExcDimensionMismatch(b.size(), 4*degree+(degree-1)*(degree-1)));
-             
+
              for (unsigned int i=0; i<(degree-1)*(degree-1); ++i)
                a.push_back(b[4*degree+i]);
            }
@@ -1261,7 +1261,7 @@ add_quad_support_points(const Triangulation<2,3>::cell_iterator &cell,
                         std::vector<Point<3> >                &a) const
 {
   std::vector<Point<3> > quad_points ((degree-1)*(degree-1));
-  
+
   cell->get_triangulation().get_boundary(cell->material_id())
     .get_intermediate_points_on_quad (cell, quad_points);
   for (unsigned int i=0; i<quad_points.size(); ++i)
@@ -1302,7 +1302,7 @@ MappingQ<dim,spacedim>::transform (
                                   // to test further
   if (!q1_data->is_mapping_q1_data)
     {
-      Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0, 
+      Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0,
              ExcInternalError());
       const InternalData &data = static_cast<const InternalData&>(mapping_data);
                                       // If we only use the
@@ -1315,7 +1315,7 @@ MappingQ<dim,spacedim>::transform (
                                   // right tensors in it and we call
                                   // the base classes transform
                                   // function
-  MappingQ1<dim,spacedim>::transform(input, output, *q1_data, mapping_type);  
+  MappingQ1<dim,spacedim>::transform(input, output, *q1_data, mapping_type);
 }
 
 
@@ -1340,7 +1340,7 @@ MappingQ<dim,spacedim>::transform (
                                   // to test further
   if (!q1_data->is_mapping_q1_data)
     {
-      Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0, 
+      Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0,
              ExcInternalError());
       const InternalData &data = static_cast<const InternalData&>(mapping_data);
                                       // If we only use the
@@ -1353,7 +1353,7 @@ MappingQ<dim,spacedim>::transform (
                                   // right tensors in it and we call
                                   // the base classes transform
                                   // function
-  MappingQ1<dim,spacedim>::transform(input, output, *q1_data, mapping_type);  
+  MappingQ1<dim,spacedim>::transform(input, output, *q1_data, mapping_type);
 }
 
 
@@ -1373,7 +1373,7 @@ transform_unit_to_real_cell (const typename Triangulation<dim,spacedim>::cell_it
   std::auto_ptr<InternalData>
     mdata (dynamic_cast<InternalData *> (
              get_data(update_transformation_values, point_quadrature)));
-  
+
   mdata->use_mapping_q1_on_current_cell = !(use_mapping_q_on_all_cells ||
                                            cell->has_boundary_lines());
 
@@ -1383,7 +1383,7 @@ transform_unit_to_real_cell (const typename Triangulation<dim,spacedim>::cell_it
                &*mdata);
 
   compute_mapping_support_points(cell, p_data->mapping_support_points);
-  
+
   return this->transform_unit_to_real_cell_internal(*p_data);
 }
 
@@ -1398,7 +1398,7 @@ transform_real_to_unit_cell (const typename Triangulation<dim,spacedim>::cell_it
                                   // first a Newton iteration based
                                   // on a Q1 mapping
   Point<dim> p_unit = MappingQ1<dim,spacedim>::transform_real_to_unit_cell(cell, p);
-  
+
                                    // then a Newton iteration based on
                                    // the full MappingQ if we need
                                    // this
@@ -1412,7 +1412,7 @@ transform_real_to_unit_cell (const typename Triangulation<dim,spacedim>::cell_it
                  get_data(update_transformation_values |
                           update_transformation_gradients,
                           point_quadrature)));
-      
+
       mdata->use_mapping_q1_on_current_cell = false;
 
       std::vector<Point<spacedim> > &points = mdata->mapping_support_points;
@@ -1420,7 +1420,7 @@ transform_real_to_unit_cell (const typename Triangulation<dim,spacedim>::cell_it
 
       this->transform_real_to_unit_cell_internal(cell, p, *mdata, p_unit);
     }
-  
+
   return p_unit;
 }
 
@@ -1442,7 +1442,7 @@ MappingQ<dim,spacedim>::clone () const
   return new MappingQ<dim,spacedim>(*this);
 }
 
-  
+
 // explicit instantiation
 template class MappingQ<deal_II_dimension>;
 
index 75f67acef67c91536bc3a4f34d17a633c23d5716..ee72c8f37b29a9ae96510efa6609d168f7084d34 100644 (file)
@@ -93,13 +93,13 @@ namespace internal
     void
     compute_shapes_virtual (const unsigned int            n_shape_functions,
                            const std::vector<Point<1> > &unit_points,
-                           typename dealii::MappingQ1<1,spacedim>::InternalData& data) 
+                           typename dealii::MappingQ1<1,spacedim>::InternalData& data)
     {
       const unsigned int n_points=unit_points.size();
       for (unsigned int k = 0 ; k < n_points ; ++k)
        {
          double x = unit_points[k](0);
-      
+
          if (data.shape_values.size()!=0)
            {
              Assert(data.shape_values.size()==n_shape_functions*n_points,
@@ -119,7 +119,7 @@ namespace internal
                                               // the following may or may not
                                               // work if dim != spacedim
              Assert (spacedim == 1, ExcNotImplemented());
-             
+
              Assert(data.shape_second_derivatives.size()==n_shape_functions*n_points,
                     ExcInternalError());
              data.second_derivative(k,0)[0][0] = 0;
@@ -133,14 +133,14 @@ namespace internal
     void
     compute_shapes_virtual (const unsigned int            n_shape_functions,
                            const std::vector<Point<2> > &unit_points,
-                           typename dealii::MappingQ1<2,spacedim>::InternalData& data) 
+                           typename dealii::MappingQ1<2,spacedim>::InternalData& data)
     {
       const unsigned int n_points=unit_points.size();
       for (unsigned int k = 0 ; k < n_points ; ++k)
        {
          double x = unit_points[k](0);
          double y = unit_points[k](1);
-      
+
          if (data.shape_values.size()!=0)
            {
              Assert(data.shape_values.size()==n_shape_functions*n_points,
@@ -168,7 +168,7 @@ namespace internal
                                               // the following may or may not
                                               // work if dim != spacedim
              Assert (spacedim == 2, ExcNotImplemented());
-             
+
              Assert(data.shape_second_derivatives.size()==n_shape_functions*n_points,
                     ExcInternalError());
              data.second_derivative(k,0)[0][0] = 0;
@@ -197,7 +197,7 @@ namespace internal
     void
     compute_shapes_virtual (const unsigned int            n_shape_functions,
                            const std::vector<Point<3> > &unit_points,
-                           typename dealii::MappingQ1<3,spacedim>::InternalData& data) 
+                           typename dealii::MappingQ1<3,spacedim>::InternalData& data)
     {
       const unsigned int n_points=unit_points.size();
       for (unsigned int k = 0 ; k < n_points ; ++k)
@@ -205,7 +205,7 @@ namespace internal
          double x = unit_points[k](0);
          double y = unit_points[k](1);
          double z = unit_points[k](2);
-      
+
          if (data.shape_values.size()!=0)
            {
              Assert(data.shape_values.size()==n_shape_functions*n_points,
@@ -253,7 +253,7 @@ namespace internal
                                               // the following may or may not
                                               // work if dim != spacedim
              Assert (spacedim == 3, ExcNotImplemented());
-             
+
              Assert(data.shape_second_derivatives.size()==n_shape_functions*n_points,
                     ExcInternalError());
              data.second_derivative(k,0)[0][0] = 0;
@@ -314,7 +314,7 @@ namespace internal
              data.second_derivative(k,5)[2][0] = (1.-y);
              data.second_derivative(k,6)[2][0] = -y;
              data.second_derivative(k,7)[2][0] = y;
-         
+
              data.second_derivative(k,0)[1][2] = (1.-x);
              data.second_derivative(k,1)[1][2] = x;
              data.second_derivative(k,2)[1][2] = -(1.-x);
@@ -418,7 +418,7 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
       if (out & (update_JxW_values
                 | update_normal_vectors))
        out |= update_boundary_forms;
-  
+
       if (out & (update_covariant_transformation
                 | update_JxW_values
                 | update_jacobians
@@ -426,7 +426,7 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
                 | update_boundary_forms
                 | update_normal_vectors))
        out |= update_contravariant_transformation;
-      
+
       if (out & (update_inverse_jacobians))
        out |= update_covariant_transformation;
 
@@ -443,7 +443,7 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
        out |= update_JxW_values;
 
     }
-  
+
   return out;
 }
 
@@ -462,7 +462,7 @@ MappingQ1<dim,spacedim>::compute_data (const UpdateFlags      update_flags,
   data.update_flags = data.update_once | data.update_each;
 
   const UpdateFlags flags(data.update_flags);
-  
+
   if (flags & update_transformation_values)
     data.shape_values.resize(data.n_shape_functions * n_q_points);
 
@@ -480,7 +480,7 @@ MappingQ1<dim,spacedim>::compute_data (const UpdateFlags      update_flags,
 
   if (flags & update_jacobian_grads)
     data.shape_second_derivatives.resize(data.n_shape_functions * n_q_points);
-  
+
   compute_shapes (q.get_points(), data);
 }
 
@@ -512,7 +512,7 @@ MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
       if (data.update_flags & update_boundary_forms)
        {
          data.aux.resize (dim-1, std::vector<Tensor<1,spacedim> > (n_original_q_points));
-      
+
                                           // Compute tangentials to the
                                           // unit cell.
          const unsigned int nfaces = GeometryInfo<dim>::faces_per_cell;
@@ -525,7 +525,7 @@ MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
              static const int tangential_orientation[4]={-1,1,1,-1};
              for (unsigned int i=0; i<nfaces; ++i)
                {
-                 Tensor<1,spacedim> tang;            
+                 Tensor<1,spacedim> tang;
                  tang[1-i/2]=tangential_orientation[i];
                  std::fill (data.unit_tangentials[i].begin(),
                             data.unit_tangentials[i].end(), tang);
@@ -539,7 +539,7 @@ MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
 
                  const unsigned int nd=
                    GeometryInfo<dim>::unit_normal_direction[i];
-             
+
                                                   // first tangential
                                                   // vector in direction
                                                   // of the (nd+1)%3 axis
@@ -564,7 +564,7 @@ MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
     }
 }
 
-  
+
 
 template<int dim, int spacedim>
 typename Mapping<dim,spacedim>::InternalDataBase *
@@ -607,7 +607,6 @@ MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim
                                       InternalData  &data,
                                       std::vector<Point<spacedim> > &quadrature_points) const
 {
-
   const UpdateFlags update_flags(data.current_update_flags());
 
                                   // if necessary, recompute the
@@ -640,7 +639,7 @@ MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim
 
       for (unsigned int point=0; point<n_q_points; ++point)
        for (unsigned int k=0; k<data.n_shape_functions; ++k)
-         quadrature_points[point] 
+         quadrature_points[point]
            += data.shape(point+data_set,k) * data.mapping_support_points[k];
     }
 
@@ -690,7 +689,7 @@ MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim
 
          else if (dim == spacedim - 1)
            {
-                                   // CODIMENSION 1 
+                                   // CODIMENSION 1
                                    // calculate left-inversion of the
                                    // contravariant matrices to obtain
                                    // covariant ones (auxiliary
@@ -701,7 +700,7 @@ MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim
              for (unsigned int point=0; point<n_q_points; ++point)
                {
                  contravariant_matrix.
-                   copy_from(data.contravariant[point],0,spacedim-1,0,dim-1);  
+                   copy_from(data.contravariant[point],0,spacedim-1,0,dim-1);
                  covariant_matrix.
                    left_invert(contravariant_matrix);
                  covariant_matrix.
@@ -745,7 +744,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
   CellSimilarity::Similarity                           &cell_similarity) const
 {
   // ensure that the following cast is really correct:
-  Assert (dynamic_cast<InternalData *>(&mapping_data) != 0, 
+  Assert (dynamic_cast<InternalData *>(&mapping_data) != 0,
          ExcInternalError());
   InternalData &data = static_cast<InternalData&>(mapping_data);
 
@@ -754,7 +753,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
   compute_fill (cell, n_q_points, DataSetDescriptor::cell (), cell_similarity,
                 data, quadrature_points);
 
-  
+
   const UpdateFlags update_flags(data.current_update_flags());
   const std::vector<double> &weights=q.get_weights();
 
@@ -766,7 +765,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
                                   // the case <2,3>
   if (update_flags & (update_normal_vectors
                      | update_JxW_values))
-    {      
+    {
       Assert (JxW_values.size() == n_q_points,
               ExcDimensionMismatch(JxW_values.size(), n_q_points));
 
@@ -776,7 +775,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
 
       if (cell_similarity != CellSimilarity::translation)
        for (unsigned int point=0; point<n_q_points; ++point)
-         { 
+         {
 
            if (dim==spacedim)
              JxW_values[point]
@@ -787,7 +786,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
                data.contravariant[point]=transpose(data.contravariant[point]);
                JxW_values[point]
                  = data.contravariant[point][0].norm()*weights[point];
-               if(update_flags & update_normal_vectors) 
+               if(update_flags & update_normal_vectors)
                {
                  normal_vectors[point][0]=-data.contravariant[point][0][1]
                      /data.contravariant[point][0].norm();
@@ -798,7 +797,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
            if ( (dim==2) && (spacedim==3) )
              {
                data.contravariant[point]=transpose(data.contravariant[point]);
-               cross_product(data.contravariant[point][2], 
+               cross_product(data.contravariant[point][2],
                              data.contravariant[point][0],
                              data.contravariant[point][1]
                              );
@@ -809,17 +808,17 @@ MappingQ1<dim,spacedim>::fill_fe_values (
                                           //is stored in the 3d
                                          //subtensor of the contravariant tensor
                data.contravariant[point][2]/=data.contravariant[point][2].norm();
-               if(update_flags & update_normal_vectors) 
+               if(update_flags & update_normal_vectors)
                  normal_vectors[point]=data.contravariant[point][2];
              }
-           
+
          }
     }
 
                                   // copy values from InternalData to vector
                                   // given by reference
   if (update_flags & update_jacobians)
-    {      
+    {
       Assert (jacobians.size() == n_q_points,
              ExcDimensionMismatch(jacobians.size(), n_q_points));
       if (cell_similarity != CellSimilarity::translation)
@@ -854,7 +853,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
                                   // copy values from InternalData to vector
                                   // given by reference
   if (update_flags & update_inverse_jacobians)
-    {      
+    {
       Assert (inverse_jacobians.size() == n_q_points,
              ExcDimensionMismatch(inverse_jacobians.size(), n_q_points));
       if (cell_similarity != CellSimilarity::translation)
@@ -875,10 +874,10 @@ MappingQ1<2,3>::compute_fill_face (const Triangulation<2,3>::cell_iterator &,
                                   const DataSetDescriptor,
                                   const std::vector<double>&,
                                   InternalData &,
-                                  std::vector<Point<2> >&,
+                                  std::vector<Point<3> >&,
                                   std::vector<double>&,
-                                  std::vector<Tensor<1,2> > &,
-                                  std::vector<Point<3> > &) const 
+                                  std::vector<Tensor<1,3> > &,
+                                  std::vector<Point<3> > &) const
 {
        Assert(false, ExcNotImplemented());
 }
@@ -897,9 +896,9 @@ MappingQ1<dim,spacedim>::compute_fill_face (
   const DataSetDescriptor,
   const std::vector<double> &,
   InternalData &,
-  std::vector<Point<dim> > &,
+  std::vector<Point<spacedim> > &,
   std::vector<double> &,
-  std::vector<Tensor<1,dim> > &,
+  std::vector<Tensor<1,spacedim> > &,
   std::vector<Point<spacedim> > &) const
 {
   Assert(false, ExcNotImplemented());
@@ -913,9 +912,9 @@ MappingQ1<dim,spacedim>::fill_fe_face_values (
   const unsigned,
   const Quadrature<dim-1>&,
   typename Mapping<dim,spacedim>::InternalDataBase&,
-  std::vector<Point<dim> >&,
+  std::vector<Point<spacedim> >&,
   std::vector<double>&,
-  std::vector<Tensor<1,dim> >&,
+  std::vector<Tensor<1,spacedim> >&,
   std::vector<Point<spacedim> >&) const
 {
   Assert(false, ExcNotImplemented());
@@ -930,9 +929,9 @@ MappingQ1<dim,spacedim>::fill_fe_subface_values (
   const unsigned,
   const Quadrature<dim-1>&,
   typename Mapping<dim,spacedim>::InternalDataBase&,
-  std::vector<Point<dim> >&,
+  std::vector<Point<spacedim> >&,
   std::vector<double>&,
-  std::vector<Tensor<1,dim> >&,
+  std::vector<Tensor<1,spacedim> >&,
   std::vector<Point<spacedim> >&) const
 {
   Assert(false, ExcNotImplemented());
@@ -952,16 +951,16 @@ MappingQ1<dim,spacedim>::compute_fill_face (
   const DataSetDescriptor          data_set,
   const std::vector<double>        &weights,
   InternalData                     &data,
-  std::vector<Point<dim> >         &quadrature_points,
+  std::vector<Point<spacedim> >         &quadrature_points,
   std::vector<double>              &JxW_values,
-  std::vector<Tensor<1,dim> > &boundary_forms,
+  std::vector<Tensor<1,spacedim> > &boundary_forms,
   std::vector<Point<spacedim> >    &normal_vectors) const
 {
   compute_fill (cell, n_q_points, data_set, CellSimilarity::none,
                data, quadrature_points);
 
   const UpdateFlags update_flags(data.current_update_flags());
-  
+
   if (update_flags & update_boundary_forms)
     {
       Assert (boundary_forms.size()==n_q_points,
@@ -971,16 +970,16 @@ MappingQ1<dim,spacedim>::compute_fill_face (
                ExcDimensionMismatch(normal_vectors.size(), n_q_points));
       if (update_flags & update_JxW_values)
        Assert (JxW_values.size() == n_q_points,
-               ExcDimensionMismatch(JxW_values.size(), n_q_points));      
-      
+               ExcDimensionMismatch(JxW_values.size(), n_q_points));
+
       transform(data.unit_tangentials[face_no], data.aux[0],
                data, mapping_contravariant);
 
-      typename std::vector<Tensor<1,dim> >::iterator
+      typename std::vector<Tensor<1,spacedim> >::iterator
        result = boundary_forms.begin();
-      const typename std::vector<Tensor<1,dim> >::iterator
+      const typename std::vector<Tensor<1,spacedim> >::iterator
        end = boundary_forms.end();
-      
+
       switch (dim)
        {
          case 2:
@@ -1009,7 +1008,7 @@ MappingQ1<dim,spacedim>::compute_fill_face (
          default:
                 Assert(false, ExcNotImplemented());
        }
-      
+
       if (update_flags & (update_normal_vectors
                          | update_JxW_values))
        for (unsigned int i=0;i<boundary_forms.size();++i)
@@ -1039,18 +1038,18 @@ MappingQ1<dim,spacedim>::fill_fe_face_values (const typename Triangulation<dim,s
                                     const unsigned int                               face_no,
                                     const Quadrature<dim-1>                          &q,
                                     typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                                    std::vector<Point<dim> >                         &quadrature_points,
+                                    std::vector<Point<spacedim> >                         &quadrature_points,
                                     std::vector<double>                              &JxW_values,
-                                    std::vector<Tensor<1,dim> >                      &boundary_forms,
+                                    std::vector<Tensor<1,spacedim> >                      &boundary_forms,
                                     std::vector<Point<spacedim> >                    &normal_vectors) const
 {
   // ensure that the following cast is really correct:
-  Assert (dynamic_cast<InternalData *>(&mapping_data) != 0, 
+  Assert (dynamic_cast<InternalData *>(&mapping_data) != 0,
          ExcInternalError());
   InternalData &data = static_cast<InternalData&>(mapping_data);
 
   const unsigned int n_q_points = q.size();
-  
+
   compute_fill_face (cell, face_no, deal_II_numbers::invalid_unsigned_int,
                     n_q_points,
                     DataSetDescriptor::face (face_no,
@@ -1074,18 +1073,18 @@ MappingQ1<dim,spacedim>::fill_fe_subface_values (const typename Triangulation<di
                                        const unsigned int       sub_no,
                                        const Quadrature<dim-1> &q,
                                        typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                                       std::vector<Point<dim> >     &quadrature_points,
+                                       std::vector<Point<spacedim> >     &quadrature_points,
                                        std::vector<double>          &JxW_values,
-                                       std::vector<Tensor<1,dim> >  &boundary_forms,
+                                       std::vector<Tensor<1,spacedim> >  &boundary_forms,
                                        std::vector<Point<spacedim> >     &normal_vectors) const
 {
   // ensure that the following cast is really correct:
-  Assert (dynamic_cast<InternalData *>(&mapping_data) != 0, 
+  Assert (dynamic_cast<InternalData *>(&mapping_data) != 0,
          ExcInternalError());
   InternalData &data = static_cast<InternalData&>(mapping_data);
 
   const unsigned int n_q_points = q.size();
-  
+
   compute_fill_face (cell, face_no, sub_no,
                     n_q_points,
                     DataSetDescriptor::subface (face_no, sub_no,
@@ -1115,19 +1114,19 @@ MappingQ1<dim,spacedim>::transform (
   const MappingType mapping_type) const
 {
   AssertDimension (input.size(), output.size());
-  Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0, 
+  Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0,
          ExcInternalError());
   const InternalData &data = static_cast<const InternalData&>(mapping_data);
 
-  Tensor<1, spacedim> auxiliary;  
-  
+  Tensor<1, spacedim> auxiliary;
+
   switch (mapping_type)
     {
       case mapping_covariant:
       {
        Assert (data.update_flags & update_covariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d=0;d<dim;++d)
@@ -1136,12 +1135,12 @@ MappingQ1<dim,spacedim>::transform (
          }
        return;
       }
-      
+
       case mapping_contravariant:
       {
        Assert (data.update_flags & update_contravariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d=0;d<dim;++d)
@@ -1150,14 +1149,14 @@ MappingQ1<dim,spacedim>::transform (
          }
        return;
       }
-      
+
       case mapping_piola:
       {
        Assert (data.update_flags & update_contravariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
        Assert (data.update_flags & update_volume_elements,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d=0;d<dim;++d)
@@ -1166,7 +1165,7 @@ MappingQ1<dim,spacedim>::transform (
          }
        return;
       }
-      
+
       default:
            Assert(false, ExcNotImplemented());
     }
@@ -1182,20 +1181,20 @@ MappingQ1<dim,spacedim>::transform (
   const MappingType mapping_type) const
 {
   AssertDimension (input.size(), output.size());
-  Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0, 
+  Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0,
          ExcInternalError());
   const InternalData &data = static_cast<const InternalData&>(mapping_data);
 
   Tensor<2, spacedim> aux1;
   Tensor<2, spacedim> aux2;
-  
+
   switch (mapping_type)
     {
       case mapping_covariant:
       {
        Assert (data.update_flags & update_covariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d1=0;d1<dim;++d1)
@@ -1205,12 +1204,12 @@ MappingQ1<dim,spacedim>::transform (
          }
        return;
       }
-      
+
       case mapping_contravariant:
       {
        Assert (data.update_flags & update_contravariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d1=0;d1<dim;++d1)
@@ -1220,45 +1219,45 @@ MappingQ1<dim,spacedim>::transform (
          }
        return;
       }
-      
+
       case mapping_covariant_gradient:
       {
        Assert (data.update_flags & update_contravariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d1=0;d1<dim;++d1)
              for (unsigned int d2=0;d2<dim;++d2)
                aux1[d1][d2] = input[i][d1][d2];
-           
+
            contract(aux2, aux1, data.covariant[i]);
            contract(output[i], data.covariant[i], aux2);
          }
-       
+
        return;
       }
-      
+
       case mapping_contravariant_gradient:
       {
        Assert (data.update_flags & update_covariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
        Assert (data.update_flags & update_contravariant_transformation,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d1=0;d1<dim;++d1)
              for (unsigned int d2=0;d2<dim;++d2)
                aux1[d1][d2] = input[i][d1][d2];
-           
+
            contract(aux2, aux1, data.covariant[i]);
            contract(output[i], data.contravariant[i], aux2);
          }
-       
+
        return;
       }
-      
+
       case mapping_piola:
       {
        Assert (data.update_flags & update_covariant_transformation,
@@ -1267,13 +1266,13 @@ MappingQ1<dim,spacedim>::transform (
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
        Assert (data.update_flags & update_volume_elements,
                typename FEValuesBase<dim>::ExcAccessToUninitializedField());
-       
+
        for (unsigned int i=0; i<output.size(); ++i)
          {
            for (unsigned int d1=0;d1<dim;++d1)
              for (unsigned int d2=0;d2<dim;++d2)
                aux1[d1][d2] = input[i][d1][d2] / data.volume_elements[i];
-           
+
            contract (aux2, aux1, data.covariant[i]);
            contract (output[i], data.contravariant[i], aux2);
          }
@@ -1319,7 +1318,7 @@ transform_unit_to_real_cell_internal (const InternalData &data) const
 {
   const unsigned int n_mapping_points=data.mapping_support_points.size();
   Assert(data.shape_values.size()==n_mapping_points, ExcInternalError());
-  
+
                                   // use now the InternalData to
                                   // compute the point in real space.
   Point<spacedim> p_real;
@@ -1365,7 +1364,7 @@ transform_real_to_unit_cell (const typename Triangulation<dim,spacedim>::cell_it
 
                                   // perform the newton iteration.
   transform_real_to_unit_cell_internal(cell, p, *mdata, p_unit);
-  
+
   return p_unit;
 }
 
@@ -1415,22 +1414,22 @@ transform_real_to_unit_cell_internal
   const unsigned int n_shapes=mdata.shape_values.size();
   Assert(n_shapes!=0, ExcInternalError());
   Assert(mdata.shape_derivatives.size()==n_shapes, ExcInternalError());
-  
+
   std::vector<Point<spacedim> > &points=mdata.mapping_support_points;
   Assert(points.size()==n_shapes, ExcInternalError());
-  
+
                                   // Newton iteration to solve
                                   // f(x)=p(x)-p=0
                                   // x_{n+1}=x_n-[f'(x)]^{-1}f(x)
-  
+
                                   // The start value is set to be the
                                   // center of the unit cell.
-  
+
                                   // The shape values and derivatives
                                   // of the mapping at this point are
                                   // previously computed.
 
-  
+
                                   // f(x)
   Point<spacedim> p_real(transform_unit_to_real_cell_internal(mdata));
   Point<spacedim> f = p_real-p;
@@ -1438,19 +1437,19 @@ transform_real_to_unit_cell_internal
   const double eps=1e-15*cell->diameter();
   unsigned int loop=0;
   while (f.square()>eps*eps && loop++<10)
-    {      
+    {
                                       // f'(x)
       Tensor<2,dim> df;
       for (unsigned int k=0; k<mdata.n_shape_functions; ++k)
        {
          const Tensor<1,dim> &grad_transform=mdata.derivative(0,k);
          const Point<spacedim> &point=points[k];
-         
+
          for (unsigned int i=0; i<dim; ++i)
            for (unsigned int j=0; j<dim; ++j)
              df[i][j]+=point[i]*grad_transform[j];
        }
-      
+
                                       // Solve  [f'(x)]d=f(x)
       Tensor<1,dim> d;
       Tensor<2,dim> df_1;
@@ -1463,7 +1462,7 @@ transform_real_to_unit_cell_internal
                                       // shape values and derivatives
                                       // at new p_unit point
       compute_shapes(std::vector<Point<dim> > (1, p_unit), mdata);
-      
+
                                       // f(x)
       p_real = transform_unit_to_real_cell_internal(mdata);
       f = p_real-p;
index c4b67767e9a77ff16a6ebd5b85bc52eabfc636eb..372d4b97c0fd6fd28465c474baf4c70a923fc606 100644 (file)
@@ -309,33 +309,33 @@ namespace hp
 // -------------------------- FEFaceValues -------------------------
 
 
-  template <int dim>
-  FEFaceValues<dim>::FEFaceValues (const hp::MappingCollection<dim> &mapping,
-                                   const hp::FECollection<dim>  &fe_collection,
+  template <int dim, int spacedim>
+  FEFaceValues<dim,spacedim>::FEFaceValues (const hp::MappingCollection<dim,spacedim> &mapping,
+                                   const hp::FECollection<dim,spacedim>  &fe_collection,
                                    const hp::QCollection<dim-1> &q_collection,
                                    const UpdateFlags         update_flags)
                   :
-                  internal::hp::FEValuesBase<dim,dim-1,dealii::FEFaceValues<dim> > (mapping,
+                  internal::hp::FEValuesBase<dim,dim-1,dealii::FEFaceValues<dim,spacedim> > (mapping,
                                                                              fe_collection,
                                                                              q_collection,
                                                                              update_flags)
   {}
 
 
-  template <int dim>
-  FEFaceValues<dim>::FEFaceValues (const hp::FECollection<dim>  &fe_collection,
+  template <int dim, int spacedim>
+  FEFaceValues<dim,spacedim>::FEFaceValues (const hp::FECollection<dim,spacedim>  &fe_collection,
                                    const hp::QCollection<dim-1> &q_collection,
                                    const UpdateFlags         update_flags)
                   :
-                  internal::hp::FEValuesBase<dim,dim-1,dealii::FEFaceValues<dim> > (fe_collection,
+                  internal::hp::FEValuesBase<dim,dim-1,dealii::FEFaceValues<dim,spacedim> > (fe_collection,
                                                                              q_collection,
                                                                              update_flags)
   {}
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FEFaceValues<dim>::reinit (const typename hp::DoFHandler<dim>::cell_iterator &cell,
+  FEFaceValues<dim,spacedim>::reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell,
                              const unsigned int face_no,
                              const unsigned int q_index,
                              const unsigned int mapping_index,
@@ -384,9 +384,9 @@ namespace hp
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FEFaceValues<dim>::reinit (const typename dealii::DoFHandler<dim>::cell_iterator &cell,
+  FEFaceValues<dim,spacedim>::reinit (const typename dealii::DoFHandler<dim,spacedim>::cell_iterator &cell,
                              const unsigned int face_no,
                              const unsigned int q_index,
                              const unsigned int mapping_index,
@@ -425,9 +425,9 @@ namespace hp
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FEFaceValues<dim>::reinit (const typename MGDoFHandler<dim>::cell_iterator &cell,
+  FEFaceValues<dim,spacedim>::reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell,
                              const unsigned int face_no,
                              const unsigned int q_index,
                              const unsigned int mapping_index,
@@ -466,9 +466,9 @@ namespace hp
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FEFaceValues<dim>::reinit (const typename Triangulation<dim>::cell_iterator &cell,
+  FEFaceValues<dim,spacedim>::reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                              const unsigned int face_no,
                              const unsigned int q_index,
                              const unsigned int mapping_index,
@@ -509,33 +509,33 @@ namespace hp
 // -------------------------- FESubfaceValues -------------------------
 
 
-  template <int dim>
-  FESubfaceValues<dim>::FESubfaceValues (const hp::MappingCollection<dim> &mapping,
-                                         const hp::FECollection<dim>  &fe_collection,
+  template <int dim, int spacedim>
+  FESubfaceValues<dim,spacedim>::FESubfaceValues (const hp::MappingCollection<dim,spacedim> &mapping,
+                                         const hp::FECollection<dim,spacedim>  &fe_collection,
                                          const hp::QCollection<dim-1> &q_collection,
                                          const UpdateFlags         update_flags)
                   :
-                  internal::hp::FEValuesBase<dim,dim-1,dealii::FESubfaceValues<dim> > (mapping,
+                  internal::hp::FEValuesBase<dim,dim-1,dealii::FESubfaceValues<dim,spacedim> > (mapping,
                                                                                 fe_collection,
                                                                                 q_collection,
                                                                                 update_flags)
   {}
 
 
-  template <int dim>
-  FESubfaceValues<dim>::FESubfaceValues (const hp::FECollection<dim>  &fe_collection,
+  template <int dim, int spacedim>
+  FESubfaceValues<dim,spacedim>::FESubfaceValues (const hp::FECollection<dim,spacedim>  &fe_collection,
                                          const hp::QCollection<dim-1> &q_collection,
                                          const UpdateFlags         update_flags)
                   :
-                  internal::hp::FEValuesBase<dim,dim-1,dealii::FESubfaceValues<dim> > (fe_collection,
+                  internal::hp::FEValuesBase<dim,dim-1,dealii::FESubfaceValues<dim,spacedim> > (fe_collection,
                                                                                 q_collection,
                                                                                 update_flags)
   {}
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FESubfaceValues<dim>::reinit (const typename hp::DoFHandler<dim>::cell_iterator &cell,
+  FESubfaceValues<dim,spacedim>::reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell,
                                 const unsigned int face_no,
                                 const unsigned int subface_no,
                                 const unsigned int q_index,
@@ -585,9 +585,9 @@ namespace hp
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FESubfaceValues<dim>::reinit (const typename dealii::DoFHandler<dim>::cell_iterator &cell,
+  FESubfaceValues<dim,spacedim>::reinit (const typename dealii::DoFHandler<dim,spacedim>::cell_iterator &cell,
                                 const unsigned int face_no,
                                 const unsigned int subface_no,
                                 const unsigned int q_index,
@@ -627,9 +627,9 @@ namespace hp
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FESubfaceValues<dim>::reinit (const typename MGDoFHandler<dim>::cell_iterator &cell,
+  FESubfaceValues<dim,spacedim>::reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell,
                                 const unsigned int face_no,
                                 const unsigned int subface_no,
                                 const unsigned int q_index,
@@ -669,9 +669,9 @@ namespace hp
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  FESubfaceValues<dim>::reinit (const typename Triangulation<dim>::cell_iterator &cell,
+  FESubfaceValues<dim,spacedim>::reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                                 const unsigned int face_no,
                                 const unsigned int subface_no,
                                 const unsigned int q_index,
@@ -746,13 +746,12 @@ namespace internal
   {
     template class FEValuesBase<deal_II_dimension,deal_II_dimension,
                                 dealii::FEValues<deal_II_dimension,deal_II_dimension+1> >;
-// not yet implemented:
-// #if deal_II_dimension == 2
-//     template class FEValuesBase<deal_II_dimension,deal_II_dimension-1,
-//                                 dealii::FEFaceValues<deal_II_dimension> >;
-//     template class FEValuesBase<deal_II_dimension,deal_II_dimension-1,
-//                                 dealii::FESubfaceValues<deal_II_dimension> >;
-// #endif
+#if deal_II_dimension == 2
+    template class FEValuesBase<deal_II_dimension,deal_II_dimension-1,
+                                dealii::FEFaceValues<deal_II_dimension,deal_II_dimension+1> >;
+    template class FEValuesBase<deal_II_dimension,deal_II_dimension-1,
+                                dealii::FESubfaceValues<deal_II_dimension,deal_II_dimension+1> >;
+#endif
   }
 }
 
@@ -760,11 +759,10 @@ namespace hp
 {
   template class FEValues<deal_II_dimension, deal_II_dimension+1>;
 
-// not yet implemented:
-//#if deal_II_dimension == 2
-//   template class FEFaceValues<deal_II_dimension, deal_II_dimension+1>;
-//   template class FESubfaceValues<deal_II_dimension, deal_II_dimension+1>;
-//#endif
+#if deal_II_dimension == 2
+  template class FEFaceValues<deal_II_dimension, deal_II_dimension+1>;
+  template class FESubfaceValues<deal_II_dimension, deal_II_dimension+1>;
+#endif
 }
 #endif
 
index c0bf0b59f7de5c7fbe59a24dfe02da4cd4c91f53..c30a9f4c51c223d5b76d8590c133443fa4496031 100644 (file)
@@ -206,6 +206,101 @@ void VectorTools::interpolate_boundary_values (
   ConstraintMatrix                    &,
   const std::vector<bool>    &);
 
+#if deal_II_dimension < 3
+template
+void VectorTools::interpolate_boundary_values (
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  std::map<unsigned int,double>       &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  std::map<unsigned int,double>       &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  std::map<unsigned int,double>       &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const unsigned char,
+  const Function<deal_II_dimension+1>   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+  const FunctionMap<deal_II_dimension+1>::type   &,
+  ConstraintMatrix                    &,
+  const std::vector<bool>    &);
+#endif
+
 template
 void VectorTools::project_boundary_values<deal_II_dimension>
 (const Mapping<deal_II_dimension>     &,
@@ -266,22 +361,6 @@ VectorTools::compute_no_normal_flux_constraints (const DoFHandler<deal_II_dimens
 
 
 
-// // Due to introducing the DoFHandler as a template parameter,
-// // the following instantiations are required in 1d
-// #if deal_II_dimension == 1
-// template
-// void VectorTools::interpolate_boundary_values<deal_II_dimension>
-// (const Mapping<1>         &,
-//  const DoFHandler<1>      &,
-//  const unsigned char,
-//  const Function<1>        &,
-//  std::map<unsigned int,double> &,
-//  const std::vector<bool>       &);
-// #endif
-
-// the following two functions are not derived from a template in 1d
-// and thus need no explicit instantiation
-//#if deal_II_dimension > 1
 template
 void VectorTools::interpolate_boundary_values
 (const Mapping<deal_II_dimension>    &,
@@ -299,7 +378,24 @@ void VectorTools::interpolate_boundary_values
  std::map<unsigned int,double>       &,
  const std::vector<bool>    &);
 
-//#endif
+#if deal_II_dimension < 3
+template
+void VectorTools::interpolate_boundary_values
+(const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+ const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+ const FunctionMap<deal_II_dimension+1>::type &,
+ std::map<unsigned int,double>       &,
+ const std::vector<bool>    &);
+
+template
+void VectorTools::interpolate_boundary_values
+(const Mapping<deal_II_dimension,deal_II_dimension+1>    &,
+ const DoFHandler<deal_II_dimension,deal_II_dimension+1> &,
+ const unsigned char,
+ const Function<deal_II_dimension+1>   &,
+ std::map<unsigned int,double>       &,
+ const std::vector<bool>    &);
+#endif
 
 
 template

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.