17 #ifndef __deal2__grid_tools_H 18 #define __deal2__grid_tools_H 21 #include <deal.II/base/config.h> 22 #include <deal.II/fe/mapping.h> 23 #include <deal.II/grid/tria.h> 24 #include <deal.II/grid/tria_accessor.h> 25 #include <deal.II/grid/tria_iterator.h> 26 #include <deal.II/fe/mapping_q1.h> 35 template <
int,
int>
class Mapping;
39 template <
int,
int>
class MappingCollection;
67 template <
int dim,
int spacedim>
80 template <
int dim,
int spacedim>
136 template <
int dim,
int spacedim>
164 template <
int dim,
int spacedim>
168 std::vector<unsigned int> &considered_vertices,
169 const double tol=1e-12);
201 template <
int dim,
typename Transformation,
int spacedim>
202 void transform (
const Transformation &transformation,
216 template <
int dim,
int spacedim>
236 void rotate (
const double angle,
285 template <
int dim,
int spacedim>
286 void scale (
const double scaling_factor,
307 template <
int dim,
int spacedim>
310 const bool keep_boundary=
true);
324 template <
int dim,
template <
int,
int>
class Container,
int spacedim>
355 template<
int dim,
template <
int,
int>
class Container,
int spacedim>
356 std::vector<typename Container<dim,spacedim>::active_cell_iterator>
358 const unsigned int vertex_index);
388 template <
int dim,
template <
int,
int>
class Container,
int spacedim>
389 typename Container<dim,spacedim>::active_cell_iterator
436 template <
int dim,
template<
int,
int>
class Container,
int spacedim>
437 std::pair<typename Container<dim,spacedim>::active_cell_iterator,
Point<dim> >
439 const Container<dim,spacedim> &container,
462 template <
int dim,
int spacedim>
463 std::pair<typename hp::DoFHandler<dim,spacedim>::active_cell_iterator,
Point<dim> >
491 template <
class Container>
492 std::vector<typename Container::active_cell_iterator>
500 template <
class Container>
503 std::vector<typename Container::active_cell_iterator> &active_neighbors);
516 template <
int dim,
int spacedim>
541 template <
int dim,
int spacedim>
622 template <
int dim,
int spacedim>
624 partition_triangulation (
const unsigned int n_partitions,
643 template <
int dim,
int spacedim>
646 std::vector<types::subdomain_id> &subdomain);
670 template <
int dim,
int spacedim>
731 template <
typename Container>
732 std::list<std::pair<
typename Container::cell_iterator,
733 typename Container::cell_iterator> >
735 const Container &mesh_2);
751 template <
int dim,
int spacedim>
769 template <
typename Container>
771 have_same_coarse_mesh (
const Container &mesh_1,
772 const Container &mesh_2);
780 template <
int dim,
int spacedim>
788 template <
int dim,
int spacedim>
841 template <
int dim,
int spacedim>
872 template <
int dim,
int spacedim>
963 template <
template <
int,
int>
class Container,
int dim,
int spacedim>
964 std::map<
typename Container<dim-1,spacedim>::cell_iterator,
965 typename Container<dim,spacedim>::face_iterator>
967 Container<dim-1,spacedim> &surface_mesh,
968 const std::set<types::boundary_id> &boundary_ids
969 = std::set<types::boundary_id>());
976 template<
typename CellIterator>
979 CellIterator cell[2];
980 unsigned int face_idx[2];
981 std::bitset<3> orientation;
1047 template<
typename FaceIterator>
1050 const FaceIterator &face1,
1051 const FaceIterator &face2,
1052 const int direction,
1053 const ::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset);
1059 template<
typename FaceIterator>
1061 orthogonal_equality (
const FaceIterator &face1,
1062 const FaceIterator &face2,
1063 const int direction,
1064 const ::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset);
1105 template<
typename CONTAINER>
1108 (
const CONTAINER &container,
1111 const int direction,
1136 template<
typename CONTAINER>
1138 collect_periodic_faces
1139 (
const CONTAINER &container,
1141 const int direction,
1151 <<
"The number of partitions you gave is " << arg1
1152 <<
", but must be greater than zero.");
1158 <<
"The subdomain id " << arg1
1159 <<
" has no cells associated with it.");
1169 <<
"The scaling factor must be positive, but is " << arg1);
1176 <<
"The point <" << arg1
1177 <<
"> could not be found inside any of the " 1178 <<
"coarse grid cells.");
1185 <<
"The point <" << arg1
1186 <<
"> could not be found inside any of the " 1187 <<
"subcells of a coarse grid cell.");
1191 <<
"The given vertex " << arg1
1192 <<
" is not used in the given triangulation");
1203 template <
int dim,
typename Predicate,
int spacedim>
1204 void transform (
const Predicate &predicate,
1210 static_cast<unsigned int>(triangulation.
begin_active()->level()+1),
1211 ExcMessage (
"Not all cells of this triangulation are at the same " 1212 "refinement level, as is required for this function."));
1214 std::vector<bool> treated_vertices (triangulation.
n_vertices(),
1225 endc = triangulation.
end ();
1226 for (; cell!=endc; ++cell)
1227 for (
unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
1228 if (treated_vertices[cell->vertex_index(v)] ==
false)
1231 cell->vertex(v) = predicate(cell->vertex(v));
1233 treated_vertices[cell->vertex_index(v)] =
true;
1240 std::vector<typename DH::active_cell_iterator>
1243 std::vector<typename DH::active_cell_iterator> child_cells;
1245 if (cell->has_children())
1247 for (
unsigned int child=0;
1248 child<cell->n_children(); ++child)
1249 if (cell->child (child)->has_children())
1251 const std::vector<typename DH::active_cell_iterator>
1252 children = get_active_child_cells<DH> (cell->child(child));
1253 child_cells.insert (child_cells.end(),
1254 children.begin(), children.end());
1257 child_cells.push_back (cell->child(child));
1265 template <
class Container>
1268 std::vector<typename Container::active_cell_iterator> &active_neighbors)
1270 active_neighbors.clear ();
1271 for (
unsigned int n=0; n<GeometryInfo<Container::dimension>::faces_per_cell; ++n)
1272 if (! cell->at_boundary(n))
1274 if (Container::dimension == 1)
1281 typename Container::cell_iterator
1282 neighbor_child = cell->neighbor(n);
1283 if (!neighbor_child->active())
1285 while (neighbor_child->has_children())
1286 neighbor_child = neighbor_child->child (n==0 ? 1 : 0);
1288 Assert (neighbor_child->neighbor(n==0 ? 1 : 0)==cell,
1291 active_neighbors.push_back (neighbor_child);
1295 if (cell->face(n)->has_children())
1299 for (
unsigned int c=0; c<cell->face(n)->number_of_children(); ++c)
1300 active_neighbors.push_back (cell->neighbor_child_on_subface(n,c));
1306 active_neighbors.push_back(cell->neighbor(n));
1319 cell_measure<3>(
const std::vector<Point<3> > &all_vertices,
1324 cell_measure<2>(
const std::vector<Point<2> > &all_vertices,
1330 DEAL_II_NAMESPACE_CLOSE
active_cell_iterator begin_active(const unsigned int level=0) const
::ExceptionBase & ExcMessage(std::string arg1)
unsigned int n_levels() const
cell_iterator end() const
#define DeclException1(Exception1, type1, outsequence)
#define Assert(cond, exc)
#define DeclException0(Exception0)
unsigned int subdomain_id
unsigned char boundary_id
::ExceptionBase & ExcInternalError()
unsigned int n_vertices() const