]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Reimplement the changes from patch 14330: Do not permute the shape functions on faces...
authorleicht <leicht@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 18 Jan 2007 10:43:46 +0000 (10:43 +0000)
committerleicht <leicht@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 18 Jan 2007 10:43:46 +0000 (10:43 +0000)
git-svn-id: https://svn.dealii.org/trunk@14339 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_accessor.templates.h
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_q.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/source/dofs/dof_accessor.cc
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_system.cc
deal.II/deal.II/source/fe/fe_values.cc
deal.II/doc/news/changes.html

index 4a732ec9149093d028b6a10f351df1d66141f9fb..124e2f9065ff46d756df08f31f44bc3510a2c994 100644 (file)
@@ -589,8 +589,13 @@ DoFObjectAccessor<3,DH>::get_dof_indices (std::vector<unsigned int> &dof_indices
     for (unsigned int d=0; d<dofs_per_line; ++d)
       *next++ = this->line(line)->dof_index(d,fe_index);
   for (unsigned int quad=0; quad<6; ++quad)
-    for (unsigned int d=0; d<dofs_per_quad; ++d)
-      *next++ = this->quad(quad)->dof_index(d,fe_index);
+    if (this->face_orientation(quad))
+      for (unsigned int d=0; d<dofs_per_quad; ++d)
+       *next++ = this->quad(quad)->dof_index(d,fe_index);
+    else
+      for (unsigned int d=0; d<dofs_per_quad; ++d)
+       *next++ = this->quad(quad)->dof_index(this->dof_handler->get_fe()[fe_index].
+                                             adjust_quad_dof_index_for_face_orientation(d),fe_index);
   for (unsigned int d=0; d<dofs_per_hex; ++d)
     *next++ = this->dof_index(d,fe_index);
 }
index 76a768dbb36da634173cfe43a804e292c782fa47..25a3752841abb86ffc2b89340b3cb432e4f0b2b9 100644 (file)
@@ -1195,6 +1195,20 @@ class FiniteElement : public Subscriptor,
     std::pair<unsigned int, unsigned int>
     face_system_to_component_index (const unsigned int index) const;
 
+                                    /**
+                                     * For faces with non-standard
+                                     * face_orientation in 3D, the dofs on
+                                     * faces (quads) have to be permuted in
+                                     * order to be combined with the correct
+                                     * shape functions. Given a local dof @p
+                                     * index on a quad, return the local index,
+                                     * if the face has non-standard
+                                     * face_orientation. In 2D and 1D there is
+                                     * no need for permutation so the identity
+                                     * is returned.
+                                     */
+    unsigned int adjust_quad_dof_index_for_face_orientation (const unsigned int index) const;
+
                                     /**
                                      * Return in which of the vector
                                      * components of this finite
@@ -1459,31 +1473,6 @@ class FiniteElement : public Subscriptor,
     unsigned int
     component_to_block_index (const unsigned int component) const;
 
-                                    /**
-                                     * For faces with non-standard
-                                     * face_orientation in 3D, the shape
-                                     * functions on faces have to be permutat
-                                     * in order to be combined with the correct
-                                     * dofs. This function returns a vector of
-                                     * integer values, that have to be added to
-                                     * the index of a shape function in order
-                                     * to result in the permuted index. Prior
-                                     * content of the vector @p shifts is
-                                     * erased. In 3D a vector of length @p
-                                     * dofs_per_quad is returned, in 2D and 1D
-                                     * there is no need for permutation and a
-                                     * vector of length 0 is returned, the same
-                                     * is true for elements which have no dofs
-                                     * on quads. The general implementation
-                                     * returns a vector of zeros, resulting in
-                                     * no permutation at all. This has to be
-                                     * overloaded by derived finite element
-                                     * classes.
-                                     */
-    virtual
-    void
-    get_face_shape_function_shifts (std::vector<int> &shifts) const;
-
                                     //@}
     
                                     /**
@@ -2078,6 +2067,29 @@ class FiniteElement : public Subscriptor,
                                      */    
     std::vector<Point<dim-1> > generalized_face_support_points;
     
