]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a new test that checks manifold clones.
authorDavid Wells <drwells@vt.edu>
Wed, 16 May 2018 22:10:12 +0000 (18:10 -0400)
committerDavid Wells <drwells@vt.edu>
Wed, 16 May 2018 22:18:45 +0000 (18:18 -0400)
We want to be sure that we can clone all of the renamed opencascade
manifolds and end up with an object that really has the original type:
put another way, we want to be sure that these are not implemented as
aliases but as their own classes.

tests/opencascade/manifold_clone_01.cc [new file with mode: 0644]
tests/opencascade/manifold_clone_01.output [new file with mode: 0644]

diff --git a/tests/opencascade/manifold_clone_01.cc b/tests/opencascade/manifold_clone_01.cc
new file mode 100644 (file)
index 0000000..5a63218
--- /dev/null
@@ -0,0 +1,128 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Check that we get the correct type of manifold when we clone it. This
+// checks that the renamed manifolds are not just aliases but, when cloned,
+// give the correct object type.
+
+#include <deal.II/opencascade/boundary_lib.h>
+#include <deal.II/opencascade/manifold_lib.h>
+
+#include <boost/core/demangle.hpp>
+
+#include <gp_Ax2.hxx>
+#include <gp_Dir.hxx>
+#include <gp_Pnt.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepFill.hxx>
+#include <GC_MakeCircle.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+
+#include <typeinfo>
+
+#include "../tests.h"
+
+int main()
+{
+  using namespace OpenCASCADE;
+
+  initlog();
+
+  {
+    // The circle passing through the vertices of the unit square
+    gp_Dir z_axis(0.,0.,1.);
+    gp_Pnt center(.5,.5,0.);
+    gp_Ax2 axis(center, z_axis);
+    Standard_Real radius(std::sqrt(2.)/2.);
+
+    GC_MakeCircle make_circle(axis, radius);
+    Handle(Geom_Circle) circle = make_circle.Value();
+    TopoDS_Edge edge = BRepBuilderAPI_MakeEdge(circle);
+
+    NormalProjectionBoundary<2, 3> manifold_b(edge);
+    std::unique_ptr<Manifold<2, 3>> clone_b = manifold_b.clone();
+    deallog << "typeid of NormalProjectionBoundary<2, 3> is "
+            << boost::core::demangle(typeid(*clone_b).name())
+            << std::endl;
+
+    NormalProjectionManifold<2, 3> manifold_m(edge);
+    std::unique_ptr<Manifold<2, 3>> clone_m = manifold_m.clone();
+    deallog << "typeid of NormalProjectionManifold<2, 3> is "
+            << boost::core::demangle(typeid(*clone_m).name())
+            << std::endl;
+  }
+
+  {
+    // Create a bspline passing through the points
+    std::vector<Point<3> > pts1, pts2;
+    pts1.push_back(Point<3>(0,0,0));
+    pts1.push_back(Point<3>(1,0,0));
+
+    pts2.push_back(Point<3>(0,1,0));
+    pts2.push_back(Point<3>(.5,1,1));
+    pts2.push_back(Point<3>(1,1,0));
+
+    TopoDS_Edge edge1 = interpolation_curve(pts1);
+    TopoDS_Edge edge2 = interpolation_curve(pts2);
+
+    TopoDS_Face face = BRepFill::Face (edge1, edge2);
+
+    DirectionalProjectionBoundary<2, 3> manifold_b(face, Point<3>(0,0,1));
+
+    std::unique_ptr<Manifold<2, 3>> clone_b = manifold_b.clone();
+    deallog << "typeid of DirectionalProjectionBoundary<2, 3> is "
+            << boost::core::demangle(typeid(*clone_b).name())
+            << std::endl;
+
+    DirectionalProjectionManifold<2, 3> manifold_m(face, Point<3>(0,0,1));
+
+    std::unique_ptr<Manifold<2, 3>> clone_m = manifold_m.clone();
+    deallog << "typeid of DirectionalProjectionManifold<2, 3> is "
+            << boost::core::demangle(typeid(*clone_m).name())
+            << std::endl;
+  }
+
+  // Create a bspline passing through the points
+  {
+    std::vector<Point<3> > pts;
+    pts.push_back(Point<3>(0,0,0));
+    pts.push_back(Point<3>(0,1,0));
+    pts.push_back(Point<3>(1,1,0));
+    pts.push_back(Point<3>(1,0,0));
+
+    TopoDS_Edge edge1 = interpolation_curve(pts);
+    for (unsigned int i=0; i<pts.size(); ++i)
+      pts[i] += Point<3>(0,0,1);
+    TopoDS_Edge edge2 = interpolation_curve(pts);
+
+    TopoDS_Face face = BRepFill::Face (edge1, edge2);
+
+    NormalToMeshProjectionBoundary<1, 3> manifold_b(face);
+
+    std::unique_ptr<Manifold<1, 3>> clone_b = manifold_b.clone();
+    deallog << "typeid of NormalProjectionBoundary<2, 3> is "
+            << boost::core::demangle(typeid(*clone_b).name())
+            << std::endl;
+
+    NormalToMeshProjectionManifold<1, 3> manifold_m(face);
+    std::unique_ptr<Manifold<1, 3>> clone_m = manifold_m.clone();
+    deallog << "typeid of NormalProjectionManifold<2, 3> is "
+            << boost::core::demangle(typeid(*clone_m).name())
+            << std::endl;
+  }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/opencascade/manifold_clone_01.output b/tests/opencascade/manifold_clone_01.output
new file mode 100644 (file)
index 0000000..7dc30c2
--- /dev/null
@@ -0,0 +1,8 @@
+
+DEAL::typeid of NormalProjectionBoundary<2, 3> is dealii::OpenCASCADE::NormalProjectionBoundary<2, 3>
+DEAL::typeid of NormalProjectionManifold<2, 3> is dealii::OpenCASCADE::NormalProjectionManifold<2, 3>
+DEAL::typeid of DirectionalProjectionBoundary<2, 3> is dealii::OpenCASCADE::DirectionalProjectionBoundary<2, 3>
+DEAL::typeid of DirectionalProjectionManifold<2, 3> is dealii::OpenCASCADE::DirectionalProjectionManifold<2, 3>
+DEAL::typeid of NormalProjectionBoundary<2, 3> is dealii::OpenCASCADE::NormalToMeshProjectionBoundary<1, 3>
+DEAL::typeid of NormalProjectionManifold<2, 3> is dealii::OpenCASCADE::NormalToMeshProjectionManifold<1, 3>
+DEAL::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.