]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move FiniteElementData::face_to_cell_index to FiniteElement. The function needs to...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 9 Aug 2013 22:37:31 +0000 (22:37 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 9 Aug 2013 22:37:31 +0000 (22:37 +0000)
Also remove FiniteElementData::face_to_equivalent_cell_index: this function had been deprecated a while back already.

git-svn-id: https://svn.dealii.org/trunk@30271 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/fe/fe.h
deal.II/include/deal.II/fe/fe_base.h
deal.II/source/fe/fe.cc
deal.II/source/fe/fe_data.cc

index c1c4818f4397f6ab4ca9cc26891bb4d55e498ba5..ad91e1e9855951466bd97d309c98d3a465761bc6 100644 (file)
@@ -23,7 +23,11 @@ modifications to application programs. We apologize for the
 inconvenience this causes.
 </p>
 
-<ol>
+<ol> Removed: The member function face_to_equivalent_cell_index() in
+FiniteElementData has been removed. It had been deprecated a while
+back already. Please use FiniteElement::face_to_cell_index() instead.
+<br>
+(Wolfgang Bangerth, 2013/08/09)
 
 </ol>
 
index 8d52d8c9430f9e71f34c2f8b440089696beecfd8..31fb8c8cb7d0fd5bfe85246c8869e39c156175ba 100644 (file)
@@ -924,6 +924,65 @@ public:
   std::pair<unsigned int, unsigned int>
   face_system_to_component_index (const unsigned int index) const;
 
+  /**
+   * Given an index in the natural ordering of indices on a face, return the
+   * index of the same degree of freedom on the cell.
+   *
+   * To explain the concept, consider the case where we would like to know
+   * whether a degree of freedom on a face, for example as part of an FESystem
+   * element, is primitive. Unfortunately, the
+   * is_primitive() function in the FiniteElement class takes a cell index, so
+   * we would need to find the cell index of the shape function that
+   * corresponds to the present face index. This function does that.
+   *
+   * Code implementing this would then look like this:
+   * @code
+   * for (i=0; i<dofs_per_face; ++i)
+   *  if (fe.is_primitive(fe.face_to_equivalent_cell_index(i, some_face_no)))
+   *   ... do whatever
+   * @endcode
+   * The function takes additional arguments that account for the fact that
+   * actual faces can be in their standard ordering with respect to the cell
+   * under consideration, or can be flipped, oriented, etc.
+   *
+   * @param face_dof_index The index of the degree of freedom on a face.
+   *   This index must be between zero and dofs_per_face.
+   * @param face The number of the face this degree of freedom lives on.
+   *   This number must be between zero and GeometryInfo::faces_per_cell.
+   * @param face_orientation One part of the description of the orientation
+   *   of the face. See @ref GlossFaceOrientation .
+   * @param face_flip One part of the description of the orientation
+   *   of the face. See @ref GlossFaceOrientation .
+   * @param face_rotation One part of the description of the orientation
+   *   of the face. See @ref GlossFaceOrientation .
+   * @return The index of this degree of freedom within the set
+   *   of degrees of freedom on the entire cell. The returned value
+   *   will be between zero and dofs_per_cell.
+   *
+   * @note This function exists in this class because that is where it
+   * was first implemented. However, it can't really work in the most
+   * general case without knowing what element we have. The reason is that
+   * when a face is flipped or rotated, we also need to know whether we
+   * need to swap the degrees of freedom on this face, or whether they
+   * are immune from this. For this, consider the situation of a $Q_3$
+   * element in 2d. If face_flip is true, then we need to consider
+   * the two degrees of freedom on the edge in reverse order. On the other
+   * hand, if the element were a $Q_1^2$, then because the two degrees of
+   * freedom on this edge belong to different vector components, they
+   * should not be considered in reverse order. What all of this shows is
+   * that the function can't work if there are more than one degree of
+   * freedom per line or quad, and that in these cases the function will
+   * throw an exception pointing out that this functionality will need
+   * to be provided by a derived class that knows what degrees of freedom
+   * actually represent.
+   */
+  virtual
+  unsigned int face_to_cell_index (const unsigned int face_dof_index,
+                                   const unsigned int face,
+                                   const bool face_orientation = true,
+                                   const bool face_flip        = false,
+                                   const bool face_rotation    = false) 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
index c317c45b323f473f0510370f91c310dd12a84832..3cd065bb95eb54a54ad9c3a3aea146f95ad91e17 100644 (file)
@@ -409,72 +409,6 @@ public:
    */
   bool conforms (const Conformity) const;
 
-  /**
-   * Given an index in the natural ordering of indices on a face, return the
-   * index of the same degree of freedom on the cell.
-   *
-   * @param face_dof_index The index of the degree of freedom on a face.
-   *   This index must be between zero and dofs_per_face.
-   * @param face The number of the face this degree of freedom lives on.
-   *   This number must be between zero and GeometryInfo::faces_per_cell.
-   * @param face_orientation One part of the description of the orientation
-   *   of the face. See @ref GlossFaceOrientation .
-   * @param face_flip One part of the description of the orientation
-   *   of the face. See @ref GlossFaceOrientation .
-   * @param face_rotation One part of the description of the orientation
-   *   of the face. See @ref GlossFaceOrientation .
-   * @return The index of this degree of freedom within the set
-   *   of degrees of freedom on the entire cell. The returned value
-   *   will be between zero and dofs_per_cell.
-   *
-   * @note This function exists in this class because that is where it
-   * was first implemented. However, it can't really work in the most
-   * general case without knowing what element we have. The reason is that
-   * when a face is flipped or rotated, we also need to know whether we
-   * need to swap the degrees of freedom on this face, or whether they
-   * are immune from this. For this, consider the situation of a $Q_3$
-   * element in 2d. If face_flip is true, then we need to consider
-   * the two degrees of freedom on the edge in reverse order. On the other
-   * hand, if the element were a $Q_1^2$, then because the two degrees of
-   * freedom on this edge belong to different vector components, they
-   * should not be considered in reverse order. What all of this shows is
-   * that the function can't work if there are more than one degree of
-   * freedom per line or quad, and that in these cases the function will
-   * throw an exception pointing out that this functionality will need
-   * to be provided by a derived class that knows what degrees of freedom
-   * actually represent.
-   */
-  virtual
-  unsigned int face_to_cell_index (const unsigned int face_dof_index,
-                                   const unsigned int face,
-                                   const bool face_orientation = true,
-                                   const bool face_flip        = false,
-                                   const bool face_rotation    = false) const;
-
-  /**
-   * @deprecated This function is just a special version of face_to_cell_index
-   * for the face zero. It is therefore of limited use in aplications and most
-   * of the time, the other function is what is required.
-   *
-   * Given an index in the natural ordering of indices on a face, return the
-   * index of an equivalent degree of freedom on the cell.
-   *
-   * To explain the concept, consider the case where we would like to know
-   * whether a degree of freedom on a face is primitive. Unfortunately, the
-   * is_primitive() function in the FiniteElement class takes a cell index, so
-   * we would need to find the cell index of the shape function that
-   * corresponds to the present face index. This function does that.
-   *
-   * Code implementing this would then look like this:
-   * @code
-   * for (i=0; i<dofs_per_face; ++i)
-   *  if (fe.is_primitive(fe.face_to_equivalent_cell_index(i)))
-   *   ... do whatever
-   * @endcode
-   *
-   */
-  unsigned int face_to_equivalent_cell_index (const unsigned int index) const DEAL_II_DEPRECATED;
-
   /**
    * Comparison operator.
    */
@@ -707,53 +641,6 @@ FiniteElementData<dim>::conforms (const Conformity space) const
 }
 
 