+                                    /**
+                                     * For faces with non-standard
+                                     * face_orientation in 3D, the dofs on
+                                     * faces (quads) have to be permuted in
+                                     * order to be combined with the correct
+                                     * shape functions. Given a local dof @p
+                                     * index on a quad, return the shift in the
+                                     * local index, if the face has
+                                     * non-standard face_orientation,
+                                     * i.e. <code>old_index + shift =
+                                     * new_index</code>. In 2D and 1D there is
+                                     * no need for permutation so the vector is
+                                     * empty. In 3D it has the size of @p
+                                     * dofs_per_quad.
+                                     *
+                                     * The standard implementation fills this
+                                     * with zeros, i.e. no permuatation at
+                                     * all. Derived finite element classes have
+                                     * to fill this vector with the correct
+                                     * values.
+                                     */
+    std::vector<int> adjust_quad_dof_index_for_face_orientation_table;
+
   private:
                                     /**
                                      * Store what
index cd57fcf2b6e63be58c1855a6d85779cc6e484042..075439ce2461d0b09773c37e9fcb8147a555ba31 100644 (file)
@@ -438,26 +438,6 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
     virtual
     FiniteElementDomination::Domination
     compare_for_face_domination (const FiniteElement<dim> &fe_other) const;
-
-                                    /**
-                                     * For faces with non-standard
-                                     * face_orientation in 3D, the shape
-                                     * functions on faces have to be permuted
-                                     * in order to be combined with the correct
-                                     * dofs. This function returns a vector of
-                                     * integer values, that have to be added to
-                                     * the index of a shape function in order
-                                     * to result in the permuted index. Prior
-                                     * content of the vector @p shifts is
-                                     * erased. In 3D a vector of length @p
-                                     * dofs_per_quad is returned, in 2D and 1D
-                                     * there is no need for permutation and a
-                                     * vector of length 0 is returned.
-                                     */
-    virtual
-    void
-    get_face_shape_function_shifts (std::vector<int> &shifts) const;
-
                                     //@}
 
                                     /**
@@ -546,6 +526,15 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * constructor.
                                      */
     void initialize_unit_face_support_points ();
+
+                                    /**
+                                     * Initialize the
+                                     * @p adjust_quad_dof_index_for_face_orientation_table field
+                                     * of the FiniteElement
+                                     * class. Called from the
+                                     * constructor.
+                                     */
+    void initialize_quad_dof_index_permutation ();
              
                                     /**
                                      * Mapping from hierarchic to
index ff0eb1d84bcc43231f6731c8d65ed11bdd7177cb..901cda21947dc7a5513c38df0504343923e53fd0 100644 (file)
@@ -539,25 +539,6 @@ class FESystem : public FiniteElement<dim>
     virtual
     FiniteElementDomination::Domination
     compare_for_face_domination (const FiniteElement<dim> &fe_other) const;
-
-                                    /**
-                                     * For faces with non-standard
-                                     * face_orientation in 3D, the shape
-                                     * functions on faces have to be permuted
-                                     * in order to be combined with the correct
-                                     * dofs. This function returns a vector of
-                                     * integer values, that have to be added to
-                                     * the index of a shape function in order
-                                     * to result in the permuted index. Prior
-                                     * content of the vector @p shifts is
-                                     * erased. In 3D a vector of length @p
-                                     * dofs_per_quad is returned, in 2D and 1D
-                                     * there is no need for permutation and a
-                                     * vector of length 0 is returned.
-                                     */
-    virtual
-    void
-    get_face_shape_function_shifts (std::vector<int> &shifts) const;
                                     //@}
     
                                     /**
@@ -743,6 +724,15 @@ class FESystem : public FiniteElement<dim>
                                      */
     void initialize_unit_face_support_points ();
 
+                                    /**
+                                     * Initialize the
+                                     * @p adjust_quad_dof_index_for_face_orientation_table field
+                                     * of the FiniteElement
+                                     * class. Called from the
+                                     * constructor.
+                                     */
+    void initialize_quad_dof_index_permutation ();
+    
                                     /**
                                      * Helper function used in the constructor:
                                      * take a @p FiniteElementData object
index a0fcde80468a2571587a099bd468441900139356..addb01b923b6affd47edbc3d1e703924ad4ce4bc 100644 (file)
@@ -263,38 +263,7 @@ class FEValuesData
                                      * non-zero components.
                                      */
     std::vector<unsigned int> shape_function_to_row_table;
