]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Impemlement get_normals_at_vertices fpr HyperShell and HalfHyperShell
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 29 Mar 2001 09:22:12 +0000 (09:22 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 29 Mar 2001 09:22:12 +0000 (09:22 +0000)
as well.

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

deal.II/deal.II/include/grid/tria_boundary_lib.h
deal.II/deal.II/source/grid/tria_boundary_lib.cc

index 92cb63666617753d1517d39a9b4df50356d0b45c..5ea463b7ad30490d9c07f3590d4d4725817f56ab 100644 (file)
@@ -248,7 +248,21 @@ class HyperShellBoundary : public StraightBoundary<dim>
                                      */
     virtual Point<dim>
     get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
-    
+
+                                    /**
+                                     * Compute the normals to the
+                                     * boundary at the vertices of
+                                     * the given face.
+                                     *
+                                     * Refer to the general
+                                     * documentation of this class
+                                     * and the documentation of the
+                                     * base class.
+                                     */
+    virtual void
+    get_normals_at_vertices (const typename Triangulation<dim>::face_iterator &face,
+                            typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const;
+
   private:
                                     /**
                                      * Store the center of the spheres.
@@ -290,6 +304,20 @@ class HalfHyperShellBoundary : public HyperShellBoundary<dim>
                                      */
     virtual Point<dim>
     get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
+
+                                    /**
+                                     * Compute the normals to the
+                                     * boundary at the vertices of
+                                     * the given face.
+                                     *
+                                     * Refer to the general
+                                     * documentation of this class
+                                     * and the documentation of the
+                                     * base class.
+                                     */
+    virtual void
+    get_normals_at_vertices (const typename Triangulation<dim>::face_iterator &face,
+                            typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const;
     
   private:
                                     /**
index db57d1d01a09ab22aaa914d34ad6d9903e408540..ac081480ed4f0dd82a4c778592b20726560a8233 100644 (file)
@@ -469,6 +469,37 @@ get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) c
 };
 
 
+
+#if deal_II_dimension == 1
+
+template <>
+void
+HyperShellBoundary<1>::
+get_normals_at_vertices (const Triangulation<1>::face_iterator &,
+                        Boundary<1>::FaceVertexNormals &) const
+{
+  Assert (false, Boundary<1>::ExcFunctionNotUseful(1));
+};
+
+#endif
+
+
+template <int dim>
+void
+HyperShellBoundary<dim>::
+get_normals_at_vertices (const typename Triangulation<dim>::face_iterator &face,
+                        typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const
+{
+                                  // this is equivalent to the case
+                                  // in the hyperball boundary. note
+                                  // that we need not normalize nor
+                                  // direct the normal
+  for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_face; ++vertex)
+    face_vertex_normals[vertex] = face->vertex(vertex)-center;
+};
+
+
+
 /* ---------------------------------------------------------------------- */
 
 
@@ -541,6 +572,35 @@ get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) c
 
 
 
+#if deal_II_dimension == 1
+
+template <>
+void
+HalfHyperShellBoundary<1>::
+get_normals_at_vertices (const Triangulation<1>::face_iterator &,
+                        Boundary<1>::FaceVertexNormals &) const
+{
+  Assert (false, Boundary<1>::ExcFunctionNotUseful(1));
+};
+
+#endif
+
+
+
+template <int dim>
+void
+HalfHyperShellBoundary<dim>::
+get_normals_at_vertices (const typename Triangulation<dim>::face_iterator &face,
+                        typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const
+{
+  if (face->center()(0) == center(0))
+    StraightBoundary<dim>::get_normals_at_vertices (face, face_vertex_normals);
+  else
+    HyperShellBoundary<dim>::get_normals_at_vertices (face, face_vertex_normals);
+};
+
+
+
 // explicit instantiations
 template class HyperBallBoundary<deal_II_dimension>;
 template class HalfHyperBallBoundary<deal_II_dimension>;

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.