From 50ebc6d8e1957989238a6a7db5354810427c64f4 Mon Sep 17 00:00:00 2001 From: third-party Date: Wed, 11 Oct 2000 13:36:31 +0000 Subject: [PATCH] Implementation of 3d boundary data in UCD format and in Triangulation<3>::create_triangulation by Michael Stadler . git-svn-id: https://svn.dealii.org/trunk@3405 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/grid/tria.h | 12 ++ deal.II/deal.II/source/grid/grid_in.cc | 53 +++++++-- deal.II/deal.II/source/grid/tria.cc | 155 ++++++++++++++++++++++++- deal.II/doc/news/2000/c-3-0.html | 18 +++ 4 files changed, 223 insertions(+), 15 deletions(-) diff --git a/deal.II/deal.II/include/grid/tria.h b/deal.II/deal.II/include/grid/tria.h index 4ac2484c39..f5e6f05e5a 100644 --- a/deal.II/deal.II/include/grid/tria.h +++ b/deal.II/deal.II/include/grid/tria.h @@ -2923,10 +2923,22 @@ class Triangulation : public TriaDimensionInfo, /** * Exception */ + DeclException4 (ExcQuadInexistant, + int, int, int, int, + << "When trying to give a boundary indicator to a quad: " + << "the quad with bounding lines " << arg1 << ", " << arg2 + << ", " << arg3 << ", " << arg4 << " does not exist."); + /** + * Exception + */ DeclException0 (ExcInteriorLineCantBeBoundary); /** * Exception */ + DeclException0 (ExcInteriorQuadCantBeBoundary); + /** + * Exception + */ DeclException0 (ExcInternalError); /** * Exception diff --git a/deal.II/deal.II/source/grid/grid_in.cc b/deal.II/deal.II/source/grid/grid_in.cc index f31ba93d30..8205cf08d2 100644 --- a/deal.II/deal.II/source/grid/grid_in.cc +++ b/deal.II/deal.II/source/grid/grid_in.cc @@ -121,7 +121,7 @@ void GridIn::read_ucd (istream &in) }; } else - if ((cell_type == "line") && (dim == 2)) + if ((cell_type == "line") && ((dim == 2) || (dim == 3))) // boundary info { subcelldata.boundary_lines.push_back (CellData<1>()); @@ -137,20 +137,51 @@ void GridIn::read_ucd (istream &in) // vertex with this index exists subcelldata.boundary_lines.back().vertices[i] = vertex_indices[subcelldata.boundary_lines.back().vertices[i]]; - else - { - // no such vertex index - AssertThrow (false, - ExcInvalidVertexIndex(cell, - subcelldata.boundary_lines.back().vertices[i])); - subcelldata.boundary_lines.back().vertices[i] = -1; - }; + else + { + // no such vertex index + AssertThrow (false, + ExcInvalidVertexIndex(cell, + subcelldata.boundary_lines.back().vertices[i])); + subcelldata.boundary_lines.back().vertices[i] = -1; + }; } else - // cannot read this - AssertThrow (false, ExcUnknownIdentifier(cell_type)); + if ((cell_type == "quad") && (dim == 3)) + // boundary info + { + subcelldata.boundary_quads.push_back (CellData<2>()); + in >> subcelldata.boundary_quads.back().vertices[0] + >> subcelldata.boundary_quads.back().vertices[1] + >> subcelldata.boundary_quads.back().vertices[2] + >> subcelldata.boundary_quads.back().vertices[3]; + + subcelldata.boundary_quads.back().material_id = material_id; + + // transform from ucd to + // consecutive numbering + for (unsigned int i=0; i<4; ++i) + if (vertex_indices.find (subcelldata.boundary_quads.back().vertices[i]) != + vertex_indices.end()) + // vertex with this index exists + subcelldata.boundary_quads.back().vertices[i] + = vertex_indices[subcelldata.boundary_quads.back().vertices[i]]; + else + { + // no such vertex index + Assert (false, + ExcInvalidVertexIndex(cell, + subcelldata.boundary_quads.back().vertices[i])); + subcelldata.boundary_quads.back().vertices[i] = -1; + }; + + } + else + // cannot read this + AssertThrow (false, ExcUnknownIdentifier(cell_type)); }; + // check that no forbidden arrays are used Assert (subcelldata.check_consistency(dim), ExcInternalError()); diff --git a/deal.II/deal.II/source/grid/tria.cc b/deal.II/deal.II/source/grid/tria.cc index 6ba3e6ad4e..fd4c96b52b 100644 --- a/deal.II/deal.II/source/grid/tria.cc +++ b/deal.II/deal.II/source/grid/tria.cc @@ -1322,10 +1322,157 @@ void Triangulation<3>::create_triangulation (const vector > &v, // now set boundary indicators // where given // - // not yet implemented - AssertThrow ((subcelldata.boundary_lines.size() == 0) && - (subcelldata.boundary_quads.size() == 0), - ExcNotImplemented()); + // first do so for lines + vector >::const_iterator boundary_line + = subcelldata.boundary_lines.begin(); + vector >::const_iterator end_boundary_line + = subcelldata.boundary_lines.end(); + for (; boundary_line!=end_boundary_line; ++boundary_line) + { + line_iterator line; + pair line_vertices(make_pair(boundary_line->vertices[0], + boundary_line->vertices[1])); + if (needed_lines.find(line_vertices) != needed_lines.end()) + // line found in this + // direction + line = needed_lines[line_vertices]; + + else + { + // look wether it exists in + // reverse direction + swap (line_vertices.first, line_vertices.second); + if (needed_lines.find(line_vertices) != needed_lines.end()) + line = needed_lines[line_vertices]; + else + { + // line does not exist + AssertThrow (false, ExcLineInexistant(line_vertices.first, + line_vertices.second)); + line = end_line(); + }; + }; + // Assert that only exterior + // lines are given a boundary + // indicator + AssertThrow (line->boundary_indicator() == 0, + ExcInteriorLineCantBeBoundary()); + + line -> set_boundary_indicator (boundary_line->material_id); + }; + + + // now go on with boundary faces + vector >::const_iterator boundary_quad + = subcelldata.boundary_quads.begin(); + vector >::const_iterator end_boundary_quad + = subcelldata.boundary_quads.end(); + for (; boundary_quad!=end_boundary_quad; ++boundary_quad) + { + quad_iterator quad; + line_iterator line[4]; + vector quad_vertices(4); + + // first find the lines that + // are made up of the given + // vertices, then build up a + // quad from these lines + // finally use the find + // function of the map template + // to find the quad + pair line_vertices[4] + = { make_pair (boundary_quad->vertices[0], boundary_quad->vertices[1]), + make_pair (boundary_quad->vertices[1], boundary_quad->vertices[2]), + make_pair (boundary_quad->vertices[3], boundary_quad->vertices[2]), + make_pair (boundary_quad->vertices[0], boundary_quad->vertices[3]) }; + + for (unsigned int i=0; i<4; ++i) + { + // check whether line + // already exists + if (needed_lines.find(line_vertices[i]) != needed_lines.end()) + line[i] = needed_lines[line_vertices[i]]; + else + // look wether it exists + // in reverse direction + { + swap (line_vertices[i].first, line_vertices[i].second); + if (needed_lines.find(line_vertices[i]) != needed_lines.end()) + line[i] = needed_lines[line_vertices[i]]; + else + { + // line does not exist + AssertThrow (false, ExcLineInexistant(line_vertices[i].first, + line_vertices[i].second)); + line[i] = end_line(); + }; + }; + }; + + + // Set up 2 quads that are + // built up fromt the lines for + // reasons of comparison to + // needed_quads. The second + // quad is the reversed version + // of the first quad in order + // find the quad regardless of + // its orientation. This is + // introduced for convenience + // and because boundary quad + // orientation does not carry + // any information. + Quad quad_compare_1(line[0]->index(), line[1]->index(), line[2]->index(), line[3]->index()); + Quad quad_compare_2(line[3]->index(), line[2]->index(), line[1]->index(), line[0]->index()); + + // try to find the quad with + // lines situated as + // constructed above. if it + // could not be found, rotate + // the boundary lines 3 times + // until it is found or it does + // not exist. + unsigned int n_rotations=0; + bool not_found_quad_1; + while ( (not_found_quad_1=(needed_quads.find(quad_compare_1) == needed_quads.end())) && + ( needed_quads.find(quad_compare_2) == needed_quads.end()) && + (n_rotations<4)) + { + rotate(line, line+1, line+4); + // update the quads with + // rotated lines + for (unsigned int i=0; i<4; ++i) + { + quad_compare_1.set_line(i, line[i]->index()); + quad_compare_2.set_line(3-i, line[i]->index()); + }; + + ++n_rotations; + }; + + if (n_rotations==4) + { + // quad does not exist + AssertThrow (false, ExcQuadInexistant(line[0]->index(), line[1]->index(), + line[2]->index(), line[3]->index())); + quad = end_quad(); + } + + // exception not thrown, + // therefore assign the quad + // appropriately + if (not_found_quad_1) + quad = needed_quads[quad_compare_2]; + else + quad = needed_quads[quad_compare_1]; + + // check whether this face is + // really an exterior one + AssertThrow(quad->boundary_indicator()==0, + ExcInteriorQuadCantBeBoundary()); + quad->set_boundary_indicator (boundary_quad->material_id); + }; + ///////////////////////////////////////// diff --git a/deal.II/doc/news/2000/c-3-0.html b/deal.II/doc/news/2000/c-3-0.html index f50c45014d..15ad0795c2 100644 --- a/deal.II/doc/news/2000/c-3-0.html +++ b/deal.II/doc/news/2000/c-3-0.html @@ -469,6 +469,24 @@ documentation, etc.

deal.II

    +

  1. + Extended: The Triangulation + class can handle boundary information in 3d as well + (i.e. lines and quads in 3d with special material IDs). +
    + (Michael + Stadler 2000/10/11) +

    + +

  2. + Extended: The GridIn + class can now read 3D UCD data, including boundary information + (i.e. lines and quads in 3d with special material IDs). +
    + (Michael + Stadler 2000/10/11) +

    +

  3. Extended: The GridIn::delete_unused_vertices function now eliminates vertices from the input that are not -- 2.39.5