-
-                                    /**
-                                     * Vector containing the permutation of
-                                     * shape functions necessary if the faces
-                                     * of a cell have the wrong
-                                     * face_orientation. This is computed only
-                                     * once. Actually, this does not contain
-                                     * the permutation itself but rather the
-                                     * shift of indices needed to calculate the
-                                     * permutation.
-                                     */
-    std::vector<int> shift_in_face_shape_functions;
-
-                                    /**
-                                     * Vector containing the permutation of
-                                     * shape functions due to faces with
-                                     * non-standard face_orientation on a given
-                                     * cell (recomputed on each cell.) If all
-                                     * faces are oriented according to the
-                                     * standard, this is the identity mapping.
-                                     */
-    std::vector<unsigned int> permuted_shape_functions;
-
-                                    /**
-                                     * Bool flag indicating the need to update
-                                     * the @p permuted_shape_functions vector
-                                     * on each cell. This is only necessary in
-                                     * 3d and if the finite element has
-                                     * shape_functions on the face.
-                                     */
-    bool update_shape_function_permutation;
-        
+    
                                      /**
                                      * Original update flags handed
                                      * to the constructor of
@@ -657,19 +626,7 @@ class FEValuesBase : protected FEValuesData<dim>,
                                    const unsigned int point_no,
                                    const unsigned int component) const;
     
-                                    /**
-                                     * If shape functions belong to a face in
-                                     * 3D, they have to be permuted, if the
-                                     * face has non-standard face
-                                     * orientation. This functuion takes an
-                                     * index of a shape function (on a standard
-                                     * cell) and returns the corresponding
-                                     * shape function on the real cell.
-                                     */
-    unsigned int
-    shift_shape_function_index (const unsigned int i) const;
-    
-    
+
                                     //@}
                                     /// @name FunctionAccess Access to values of global finite element functions
                                     //@{
@@ -1666,12 +1623,6 @@ class FEValuesBase : protected FEValuesData<dim>,
                                      */
     UpdateFlags compute_update_flags (const UpdateFlags update_flags) const;
 
