Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
tria_boundary.h
1 // ---------------------------------------------------------------------
2 // @f$Id: tria_boundary.h 30450 2013-08-23 15:48:29Z kronbichler @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__tria_boundary_h
18 #define __deal2__tria_boundary_h
19 
20 
21 /*---------------------------- boundary-function.h ---------------------------*/
22 
23 #include <deal.II/base/config.h>
24 #include <deal.II/base/subscriptor.h>
25 #include <deal.II/base/quadrature_lib.h>
26 #include <deal.II/base/thread_management.h>
27 #include <deal.II/base/point.h>
28 #include <deal.II/grid/tria.h>
29 
31 
32 template <int dim, int space_dim> class Triangulation;
33 
34 
35 
84 template <int dim, int spacedim=dim>
85 class Boundary : public Subscriptor
86 {
87 public:
88 
101 
106  virtual ~Boundary ();
107 
114  virtual
116  get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const = 0;
117 
132  virtual
134  get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
135 
144 
160  virtual
161  void
162  get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
163  std::vector<Point<spacedim> > &points) const;
164 
184  virtual
185  void
186  get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
187  std::vector<Point<spacedim> > &points) const;
188 
195  void
197  std::vector<Point<spacedim> > &points) const;
198 
221  virtual
224  const Point<spacedim> &p) const;
225 
239  virtual
240  void
242  FaceVertexNormals &face_vertex_normals) const;
243 
260  virtual
262  project_to_surface (const typename Triangulation<dim,spacedim>::line_iterator &line,
263  const Point<spacedim> &candidate) const;
264 
273  virtual
275  project_to_surface (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
276  const Point<spacedim> &candidate) const;
277 
286  virtual
288  project_to_surface (const typename Triangulation<dim,spacedim>::hex_iterator &hex,
289  const Point<spacedim> &candidate) const;
290 
291 protected:
300  const std::vector<Point<1> > &
301  get_line_support_points (const unsigned int n_intermediate_points) const;
302 
303 private:
307  mutable std::vector<std_cxx1x::shared_ptr<QGaussLobatto<1> > > points;
308 
313 };
314 
315 
316 
332 template <int dim, int spacedim=dim>
333 class StraightBoundary : public Boundary<dim,spacedim>
334 {
335 public:
339  StraightBoundary ();
340 
347  virtual Point<spacedim>
348  get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
349 
358  virtual
360  get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
361 
369  virtual
370  void
371  get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
372  std::vector<Point<spacedim> > &points) const;
373 
381  virtual
382  void
383  get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
384  std::vector<Point<spacedim> > &points) const;
385 
392  virtual
395  const Point<spacedim> &p) const;
396 
403  virtual
404  void
406  typename Boundary<dim,spacedim>::FaceVertexNormals &face_vertex_normals) const;
407 
422  virtual
424  project_to_surface (const typename Triangulation<dim,spacedim>::line_iterator &line,
425  const Point<spacedim> &candidate) const;
426 
438  virtual
440  project_to_surface (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
441  const Point<spacedim> &candidate) const;
442 
455  virtual
457  project_to_surface (const typename Triangulation<dim,spacedim>::hex_iterator &hex,
458  const Point<spacedim> &candidate) const;
459 };
460 
461 
462 
463 /* -------------- declaration of explicit specializations ------------- */
464 
465 #ifndef DOXYGEN
466 
467 template <>
468 Point<1>
471 
472 template <>
473 void
476  std::vector<Point<1> > &) const;
477 
478 template <>
479 Point<2>
482 
483 template <>
484 void
487  std::vector<Point<2> > &) const;
488 
489 
490 
491 template <>
492 Point<3>
495 
496 template <>
497 void
500  std::vector<Point<3> > &) const;
501 template <>
502 void
506 template <>
507 void
510  Boundary<2,2>::FaceVertexNormals &face_vertex_normals) const;
511 template <>
512 void
515  Boundary<3,3>::FaceVertexNormals &face_vertex_normals) const;
516 
517 template <>
518 Point<3>
520 get_new_point_on_quad (const Triangulation<3,3>::quad_iterator &quad) const;
521 
522 template <>
523 void
525 get_intermediate_points_on_line (const Triangulation<1,1>::line_iterator &,
526  std::vector<Point<1> > &) const;
527 
528 template <>
529 void
531 get_intermediate_points_on_quad (const Triangulation<3,3>::quad_iterator &quad,
532  std::vector<Point<3> > &points) const;
533 
534 template <>
535 Point<3>
537 project_to_surface (const Triangulation<1, 3>::quad_iterator &quad,
538  const Point<3> &y) const;
539 
540 
541 #endif // DOXYGEN
542 
543 DEAL_II_NAMESPACE_CLOSE
544 
545 #endif
virtual Point< spacedim > get_new_point_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line) const
virtual void get_normals_at_vertices(const typename Triangulation< dim, spacedim >::face_iterator &face, FaceVertexNormals &face_vertex_normals) const
virtual Tensor< 1, spacedim > normal_vector(const typename Triangulation< dim, spacedim >::face_iterator &face, const Point< spacedim > &p) const
virtual Point< spacedim > get_new_point_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line) const =0
virtual Point< spacedim > project_to_surface(const typename Triangulation< dim, spacedim >::line_iterator &line, const Point< spacedim > &candidate) const
Tensor< 1, spacedim > FaceVertexNormals[GeometryInfo< dim >::vertices_per_face]
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad, std::vector< Point< spacedim > > &points) const
std::vector< std_cxx1x::shared_ptr< QGaussLobatto< 1 > > > points
Point< spacedim > get_new_point_on_face(const typename Triangulation< dim, spacedim >::face_iterator &face) const
virtual void get_normals_at_vertices(const typename Triangulation< dim, spacedim >::face_iterator &face, typename Boundary< dim, spacedim >::FaceVertexNormals &face_vertex_normals) const
virtual void get_intermediate_points_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line, std::vector< Point< spacedim > > &points) const
virtual Point< spacedim > get_new_point_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad) const
virtual Tensor< 1, spacedim > normal_vector(const typename Triangulation< dim, spacedim >::face_iterator &face, const Point< spacedim > &p) const
void get_intermediate_points_on_face(const typename Triangulation< dim, spacedim >::face_iterator &face, std::vector< Point< spacedim > > &points) const
virtual Point< spacedim > project_to_surface(const typename Triangulation< dim, spacedim >::line_iterator &line, const Point< spacedim > &candidate) const
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad, std::vector< Point< spacedim > > &points) const
Threads::Mutex mutex
const std::vector< Point< 1 > > & get_line_support_points(const unsigned int n_intermediate_points) const
virtual ~Boundary()
virtual Point< spacedim > get_new_point_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad) const
virtual void get_intermediate_points_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line, std::vector< Point< spacedim > > &points) const