-
-template <>
-inline
-unsigned int
-FiniteElementData<1>::
-face_to_equivalent_cell_index (const unsigned int index) const
-{
-  Assert (index < dofs_per_face,
-          ExcIndexRange (index, 0, dofs_per_face));
-  return index;
-}
-
-
-
-template <>
-inline
-unsigned int
-FiniteElementData<2>::
-face_to_equivalent_cell_index (const unsigned int index) const
-{
-  Assert (index < dofs_per_face,
-          ExcIndexRange (index, 0, dofs_per_face));
-
-  // on face 0, the vertices are 0 and 2
-  if (index < this->dofs_per_vertex)
-    return index;
-  else if (index < 2*this->dofs_per_vertex)
-    return index + this->dofs_per_vertex;
-  else
-    // this is a dof on line 0, so on the cell there are now 4 vertices
-    // instead of only 2 ahead of this one
-    return index + 2*this->dofs_per_vertex;
-}
-
-
-
-
-template <>
-inline
-unsigned int
-FiniteElementData<3>::
-face_to_equivalent_cell_index (const unsigned int index) const
-{
-  // this case is just way too complicated. fall back to face_to_cell_index
-  return face_to_cell_index(index, 0, true);
-}
-
 #endif // DOXYGEN
 
 
index cdb97b5104b3bbd2dfa062be453fc04490006ee0..f5332eb778801ff6e158a983eb468695e1a5edf0 100644 (file)
@@ -544,6 +544,68 @@ block_mask (const ComponentMask &component_mask) const
 
 
 