-                                    /**
-                                     * Reinit the permutation of (face) shape
-                                     * functions to match the present cell.
-                                     */
-    void reinit();
-
   private:
                                      /**
                                       * Copy constructor. Since
@@ -2371,17 +2322,15 @@ double
 FEValuesBase<dim>::shape_value (const unsigned int i,
                                const unsigned int j) const
 {
-  const unsigned int I=shift_shape_function_index(i);
-  
   Assert (this->update_flags & update_values,
          ExcAccessToUninitializedField());
-  Assert (fe->is_primitive (I),
-         ExcShapeFunctionNotPrimitive(I));
+  Assert (fe->is_primitive (i),
+         ExcShapeFunctionNotPrimitive(i));
 
                                   // if the entire FE is primitive,
                                   // then we can take a short-cut:
   if (fe->is_primitive())
-    return this->shape_values(I,j);
+    return this->shape_values(i,j);
   else
                                     // otherwise, use the mapping
                                     // between shape function numbers
@@ -2392,7 +2341,7 @@ FEValuesBase<dim>::shape_value (const unsigned int i,
                                     // question to which vector
                                     // component the call of this
                                     // function refers
-    return this->shape_values(this->shape_function_to_row_table[I], j);
+    return this->shape_values(this->shape_function_to_row_table[i], j);
 }
 
 
@@ -2404,8 +2353,6 @@ FEValuesBase<dim>::shape_value_component (const unsigned int i,
                                          const unsigned int j,
                                          const unsigned int component) const
 {
-  const unsigned int I=shift_shape_function_index(i);
-
   Assert (this->update_flags & update_values,
          ExcAccessToUninitializedField());
   Assert (component < fe->n_components(),
@@ -2420,10 +2367,10 @@ FEValuesBase<dim>::shape_value_component (const unsigned int i,
                                   // system_to_component_table only
                                   // works if the shape function is
                                   // primitive):
-  if (fe->is_primitive(I))
+  if (fe->is_primitive(i))
     {
-      if (component == fe->system_to_component_index(I).first)
-       return this->shape_values(this->shape_function_to_row_table[I],j);
+      if (component == fe->system_to_component_index(i).first)
+       return this->shape_values(this->shape_function_to_row_table[i],j);
       else
        return 0;
     }
@@ -2438,7 +2385,7 @@ FEValuesBase<dim>::shape_value_component (const unsigned int i,
                                       // whether the shape function
                                       // is non-zero at all within
                                       // this component:
-      if (fe->get_nonzero_components(I)[component] == false)
+      if (fe->get_nonzero_components(i)[component] == false)
        return 0.;
 
                                       // count how many non-zero
@@ -2450,10 +2397,10 @@ FEValuesBase<dim>::shape_value_component (const unsigned int i,
                                       // shape function in the arrays
                                       // we index presently:
       const unsigned int
-       row = (this->shape_function_to_row_table[I]
+       row = (this->shape_function_to_row_table[i]
               +
-              std::count (fe->get_nonzero_components(I).begin(),
-                          fe->get_nonzero_components(I).begin()+component,
+              std::count (fe->get_nonzero_components(i).begin(),
+                          fe->get_nonzero_components(i).begin()+component,
                           true));
       return this->shape_values(row, j);
     };
@@ -2467,21 +2414,19 @@ const Tensor<1,dim> &
 FEValuesBase<dim>::shape_grad (const unsigned int i,
                               const unsigned int j) const
 {
-  const unsigned int I=shift_shape_function_index(i);
-
   Assert (this->update_flags & update_gradients,
          ExcAccessToUninitializedField());
-  Assert (fe->is_primitive (I),
-         ExcShapeFunctionNotPrimitive(I));
+  Assert (fe->is_primitive (i),
+         ExcShapeFunctionNotPrimitive(i));
   Assert (i<this->shape_gradients.size(),
-         ExcIndexRange (I, 0, this->shape_gradients.size()));
+         ExcIndexRange (i, 0, this->shape_gradients.size()));
   Assert (j<this->shape_gradients[0].size(),
          ExcIndexRange (j, 0, this->shape_gradients[0].size()));
 
                                   // if the entire FE is primitive,
                                   // then we can take a short-cut:
   if (fe->is_primitive())
-    return this->shape_gradients[I][j];
+    return this->shape_gradients[i][j];
   else
                                     // otherwise, use the mapping
                                     // between shape function numbers
@@ -2492,7 +2437,7 @@ FEValuesBase<dim>::shape_grad (const unsigned int i,
                                     // question to which vector
                                     // component the call of this
                                     // function refers
-    return this->shape_gradients[this->shape_function_to_row_table[I]][j];
+    return this->shape_gradients[this->shape_function_to_row_table[i]][j];
 }
 
 
@@ -2504,8 +2449,6 @@ FEValuesBase<dim>::shape_grad_component (const unsigned int i,
                                         const unsigned int j,
                                         const unsigned int component) const
 {
-  const unsigned int I=shift_shape_function_index(i);
-
   Assert (this->update_flags & update_gradients,
          ExcAccessToUninitializedField());
   Assert (component < fe->n_components(),
@@ -2520,10 +2463,10 @@ FEValuesBase<dim>::shape_grad_component (const unsigned int i,
                                   // system_to_component_table only
                                   // works if the shape function is
                                   // primitive):
-  if (fe->is_primitive(I))
+  if (fe->is_primitive(i))
     {
-      if (component == fe->system_to_component_index(I).first)
-       return this->shape_gradients[this->shape_function_to_row_table[I]][j];
+      if (component == fe->system_to_component_index(i).first)
+       return this->shape_gradients[this->shape_function_to_row_table[i]][j];
       else
        return Tensor<1,dim>();
     }
@@ -2538,7 +2481,7 @@ FEValuesBase<dim>::shape_grad_component (const unsigned int i,
                                       // whether the shape function
                                       // is non-zero at all within
                                       // this component:
-      if (fe->get_nonzero_components(I)[component] == false)
+      if (fe->get_nonzero_components(i)[component] == false)
        return Tensor<1,dim>();
 
                                       // count how many non-zero
@@ -2550,10 +2493,10 @@ FEValuesBase<dim>::shape_grad_component (const unsigned int i,
                                       // shape function in the arrays
                                       // we index presently:
       const unsigned int
-       row = (this->shape_function_to_row_table[I]
+       row = (this->shape_function_to_row_table[i]
               +
-              std::count (fe->get_nonzero_components(I).begin(),
-                          fe->get_nonzero_components(I).begin()+component,
+              std::count (fe->get_nonzero_components(i).begin(),
+                          fe->get_nonzero_components(i).begin()+component,
                           true));
       return this->shape_gradients[row][j];
     };
@@ -2567,21 +2510,19 @@ const Tensor<2,dim> &
 FEValuesBase<dim>::shape_2nd_derivative (const unsigned int i,
                                         const unsigned int j) const
 {
-  const unsigned int I=shift_shape_function_index(i);
-
   Assert (this->update_flags & update_second_derivatives,
          ExcAccessToUninitializedField());
-  Assert (fe->is_primitive (I),
+  Assert (fe->is_primitive (i),
          ExcShapeFunctionNotPrimitive(i));
-  Assert (I<this->shape_2nd_derivatives.size(),
-         ExcIndexRange (I, 0, this->shape_2nd_derivatives.size()));
+  Assert (i<this->shape_2nd_derivatives.size(),
+         ExcIndexRange (i, 0, this->shape_2nd_derivatives.size()));
   Assert (j<this->shape_2nd_derivatives[0].size(),
          ExcIndexRange (j, 0, this->shape_2nd_derivatives[0].size()));
 
                                   // if the entire FE is primitive,
                                   // then we can take a short-cut:
   if (fe->is_primitive())
-    return this->shape_2nd_derivatives[I][j];
+    return this->shape_2nd_derivatives[i][j];
   else
                                     // otherwise, use the mapping
                                     // between shape function numbers
@@ -2592,7 +2533,7 @@ FEValuesBase<dim>::shape_2nd_derivative (const unsigned int i,
                                     // question to which vector
                                     // component the call of this
                                     // function refers
-    return this->shape_2nd_derivatives[this->shape_function_to_row_table[I]][j];
+    return this->shape_2nd_derivatives[this->shape_function_to_row_table[i]][j];
 }
 
 
@@ -2604,8 +2545,6 @@ FEValuesBase<dim>::shape_2nd_derivative_component (const unsigned int i,
                                                   const unsigned int j,
                                                   const unsigned int component) const
 {
-  const unsigned int I=shift_shape_function_index(i);
-
   Assert (this->update_flags & update_second_derivatives,
          ExcAccessToUninitializedField());
   Assert (component < fe->n_components(),
@@ -2620,10 +2559,10 @@ FEValuesBase<dim>::shape_2nd_derivative_component (const unsigned int i,
                                   // system_to_component_table only
                                   // works if the shape function is
                                   // primitive):
-  if (fe->is_primitive(I))
+  if (fe->is_primitive(i))
     {
-      if (component == fe->system_to_component_index(I).first)
-       return this->shape_2nd_derivatives[this->shape_function_to_row_table[I]][j];
+      if (component == fe->system_to_component_index(i).first)
+       return this->shape_2nd_derivatives[this->shape_function_to_row_table[i]][j];
       else
        return Tensor<2,dim>();
     }
@@ -2650,10 +2589,10 @@ FEValuesBase<dim>::shape_2nd_derivative_component (const unsigned int i,
                                       // shape function in the arrays
                                       // we index presently:
       const unsigned int
-       row = (this->shape_function_to_row_table[I]
+       row = (this->shape_function_to_row_table[i]
               +
-              std::count (fe->get_nonzero_components(I).begin(),
-                          fe->get_nonzero_components(I).begin()+component,
+              std::count (fe->get_nonzero_components(i).begin(),
+                          fe->get_nonzero_components(i).begin()+component,
                           true));
       return this->shape_2nd_derivatives[row][j];
     };
@@ -2737,28 +2676,6 @@ FEValuesBase<dim>::JxW (const unsigned int i) const
 }
 
 
-
-template <int dim>
-inline
-unsigned int
-FEValuesBase<dim>::shift_shape_function_index (const unsigned int i) const
-{
-                                  // standard implementation for 1D and 2D
-  Assert(i<fe->dofs_per_cell, ExcInternalError());
-  return i;
-}
-
-template <>
-inline
-unsigned int
-FEValuesBase<3>::shift_shape_function_index (const unsigned int i) const
-{
-  Assert(i<fe->dofs_per_cell, ExcInternalError());
-  return this->permuted_shape_functions[i];
-}
-
-
-
 /*------------------------ Inline functions: FEValues ----------------------------*/
 
 
