Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
tria_boundary_lib.h
1 // ---------------------------------------------------------------------
2 // @f$Id: tria_boundary_lib.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 1999 - 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_lib_h
18 #define __deal2__tria_boundary_lib_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/grid/tria_boundary.h>
23 
25 
50 template <int dim, int spacedim = dim>
51 class CylinderBoundary : public StraightBoundary<dim,spacedim>
52 {
53 public:
59  CylinderBoundary (const double radius = 1.0,
60  const unsigned int axis = 0);
61 
74  CylinderBoundary (const double radius,
77 
83  virtual Point<spacedim>
84  get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
85 
91  virtual Point<spacedim>
92  get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
93 
103  virtual void
104  get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
105  std::vector<Point<spacedim> > &points) const;
106 
116  virtual void
117  get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
118  std::vector<Point<spacedim> > &points) const;
119 
130  virtual void
132  typename Boundary<dim,spacedim>::FaceVertexNormals &face_vertex_normals) const;
133 
137  double get_radius () const;
138 
145  DeclException0 (ExcRadiusNotSet);
146 
147 
148 protected:
152  const double radius;
153 
158 
163 
164 private:
165 
179  std::vector<Point<spacedim> > &points) const;
180 
186  static Point<spacedim> get_axis_vector (const unsigned int axis);
187 };
188 
189 
190 
224 template <int dim>
225 class ConeBoundary : public StraightBoundary<dim>
226 {
227 public:
241  ConeBoundary (const double radius_0,
242  const double radius_1,
243  const Point<dim> x_0,
244  const Point<dim> x_1);
245 
251  double get_radius (const Point<dim> x) const;
252 
259  virtual
260  Point<dim>
261  get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const;
262 
269  virtual
270  Point<dim>
271  get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
272 
282  virtual
283  void
284  get_intermediate_points_on_line (const typename Triangulation<dim>::line_iterator &line,
285  std::vector<Point<dim> > &points) const;
286 
297  virtual
298  void
299  get_intermediate_points_on_quad (const typename Triangulation<dim>::quad_iterator &quad,
300  std::vector<Point<dim> > &points) const;
301 
312  virtual
313  void
315  typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const;
316 
317 protected:
322  const double radius_0;
323 
328  const double radius_1;
329 
334 
339 
340 private:
353  void
355  const Point<dim> &p1,
356  std::vector<Point<dim> > &points) const;
357 };
358 
359 
360 
378 template <int dim, int spacedim=dim>
379 class HyperBallBoundary : public StraightBoundary<dim,spacedim>
380 {
381 public:
386  const double radius = 1.0);
387 
393  virtual
395  get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
396 
402  virtual
404  get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
405 
415  virtual
416  void
417  get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
418  std::vector<Point<spacedim> > &points) const;
419 
429  virtual
430  void
431  get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
432  std::vector<Point<spacedim> > &points) const;
433 
443  virtual
446  const Point<spacedim> &p) const;
447 
458  virtual
459  void
461  typename Boundary<dim,spacedim>::FaceVertexNormals &face_vertex_normals) const;
462 
467  get_center () const;
468 
472  double
473  get_radius () const;
474 
481  DeclException0 (ExcRadiusNotSet);
482 
483 
484 protected:
485 
490 
494  const double radius;
495 
514 
515 private:
516 
530  std::vector<Point<spacedim> > &points) const;
531 };
532 
533 
534 
547 template <int dim>
549 {
550 public:
555  const double radius = 1.0);
556 
562  virtual Point<dim>
563  get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const;
564 
570  virtual Point<dim>
571  get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
572 
582  virtual void
583  get_intermediate_points_on_line (const typename Triangulation<dim>::line_iterator &line,
584  std::vector<Point<dim> > &points) const;
585 
595  virtual void
596  get_intermediate_points_on_quad (const typename Triangulation<dim>::quad_iterator &quad,
597  std::vector<Point<dim> > &points) const;
598 
609  virtual void
611  typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const;
612 };
613 
614 
615 
626 template <int dim>
628 {
629 public:
642 };
643 
644 
645 
658 template <int dim>
660 {
661 public:
675  const double inner_radius = -1,
676  const double outer_radius = -1);
677 
681  virtual Point<dim>
682  get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const;
683 
687  virtual Point<dim>
688  get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
689 
699  virtual void
700  get_intermediate_points_on_line (const typename Triangulation<dim>::line_iterator &line,
701  std::vector<Point<dim> > &points) const;
702 
712  virtual void
713  get_intermediate_points_on_quad (const typename Triangulation<dim>::quad_iterator &quad,
714  std::vector<Point<dim> > &points) const;
715 
726  virtual void
728  typename Boundary<dim>::FaceVertexNormals &face_vertex_normals) const;
729 
730 private:
734  const double inner_radius;
735  const double outer_radius;
736 };
737 
738 
747 template <int dim, int spacedim>
748 class TorusBoundary : public Boundary<dim,spacedim>
749 {
750 public:
755  TorusBoundary (const double R, const double r);
756 
757 //Boundary Refinenment Functions
761  virtual Point<spacedim>
762  get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
763 
767  virtual Point<spacedim>
768  get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
769 
773  virtual void get_intermediate_points_on_line (
774  const typename Triangulation< dim, spacedim >::line_iterator &line,
775  std::vector< Point< spacedim > > &points) const ;
776 
780  virtual void get_intermediate_points_on_quad (
781  const typename Triangulation< dim, spacedim >::quad_iterator &quad,
782  std::vector< Point< spacedim > > &points ) const ;
783 
789  virtual void get_normals_at_vertices (
791  typename Boundary<dim,spacedim>::FaceVertexNormals &face_vertex_normals) const ;
792 
793 private:
794  //Handy functions
804  double get_correct_angle(const double angle,const double x,const double y) const;
805 
811  Point<spacedim> get_real_coord(const Point<dim> &surfP) const ;
812 
818  Point<dim> get_surf_coord(const Point<spacedim> &p) const ;
819 
825  Point<spacedim> get_surf_norm_from_sp(const Point<dim> &surfP) const ;
826 
833 
837  const double R;
838  const double r;
839 };
840 
841 
842 
843 /* -------------- declaration of explicit specializations ------------- */
844 
845 #ifndef DOXYGEN
846 
847 template <>
848 Point<1>
850 get_new_point_on_quad (const Triangulation<1>::quad_iterator &) const;
851 template <>
852 void
854  const Triangulation<1>::line_iterator &,
855  std::vector<Point<1> > &) const;
856 template <>
857 void
859  const Triangulation<3>::quad_iterator &quad,
860  std::vector<Point<3> > &points) const;
861 template <>
862 void
866 template <>
867 Point<1>
869 get_new_point_on_quad (const Triangulation<1>::quad_iterator &) const;
870 template <>
871 void
873 get_intermediate_points_on_quad (const Triangulation<1>::quad_iterator &,
874  std::vector<Point<1> > &) const;
875 template <>
876 void
880 template <>
881 Point<1>
883 get_new_point_on_quad (const Triangulation<1>::quad_iterator &) const;
884 template <>
885 void
887 get_intermediate_points_on_quad (const Triangulation<1>::quad_iterator &,
888  std::vector<Point<1> > &) const;
889 template <>
890 void
894 
895 
896 #endif // DOXYGEN
897 
898 DEAL_II_NAMESPACE_CLOSE
899 
900 #endif
virtual void get_intermediate_points_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line, std::vector< Point< spacedim > > &points) const
virtual void get_intermediate_points_on_line(const typename Triangulation< dim >::line_iterator &line, std::vector< Point< dim > > &points) const
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad, std::vector< Point< spacedim > > &points) 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_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad, std::vector< Point< spacedim > > &points) const
virtual Point< dim > get_new_point_on_line(const typename Triangulation< dim >::line_iterator &line) const
virtual void get_normals_at_vertices(const typename Triangulation< dim >::face_iterator &face, typename Boundary< dim >::FaceVertexNormals &face_vertex_normals) const
virtual Point< dim > get_new_point_on_quad(const typename Triangulation< dim >::quad_iterator &quad) const
virtual Point< spacedim > get_new_point_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad) const
Point< spacedim > get_surf_norm_from_sp(const Point< dim > &surfP) const
const Point< spacedim > point_on_axis
void get_intermediate_points_between_points(const Point< spacedim > &p0, const Point< spacedim > &p1, std::vector< Point< spacedim > > &points) const
const Point< spacedim > direction
const double R
double get_radius(const Point< dim > x) const
const Point< dim > x_1
std::vector< std_cxx1x::shared_ptr< QGaussLobatto< 1 > > > points
virtual void get_normals_at_vertices(const typename Triangulation< dim >::face_iterator &face, typename Boundary< dim >::FaceVertexNormals &face_vertex_normals) const
Point< dim > get_surf_coord(const Point< spacedim > &p) 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 Point< spacedim > get_new_point_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad) const
double get_radius() const
virtual Point< dim > get_new_point_on_line(const typename Triangulation< dim >::line_iterator &line) const
HalfHyperShellBoundary(const Point< dim > &center=Point< dim >(), const double inner_radius=-1, const double outer_radius=-1)
Point< spacedim > get_surf_norm(const Point< spacedim > &p) const
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, typename Boundary< dim, spacedim >::FaceVertexNormals &face_vertex_normals) const
HalfHyperBallBoundary(const Point< dim > p=Point< dim >(), const double radius=1.0)
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim >::quad_iterator &quad, std::vector< Point< dim > > &points) const
virtual Point< spacedim > get_new_point_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line) const
Point< spacedim > get_real_coord(const Point< dim > &surfP) const
double get_correct_angle(const double angle, const double x, const double y) const
virtual Point< dim > get_new_point_on_line(const typename Triangulation< dim >::line_iterator &line) const
virtual void get_intermediate_points_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line, std::vector< Point< spacedim > > &points) const
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim >::quad_iterator &quad, std::vector< Point< dim > > &points) const
virtual Tensor< 1, spacedim > normal_vector(const typename Triangulation< dim, spacedim >::face_iterator &face, const Point< spacedim > &p) const
virtual Point< dim > get_new_point_on_quad(const typename Triangulation< dim >::quad_iterator &quad) const
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim >::quad_iterator &quad, std::vector< Point< dim > > &points) const
void get_intermediate_points_between_points(const Point< dim > &p0, const Point< dim > &p1, std::vector< Point< dim > > &points) const
Point< spacedim > get_center() const
void get_intermediate_points_between_points(const Point< spacedim > &p0, const Point< spacedim > &p1, std::vector< Point< spacedim > > &points) const
ConeBoundary(const double radius_0, const double radius_1, const Point< dim > x_0, const Point< dim > x_1)
const Point< spacedim > center
virtual void get_intermediate_points_on_line(const typename Triangulation< dim >::line_iterator &line, std::vector< Point< dim > > &points) const
const double radius_1
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 >::line_iterator &line, std::vector< Point< dim > > &points) const
virtual void get_intermediate_points_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line, std::vector< Point< spacedim > > &points) const
virtual Point< dim > get_new_point_on_quad(const typename Triangulation< dim >::quad_iterator &quad) const
DeclException0(ExcRadiusNotSet)
virtual void get_intermediate_points_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad, std::vector< Point< spacedim > > &points) const
HyperShellBoundary(const Point< dim > &center=Point< dim >())
TorusBoundary(const double R, const double r)
const Point< dim > x_0
virtual void get_normals_at_vertices(const typename Triangulation< dim >::face_iterator &face, typename Boundary< dim >::FaceVertexNormals &face_vertex_normals) const
HyperBallBoundary(const Point< spacedim > p=Point< spacedim >(), const double radius=1.0)
CylinderBoundary(const double radius=1.0, const unsigned int axis=0)
const double radius_0
virtual Point< spacedim > get_new_point_on_line(const typename Triangulation< dim, spacedim >::line_iterator &line) const
DeclException0(ExcRadiusNotSet)
double get_radius() const
static Point< spacedim > get_axis_vector(const unsigned int axis)