16 #ifndef dealii__tria_manifold_h 17 #define dealii__tria_manifold_h 22 #include <deal.II/base/config.h> 23 #include <deal.II/base/subscriptor.h> 24 #include <deal.II/base/quadrature_lib.h> 25 #include <deal.II/base/thread_management.h> 26 #include <deal.II/base/point.h> 27 #include <deal.II/base/derivative_form.h> 28 #include <deal.II/grid/tria.h> 30 DEAL_II_NAMESPACE_OPEN
33 template <
int,
typename>
class Table;
82 template <
typename MeshIteratorType>
85 const bool with_laplace =
false) DEAL_II_DEPRECATED;
128 template <
typename MeshIteratorType>
129 std::pair<std::vector<Point<MeshIteratorType::AccessorType::space_dimension> >,
130 std::vector<double> >
132 const bool with_laplace =
false);
299 template <
int dim,
int spacedim=dim>
305 #ifdef DEAL_II_WITH_CXX11 306 static_assert (dim<=spacedim,
307 "The dimension <dim> of a Manifold must be less than or " 308 "equal to the space dimension <spacedim> in which it lives.");
358 const double w)
const;
399 get_new_point (
const std::vector<
Point<spacedim> > &surrounding_points,
400 const std::vector<double> &weights)
const;
427 add_new_points (
const std::vector<
Point<spacedim> > &surrounding_points,
463 get_new_point_on_line (
const typename Triangulation<dim,spacedim>::line_iterator &line)
const;
484 get_new_point_on_quad (
const typename Triangulation<dim,spacedim>::quad_iterator &quad)
const;
506 get_new_point_on_hex (
const typename Triangulation<dim,spacedim>::hex_iterator &hex)
const;
656 FaceVertexNormals &face_vertex_normals)
const;
673 template <
int dim,
int spacedim=dim>
705 const double tolerance=1e-10);
756 const std::vector<double> &weights)
const;
773 add_new_points (
const std::vector<
Point<spacedim> > &surrounding_points,
837 <<
"The component number " << arg1 <<
" of the point [ " << arg2
838 <<
" ] is not in the interval [ 0, " << arg3 <<
"), bailing out.");
936 template <
int dim,
int spacedim=dim,
int chartdim=dim>
941 #ifdef DEAL_II_WITH_CXX11 942 static_assert (dim<=spacedim,
943 "The dimension <dim> of a ChartManifold must be less than or " 944 "equal to the space dimension <spacedim> in which it lives.");
985 const std::vector<double> &weights)
const;
1013 add_new_points (
const std::vector<
Point<spacedim> > &surrounding_points,
1187 template <
typename MeshIteratorType>
1190 const bool with_laplace)
1192 const std::pair<std::vector<Point<MeshIteratorType::AccessorType::space_dimension> >,
1196 points_and_weights.second);
1199 template <
typename MeshIteratorType>
1200 std::pair<std::vector<Point<MeshIteratorType::AccessorType::space_dimension> >,
1201 std::vector<double> >
1203 const bool with_laplace)
1205 const int spacedim = MeshIteratorType::AccessorType::space_dimension;
1206 const int dim = MeshIteratorType::AccessorType::structure_dimension;
1208 std::pair<std::vector<Point<spacedim> >,
1209 std::vector<double> > points_weights;
1219 points_weights.first.resize(2);
1220 points_weights.second.resize(2);
1221 points_weights.first[0] = iterator->vertex(0);
1222 points_weights.second[0] = .5;
1223 points_weights.first[1] = iterator->vertex(1);
1224 points_weights.second[1] = .5;
1227 points_weights.first.resize(8);
1228 points_weights.second.resize(8);
1230 for (
unsigned int i=0; i<4; ++i)
1232 points_weights.first[i] = iterator->vertex(i);
1233 points_weights.first[4+i] = ( iterator->line(i)->has_children() ?
1234 iterator->line(i)->child(0)->vertex(1) :
1235 iterator->line(i)->get_manifold().get_new_point_on_line(iterator->line(i)) );
1240 std::fill(points_weights.second.begin(), points_weights.second.begin()+4, 1.0/16.0);
1241 std::fill(points_weights.second.begin()+4, points_weights.second.end(), 3.0/16.0);
1244 std::fill(points_weights.second.begin(), points_weights.second.end(), 1.0/8.0);
1250 const unsigned int np =
1254 points_weights.first.resize(np);
1255 points_weights.second.resize(np);
1256 std::vector<Point<3> > *sp3 =
reinterpret_cast<std::vector<Point<3>
> *>(&points_weights.first);
1264 for (
unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i, ++j)
1266 (*sp3)[j] = hex->vertex(i);
1267 points_weights.second[j] = 1.0/128.0;
1269 for (
unsigned int i=0; i<GeometryInfo<dim>::lines_per_cell; ++i, ++j)
1271 (*sp3)[j] = (hex->line(i)->has_children() ?
1272 hex->line(i)->child(0)->vertex(1) :
1273 hex->line(i)->get_manifold().get_new_point_on_line(hex->line(i)));
1274 points_weights.second[j] = 7.0/192.0;
1276 for (
unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i, ++j)
1278 (*sp3)[j] = (hex->quad(i)->has_children() ?
1279 hex->quad(i)->isotropic_child(0)->vertex(3) :
1280 hex->quad(i)->get_manifold().get_new_point_on_quad(hex->quad(i)));
1281 points_weights.second[j] = 1.0/12.0;
1285 if (with_laplace ==
false)
1286 std::fill(points_weights.second.begin(), points_weights.second.end(), 1.0/np);
1293 return points_weights;
1299 DEAL_II_NAMESPACE_CLOSE
const Tensor< 1, spacedim > periodicity
const FlatManifold< chartdim, chartdim > sub_manifold
std::pair< std::vector< Point< MeshIteratorType::AccessorType::space_dimension > >, std::vector< double > > get_default_points_and_weights(const MeshIteratorType &iterator, const bool with_laplace=false)
Quadrature< MeshIteratorType::AccessorType::space_dimension > get_default_quadrature(const MeshIteratorType &iterator, const bool with_laplace=false) 1
Point< spacedim > get_new_point_on_face(const typename Triangulation< dim, spacedim >::face_iterator &face) const
#define Assert(cond, exc)
virtual Point< spacedim > get_new_point_on_quad(const typename Triangulation< dim, spacedim >::quad_iterator &quad) const
virtual Point< spacedim > get_new_point_on_hex(const typename Triangulation< dim, spacedim >::hex_iterator &hex) const
#define DeclException3(Exception3, type1, type2, type3, outsequence)
static::ExceptionBase & ExcInternalError()