@@ -2844,8 +2761,6 @@ FEFaceValuesBase<dim>::boundary_form (const unsigned int i) const
   return this->boundary_forms[i];
 }
 
-
-
 #endif // DOXYGEN
 
 DEAL_II_NAMESPACE_CLOSE
index b0f24874c5d6f613772b67bd33e3e8da66431320..13043e4642fe1c4eac924628a0abc08cd6ad4381 100644 (file)
@@ -253,8 +253,13 @@ DoFCellAccessor<DoFHandler<3> >::update_cell_dof_indices_cache () const
     for (unsigned int d=0; d<dofs_per_line; ++d)
       *next++ = this->line(line)->dof_index(d);
   for (unsigned int quad=0; quad<6; ++quad)
-    for (unsigned int d=0; d<dofs_per_quad; ++d)
-      *next++ = this->quad(quad)->dof_index(d);
+    if (this->face_orientation(quad))
+      for (unsigned int d=0; d<dofs_per_quad; ++d)
+       *next++ = this->quad(quad)->dof_index(d);
+    else
+      for (unsigned int d=0; d<dofs_per_quad; ++d)
+       *next++ = this->quad(quad)->dof_index(this->dof_handler->get_fe().
+                                             adjust_quad_dof_index_for_face_orientation(d));
   for (unsigned int d=0; d<dofs_per_hex; ++d)
     *next++ = dof_index(d);
 }
