]> https://gitweb.dealii.org/ - dealii.git/commitdiff
remove cell_normal and move normal into FEValuesBase
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 24 Feb 2010 18:12:25 +0000 (18:12 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 24 Feb 2010 18:12:25 +0000 (18:12 +0000)
git-svn-id: https://svn.dealii.org/trunk@20688 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_update_flags.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/source/fe/fe_values.cc
deal.II/deal.II/source/fe/mapping_q.cc
deal.II/deal.II/source/fe/mapping_q1.cc
deal.II/deal.II/source/fe/mapping_q1_eulerian.cc
deal.II/deal.II/source/fe/mapping_q_eulerian.cc

index 33fe2b8384ded798f311dfc5847a5030da49c1ca..c9438f008d5a5da00b8e5f220c1f306c9c80ab0b 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -113,7 +113,7 @@ enum UpdateFlags
                                        * the real cell.
                                        */
       update_hessians                     = 0x0004,
-                                      //! Outter normal vector, not normalized
+                                      //! Outer normal vector, not normalized
                                       /**
                                        * Vector product of tangential
                                        * vectors, yielding a normal
@@ -143,21 +143,24 @@ enum UpdateFlags
                                        * reference to realcell.
                                        */
       update_JxW_values                   = 0x0020,
-                                      //! Normal vectors to the faces
+                                       //! Normal vectors
+                                      /**
+                                       * Compute the normal vectors,
+                                       * either for a face or for a
+                                       * cell of codimension
+                                       * one. Setting this flag for
+                                       * any other object will raise
+                                       * an error.
+                                       */
+      update_normal_vectors               = 0x0040,
                                       /**
-                                       * Compute the unit outer
-                                       * normal vector on the face of
-                                       * a cell.
+                                       * @deprecated Use #update_normal_vectors
                                        */
-      update_face_normal_vectors          = 0x0040,
-                                      //! Normal vectors to the cells
+      update_face_normal_vectors          = update_normal_vectors,
                                       /**
-                                       * Compute the unit outer
-                                       * normal vector on the cell
-                                       * itself. Only possible if
-                                       * dim=spacedim-1
+                                       * @deprecated Use #update_normal_vectors
                                        */      
-      update_cell_normal_vectors          = 0x20000,
+      update_cell_normal_vectors          = update_normal_vectors,
                                       //! Volume element
                                       /**
                                        * Compute the Jacobian of the
@@ -255,13 +258,6 @@ enum UpdateFlags
                                        * derivatives.
                                        */
       update_second_derivatives = update_hessians,
-                                       //! Normal vectors
-                                      /**
-                                       * @deprecated Update normal
-                                       * vectors. Use
-                                       * update_face_normal_vectors
-                                       */
-      update_normal_vectors               = update_face_normal_vectors,
                                       //! Values needed for Piola transform
                                       /**
                                        * Combination of the flags
index 79501b6fe3fc5b60c7f77ebebea511d31950a6c4..e4d8148b758e2081f5159fb2c90eb01f40a2a0a4 100644 (file)
@@ -1282,14 +1282,6 @@ class FEValuesData
                                      */
     std::vector<Point<spacedim> >  normal_vectors;
 
-                                    /**
-                                     * List of outward vectors normal to the cell
-                                     * surface (line) at the quadrature points
-                                     * for the codimension 1 case,
-                                     * when spacedim=3 (=2).
-                                     */
-    std::vector<Point<spacedim> >  cell_normal_vectors;
-
                                      /**
                                      * List of boundary forms at the
                                      * quadrature points. This field is filled
@@ -2550,15 +2542,38 @@ class FEValuesBase : protected FEValuesData<dim,spacedim>,
                                      */
     const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
 
-                                    /**
-                                     * Returns the vectors normal to
-                                     * the cell in each of the
-                                     * quadrature points.
+                                    /**
+                                     * For a face, return the outward
+                                     * normal vector to the cell at
+                                     * the <tt>i</tt>th quadrature
+                                     * point.
+                                     *
+                                     * For a cell of codimension one,
+                                     * return the normal vector, as
+                                     * it is specified by the
+                                     * numbering of the vertices.
+                                     *
+                                     * The length of the vector
+                                     * is normalized to one.
                                      */
+    const Point<spacedim> & normal_vector (const unsigned int i) const;
 
-    const std::vector<Point<spacedim> > & get_cell_normal_vectors () const;
+                                    /**
+                                     * Return the normal vectors at
+                                     * the quadrature points. For a
+                                     * face, these are the outward
+                                     * normal vectors to the
+                                     * cell. For a cell of
+                                     * codimension one, the
+                                     * orientation is given by the
+                                     * numbering of vertices.
+                                     */
+    const std::vector<Point<spacedim> > & get_normal_vectors () const;
 
                                     /**
+                                     * @deprecated Use
+                                     * normal_vector() instead.
+                                     *
                                      * Return the outward normal vector to
                                      * the cell at the <tt>i</tt>th quadrature
                                      * point. The length of the vector
@@ -2566,6 +2581,16 @@ class FEValuesBase : protected FEValuesData<dim,spacedim>,
                                      */
     const Point<spacedim> & cell_normal_vector (const unsigned int i) const;
 
+                                    /**
+                                     * @deprecated Use
+                                     * get_normal_vectors() instead.
+                                     *
+                                     * Returns the vectors normal to
+                                     * the cell in each of the
+                                     * quadrature points.
+                                     */
+    const std::vector<Point<spacedim> > & get_cell_normal_vectors () const;
+
                                     /**
                                      * Return the relation of the current
                                      * cell to the previous cell. This
@@ -2998,7 +3023,7 @@ class FEValues : public FEValuesBase<dim,spacedim>
  * when evaluating something on the surface of a cell. All the data
  * that is available in the interior of cells is also available here.
  *
- * On surfaces of mesh cells, normal vectors and boundary forms are
+ * On surfaces of mesh cells, boundary forms are
  * additional values that can be computed. This class provides the
  * interface to access those. Implementations are in derived classes
  * FEFaceValues and FESubfaceValues.
@@ -3050,14 +3075,6 @@ class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
                      const FiniteElement<dim,spacedim> &fe,
                      const Quadrature<dim-1>&           quadrature);
 
-                                    /**
-                                     * Return the outward normal vector to
-                                     * the cell at the <tt>i</tt>th quadrature
-                                     * point. The length of the vector
-                                     * is normalized to one.
-                                     */
-    const Point<dim> & normal_vector (const unsigned int i) const;
-
                                     /**
                                      * Boundary form of the
                                      * transformation of the cell at
@@ -3077,13 +3094,6 @@ class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
                                      */
     const Tensor<1,spacedim> & boundary_form (const unsigned int i) const;
 
-                                    /**
-                                     * Return the list of outward normal
-                                     * vectors to the cell at the
-                                     * quadrature points.
-                                     */
-    const std::vector<Point<dim> > & get_normal_vectors () const;
-
                                     /**
                                      * Return the list of outward
                                      * normal vectors times the
@@ -4497,15 +4507,25 @@ get_function_2nd_derivatives (const InputVector                         &fe_func
 template <int dim, int spacedim>
 inline
 const Point<spacedim> &
-FEValuesBase<dim,spacedim>::cell_normal_vector (const unsigned int i) const
+FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
 {
   typedef FEValuesBase<dim,spacedim> FVB;
-  Assert (i<this->cell_normal_vectors.size(),
-         ExcIndexRange(i, 0, this->cell_normal_vectors.size()));
-  Assert (this->update_flags & update_cell_normal_vectors,
+  Assert (this->update_flags & update_normal_vectors,
          typename FVB::ExcAccessToUninitializedField());
+  Assert (i<this->normal_vectors.size(),
+         ExcIndexRange(i, 0, this->normal_vectors.size()));
+
+  return this->normal_vectors[i];
+}
 
-  return this->cell_normal_vectors[i];
+
+
+template <int dim, int spacedim>
+inline
+const Point<spacedim> &
+FEValuesBase<dim,spacedim>::cell_normal_vector (const unsigned int i) const
+{
+  return this->normal_vector(i);
 }
 
 
@@ -4536,21 +4556,6 @@ FEValues<dim,spacedim>::get_present_fe_values () const
 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
 
 
-template <int dim, int spacedim>
-inline
-const Point<dim> &
-FEFaceValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
-{
-  typedef FEValuesBase<dim,spacedim> FVB;
-  Assert (i<this->normal_vectors.size(),
-         ExcIndexRange(i, 0, this->normal_vectors.size()));
-  Assert (this->update_flags & update_normal_vectors,
-         typename FVB::ExcAccessToUninitializedField());
-
-  return this->normal_vectors[i];
-}
-
-
 template <int dim, int spacedim>
 inline
 unsigned int
index 983102bcaf51a74e85082e546cfb39ad13e54932..951b5b5c414cb264202d745ffaab6a01f5192f78 100644 (file)
@@ -1463,9 +1463,6 @@ FEValuesData<dim,spacedim>::initialize (const unsigned int        n_quadrature_p
 
   if (flags & update_normal_vectors)
     this->normal_vectors.resize(n_quadrature_points);
-
-  if (flags & update_cell_normal_vectors)
-    this->cell_normal_vectors.resize(n_quadrature_points);
 }
 
 
@@ -2959,12 +2956,21 @@ void FEValuesBase<dim,spacedim>::get_function_laplacians (
 
 template <int dim, int spacedim>
 const std::vector<Point<spacedim> > &
-FEValuesBase<dim,spacedim>::get_cell_normal_vectors () const
+FEValuesBase<dim,spacedim>::get_normal_vectors () const
 {
   typedef FEValuesBase<dim,spacedim> FEVB;
-  Assert (this->update_flags & update_cell_normal_vectors,
+  Assert (this->update_flags & update_normal_vectors,
          typename FEVB::ExcAccessToUninitializedField());
-  return this->cell_normal_vectors;
+  return this->normal_vectors;
+}
+
+
+
+template <int dim, int spacedim>
+const std::vector<Point<spacedim> > &
+FEValuesBase<dim,spacedim>::get_cell_normal_vectors () const
+{
+  return this->get_normal_vectors ();
 }
 
 
@@ -3114,20 +3120,14 @@ template <int dim, int spacedim>
 void
 FEValues<dim,spacedim>::initialize (const UpdateFlags update_flags)
 {
-                                  // you can't compute normal vectors
-                                  // on cells, only on faces
-  typedef FEValuesBase<dim,spacedim> FEVB;
-  Assert ((update_flags & update_normal_vectors) == false,
-         typename FEVB::ExcInvalidUpdateFlag());
-
                                   // You can compute normal vectors
                                   // to the cells only in the
                                   // codimension one case.
-  if(dim == spacedim)
-    Assert ( (update_flags & update_cell_normal_vectors) == false,
-            typename FEVB::ExcInvalidUpdateFlag());
-
-
+  typedef FEValuesBase<dim,spacedim> FEVB;
+  if (dim != spacedim-1)
+    Assert ((update_flags & update_normal_vectors) == false,
+           typename FEVB::ExcInvalidUpdateFlag());
+  
   const UpdateFlags flags = this->compute_update_flags (update_flags);
 
                                   // then get objects into which the
@@ -3304,7 +3304,7 @@ void FEValues<dim,spacedim>::do_reinit ()
                                     this->jacobians,
                                     this->jacobian_grads,
                                     this->inverse_jacobians,
-                                    this->cell_normal_vectors,
+                                    this->normal_vectors,
                                     this->cell_similarity);
 
   this->get_fe().fill_fe_values(this->get_mapping(),
@@ -3351,18 +3351,6 @@ FEFaceValuesBase<dim,spacedim>::FEFaceValuesBase (const unsigned int n_q_points,
 
 
 
-template <int dim, int spacedim>
-const std::vector<Point<dim> > &
-FEFaceValuesBase<dim,spacedim>::get_normal_vectors () const
-{
-  typedef FEValuesBase<dim,spacedim> FEVB;
-  Assert (this->update_flags & update_normal_vectors,
-         typename FEVB::ExcAccessToUninitializedField());
-  return this->normal_vectors;
-}
-
-
-
 template <int dim, int spacedim>
 const std::vector<Tensor<1,spacedim> > &
 FEFaceValuesBase<dim,spacedim>::get_boundary_forms () const
index a77a6b8c84e08013c80445178d7f5128d3b7b386..ce211ca58b7507ac716602015908c0ad94cd3785 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -313,7 +313,7 @@ MappingQ<dim,spacedim>::fill_fe_values (
   std::vector<Tensor<2,spacedim> >                          &jacobians,
   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
-  std::vector<Point<spacedim> >                             &cell_normal_vectors,
+  std::vector<Point<spacedim> >                             &normal_vectors,
   CellSimilarity::Similarity                           &cell_similarity) const
 {
                                   // convert data object to internal
@@ -358,7 +358,7 @@ MappingQ<dim,spacedim>::fill_fe_values (
   MappingQ1<dim,spacedim>::fill_fe_values(cell, q, *p_data,
                                          quadrature_points, JxW_values,
                                          jacobians, jacobian_grads, inverse_jacobians,
-                                         cell_normal_vectors, cell_similarity);
+                                         normal_vectors, cell_similarity);
 }
 
 
index ca0f2c818a5fb99a3443398408adba2770b4c9e8..75f67acef67c91536bc3a4f34d17a633c23d5716 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -366,7 +366,6 @@ MappingQ1<dim,spacedim>::update_once (const UpdateFlags in) const
   if (in & (update_covariant_transformation
            | update_contravariant_transformation
            | update_JxW_values
-           | update_cell_normal_vectors
            | update_boundary_forms
            | update_normal_vectors
            | update_jacobians
@@ -389,7 +388,6 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
                                      | update_covariant_transformation
                                      | update_contravariant_transformation
                                      | update_JxW_values
-                                     | update_cell_normal_vectors
                                      | update_boundary_forms
                                      | update_normal_vectors
                                      | update_volume_elements
@@ -423,7 +421,6 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
   
       if (out & (update_covariant_transformation
                 | update_JxW_values
-                | update_cell_normal_vectors
                 | update_jacobians
                 | update_jacobian_grads
                 | update_boundary_forms
@@ -442,8 +439,8 @@ MappingQ1<dim,spacedim>::update_each (const UpdateFlags in) const
       if (out & update_contravariant_transformation)
        out |= update_JxW_values;
 
-      if (out & update_cell_normal_vectors)
-       out |= update_JxW_values | update_cell_normal_vectors;
+      if (out & update_normal_vectors)
+       out |= update_JxW_values;
 
     }
   
@@ -744,7 +741,7 @@ MappingQ1<dim,spacedim>::fill_fe_values (
   std::vector<Tensor<2,spacedim> >                          &jacobians,
   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
-  std::vector<Point<spacedim> >                             &cell_normal_vectors,
+  std::vector<Point<spacedim> >                             &normal_vectors,
   CellSimilarity::Similarity                           &cell_similarity) const
 {
   // ensure that the following cast is really correct:
@@ -767,15 +764,15 @@ MappingQ1<dim,spacedim>::fill_fe_values (
                                   // columns of the Jacobian and
                                   // store the cell normal vectors in
                                   // the case <2,3>
-  if (update_flags & (update_cell_normal_vectors
+  if (update_flags & (update_normal_vectors
                      | update_JxW_values))
     {      
       Assert (JxW_values.size() == n_q_points,
               ExcDimensionMismatch(JxW_values.size(), n_q_points));
 
-      Assert( !(update_flags & update_cell_normal_vectors ) ||
-             (cell_normal_vectors.size() == n_q_points),
-            ExcDimensionMismatch(cell_normal_vectors.size(), n_q_points));
+      Assert( !(update_flags & update_normal_vectors ) ||
+             (normal_vectors.size() == n_q_points),
+            ExcDimensionMismatch(normal_vectors.size(), n_q_points));
 
       if (cell_similarity != CellSimilarity::translation)
        for (unsigned int point=0; point<n_q_points; ++point)
@@ -790,11 +787,11 @@ 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_cell_normal_vectors) 
+               if(update_flags & update_normal_vectors) 
                {
-                 cell_normal_vectors[point][0]=-data.contravariant[point][0][1]
+                 normal_vectors[point][0]=-data.contravariant[point][0][1]
                      /data.contravariant[point][0].norm();
-                 cell_normal_vectors[point][1]=data.contravariant[point][0][0]
+                 normal_vectors[point][1]=data.contravariant[point][0][0]
                      /data.contravariant[point][0].norm();
                }
              }
@@ -812,8 +809,8 @@ 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_cell_normal_vectors) 
-                 cell_normal_vectors[point]=data.contravariant[point][2];
+               if(update_flags & update_normal_vectors) 
+                 normal_vectors[point]=data.contravariant[point][2];
              }
            
          }
index fa37567b27980ef15330aaebc03bf8fe68335100..1f85b2c50c7487012132e105055c0f2189055ab7 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2009 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -126,7 +126,7 @@ MappingQ1Eulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
   std::vector<Tensor<2,spacedim> >                          &jacobians,
   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
-  std::vector<Point<spacedim> >                             &cell_normal_vectors,
+  std::vector<Point<spacedim> >                             &normal_vectors,
   CellSimilarity::Similarity                           &cell_similarity) const
 {
                                   // disable any previously detected
@@ -136,7 +136,7 @@ MappingQ1Eulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
   MappingQ1<dim,spacedim>::fill_fe_values (cell, q, mapping_data,
                                           quadrature_points, JxW_values, jacobians,
                                           jacobian_grads, inverse_jacobians,
-                                          cell_normal_vectors, cell_similarity);
+                                          normal_vectors, cell_similarity);
 }
 
 
index e810fb4c4153ee8290d5897bc17813b939bd9a6b..4c85e5eb4a47a347f9ecc8fd682129db76bb3e16 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -193,7 +193,7 @@ MappingQEulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
   std::vector<Tensor<2,spacedim> >                          &jacobians,
   std::vector<Tensor<3,spacedim> >                          &jacobian_grads,
   std::vector<Tensor<2,spacedim> >                          &inverse_jacobians,
-  std::vector<Point<spacedim> >                             &cell_normal_vectors,
+  std::vector<Point<spacedim> >                             &normal_vectors,
   CellSimilarity::Similarity                           &cell_similarity) const
 {
                                   // disable any previously detected
@@ -203,7 +203,7 @@ MappingQEulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
   MappingQ<dim,spacedim>::fill_fe_values (cell, q, mapping_data,
                                          quadrature_points, JxW_values, jacobians,
                                          jacobian_grads, inverse_jacobians,
-                                         cell_normal_vectors, cell_similarity);
+                                         normal_vectors, cell_similarity);
 }
 
 

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.