Compute normal vectors in no_normal_flux by querying the Boundary object, rather...

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

diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index d8a307474c039818c84c1c4fab5b6c6bd7d7834f..7c28b89c31fee8c6014f8a498a157bf07fdb437b 100644 (file)
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -281,6 +281,23 @@ and DoF handlers embedded in higher dimensional space have been added.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Changed: The function VectorTools::compute_no_normal_flux_constraints
+used to compute its constraints by evaluating the normal vector to the
+surface as described by the mapping, rather than using the normal to
+the surface described by the Boundary object associated with this face.
+(Note that the Mapping computes its approximation by polynomial interpolation
+of the surface described by the Boundary object.) This has now been changed:
+the normal vector is now obtained from the Boundary object directly, at
+points computed by the Mapping.
+<br>
+(Wolfgang Bangerth, 2011/08/25)
+
+<li> New: The Boundary base class now has a function Boundary::normal_vector
+that returns the normal vector to the surface at a given location. Derived
+classes need to implement it, of course, if they want it to be used.
+<br>
+(Wolfgang Bangerth, 2011/08/25)
+
 <li> Fixed: The function VectorTools::compute_no_normal_flux_constraints had
 a problem that led to extremely difficult to pin down bugs when running
 with sufficiently many processors. Basically, the constraints computed

diff --git a/deal.II/include/deal.II/numerics/vectors.h b/deal.II/include/deal.II/numerics/vectors.h
index 3f40890cff3ab86769216e118cc35741b80333a0..3b01ca30a385c8689bb8326f876225597b7ee086 100644 (file)
--- a/deal.II/include/deal.II/numerics/vectors.h
+++ b/deal.II/include/deal.II/numerics/vectors.h
@@ -1267,9 +1267,10 @@ class VectorTools
                                      * function onces with the whole set of
                                      * boundary indicators at once.
                                      *
-                                     * The last argument is denoted to
-                                     * compute the normal vector $\vec n$ at
-                                     * the boundary points.
+                                     * The mapping argument is used to
+                                     * compute the boundary points where the function
+                                     * needs to request the normal vector $\vec n$
+                                     * from the boundary description.
                                      *
                                      * @note When combining adaptively
                                      * refined meshes with hanging node

diff --git a/deal.II/include/deal.II/numerics/vectors.templates.h b/deal.II/include/deal.II/numerics/vectors.templates.h
index 9763b62fc97b19e886a80a03bea30abb0f22cd88..78a23c855a89bf8ea3f7889694587c21d194d5bd 100644 (file)
--- a/deal.II/include/deal.II/numerics/vectors.templates.h
+++ b/deal.II/include/deal.II/numerics/vectors.templates.h
@@ -30,6 +30,7 @@
 #include <deal.II/lac/constraint_matrix.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_boundary.h>
 #include <deal.II/grid/grid_tools.h>
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/dofs/dof_handler.h>
@@ -4204,7 +4205,7 @@ compute_no_normal_flux_constraints (const DH<dim,spacedim>         &dof_handler,
   hp::FEFaceValues<dim,spacedim> x_fe_face_values (mapping_collection,
                                                   fe_collection,
                                                   face_quadrature_collection,
-                                                  update_normal_vectors);
+                                                  update_q_points);
 
                                   // have a map that stores normal
                                   // vectors for each vector-dof
@@ -4345,7 +4346,10 @@ compute_no_normal_flux_constraints (const DH<dim,spacedim>         &dof_handler,
                                                   // length of the vector
                                                   // back to one, just to be
                                                   // on the safe side
-                 Point<dim> normal_vector = fe_values.normal_vector(i);
+                 Point<dim> normal_vector 
+                 = (cell->face(face_no)->get_boundary()
+                    .normal_vector (cell->face(face_no),
+                                    fe_values.quadrature_point(i)));
                  Assert (std::fabs(normal_vector.norm() - 1) < 1e-14,
                          ExcInternalError());
                  for (unsigned int d=0; d<dim; ++d)