index d08e3ff599ad65e95f54fbbe440936f782878baa..b42264f063f64f792d7bc57395948077c1c6ee56 100644 (file)
@@ -114,6 +114,7 @@ FiniteElement<dim>::FiniteElement (
                :
                FiniteElementData<dim> (fe_data),
                cached_primitivity(false),
+               adjust_quad_dof_index_for_face_orientation_table (this->dofs_per_quad, 0),
                 system_to_base_table(this->dofs_per_cell),
                 face_system_to_base_table(this->dofs_per_face),                
                 component_to_base_table (this->components,
@@ -333,25 +334,25 @@ FiniteElement<dim>::component_to_block_index (const unsigned int index) const
 }
 
 
-#if deal_II_dimension < 3
-
 template <int dim>
-void
-FiniteElement<dim>::get_face_shape_function_shifts (std::vector<int> &shifts) const
+unsigned int
+FiniteElement<dim>::adjust_quad_dof_index_for_face_orientation (const unsigned int) const
 {
-                                  // general template for 1D and 2D, return an
-                                  // empty vector
-  shifts.clear();
+                                  // general template for 1D and 2D: not implemented
+  Assert (false, ExcNotImplemented());
+  return deal_II_numbers::invalid_unsigned_int;
 }
 
-#else
+#if deal_II_dimension == 3
 
 template <>
-void
-FiniteElement<3>::get_face_shape_function_shifts (std::vector<int> &shifts) const
+unsigned int
+FiniteElement<3>::adjust_quad_dof_index_for_face_orientation (const unsigned int index) const
 {
-  shifts.clear();
-  shifts.resize(this->dofs_per_quad,0);
+  Assert (index<this->dofs_per_quad, ExcIndexRange(index,0,this->dofs_per_quad));
+  Assert (adjust_quad_dof_index_for_face_orientation_table.size()==this->dofs_per_quad,
+         ExcInternalError());
+  return index+adjust_quad_dof_index_for_face_orientation_table[index];
 }
 
 #endif
index b584801c386bf89fc5b9356d198226fcefbb42b8..199ecb5e5c9bc6b51772bf0ecf72622ac21f26b4 100644 (file)
@@ -210,6 +210,8 @@ FE_Q<dim>::FE_Q (const unsigned int degree)
   initialize_constraints ();
   initialize_embedding ();
   initialize_restriction ();
+
+  initialize_quad_dof_index_permutation();
 }
 
 
@@ -709,38 +711,6 @@ compare_for_face_domination (const FiniteElement<dim> &fe_other) const
 }
 
 
-
-template <int dim>
-void
-FE_Q<dim>::get_face_shape_function_shifts (std::vector<int> &shifts) const
-{
-                                  // general template for 1D and 2D, return an
-                                  // empty vector
-  shifts.clear();
-}
-
-
-
-#if deal_II_dimension == 3
-
-template <>
-void
-FE_Q<3>::get_face_shape_function_shifts (std::vector<int> &shifts) const
-{
-  shifts.resize(this->dofs_per_quad);
-
-  unsigned int points=this->degree-1;
-  Assert(points*points==this->dofs_per_quad, ExcInternalError());
-               
-  for (unsigned int local=0; local<this->dofs_per_quad; ++local)
-                                    // face support points are in lexicographic
-                                    // ordering with x running fastest. invert
-                                    // that (y running fastest)
-    shifts[local] = (local%points)*points + local/points - local;
-}
-
-#endif
-
 //---------------------------------------------------------------------------
 // Auxiliary functions
 //---------------------------------------------------------------------------
@@ -825,6 +795,40 @@ void FE_Q<dim>::initialize_unit_face_support_points ()
 
 
 
+template <int dim>
+void
+FE_Q<dim>::initialize_quad_dof_index_permutation ()
+{
+                                  // general template for 1D and 2D, do nothing
+}
+
+
+
+#if deal_II_dimension == 3
+
+template <>
+void
+FE_Q<3>::initialize_quad_dof_index_permutation ()
+{
+
+  Assert (adjust_quad_dof_index_for_face_orientation_table.size()==this->dofs_per_quad,
+         ExcInternalError());
+
+  unsigned int points=this->degree-1;
+  Assert(points*points==this->dofs_per_quad, ExcInternalError());
+               
+  for (unsigned int local=0; local<this->dofs_per_quad; ++local)
+                                    // face support points are in lexicographic
+                                    // ordering with x running fastest. invert
+                                    // that (y running fastest)
+    this->adjust_quad_dof_index_for_face_orientation_table[local]
+      =(local%points)*points + local/points - local;
+}
+
+#endif
+
+
+
 template <int dim>
 std::vector<unsigned int>
 FE_Q<dim>::get_dpo_vector(const unsigned int deg)
