--- /dev/null
+New: GridGenerator_merge_triangulations gained an option to copy manifold ids.
+<br>
+(Daniel Arndt, 2018/07/08)
* refined cells.
*
* @note The function copies the material ids of the cells of the two input
- * triangulations into the output triangulation but it currently makes no
- * attempt to do the same for boundary ids. In other words, if the two
- * coarse meshes have anything but the default boundary indicators, then you
- * will currently have to set boundary indicators again by hand in the
- * output triangulation.
+ * triangulations into the output triangulation. If @p copy_manifold_ids is
+ * activated, manifold ids will be copied. Currently, boundary_ids are never
+ * copied. In other words, if the two coarse meshes have anything but the
+ * default boundary indicators, then you will currently have to set boundary
+ * indicators again by hand in the output triangulation.
*
* @note Unlike most GridGenerator functions, this function does not attach
* any manifolds to @p result, nor does it set any manifold ids.
merge_triangulations(const Triangulation<dim, spacedim> &triangulation_1,
const Triangulation<dim, spacedim> &triangulation_2,
Triangulation<dim, spacedim> & result,
- const double duplicated_vertex_tolerance = 1.0e-12);
+ const double duplicated_vertex_tolerance = 1.0e-12,
+ const bool copy_manifold_ids = false);
/**
* Given the two triangulations specified as the first two arguments, create
merge_triangulations(const Triangulation<dim, spacedim> &triangulation_1,
const Triangulation<dim, spacedim> &triangulation_2,
Triangulation<dim, spacedim> & result,
- const double duplicated_vertex_tolerance)
+ const double duplicated_vertex_tolerance,
+ const bool copy_manifold_ids)
{
Assert(triangulation_1.n_levels() == 1,
ExcMessage("The input triangulations must be coarse meshes."));
// now form the union of the set of cells
std::vector<CellData<dim>> cells;
cells.reserve(triangulation_1.n_cells() + triangulation_2.n_cells());
- for (typename Triangulation<dim, spacedim>::cell_iterator cell =
- triangulation_1.begin();
- cell != triangulation_1.end();
- ++cell)
+ for (const auto &cell : triangulation_1.cell_iterators())
{
CellData<dim> this_cell;
for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
this_cell.vertices[v] = cell->vertex_index(v);
this_cell.material_id = cell->material_id();
- cells.push_back(this_cell);
+ this_cell.manifold_id = cell->manifold_id();
+ cells.push_back(std::move(this_cell));
}
// now do the same for the other other mesh. note that we have to
// translate the vertex indices
- for (typename Triangulation<dim, spacedim>::cell_iterator cell =
- triangulation_2.begin();
- cell != triangulation_2.end();
- ++cell)
+ for (const auto &cell : triangulation_2.cell_iterators())
{
CellData<dim> this_cell;
for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
this_cell.vertices[v] =
cell->vertex_index(v) + triangulation_1.n_vertices();
this_cell.material_id = cell->material_id();
- cells.push_back(this_cell);
+ this_cell.manifold_id = cell->manifold_id();
+ cells.push_back(std::move(this_cell));
+ }
+
+ // Now for the SubCellData
+ SubCellData subcell_data;
+ if (copy_manifold_ids)
+ {
+ if (dim == 2)
+ {
+ subcell_data.boundary_lines.reserve(triangulation_1.n_lines() +
+ triangulation_2.n_lines());
+
+ auto copy_line_manifold_ids =
+ [&](const Triangulation<dim, spacedim> &tria,
+ const unsigned int offset) {
+ for (typename Triangulation<dim, spacedim>::line_iterator line =
+ tria.begin_face();
+ line != tria.end_face();
+ ++line)
+ if (line->manifold_id() != numbers::flat_manifold_id)
+ {
+ CellData<1> boundary_line;
+ boundary_line.vertices[0] =
+ line->vertex_index(0) + offset;
+ boundary_line.vertices[1] =
+ line->vertex_index(1) + offset;
+ boundary_line.manifold_id = line->manifold_id();
+ subcell_data.boundary_lines.push_back(
+ std::move(boundary_line));
+ }
+ };
+
+ copy_line_manifold_ids(triangulation_1, 0);
+ copy_line_manifold_ids(triangulation_2,
+ triangulation_1.n_vertices());
+ }
+
+ if (dim == 3)
+ {
+ subcell_data.boundary_quads.reserve(triangulation_1.n_quads() +
+ triangulation_2.n_quads());
+ // we can't do better here than to loop over all the lines bounding
+ // a face. For regular meshes an (interior) line in 3D is part of
+ // four cells. So this should be an appropriate guess.
+ subcell_data.boundary_lines.reserve(triangulation_1.n_cells() * 4 +
+ triangulation_2.n_cells() * 4);
+
+ auto copy_face_and_line_manifold_ids =
+ [&](const Triangulation<dim, spacedim> &tria,
+ const unsigned int offset) {
+ for (typename Triangulation<dim, spacedim>::face_iterator face =
+ tria.begin_face();
+ face != tria.end_face();
+ ++face)
+ if (face->manifold_id() != numbers::flat_manifold_id)
+ {
+ CellData<2> boundary_quad;
+ boundary_quad.vertices[0] =
+ face->vertex_index(0) + offset;
+ boundary_quad.vertices[1] =
+ face->vertex_index(1) + offset;
+ boundary_quad.vertices[2] =
+ face->vertex_index(2) + offset;
+ boundary_quad.vertices[3] =
+ face->vertex_index(3) + offset;
+
+ boundary_quad.manifold_id = face->manifold_id();
+ subcell_data.boundary_quads.push_back(
+ std::move(boundary_quad));
+ }
+ for (const auto &cell : tria.cell_iterators())
+ for (unsigned int l = 0; l < 12; ++l)
+ if (cell->line(l)->manifold_id() !=
+ numbers::flat_manifold_id)
+ {
+ CellData<1> boundary_line;
+ boundary_line.vertices[0] =
+ cell->line(l)->vertex_index(0) + offset;
+ boundary_line.vertices[1] =
+ cell->line(l)->vertex_index(1) + offset;
+ boundary_line.manifold_id =
+ cell->line(l)->manifold_id();
+ subcell_data.boundary_lines.push_back(
+ std::move(boundary_line));
+ }
+ };
+
+ copy_face_and_line_manifold_ids(triangulation_1, 0);
+ copy_face_and_line_manifold_ids(triangulation_2,
+ triangulation_1.n_vertices());
+ }
}
// throw out duplicated vertices from the two meshes, reorder vertices as
// necessary and create the triangulation
- SubCellData subcell_data;
std::vector<unsigned int> considered_vertices;
GridTools::delete_duplicated_vertices(vertices,
cells,
const Triangulation<deal_II_dimension, deal_II_space_dimension>
&triangulation_2,
Triangulation<deal_II_dimension, deal_II_space_dimension> &result,
- const double duplicated_vertex_tolerance);
+ const double duplicated_vertex_tolerance,
+ const bool copy_manifold_ids);
template void
create_union_triangulation(
}
}
- // now we got a renumbering list. simply
- // renumber all vertices (non-duplicate
- // vertices get renumbered to themselves, so
- // nothing bad happens). after that, the
- // duplicate vertices will be unused, so call
- // delete_unused_vertices() to do that part
- // of the job.
- for (unsigned int c = 0; c < cells.size(); ++c)
- for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
- cells[c].vertices[v] = new_vertex_numbers[cells[c].vertices[v]];
+ // now we got a renumbering list. simply renumber all vertices
+ // (non-duplicate vertices get renumbered to themselves, so nothing bad
+ // happens). after that, the duplicate vertices will be unused, so call
+ // delete_unused_vertices() to do that part of the job.
+ for (auto &cell : cells)
+ for (auto &vertex_index : cell.vertices)
+ vertex_index = new_vertex_numbers[vertex_index];
+ for (auto &quad : subcelldata.boundary_quads)
+ for (auto &vertex_index : quad.vertices)
+ vertex_index = new_vertex_numbers[vertex_index];
+ for (auto &line : subcelldata.boundary_lines)
+ for (auto &vertex_index : line.vertices)
+ vertex_index = new_vertex_numbers[vertex_index];
delete_unused_vertices(vertices, cells, subcelldata);
}
ExcLineInexistant(line_vertices.first,
line_vertices.second));
}
- // Assert that only exterior
- // lines are given a boundary
- // indicator
- AssertThrow(
- line->at_boundary(),
- ExcInteriorLineCantBeBoundary(line->vertex_index(0),
- line->vertex_index(1),
- boundary_line->boundary_id));
-
- // and make sure that we don't
- // attempt to reset the
- // boundary indicator to a
- // different than the
- // previously set value
- if (line->boundary_id() != 0)
- AssertThrow(line->boundary_id() == boundary_line->boundary_id,
- ExcMessage(
- "Duplicate boundary lines are only allowed "
- "if they carry the same boundary indicator."));
+ // Only exterior lines can be given a boundary indicator
+ if (line->at_boundary())
+ {
+ // make sure that we don't attempt to reset the boundary
+ // indicator to a different than the previously set value
+ if (line->boundary_id() != 0)
+ AssertThrow(line->boundary_id() == boundary_line->boundary_id,
+ ExcMessage(
+ "Duplicate boundary lines are only allowed "
+ "if they carry the same boundary indicator."));
- line->set_boundary_id_internal(boundary_line->boundary_id);
+ line->set_boundary_id_internal(boundary_line->boundary_id);
+ }
// Set manifold id if given
+ if (line->manifold_id() != numbers::flat_manifold_id)
+ AssertThrow(line->manifold_id() == boundary_line->manifold_id,
+ ExcMessage(
+ "Duplicate boundary lines are only allowed "
+ "if they carry the same manifold indicator."));
line->set_manifold_id(boundary_line->manifold_id);
}
// check whether this face is
// really an exterior one
- AssertThrow(
- quad->at_boundary(),
- ExcInteriorQuadCantBeBoundary(quad->vertex_index(0),
- quad->vertex_index(1),
- quad->vertex_index(2),
- quad->vertex_index(3),
- boundary_quad->boundary_id));
-
- // and make sure that we don't
- // attempt to reset the
- // boundary indicator to a
- // different than the
- // previously set value
- if (quad->boundary_id() != 0)
- AssertThrow(quad->boundary_id() == boundary_quad->boundary_id,
+ if (quad->at_boundary())
+ {
+ // and make sure that we don't attempt to reset the boundary
+ // indicator to a different than the previously set value
+ if (quad->boundary_id() != 0)
+ AssertThrow(quad->boundary_id() == boundary_quad->boundary_id,
+ ExcMessage(
+ "Duplicate boundary quads are only allowed "
+ "if they carry the same boundary indicator."));
+
+ quad->set_boundary_id_internal(boundary_quad->boundary_id);
+ }
+ // Set manifold id if given
+ if (quad->manifold_id() != numbers::flat_manifold_id)
+ AssertThrow(quad->manifold_id() == boundary_quad->manifold_id,
ExcMessage(
"Duplicate boundary quads are only allowed "
- "if they carry the same boundary indicator."));
+ "if they carry the same manifold indicator."));
- quad->set_boundary_id_internal(boundary_quad->boundary_id);
quad->set_manifold_id(boundary_quad->manifold_id);
}
projects / dealii.git / commitdiff