#include <deal.II/grid/tria.h>
#include <deal.II/base/point.h>
+#include <deal.II/fe/mapping_q1.h>
#include <string>
* TopoDS_Shapes of one type or another, and provide interfaces to common
* deal.II objects, like Triangulation, Manifold, and so on.
*
- * Notice that these tools are only useful when spacedim is equal to three,
- * since OpenCASCADE only operates in three-dimensional mode.
+ * Notice that most of these tools are only useful when spacedim is equal to
+ * three, since OpenCASCADE only operates in three-dimensional mode. In some
+ * cases they can be used in two dimensions as well, and the third dimension
+ * will be set to zero.
*
* @author Luca Heltai, Andrea Mola, 2011--2014.
*/
void create_triangulation(const TopoDS_Face &face,
Triangulation<2,3> &tria);
+
+ /**
+ * Given a Triangulation and an optional Mapping, create a vector of smooth
+ * curves that interpolate the connected parts of the boundary vertices of
+ * the Triangulation and return them as a vector of TopoDS_Edge.
+ *
+ * The returned curves are ordered with respect to the indices of the faces
+ * that make up the triangulation boundary, i.e., the first curve is the one
+ * extracted starting from the face with the lowest index, and so on.
+ *
+ * @param[input] triangulation Input triangulation
+ * @param[input] mapping Optional input mapping
+ * @return An std::vector of TopoDS_Edge objects representing the smooth
+ * interpolation of the boundary of the `triangulation`
+ *
+ * @author Dirk Peschka, Luca Heltai, 2017.
+ */
+ template <int spacedim>
+ std::vector<TopoDS_Edge> create_curves_from_triangulation_boundary
+ (const Triangulation<2,spacedim> &triangulation,
+ const Mapping<2,spacedim> &mapping=StaticMappingQ1<2,spacedim>::mapping);
+
/**
* Extract all compound shapes from a TopoDS_Shape, and store the results
* into standard containers. If the shape does not contain a certain type of
return out_shape;
}
+
+
+ template<int spacedim>
+ std::vector<TopoDS_Edge> create_curves_from_triangulation_boundary
+ (const Triangulation<2, spacedim> &triangulation,
+ const Mapping<2, spacedim> &mapping)
+
+ {
+ // store maps from global vertex index to pairs of global face indices
+ // and from global face index to pairs of global vertex indices
+ std::map<unsigned int, std::pair<unsigned int,unsigned int> > vert_to_faces;
+ std::map<unsigned int, std::pair<unsigned int,unsigned int> > face_to_verts;
+ std::map<unsigned int, bool> visited_faces;
+ std::map<unsigned int, Point<spacedim> > vert_to_point;
+
+ unsigned int face_index;
+
+ for (auto cell: triangulation.active_cell_iterators())
+ for (unsigned int f=0; f<GeometryInfo<2>::faces_per_cell; ++f)
+ if (cell->face(f)->at_boundary())
+ {
+ // get global face and vertex indices
+ face_index = cell->face(f)->index();
+ unsigned int v0 = cell->face(f)->vertex_index(0);
+ unsigned int v1 = cell->face(f)->vertex_index(1);
+ face_to_verts[face_index].first = v0;
+ face_to_verts[face_index].second = v1;
+ visited_faces[face_index] = false;
+
+ // extract mapped vertex locations
+ std::array<Point<spacedim>,GeometryInfo<2>::vertices_per_cell> verts = mapping.get_vertices(cell);
+ vert_to_point[v0]=verts[GeometryInfo<2>::face_to_cell_vertices(f,0,true,false,false)];
+ vert_to_point[v1]=verts[GeometryInfo<2>::face_to_cell_vertices(f,1,true,false,false)];
+
+ // distribute indices into maps
+ if (vert_to_faces.find(v0) == vert_to_faces.end())
+ {
+ vert_to_faces[v0].first = face_index;
+ }
+ else
+ {
+ vert_to_faces[v0].second = face_index;
+ }
+ if (vert_to_faces.find(v1) == vert_to_faces.end())
+ {
+ vert_to_faces[v1].first = face_index;
+ }
+ else
+ {
+ vert_to_faces[v1].second = face_index;
+ }
+ }
+
+ // run through maps in an orderly fashion, i.e., through the
+ // boundary in one cycle and add points to pointlist.
+ std::vector<TopoDS_Edge> interpolation_curves;
+ bool finished = (face_to_verts.size() == 0);
+ face_index = finished ? 0 : face_to_verts.begin()->first;
+
+ while (finished == false)
+ {
+ const unsigned int start_point_index = face_to_verts[face_index].first;
+ unsigned int point_index = start_point_index;
+
+ // point_index and face_index always run together
+ std::vector<Point<3> > pointlist;
+ do
+ {
+ visited_faces[face_index] = true;
+ auto current_point = vert_to_point[point_index];
+ if (spacedim==2)
+ pointlist.push_back(Point<3>(current_point[0],current_point[1],0));
+ else
+ pointlist.push_back(Point<3>(current_point[0],current_point[1],current_point[2]));
+
+ // Get next point
+ if (face_to_verts[face_index].first != point_index )
+ point_index = face_to_verts[face_index].first;
+ else
+ point_index = face_to_verts[face_index].second;
+
+ // Get next face
+ if (vert_to_faces[point_index].first != face_index)
+ face_index = vert_to_faces[point_index].first;
+ else
+ face_index = vert_to_faces[point_index].second;
+ }
+ while (point_index != start_point_index);
+
+ interpolation_curves.push_back(interpolation_curve(pointlist, Tensor<1,3>(), true));
+ pointlist.clear();
+
+ finished = true;
+ for (const auto &f : visited_faces)
+ if (f.second == false)
+ {
+ face_index = f.first;
+ finished = false;
+ break;
+ }
+ }
+ return interpolation_curves;
+ }
+
+
+
std::tuple<Point<3>, TopoDS_Shape, double, double>
project_point_and_pull_back(const TopoDS_Shape &in_shape,
const Point<3> &origin,
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 - 2015 by the deal.II authors
+//
+// This file is subject to LGPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//-----------------------------------------------------------
+
+// Create a Triangulation, interpolate its boundary points to a close
+// smooth BSpline, and use that as a Boundary Descriptor.
+
+#include "../tests.h"
+
+#include <deal.II/opencascade/utilities.h>
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_tools.h>
+
+using namespace OpenCASCADE;
+
+void remove_iges_header(const std::string &in_filename,
+ const std::string &out_filename)
+{
+ std::ifstream in(in_filename);
+ std::ofstream out(out_filename);
+ std::string line;
+ unsigned int counter = 5;
+ while (counter--) std::getline(in,line);
+ while (std::getline(in,line))
+ out << line << std::endl;
+ in.close();
+ out.close();
+}
+
+template<int spacedim>
+void test()
+{
+ deallog << "Testing <2," << spacedim << ">" << std::endl;
+
+ Triangulation<2,spacedim> tria1;
+ Triangulation<2,spacedim> tria2;
+ Triangulation<2,spacedim> tria3;
+
+ GridGenerator::hyper_cube(tria1, 0, 1);
+ GridGenerator::hyper_cube(tria2, 2, 3);
+
+ GridGenerator::merge_triangulations(tria1, tria2, tria3);
+
+ tria3.refine_global(1);
+
+ auto v1 = create_curves_from_triangulation_boundary(tria1);
+ auto v3 = create_curves_from_triangulation_boundary(tria3);
+
+ if (spacedim == 3)
+ {
+ Triangulation<2,3> tria4;
+ GridGenerator::hyper_sphere(tria4);
+ auto v4 = create_curves_from_triangulation_boundary(tria4);
+
+ if (v4.size() != 0)
+ deallog << "Not OK!: size is " << v4.size() << std::endl;
+ else
+ deallog << "OK" << std::endl;
+ }
+
+ if (v1.size() != 1)
+ deallog << "Not OK: size is " << v1.size() << std::endl;
+ else
+ deallog << "OK" << std::endl;
+
+ if (v3.size() != 2)
+ deallog << "Not OK!: size is " << v3.size() << std::endl;
+ else
+ deallog << "OK" << std::endl;
+
+ write_IGES(v1[0], "edge1.iges");
+ write_IGES(v3[0], "edge30.iges");
+ write_IGES(v3[1], "edge31.iges");
+
+ remove_iges_header("edge1.iges", "edge1_noheader.iges");
+ remove_iges_header("edge30.iges", "edge30_noheader.iges");
+ remove_iges_header("edge31.iges", "edge31_noheader.iges");
+
+ cat_file("edge1_noheader.iges");
+ cat_file("edge30_noheader.iges");
+ cat_file("edge31_noheader.iges");
+}
+
+int main ()
+{
+ initlog();
+
+ test<2>();
+ test<3>();
+
+ return 0;
+}
+
--- /dev/null
+
+DEAL::Testing <2,2>
+DEAL::OK
+DEAL::OK
+ 126 1 0 0 0 0 0 000000000D0000001
+ 126 0 0 4 0 0D0000002
+126,6,3,1,1,1,1,-2.,-1.,0.,0.,1.,2.,3.,4.,4.,5.,6.,1.,1.,1.,1., 0000001P0000001
+1.,1.,1.,-0.5,0.388888889,0.,0.5,-0.388888889,0.,1.142857143, 0000001P0000002
+-0.19047619,0.,1.25,1.25,0.,-0.142857143,1.19047619,0.,-0.5, 0000001P0000003
+0.388888889,0.,0.5,-0.388888889,0.,0.,4.,0.,0.,1.; 0000001P0000004
+S 1G 4D 2P 4 T0000001
+
+ 126 1 0 0 0 0 0 000000000D0000001
+ 126 0 0 8 0 0D0000002
+126,10,3,1,1,1,1,-1.,-0.5,0.,0.,0.5,1.,1.5,2.,2.5,3.,3.5,4.,4., 0000001P0000001
+4.5,5.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,-0.111111111, 0000001P0000002
+-5.555555556E-02,0.,0.111111111,5.555555556E-02,0.,0.481590574, 0000001P0000003
+1.251840943E-02,0.,1.147766323,-0.127147766,0.,0.927344134, 0000001P0000004
+0.496072656,0.,1.142857143,1.142857143,0.,0.501227295, 0000001P0000005
+0.932498773,0.,-0.147766323,1.127147766,0.,8.983799705E-02, 0000001P0000006
+0.558910162,0.,-0.111111111,-5.555555556E-02,0.,0.111111111, 0000001P0000007
+5.555555556E-02,0.,0.,4.,0.,0.,1.; 0000001P0000008
+S 1G 4D 2P 8 T0000001
+
+ 126 1 0 0 0 0 0 000000000D0000001
+ 126 0 0 8 0 0D0000002
+126,10,3,1,1,1,1,-1.,-0.5,0.,0.,0.5,1.,1.5,2.,2.5,3.,3.5,4.,4., 0000001P0000001
+4.5,5.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.888888889,1.944444444, 0000001P0000002
+0.,2.111111111,2.055555556,0.,2.481590574,2.012518409,0., 0000001P0000003
+3.147766323,1.872852234,0.,2.927344134,2.496072656,0., 0000001P0000004
+3.142857143,3.142857143,0.,2.501227295,2.932498773,0., 0000001P0000005
+1.852233677,3.127147766,0.,2.089837997,2.558910162,0., 0000001P0000006
+1.888888889,1.944444444,0.,2.111111111,2.055555556,0.,0.,4.,0., 0000001P0000007
+0.,1.; 0000001P0000008
+S 1G 4D 2P 8 T0000001
+
+DEAL::Testing <2,3>
+DEAL::OK
+DEAL::OK
+DEAL::OK
+ 126 1 0 0 0 0 0 000000000D0000001
+ 126 0 0 4 0 0D0000002
+126,6,3,1,1,1,1,-2.,-1.,0.,0.,1.,2.,3.,4.,4.,5.,6.,1.,1.,1.,1., 0000001P0000001
+1.,1.,1.,-0.5,0.388888889,0.,0.5,-0.388888889,0.,1.142857143, 0000001P0000002
+-0.19047619,0.,1.25,1.25,0.,-0.142857143,1.19047619,0.,-0.5, 0000001P0000003
+0.388888889,0.,0.5,-0.388888889,0.,0.,4.,0.,0.,1.; 0000001P0000004
+S 1G 4D 2P 4 T0000001
+
+ 126 1 0 0 0 0 0 000000000D0000001
+ 126 0 0 8 0 0D0000002
+126,10,3,1,1,1,1,-1.,-0.5,0.,0.,0.5,1.,1.5,2.,2.5,3.,3.5,4.,4., 0000001P0000001
+4.5,5.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,-0.111111111, 0000001P0000002
+-5.555555556E-02,0.,0.111111111,5.555555556E-02,0.,0.481590574, 0000001P0000003
+1.251840943E-02,0.,1.147766323,-0.127147766,0.,0.927344134, 0000001P0000004
+0.496072656,0.,1.142857143,1.142857143,0.,0.501227295, 0000001P0000005
+0.932498773,0.,-0.147766323,1.127147766,0.,8.983799705E-02, 0000001P0000006
+0.558910162,0.,-0.111111111,-5.555555556E-02,0.,0.111111111, 0000001P0000007
+5.555555556E-02,0.,0.,4.,0.,0.,1.; 0000001P0000008
+S 1G 4D 2P 8 T0000001
+
+ 126 1 0 0 0 0 0 000000000D0000001
+ 126 0 0 8 0 0D0000002
+126,10,3,1,1,1,1,-1.,-0.5,0.,0.,0.5,1.,1.5,2.,2.5,3.,3.5,4.,4., 0000001P0000001
+4.5,5.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.888888889,1.944444444, 0000001P0000002
+0.,2.111111111,2.055555556,0.,2.481590574,2.012518409,0., 0000001P0000003
+3.147766323,1.872852234,0.,2.927344134,2.496072656,0., 0000001P0000004
+3.142857143,3.142857143,0.,2.501227295,2.932498773,0., 0000001P0000005
+1.852233677,3.127147766,0.,2.089837997,2.558910162,0., 0000001P0000006
+1.888888889,1.944444444,0.,2.111111111,2.055555556,0.,0.,4.,0., 0000001P0000007
+0.,1.; 0000001P0000008
+S 1G 4D 2P 8 T0000001
+
projects / dealii.git / commitdiff