index 48f2e4b7a8ec00653af28bce43b289f562910e76..b6db8a9000a3a54c9b18b35b37c7cd34408ee491 100644 (file)
@@ -1729,6 +1729,8 @@ void FESystem<dim>::initialize ()
                                   // on cell and face
   initialize_unit_support_points ();
   initialize_unit_face_support_points ();
+
+  initialize_quad_dof_index_permutation ();
 }
 
 
@@ -1847,6 +1849,40 @@ initialize_unit_face_support_points ()
 
 
 
+template <int dim>
+void
+FESystem<dim>::initialize_quad_dof_index_permutation ()
+{
+                                  // general template for 1D and 2D, do nothing
+}
+
+
+
+#if deal_II_dimension == 3
+
+template <>
+void
+FESystem<3>::initialize_quad_dof_index_permutation ()
+{
+                                  // to obtain the shifts for this composed
+                                  // element, concatenate the shift vectors of
+                                  // the base elements
+  for (unsigned int b=0; b<n_base_elements();++b)
+    {
+      const std::vector<int> &temp=this->base_element(b).adjust_quad_dof_index_for_face_orientation_table;
+      for (unsigned int c=0; c<element_multiplicity(b); ++c)
+       adjust_quad_dof_index_for_face_orientation_table.insert
+         (adjust_quad_dof_index_for_face_orientation_table.begin(),
+          temp.begin(),temp.end());
+    }
+  Assert (adjust_quad_dof_index_for_face_orientation_table.size()==this->dofs_per_quad,
+         ExcInternalError());
+}
+
+#endif
+
+
+
 template <int dim>
 bool
 FESystem<dim>::
@@ -2926,41 +2962,6 @@ FESystem<dim>::unit_face_support_point (const unsigned index) const
 
 
 
-template <int dim>
-void
-FESystem<dim>::get_face_shape_function_shifts (std::vector<int> &shifts) const
-{
-                                  // general template for 1D and 2D, return an
-                                  // empty vector
-  shifts.clear();
-}
-
-
-
-#if deal_II_dimension == 3
-
-template <>
-void
-FESystem<3>::get_face_shape_function_shifts (std::vector<int> &shifts) const
-{
-  shifts.clear();
-  std::vector<int> temp;
-                                  // to obtain the shifts for this composed
-                                  // element, concatenate the shift vectors of
-                                  // the base elements
-  for (unsigned int b=0; b<n_base_elements();++b)
-    {
-      this->base_element(b).get_face_shape_function_shifts(temp);
-      for (unsigned int c=0; c<element_multiplicity(b); ++c)
-       shifts.insert(shifts.begin(),temp.begin(),temp.end());
-    }
-  Assert (shifts.size()==this->dofs_per_quad, ExcInternalError());
-}
-
-#endif
-
-
-
 template <int dim>
 unsigned int
 FESystem<dim>::memory_consumption () const 
index 7d5e1eadb67118ec43666d087abc552f8a3b8a9e..84aeac4bf5b0340d061eca290417af9cd33dfe69 100644 (file)
@@ -316,34 +316,6 @@ FEValuesData<dim>::initialize (const unsigned int        n_quadrature_points,
 
   if (flags & update_cell_JxW_values)
     this->cell_JxW_values.resize(n_quadrature_points);
-
-                                  // initialize the permutation fields, if they
-                                  // are needed
-  if (dim==3)
-    {
-      permuted_shape_functions.resize(fe.dofs_per_cell);
-                                      // initialize cell permutation mapping
-                                      // with identity
-      for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
-       this->permuted_shape_functions[i]=i;
-
-      if (fe.dofs_per_quad>0 &&
-         (flags & update_values ||
-          flags & update_gradients ||
-          flags & update_second_derivatives))
-       {
-                                          // ask the fe to fill the vector of
-                                          // shifts
-         fe.get_face_shape_function_shifts(shift_in_face_shape_functions);
-         Assert (shift_in_face_shape_functions.size()==fe.dofs_per_quad,
-                 ExcInternalError());
-         update_shape_function_permutation=true;
-       }
-      else
-       update_shape_function_permutation=false;
-    }
-  else
-    update_shape_function_permutation=false;
 }
 
 