+template <int dim, int spacedim>
+unsigned int
+FiniteElement<dim,spacedim>::
+face_to_cell_index (const unsigned int face_index,
+                    const unsigned int face,
+                    const bool face_orientation,
+                    const bool face_flip,
+                    const bool face_rotation) const
+{
+  Assert (face_index < this->dofs_per_face,
+          ExcIndexRange(face_index, 0, this->dofs_per_face));
+  Assert (face < GeometryInfo<dim>::faces_per_cell,
+          ExcIndexRange(face, 0, GeometryInfo<dim>::faces_per_cell));
+
+  // see the function's documentation for an explanation of this
+  // assertion -- in essence, derived classes have to implement
+  // an overloaded version of this function if we are to use any
+  // other than standard orientation
+  if ((face_orientation != true) || (face_flip != false) || (face_rotation != false))
+    Assert ((this->dofs_per_line <= 1) && (this->dofs_per_quad <= 1),
+            ExcMessage ("The function in this base class can not handle this case. "
+                        "Rather, the derived class you are using must provide "
+                        "an overloaded version but apparently hasn't done so. See "
+                        "the documentation of this function for more information."));
+
+  // DoF on a vertex
+  if (face_index < this->first_face_line_index)
+    {
+      // Vertex number on the face
+      const unsigned int face_vertex = face_index / this->dofs_per_vertex;
+      return face_index % this->dofs_per_vertex
+             + GeometryInfo<dim>::face_to_cell_vertices(face, face_vertex,
+                                                        face_orientation,
+                                                        face_flip,
+                                                        face_rotation)
+             * this->dofs_per_vertex;
+    }
+  // Else, DoF on a line?
+  if (face_index < this->first_face_quad_index)
+    {
+      // Ignore vertex dofs
+      const unsigned int index = face_index - this->first_face_line_index;
+      // Line number on the face
+      const unsigned int face_line = index / this->dofs_per_line;
+      return this->first_line_index + index % this->dofs_per_line
+             + GeometryInfo<dim>::face_to_cell_lines(face, face_line,
+                                                     face_orientation,
+                                                     face_flip,
+                                                     face_rotation)
+             * this->dofs_per_line;
+    }
+  // Else, DoF is on a quad
+
+  // Ignore vertex and line dofs
+  const unsigned int index = face_index - this->first_face_quad_index;
+  return this->first_quad_index + index
+         + face * this->dofs_per_quad;
+}
+
+
+
+
 template <int dim, int spacedim>
 unsigned int
 FiniteElement<dim,spacedim>::adjust_quad_dof_index_for_face_orientation (const unsigned int index,
index 3fec91c5a9020447fca4cc3b8c389fa4b34dfcb4..894e3c3f1859fbeb10eb03059b481d6a415bb9b3 100644 (file)
@@ -99,65 +99,6 @@ bool FiniteElementData<dim>::operator== (const FiniteElementData<dim> &f) const
           (conforming_space == f.conforming_space));
 }
 
-template<int dim>
-unsigned int
-FiniteElementData<dim>::
-face_to_cell_index (const unsigned int face_index,
-                    const unsigned int face,
-                    const bool face_orientation,
-                    const bool face_flip,
-                    const bool face_rotation) const
-{
-  Assert (face_index < this->dofs_per_face,
-          ExcIndexRange(face_index, 0, this->dofs_per_face));
-  Assert (face < GeometryInfo<dim>::faces_per_cell,
-          ExcIndexRange(face, 0, GeometryInfo<dim>::faces_per_cell));
-
-  // see the function's documentation for an explanation of this
-  // assertion -- in essence, derived classes have to implement
-  // an overloaded version of this function if we are to use any
-  // other than standard orientation
-  if ((face_orientation != true) || (face_flip != false) || (face_rotation != false))
-    Assert ((this->dofs_per_line <= 1) && (this->dofs_per_quad <= 1),
-            ExcMessage ("The function in this base class can not handle this case. "
-                        "Rather, the derived class you are using must provide "
-                        "an overloaded version but apparently hasn't done so. See "
-                        "the documentation of this function for more information."));
-
-  // DoF on a vertex
-  if (face_index < this->first_face_line_index)
-    {
-      // Vertex number on the face
-      const unsigned int face_vertex = face_index / this->dofs_per_vertex;
-      return face_index % this->dofs_per_vertex
-             + GeometryInfo<dim>::face_to_cell_vertices(face, face_vertex,
-                                                        face_orientation,
-                                                        face_flip,
-                                                        face_rotation)
-             * this->dofs_per_vertex;
-    }
-  // Else, DoF on a line?
-  if (face_index < this->first_face_quad_index)
-    {
-      // Ignore vertex dofs
-      const unsigned int index = face_index - this->first_face_line_index;
-      // Line number on the face
-      const unsigned int face_line = index / this->dofs_per_line;
-      return this->first_line_index + index % this->dofs_per_line
-             + GeometryInfo<dim>::face_to_cell_lines(face, face_line,
-                                                     face_orientation,
-                                                     face_flip,
-                                                     face_rotation)
-             * this->dofs_per_line;
-    }
-  // Else, DoF is on a quad
-
-  // Ignore vertex and line dofs
-  const unsigned int index = face_index - this->first_face_quad_index;
-  return this->first_quad_index + index
-         + face * this->dofs_per_quad;
-}
-
 
 template class FiniteElementData<1>;
 template class FiniteElementData<2>;

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.