@@ -393,43 +365,6 @@ FEValuesBase<dim>::~FEValuesBase ()
 }
 
 
-#if deal_II_dimension <3
-template <int dim>
-void FEValuesBase<dim>::reinit ()
-{
-                                  // do nothing in 1D and 2D
-}
-#else
-
-template <>
-void FEValuesBase<3>::reinit ()
-{
-                                  // in 3D reinit the permutation of face shape
-                                  // functions, if necessary
-  if (update_shape_function_permutation)
-    {
-      Assert (this->shift_in_face_shape_functions.size()==this->fe->dofs_per_quad,
-             ExcInternalError());
-      Triangulation<3>::cell_iterator thiscell = *(this->present_cell);
-      unsigned int offset=fe->first_quad_index;
-      for (unsigned int face_no=0; face_no<GeometryInfo<3>::faces_per_cell; ++face_no)
-       {
-         if (thiscell->face_orientation(face_no))
-           for (unsigned int i=0; i<fe->dofs_per_quad; ++i)
-             this->permuted_shape_functions[offset+i]=offset+i;
-         else
-           for (unsigned int i=0; i<fe->dofs_per_quad; ++i)
-             this->permuted_shape_functions[offset+i]=offset+i+this->shift_in_face_shape_functions[i];
-
-         offset+=this->fe->dofs_per_quad;
-       }
-      Assert (offset-fe->first_quad_index==GeometryInfo<3>::faces_per_cell*this->fe->dofs_per_quad,
-             ExcInternalError());
-    }
-}
-#endif
-
-
 
 template <int dim>
 template <class InputVector, typename number>
@@ -1076,15 +1011,12 @@ FEValuesBase<dim>::memory_consumption () const
   return (MemoryConsumption::memory_consumption (this->shape_values) +
          MemoryConsumption::memory_consumption (this->shape_gradients) +
          MemoryConsumption::memory_consumption (this->shape_2nd_derivatives) +
-         MemoryConsumption::memory_consumption (this->permuted_shape_functions) +
-         MemoryConsumption::memory_consumption (this->shift_in_face_shape_functions) +
          MemoryConsumption::memory_consumption (this->JxW_values) +
          MemoryConsumption::memory_consumption (this->quadrature_points) +
          MemoryConsumption::memory_consumption (this->normal_vectors) +
          MemoryConsumption::memory_consumption (this->boundary_forms) +
          MemoryConsumption::memory_consumption (this->cell_JxW_values) +
          sizeof(this->update_flags) +
-         sizeof(this->update_shape_function_permutation) +
          MemoryConsumption::memory_consumption (n_quadrature_points) +
          MemoryConsumption::memory_consumption (dofs_per_cell) +
          MemoryConsumption::memory_consumption (mapping) +
@@ -1328,8 +1260,6 @@ void FEValues<dim>::do_reinit ()
 
   this->fe_data->clear_first_cell ();
   this->mapping_data->clear_first_cell ();
-  
-  FEValuesBase<dim>::reinit();
 }
 
 
@@ -1607,7 +1537,6 @@ void FEFaceValues<dim>::do_reinit (const unsigned int face_no)
 
   this->fe_data->clear_first_cell ();
   this->mapping_data->clear_first_cell ();
-  FEValuesBase<dim>::reinit();
 }
 
 
@@ -1860,7 +1789,6 @@ void FESubfaceValues<dim>::do_reinit (const unsigned int face_no,
 
   this->fe_data->clear_first_cell ();
   this->mapping_data->clear_first_cell ();
-  FEValuesBase<dim>::reinit();
 }
 
 
index 77293dbeaac517d723199e4683ad650c7a01e60f..cce5704b4ff455a268e831bb31df3123e18b16d0 100644 (file)
@@ -959,9 +959,9 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
-  <li> <p> Fixed: On faces with wrong <code>face_orientation</code> the shape
-       functions have to reordered in order to be combined with the correct
-       dofs. This is only relevant for continuous elements in 3D. At least for
+  <li> <p> Fixed: On faces with wrong <code>face_orientation</code> the dofs
+       have to reordered in order to be combined with the correct shape
+       functions. This is only relevant for continuous elements in 3D. At least for
        <code class="class">FE_Q</code> and systems of <code
        class="class">FE_Q</code> this works now, for other finite elements the
        reordering vector still has to be implemented.

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.