PERL = perl
-USE_CONTRIB_PETSC = @MAKEFILE_PETSC@
-USE_CONTRIB_SLEPC = @MAKEFILE_SLEPC@
-USE_CONTRIB_TRILINOS = @MAKEFILE_TRILINOS@
-USE_CONTRIB_HDF5 = @MAKEFILE_HDF5@
-USE_CONTRIB_BLAS = @MAKEFILE_BLAS@
-USE_CONTRIB_LAPACK = @MAKEFILE_LAPACK@
-USE_CONTRIB_MUMPS = @MAKEFILE_MUMPS@
-USE_CONTRIB_ARPACK = @MAKEFILE_ARPACK@
-USE_CONTRIB_METIS = @MAKEFILE_METIS@
-USE_CONTRIB_UMFPACK = @MAKEFILE_UMFPACK@
-USE_CONTRIB_P4EST = @MAKEFILE_P4EST@
-USE_CONTRIB_OPENCASCADE = @MAKEFILE_OPENCASCADE@
-DEAL_II_USE_MPI = @MAKEFILE_MPI@
+USE_CONTRIB_PETSC = @MAKEFILE_PETSC@
+USE_CONTRIB_SLEPC = @MAKEFILE_SLEPC@
+USE_CONTRIB_TRILINOS = @MAKEFILE_TRILINOS@
+USE_CONTRIB_HDF5 = @MAKEFILE_HDF5@
+USE_CONTRIB_BLAS = @MAKEFILE_BLAS@
+USE_CONTRIB_LAPACK = @MAKEFILE_LAPACK@
+USE_CONTRIB_MUMPS = @MAKEFILE_MUMPS@
+USE_CONTRIB_ARPACK = @MAKEFILE_ARPACK@
+USE_CONTRIB_METIS = @MAKEFILE_METIS@
+USE_CONTRIB_UMFPACK = @MAKEFILE_UMFPACK@
+USE_CONTRIB_P4EST = @MAKEFILE_P4EST@
+USE_CONTRIB_OPENCASCADE = @MAKEFILE_OPENCASCADE@
+DEAL_II_USE_MPI = @MAKEFILE_MPI@
OBJEXT = o
EXEEXT =
ENDIF()
# These seem to be pretty much the only required ones.
-SET(OPENCASCADE_LIBRARIES
+SET(_opencascade_libraries
TKFillet
TKMesh
TKernel
)
SET(_libraries "")
-FOREACH(_library ${OPENCASCADE_LIBRARIES})
- LIST(APPEND _libraries ${_library})
- DEAL_II_FIND_LIBRARY(${_library}
+FOREACH(_library ${_opencascade_libraries})
+ LIST(APPEND _libraries OPENCASCADE_${_library})
+ DEAL_II_FIND_LIBRARY(OPENCASCADE_${_library}
NAMES ${_library}
HINTS ${OPENCASCADE_DIR}
PATH_SUFFIXES lib${LIB_SUFFIX} lib64 lib
USER_INCLUDE_DIRS
REQUIRED OPENCASCADE_INCLUDE_DIR
CLEAR
- OPENCASCADE_LIBRARIES ${_libraries}
- TRILINOS_SUPPORTS_CPP11 TRILINOS_HAS_C99_TR1_WORKAROUND
+ _opencascade_libraries ${_libraries}
)
#include <deal.II/grid/occ_utilities.h>
#include <deal.II/grid/tria_boundary.h>
+#include <BRepAdaptor_Curve.hxx>
+#include <Adaptor3d_Curve.hxx>
+
DEAL_II_NAMESPACE_OPEN
-namespace OpenCASCADE
+/**
+ * @addtogroup OpenCASCADE
+ * @{
+ */
+
+namespace OpenCASCADE
{
/**
- * @addtogroup OpenCASCADE
- * @{
- *
- * A Boundary object based on OpenCASCADE TopoDS_Shape where new
- * points are first computed using the FlatManifold class, and then
- * projected in the normal direction using OpenCASCADE utilities.
+ * A Boundary object based on OpenCASCADE TopoDS_Shape where where
+ * new points are first computed by averaging the surrounding points
+ * in the same way as FlatManifold does, and then projecting in the
+ * normal direction using OpenCASCADE utilities.
*
* This class makes no assumptions on the shape you pass to it, and
* the topological dimension of the Manifold is inferred from the
* This could happen, for example, if you are trying to use a shape
* of type TopoDS_Edge when projecting on a face. In this case, the
* vertices of the face would be collapsed to the edge, and your
- * surrounding points would not be on lying on the given shape,
- * raising an exception.
- *
+ * surrounding points would not be lying on the given shape, raising
+ * an exception.
+ *
* @author Luca Heltai, Andrea Mola, 2011--2014.
*/
template <int dim, int spacedim>
- class NormalProjectionBoundary : public Boundary<dim,spacedim> {
- public:
- NormalProjectionBoundary(const TopoDS_Shape sh,
- const double tolerance=1e-7);
-
+ class NormalProjectionBoundary : public Boundary<dim,spacedim>
+ {
+ public:
+
+ /**
+ * The standard constructor takes a generic TopoDS_Shape @p sh,
+ * and a tolerance used to compute distances internally.
+ *
+ * The TopoDS_Shape can be of arbitrary, i.e., a collection of
+ * shapes, faces, edges or a single face or edge.
+ */
+ NormalProjectionBoundary(const TopoDS_Shape &sh,
+ const double tolerance=1e-7);
+
/**
* Perform the actual projection onto the manifold. This function,
- * in debug mode, checks that each of the #surrounding_points is
+ * in debug mode, checks that each of the @p surrounding_points is
* within tolerance from the given TopoDS_Shape. If this is not
* the case, an exception is thrown.
*
*/
virtual Point<spacedim>
project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
- const Point<spacedim> &candidate) const;
+ const Point<spacedim> &candidate) const;
private:
/**
* The topological shape which is used internally to project
- * points. You can construct one such a shape by calling the
+ * points. You can construct such a shape by calling the
* OpenCASCADE::read_IGES() function, which will create a
* TopoDS_Shape with the geometry contained in the IGES file.
- */
+ */
const TopoDS_Shape sh;
/**
/**
* A Boundary object based on OpenCASCADE TopoDS_Shape where new
- * points are first computed using the FlatManifold class, and then
- * projected along the direction given at construction time, using
- * OpenCASCADE utilities.
+ * points are first computed by averaging the surrounding points in
+ * the same way as FlatManifold does, and then projecting in onto
+ * the manifold along the direction specified at construction time
+ * using OpenCASCADE utilities.
*
* This class makes no assumptions on the shape you pass to it, and
* the topological dimension of the Manifold is inferred from the
* calling OpenCASCADE::closest_point() on those points leaves them
* untouched. If this is not the case, an ExcPointNotOnManifold is
* thrown.
- *
+ *
+ * Notice that this type of Boundary descriptor may fail to give
+ * results if the triangulation to be refined is close to the
+ * boundary of the given TopoDS_Shape, or when the direction you use
+ * at construction time does not intersect the shape. An exception
+ * is thrown when this appens.
+ *
* @author Luca Heltai, Andrea Mola, 2011--2014.
*/
template <int dim, int spacedim>
- class AxisProjectionBoundary : public Boundary<dim,spacedim> {
- public:
- AxisProjectionBoundary(const TopoDS_Shape sh,
- const Point<3> direction,
- const double tolerance=1e-7);
-
+ class DirectionalProjectionBoundary : public Boundary<dim,spacedim>
+ {
+ public:
+ /**
+ * Construct a Boundary object which will project points on the
+ * TopoDS_Shape @p sh, along the given @p direction.
+ */
+ DirectionalProjectionBoundary(const TopoDS_Shape &sh,
+ const Tensor<1,spacedim> &direction,
+ const double tolerance=1e-7);
+
/**
* Perform the actual projection onto the manifold. This function,
- * in debug mode, checks that each of the #surrounding_points is
+ * in debug mode, checks that each of the @p surrounding_points is
* within tolerance from the given TopoDS_Shape. If this is not
* the case, an exception is thrown.
*
- * The projected point is computed using OpenCASCADE normal
+ * The projected point is computed using OpenCASCADE directional
* projection algorithms.
*/
virtual Point<spacedim>
project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
- const Point<spacedim> &candidate) const;
+ const Point<spacedim> &candidate) const;
private:
/**
* The topological shape which is used internally to project
- * points. You can construct one such a shape by calling the
+ * points. You can construct such a shape by calling the
* OpenCASCADE::read_IGES() function, which will create a
* TopoDS_Shape with the geometry contained in the IGES file.
- */
+ */
const TopoDS_Shape sh;
/**
- * Direction used to project new points on the shape.
+ * Direction used to project new points on the shape.
*/
const Point<3> direction;
const double tolerance;
};
+ /**
+ * A Boundary object based on OpenCASCADE TopoDS_Shape objects which
+ * have topological dimension equal to one (TopoDS_Edge or
+ * TopoDS_Wire) where new points are located at the arclength
+ * average of the surrounding points. If the given TopoDS_Shape can
+ * be casted to a periodic (closed) curve, then this information is
+ * used internally to set the periodicity of the base ChartManifold
+ * class.
+ *
+ * This class can only work on TopoDS_Edge or TopoDS_Wire objects,
+ * and it only makes sense when spacedim is three. If you use an
+ * object of topological dimension different from one, an exception
+ * is throw.
+ *
+ * In debug mode there is an additional sanity check to make sure
+ * that the surrounding points actually live on the Manifold, i.e.,
+ * calling OpenCASCADE::closest_point() on those points leaves them
+ * untouched. If this is not the case, an ExcPointNotOnManifold is
+ * thrown.
+ *
+ * @author Luca Heltai, Andrea Mola, 2011--2014.
+ */
+ template <int dim, int spacedim>
+ class ArclengthProjectionLineManifold : public ChartManifold<dim,spacedim,1>
+ {
+ public:
+ /**
+ * Default constructor with a TopoDS_Edge.
+ */
+ ArclengthProjectionLineManifold(const TopoDS_Shape &sh,
+ const double tolerance=1e-7);
+
+ /**
+ * Given a point on real space, find its arclength
+ * parameter. Throws an error in debug mode, if the point is not
+ * on the TopoDS_Edge given at construction time.
+ */
+ virtual Point<1>
+ pull_back(const Point<spacedim> &space_point) const;
+
+ /**
+ * Given an arclength parameter, find its image in real space.
+ */
+ virtual Point<spacedim>
+ push_forward(const Point<1> &chart_point) const;
+
+ private:
+ /**
+ * A Curve adaptor. This is the one which is used in the
+ * computations, and it points to the right one above.
+ */
+ Handle_Adaptor3d_HCurve curve;
+
+ /**
+ * Relative tolerance used in all internal computations.
+ */
+ const double tolerance;
+
+ /**
+ * The total length of the curve. This is also used as a period if
+ * the edge is periodic.
+ */
+ const double length;
+ };
}
/*@}*/
#ifdef DEAL_II_WITH_OPENCASCADE
+#include <deal.II/grid/tria.h>
+
#include <string>
#include <TopoDS_Shape.hxx>
+#include <TopoDS_Face.hxx>
#include <TopoDS_Edge.hxx>
+#include <TopoDS_Vertex.hxx>
+#include <TopoDS_Compound.hxx>
+#include <TopoDS_CompSolid.hxx>
+#include <TopoDS_Solid.hxx>
+#include <TopoDS_Shell.hxx>
+#include <TopoDS_Wire.hxx>
#include <IFSelect_ReturnStatus.hxx>
-#include <Geom_Plane.hxx>
-#include <Geom_Curve.hxx>
#include <gp_Pnt.hxx>
#include <deal.II/base/point.h>
+
DEAL_II_NAMESPACE_OPEN
/**
* contains 3 fields: location, orientation and a myTShape handle (of
* the TopoDS_TShape type). According to OpenCASCADE documentation,
* myTShape and Location are used to share data between various shapes
- * and thus save huge amounts of memory. For example, an edge
- * belonging to two faces has equal Locations and myTShape fields but
- * different Orientations (Forward in context of one face and Reversed
- * in one of the other).
+ * to save memory. For example, an edge belonging to two faces has
+ * equal Locations and myTShape fields but different Orientations
+ * (Forward in context of one face and Reversed in one of the other).
*
* Valid shapes include collection of other shapes, solids, faces,
* edges, vertices, etc.
*
* @author Luca Heltai, Andrea Mola, 2011--2014.
*/
-namespace OpenCASCADE
+namespace OpenCASCADE
{
/**
- * Count the subobjects of a shape. This function just outputs some
- * information about the TopoDS_Shape passed as argument. It counts
- * the number of faces, edges and vertices (the only topological
- * entities associated with actual geometries) which are contained
- * in the given shape.
+ * Count the subobjects of a shape. This function is useful to
+ * gather information about the TopoDS_Shape passed as argument. It
+ * counts the number of faces, edges and vertices (the only
+ * topological entities associated with actual geometries) which are
+ * contained in the given shape.
*/
void count_elements(const TopoDS_Shape &shape,
- unsigned int &n_faces,
- unsigned int &n_edges,
- unsigned int &n_vertices);
-
+ unsigned int &n_faces,
+ unsigned int &n_edges,
+ unsigned int &n_vertices);
+
/**
* Read IGES files and translate their content into openCascade
* topological entities. The option scale_factor is used to
* compensate for different units being used in the IGES files and
* in the target application. The standard unit for IGES files is
* millimiters. The return object is a TopoDS_Shape which contains
- * all objects from the file.
- */
- TopoDS_Shape read_IGES(const std::string &filename,
- const double scale_factor=1e-3);
-
- /**
- * Write the given topological shape into an IGES file.
+ * all objects from the file.
*/
- void write_IGES(const TopoDS_Shape &shape,
- const std::string &filename);
-
-
+ TopoDS_Shape read_IGES(const std::string &filename,
+ const double scale_factor=1e-3);
+
+ /**
+ * Write the given topological shape into an IGES file.
+ */
+ void write_IGES(const TopoDS_Shape &shape,
+ const std::string &filename);
+
+
/**
* Perform the intersection of the given topological shape with the
* thrown if the intersection produces an empty shape.
*/
TopoDS_Shape intersect_plane(const TopoDS_Shape &in_shape,
- const double c_x,
- const double c_y,
- const double c_z,
- const double c,
- const double tolerance=1e-7);
+ const double c_x,
+ const double c_y,
+ const double c_z,
+ const double c,
+ const double tolerance=1e-7);
+
+ /**
+ * Try to join all edges contained in the given TopoDS_Shape into a
+ * single TopoDS_Edge, containing as few BSPlines as possible. If
+ * the input shape contains faces, they will be ignored by this
+ * function. If the contained edges cannot be joined into a single
+ * one, i.e., they form disconnected curves, an exception will be
+ * thrown.
+ */
+ TopoDS_Edge join_edges(const TopoDS_Shape &in_shape,
+ const double tolerance=1e-7);
-
/**
* Creates a 3D smooth BSpline curve passing through the points in
* the assigned vector, and store it in the returned TopoDS_Shape
*
* This class is used to interpolate a BsplineCurve passing through
* an array of points, with a C2 Continuity. If the optional
- * parameter #closed is set to true, then the curve will be C2 at
+ * parameter @p closed is set to true, then the curve will be C2 at
* all points execpt the first (where only C1 continuity will be
* given), and it will be a closed curve.
*
- * The curve is garanteed to be at distance #tolerance from the
+ * The curve is garanteed to be at distance @p tolerance from the
* input points. If the algorithm fails in generating such a curve,
* an exception is thrown.
*/
TopoDS_Edge interpolation_curve(std::vector<dealii::Point<3> > &curve_points,
- const dealii::Point<3> direction=dealii::Point<3>(),
- const bool closed=false,
- const double tolerance=1e-7);
+ const dealii::Point<3> direction=dealii::Point<3>(),
+ const bool closed=false,
+ const double tolerance=1e-7);
+
+ /**
+ * Extract all subshapes from a TopoDS_Shape, and store the results
+ * into standard containers. If the shape does not contain a certain
+ * type of shape, the respective container will be empty.
+ */
+ void extract_geometrical_shapes(const TopoDS_Shape &shape,
+ std::vector<TopoDS_Face> &faces,
+ std::vector<TopoDS_Edge> &edges,
+ std::vector<TopoDS_Vertex> &vertices);
+
+ /**
+ * Create a triangulation from a single face. This class extract the
+ * first u and v parameter of the parametric surface making up this
+ * face, and creates a Triangulation<2,3> containing a single coarse
+ * cell reflecting this face. If the surface is not a trimmed
+ * surface, the vertices of this cell will coincide with the
+ * TopoDS_Vertex vertices of the original TopoDS_Face. This,
+ * however, is often not the case, and the user should be careful on
+ * how this mesh is used.
+ */
+ void create_triangulation(const TopoDS_Face &face,
+ Triangulation<2,3> &tria);
+
+ /**
+ * 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 compound, the respective container will be empty.
+ */
+ void extract_compound_shapes(const TopoDS_Shape &shape,
+ std::vector<TopoDS_Compound> &compounds,
+ std::vector<TopoDS_CompSolid> &compsolids,
+ std::vector<TopoDS_Solid> &solids,
+ std::vector<TopoDS_Shell> &shells,
+ std::vector<TopoDS_Wire> &wires);
/**
* coordinate is filled with sensible information, and the v
* coordinate is set to zero.
*/
- Point<3> closest_point(const TopoDS_Shape in_shape,
- const Point<3> origin,
- TopoDS_Shape &out_shape,
- double &u,
- double &v,
- const double tolerance=1e-7);
+ Point<3> closest_point(const TopoDS_Shape in_shape,
+ const Point<3> origin,
+ TopoDS_Shape &out_shape,
+ double &u,
+ double &v,
+ const double tolerance=1e-7);
/**
- * Intersect a line passing through the given #origin point along
- * #direction and the given topological shape. If there is more than
+ * Intersect a line passing through the given @p origin point along
+ * @p direction and the given topological shape. If there is more than
* one intersection, it will return the closest one.
*
- * The optional #tolerance parameter is used to compute distances.
+ * The optional @p tolerance parameter is used to compute distances.
*/
- Point<3> axis_intersection(const TopoDS_Shape in_shape,
- const Point<3> origin,
- const Point<3> direction,
- const double tolerance=1e-7);
-
+ Point<3> axis_intersection(const TopoDS_Shape in_shape,
+ const Point<3> origin,
+ const Point<3> direction,
+ const double tolerance=1e-7);
+
/**
* Convert OpenCASCADE point into a Point<3>.
*/
- inline Point<3> Pnt(const gp_Pnt &p);
+ Point<3> Pnt(const gp_Pnt &p);
/**
* Convert Point<3> into OpenCASCADE point.
*/
- inline gp_Pnt Pnt(const Point<3> &p);
+ gp_Pnt Pnt(const Point<3> &p);
+
-
/**
* Sort two points according to their scalar product with
* direction. If the norm of the direction is zero, then use
* lexycographical ordering. The optional parameter is used as a
* relative tolerance when comparing objects.
*/
- inline bool point_compare(const dealii::Point<3> &p1, const dealii::Point<3> &p2,
- const dealii::Point<3> direction=Point<3>(),
- const double tolerance=1e-10);
+ bool point_compare(const dealii::Point<3> &p1, const dealii::Point<3> &p2,
+ const dealii::Point<3> direction=Point<3>(),
+ const double tolerance=1e-10);
/**
* Exception thrown when the point specified as argument does not
- * lie between #tolerance from the given TopoDS_Shape.
+ * lie between @p tolerance from the given TopoDS_Shape.
*/
DeclException1 (ExcPointNotOnManifold,
- Point<3>,
- <<"The point [ "<<arg1<<" ] is not on the manifold.");
+ Point<3>,
+ <<"The point [ "<<arg1<<" ] is not on the manifold.");
/**
* Exception thrown when the point specified as argument cannot be
* projected to the manifold.
- */
- DeclException1 (ExcProjectionFailed,
- Point<3>,
- <<"Projection of point [ "<< arg1
- << " ] failed.");
+ */
+ DeclException1 (ExcProjectionFailed,
+ Point<3>,
+ <<"Projection of point [ "<< arg1
+ << " ] failed.");
/**
* Thrown when internal OpenCASCADE utilities fail to return the OK
* status.
*/
- DeclException1 (ExcOCCError,
- IFSelect_ReturnStatus,
- <<"An OpenCASCADE routine failed with return status "
- <<arg1);
-}
+ DeclException1 (ExcOCCError,
+ IFSelect_ReturnStatus,
+ <<"An OpenCASCADE routine failed with return status "
+ <<arg1);
+
+ /**
+ * Trying to make curve operations on a degenerate edge.
+ */
+ DeclException0(ExcEdgeIsDegenerate);
+
+ /**
+ * Trying to make operations on the wrong type of shapes.
+ */
+ DeclException0(ExcUnsupportedShape);
+
+}
/*@}*/
DEAL_II_NAMESPACE_CLOSE
#ifdef DEAL_II_WITH_OPENCASCADE
+#include <GCPnts_AbscissaPoint.hxx>
+#include <BRepAdaptor_Curve.hxx>
+#include <BRepAdaptor_CompCurve.hxx>
+#include <BRepAdaptor_HCurve.hxx>
+#include <BRepAdaptor_HCompCurve.hxx>
+#include <GCPnts_AbscissaPoint.hxx>
+#include <ShapeAnalysis_Curve.hxx>
+#include <BRep_Tool.hxx>
+#include <TopoDS.hxx>
+#include <Adaptor3d_HCurve.hxx>
+#include <Handle_Adaptor3d_HCurve.hxx>
+
DEAL_II_NAMESPACE_OPEN
+
namespace OpenCASCADE
{
+
+ namespace
+ {
+ /**
+ * Return a Geometrical curve representation for the given
+ * TopoDS_Shape. This function will fail when the given shape is
+ * not of topological dimension one.
+ */
+ Handle_Adaptor3d_HCurve curve_adaptor(const TopoDS_Shape &shape)
+ {
+ Assert( (shape.ShapeType() == TopAbs_WIRE) ||
+ (shape.ShapeType() == TopAbs_EDGE),
+ ExcUnsupportedShape());
+ if (shape.ShapeType() == TopAbs_WIRE)
+ return (Handle(BRepAdaptor_HCompCurve(new BRepAdaptor_HCompCurve(TopoDS::Wire(shape)))));
+ else if (shape.ShapeType() == TopAbs_EDGE)
+ return (Handle(BRepAdaptor_HCurve(new BRepAdaptor_HCurve(TopoDS::Edge(shape)))));
+
+ Assert(false, ExcInternalError());
+ return Handle(BRepAdaptor_HCurve(new BRepAdaptor_HCurve()));
+ }
+
+
+
+ // Helper internal functions.
+ double shape_length(const TopoDS_Shape &sh)
+ {
+ Handle_Adaptor3d_HCurve adapt = curve_adaptor(sh);
+ return GCPnts_AbscissaPoint::Length(adapt->GetCurve());
+ }
+ }
+
/*============================== NormalProjectionBoundary ==============================*/
template <int dim, int spacedim>
- NormalProjectionBoundary<dim,spacedim>::NormalProjectionBoundary(const TopoDS_Shape sh,
- const double tolerance) :
+ NormalProjectionBoundary<dim,spacedim>::NormalProjectionBoundary(const TopoDS_Shape &sh,
+ const double tolerance) :
sh(sh),
- tolerance(tolerance)
+ tolerance(tolerance)
{
- Assert(spacedim == 3, ExcImpossibleInDim(spacedim));
+ Assert(spacedim == 3, ExcNotImplemented());
}
-
-
+
+
template <int dim, int spacedim>
Point<spacedim> NormalProjectionBoundary<dim,spacedim>::
project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
- const Point<spacedim> &candidate) const {
+ const Point<spacedim> &candidate) const
+ {
TopoDS_Shape out_shape;
double u=0, v=0;
- for(unsigned int i=0; i<surrounding_points.size(); ++i)
+ for (unsigned int i=0; i<surrounding_points.size(); ++i)
Assert(closest_point(sh, surrounding_points[i], out_shape, u, v)
- .distance(surrounding_points[i]) < (surrounding_points[i].norm()>0 ?
- tolerance*surrounding_points[i].norm() :
- tolerance),
- // ExcPointNotOnManifold(surrounding_points[i]));
- ExcMessage("Points are not on manifold"));
-
+ .distance(surrounding_points[i]) <
+ std::max(tolerance*surrounding_points[i].norm(), tolerance),
+ ExcPointNotOnManifold(surrounding_points[i]));
+
return closest_point(sh, candidate, out_shape, u, v);
}
- /*============================== AxisProjectionBoundary ==============================*/
+ /*============================== DirectionalProjectionBoundary ==============================*/
template <int dim, int spacedim>
- AxisProjectionBoundary<dim,spacedim>::AxisProjectionBoundary(const TopoDS_Shape sh,
- const Point<3> direction,
- const double tolerance) :
+ DirectionalProjectionBoundary<dim,spacedim>::DirectionalProjectionBoundary(const TopoDS_Shape &sh,
+ const Tensor<1,spacedim> &direction,
+ const double tolerance) :
sh(sh),
direction(direction),
- tolerance(tolerance)
+ tolerance(tolerance)
{
- Assert(spacedim == 3, ExcImpossibleInDim(spacedim));
+ Assert(spacedim == 3, ExcNotImplemented());
}
-
-
+
+
template <int dim, int spacedim>
- Point<spacedim> AxisProjectionBoundary<dim,spacedim>::
+ Point<spacedim> DirectionalProjectionBoundary<dim,spacedim>::
project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
- const Point<spacedim> &candidate) const {
+ const Point<spacedim> &candidate) const
+ {
TopoDS_Shape out_shape;
double u=0, v=0;
- for(unsigned int i=0; i<surrounding_points.size(); ++i)
+ for (unsigned int i=0; i<surrounding_points.size(); ++i)
Assert(closest_point(sh, surrounding_points[i], out_shape, u, v)
- .distance(surrounding_points[i]) < (surrounding_points[i].norm()>0 ?
- tolerance*surrounding_points[i].norm() :
- tolerance),
- ExcPointNotOnManifold(surrounding_points[i]));
-
+ .distance(surrounding_points[i]) < (surrounding_points[i].norm()>0 ?
+ tolerance*surrounding_points[i].norm() :
+ tolerance),
+ ExcPointNotOnManifold(surrounding_points[i]));
+
return axis_intersection(sh, candidate, direction, tolerance);
}
-
+ /*============================== ArclengthProjectionLineManifold ==============================*/
+
+ template <int dim, int spacedim>
+ ArclengthProjectionLineManifold<dim,spacedim>::ArclengthProjectionLineManifold(const TopoDS_Shape &sh,
+ const double tolerance):
+
+ ChartManifold<dim,spacedim,1>(sh.Closed() ?
+ Point<1>(shape_length(sh)) :
+ Point<1>()),
+ curve(curve_adaptor(sh)),
+ tolerance(tolerance),
+ length(shape_length(sh))
+ {
+ Assert(spacedim == 3, ExcNotImplemented());
+ }
+
+
+ template <int dim, int spacedim>
+ Point<1>
+ ArclengthProjectionLineManifold<dim,spacedim>::pull_back(const Point<spacedim> &space_point) const
+ {
+ ShapeAnalysis_Curve curve_analysis;
+ gp_Pnt proj;
+ double t;
+ double dist = curve_analysis.Project(curve->GetCurve(), Pnt(space_point), tolerance, proj, t, true);
+ Assert(dist < tolerance*length, ExcPointNotOnManifold(space_point));
+ return Point<1>(GCPnts_AbscissaPoint::Length(curve->GetCurve(),curve->GetCurve().FirstParameter(),t));
+ }
+
+
+
+ template <int dim, int spacedim>
+ Point<spacedim>
+ ArclengthProjectionLineManifold<dim,spacedim>::push_forward(const Point<1> &chart_point) const
+ {
+ GCPnts_AbscissaPoint AP(curve->GetCurve(), chart_point[0], curve->GetCurve().FirstParameter());
+ gp_Pnt P = curve->GetCurve().Value(AP.Parameter());
+ return Pnt(P);
+ }
+
+
// Explicit instantiations
-#include "occ_boundary_lib.inst"
-
+#include "occ_boundary_lib.inst"
+
} // end namespace OpenCASCADE
DEAL_II_NAMESPACE_CLOSE
for (deal_II_dimension : DIMENSIONS)
{
template class NormalProjectionBoundary<deal_II_dimension, 3>;
- template class AxisProjectionBoundary<deal_II_dimension, 3>;
+ template class DirectionalProjectionBoundary<deal_II_dimension, 3>;
+ template class ArclengthProjectionLineManifold<deal_II_dimension, 3>;
}
#ifdef DEAL_II_WITH_OPENCASCADE
#include <deal.II/base/point.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/exceptions.h>
+
+#include <boost/bind.hpp>
#include <stdio.h>
#include <stdlib.h>
#include <TopoDS.hxx>
#include <TopoDS_Shape.hxx>
-// #include <TopoDS_Shell.hxx>
#include <TopoDS_Face.hxx>
#include <TopoDS_Edge.hxx>
-#include <BRepTools.hxx>
+#include <TopExp_Explorer.hxx>
-// #include <XSControl_Reader.hxx>
-// #include <TopTools_SequenceOfShape.hxx>
#include <Handle_Standard_Transient.hxx>
+
#include <TColStd_SequenceOfTransient.hxx>
#include <TColStd_HSequenceOfTransient.hxx>
-#include <TopExp_Explorer.hxx>
+#include <TColgp_HArray1OfPnt.hxx>
+
#include <gp_Pnt.hxx>
#include <gp_Lin.hxx>
#include <gp_Vec.hxx>
#include <GeomAPI_ProjectPointOnSurf.hxx>
#include <GeomAPI_ProjectPointOnCurve.hxx>
#include <IntCurvesFace_ShapeIntersector.hxx>
-// #include <Standard_Real.hxx>
-// #include <Standard_Integer.hxx>
-// #include <BRep_Tool.hxx>
-// #include <Geom_Surface.hxx>
-#include <Geom_Plane.hxx>
-// #include <Prs3d_ShapeTool.hxx>
-// #include <GeomAPI_IntSS.hxx>
-// #include <Bnd_Box.hxx>
-// #include <gp_Trsf.hxx>
-// #include <gp_Ax3.hxx>
-// #include <gp_Pln.hxx>
+
+#include <BRepTools.hxx>
+#include <BRepAdaptor_Curve.hxx>
+#include <BRepAdaptor_HCurve.hxx>
+#include <BRepAdaptor_HCompCurve.hxx>
+#include <BRepAdaptor_Surface.hxx>
+#include <BRep_Builder.hxx>
#include <BRepBuilderAPI_Transform.hxx>
-#include <GeomConvert_CompCurveToBSplineCurve.hxx>
#include <BRepBuilderAPI_MakeEdge.hxx>
-// #include <TColGeom_Array1OfCurve.hxx>
-#include <TColgp_HArray1OfPnt.hxx>
-// #include <Geom_Curve.hxx>
-#include <Geom_BoundedCurve.hxx>
-// #include <Geom_TrimmedCurve.hxx>
-// #include <Geom_BSplineCurve.hxx>
-#include <GeomAPI_Interpolate.hxx>
#include <BRepAlgo_Section.hxx>
-// #include <GeomLib_Tool.hxx>
-// #include <TColGeom_Array2OfBezierSurface.hxx>
-// #include <ProjLib_ProjectOnPlane.hxx>
-// #include <Adaptor3d_HCurve.hxx>
-// #include <GeomAdaptor_HCurve.hxx>
-#include <ShapeAnalysis_Surface.hxx>
-// #include <GeomLProp_SLProps.hxx>
-// #include <BRepExtrema_DistShapeShape.hxx>
-// #include <BRepBuilderAPI_MakeVertex.hxx>
-// #include <GCPnts_AbscissaPoint.hxx>
-// #include <GeomLProp_CLProps.hxx>
-#include <deal.II/base/utilities.h>
-#include <deal.II/base/exceptions.h>
-
-#include <boost/bind.hpp>
-
-#include <BRep_Builder.hxx>
+#include <Geom_Plane.hxx>
+#include <Geom_BoundedCurve.hxx>
+#include <GeomAPI_Interpolate.hxx>
+#include <GeomConvert_CompCurveToBSplineCurve.hxx>
+#include <GCPnts_AbscissaPoint.hxx>
+#include <ShapeAnalysis_Surface.hxx>
#include <vector>
#include <algorithm>
namespace OpenCASCADE
{
void count_elements(const TopoDS_Shape &shape,
- unsigned int &n_faces,
- unsigned int &n_edges,
- unsigned int &n_vertices) {
+ unsigned int &n_faces,
+ unsigned int &n_edges,
+ unsigned int &n_vertices)
+ {
TopExp_Explorer exp;
n_faces=0, n_edges=0, n_vertices=0;
- for(exp.Init(shape, TopAbs_FACE); exp.More(); exp.Next(), ++n_faces) {
- }
- for(exp.Init(shape, TopAbs_EDGE); exp.More(); exp.Next(), ++n_edges) {
- }
- for(exp.Init(shape, TopAbs_VERTEX); exp.More(); exp.Next(), ++n_vertices) {
- }
+ for (exp.Init(shape, TopAbs_FACE);
+ exp.More(); exp.Next(), ++n_faces)
+ {}
+ for (exp.Init(shape, TopAbs_EDGE);
+ exp.More(); exp.Next(), ++n_edges)
+ {}
+ for (exp.Init(shape, TopAbs_VERTEX);
+ exp.More(); exp.Next(), ++n_vertices)
+ {}
+ }
+
+ void extract_geometrical_shapes(const TopoDS_Shape &shape,
+ std::vector<TopoDS_Face> &faces,
+ std::vector<TopoDS_Edge> &edges,
+ std::vector<TopoDS_Vertex> &vertices)
+ {
+ faces.resize(0);
+ edges.resize(0);
+ vertices.resize(0);
+
+ TopExp_Explorer exp;
+ for (exp.Init(shape, TopAbs_FACE); exp.More(); exp.Next())
+ {
+ faces.push_back(TopoDS::Face(exp.Current()));
+ }
+ for (exp.Init(shape, TopAbs_EDGE); exp.More(); exp.Next())
+ {
+ edges.push_back(TopoDS::Edge(exp.Current()));
+ }
+ for (exp.Init(shape, TopAbs_VERTEX); exp.More(); exp.Next())
+ {
+ vertices.push_back(TopoDS::Vertex(exp.Current()));
+ }
+ }
+
+
+ void extract_compound_shapes(const TopoDS_Shape &shape,
+ std::vector<TopoDS_Compound> &compounds,
+ std::vector<TopoDS_CompSolid> &compsolids,
+ std::vector<TopoDS_Solid> &solids,
+ std::vector<TopoDS_Shell> &shells,
+ std::vector<TopoDS_Wire> &wires)
+ {
+ compounds.resize(0);
+ compsolids.resize(0);
+ solids.resize(0);
+ shells.resize(0);
+ wires.resize(0);
+
+ TopExp_Explorer exp;
+ for (exp.Init(shape, TopAbs_COMPOUND); exp.More(); exp.Next())
+ {
+ compounds.push_back(TopoDS::Compound(exp.Current()));
+ }
+ for (exp.Init(shape, TopAbs_COMPSOLID); exp.More(); exp.Next())
+ {
+ compsolids.push_back(TopoDS::CompSolid(exp.Current()));
+ }
+ for (exp.Init(shape, TopAbs_SOLID); exp.More(); exp.Next())
+ {
+ solids.push_back(TopoDS::Solid(exp.Current()));
+ }
+ for (exp.Init(shape, TopAbs_SHELL); exp.More(); exp.Next())
+ {
+ shells.push_back(TopoDS::Shell(exp.Current()));
+ }
+ for (exp.Init(shape, TopAbs_WIRE); exp.More(); exp.Next())
+ {
+ wires.push_back(TopoDS::Wire(exp.Current()));
+ }
}
gp_Pnt Pnt(const Point<3> &p)
return gp_Pnt(p(0), p(1), p(2));
}
-
+
Point<3> Pnt(const gp_Pnt &p)
{
return Point<3>(p.X(), p.Y(), p.Z());
}
-
- inline bool point_compare(const dealii::Point<3> &p1, const dealii::Point<3> &p2,
- const dealii::Point<3> direction,
- const double tolerance)
+
+ bool point_compare(const dealii::Point<3> &p1, const dealii::Point<3> &p2,
+ const dealii::Point<3> direction,
+ const double tolerance)
{
const double rel_tol=std::max(tolerance, std::max(p1.norm(), p2.norm())*tolerance);
- if(direction.norm() > 0.0)
+ if (direction.norm() > 0.0)
return (p1*direction < p2*direction-rel_tol);
- else
- for(int d=2; d>=0; --d)
- if(p1[d] < p2[d]-rel_tol)
- return true;
- else if(p2[d] < p1[d]-rel_tol)
- return false;
-
+ else
+ for (int d=2; d>=0; --d)
+ if (p1[d] < p2[d]-rel_tol)
+ return true;
+ else if (p2[d] < p1[d]-rel_tol)
+ return false;
+
// If we got here, for all d, none of the conditions above was
// satisfied. The two points are equal up to tolerance
return false;
}
-
- TopoDS_Shape read_IGES(const std::string &filename,
- const double scale_factor)
+
+ TopoDS_Shape read_IGES(const std::string &filename,
+ const double scale_factor)
{
IGESControl_Reader reader;
IFSelect_ReturnStatus stat;
- stat = reader.ReadFile(filename.c_str());
- Assert(stat == IFSelect_RetDone,
- ExcMessage("Error in reading file!"));
-
+ stat = reader.ReadFile(filename.c_str());
+ Assert(stat == IFSelect_RetDone,
+ ExcMessage("Error in reading file!"));
+
Standard_Boolean failsonly = Standard_False;
IFSelect_PrintCount mode = IFSelect_ItemsByEntity;
reader.PrintCheckLoad (failsonly, mode);
Handle(TColStd_HSequenceOfTransient) myList = reader.GiveList("iges-faces");
- //selects all IGES faces in the
- //file and puts them into a list
- //called MyList,
- Standard_Integer // nIgesFaces = myList->Length(),
- nTransFaces = reader.TransferList(myList);
+ //selects all IGES faces in the
+ //file and puts them into a list
+ //called MyList,
+ Standard_Integer nTransFaces = reader.TransferList(myList);
- AssertThrow(nTransFaces > 0,
- ExcMessage("Read nothing from file."));
+ AssertThrow(nTransFaces > 0,
+ ExcMessage("Read nothing from file."));
- // Handle IGES Scale here.
+ // Handle IGES Scale here.
gp_Pnt Origin;
gp_Trsf scale;
scale.SetScale (Origin, scale_factor);
TopoDS_Shape sh = reader.OneShape();
BRepBuilderAPI_Transform trans(sh, scale);
-
+
return trans.Shape(); // this is the actual translation
}
-
- void write_IGES(const TopoDS_Shape &shape,
- const std::string &filename) {
-
+
+ void write_IGES(const TopoDS_Shape &shape,
+ const std::string &filename)
+ {
+
IGESControl_Controller::Init();
IGESControl_Writer ICW ("MM", 0);
Standard_Boolean ok = ICW.AddShape (shape);
Standard_Boolean OK = ICW.Write (filename.c_str());
Assert(OK, ExcMessage("Failed to write IGES file."));
}
-
-
- TopoDS_Shape intersect_plane(const TopoDS_Shape &in_shape,
- const double c_x,
- const double c_y,
- const double c_z,
- const double c,
- const double tolerance)
- {
- TopoDS_Shape out_shape;
- Handle(Geom_Plane) plane = new Geom_Plane(c_x,c_y,c_z,c);
-
- // This is to loop on the faces,
- // that extracts all
- // intersections.
-
- //TopExp_Explorer faceExplorer(in_shape , TopAbs_FACE);
-
- std::vector<Handle_Geom_BoundedCurve> intersections;
-
- BRepAlgo_Section section(in_shape, plane);
- TopoDS_Shape edges = section.Shape();
-
- TopoDS_Edge edge;
- TopLoc_Location L;
- Standard_Real First;
- Standard_Real Last;
- gp_Pnt PIn(0.0,0.0,0.0);
- gp_Pnt PFin(0.0,0.0,0.0);
- gp_Pnt PMid(0.0,0.0,0.0);
- TopExp_Explorer edgeExplorer(edges , TopAbs_EDGE);
- while (edgeExplorer.More())
- {
- edge = TopoDS::Edge(edgeExplorer.Current());
- Handle(Geom_Curve) curve = BRep_Tool::Curve(edge,L,First,Last);
- intersections.push_back(Handle(Geom_BoundedCurve)::DownCast(curve));
- edgeExplorer.Next();
- }
-
- // Now we build a single bspline out of all the geometrical
- // curves
- unsigned int numIntersEdges = intersections.size();
- for (unsigned int i = 0;i<intersections.size();++i)
- {
- if (intersections[i]->Value(intersections[i]->FirstParameter()).X() >
- intersections[i]->Value(intersections[i]->LastParameter()).X() )
- intersections[i]->Reverse();
- }
-
- GeomConvert_CompCurveToBSplineCurve
- convert_bspline(intersections[0], Convert_TgtThetaOver2);
- bool check = false, one_added = true, one_failed=true;
- std::vector<bool> added(numIntersEdges, false);
- added[0] = true;
- while(one_added == true)
- {
- one_added = false;
- one_failed = false;
- for (unsigned int i=1; i<numIntersEdges; ++i)
- if(added[i] == false)
- {
- Handle(Geom_Curve) curve = intersections[i];
- Handle(Geom_BoundedCurve) bcurve = Handle(Geom_BoundedCurve)::DownCast(curve);
- check = convert_bspline.Add(bcurve,tolerance,0,1,0);
- one_failed = one_failed || (check == false);
- one_added = one_added || (check == true);
- added[i] = check;
- //cout<<i<<" --> "<<added[i]<<" "<<false<<endl;
- }
- }
-
- Assert(one_failed == false,
- ExcMessage("Bspline convertion of intersection with plane has failed."));
-
- Handle(Geom_Curve) bspline = convert_bspline.BSplineCurve();
-
- out_shape = BRepBuilderAPI_MakeEdge(bspline);
-
- if (bspline->IsCN(1))
- cout<<"Intersection with plane is at least a C1 curve"<<endl;
- else
- cout<<"Intersection with plane is not a C1 curve "<<endl;
- return out_shape;
- }
-
-
- Point<3> axis_intersection(const TopoDS_Shape in_shape,
- const Point<3> origin,
- const Point<3> direction,
- const double tolerance)
+
+
+ TopoDS_Shape intersect_plane(const TopoDS_Shape &in_shape,
+ const double c_x,
+ const double c_y,
+ const double c_z,
+ const double c,
+ const double tolerance)
+ {
+ Handle(Geom_Plane) plane = new Geom_Plane(c_x,c_y,c_z,c);
+ BRepAlgo_Section section(in_shape, plane);
+ TopoDS_Shape edges = section.Shape();
+ return edges;
+ }
+
+ TopoDS_Edge join_edges(const TopoDS_Shape &in_shape,
+ const double tolerance)
+ {
+ TopoDS_Edge out_shape;
+ TopoDS_Shape edges = in_shape;
+ std::vector<Handle_Geom_BoundedCurve> intersections;
+ TopLoc_Location L;
+ Standard_Real First;
+ Standard_Real Last;
+ gp_Pnt PIn(0.0,0.0,0.0);
+ gp_Pnt PFin(0.0,0.0,0.0);
+ gp_Pnt PMid(0.0,0.0,0.0);
+ TopExp_Explorer edgeExplorer(edges , TopAbs_EDGE);
+ TopoDS_Edge edge;
+ while (edgeExplorer.More())
+ {
+ edge = TopoDS::Edge(edgeExplorer.Current());
+ Handle(Geom_Curve) curve = BRep_Tool::Curve(edge,L,First,Last);
+ intersections.push_back(Handle(Geom_BoundedCurve)::DownCast(curve));
+ edgeExplorer.Next();
+ }
+
+ // Now we build a single bspline out of all the geometrical
+ // curves, in Lexycographical order
+ unsigned int numIntersEdges = intersections.size();
+ Assert(numIntersEdges>0, ExcMessage("No curves to process!"));
+
+ GeomConvert_CompCurveToBSplineCurve convert_bspline(intersections[0]);
+
+ bool check = false, one_added = true, one_failed=true;
+ std::vector<bool> added(numIntersEdges, false);
+ added[0] = true;
+ while (one_added == true)
+ {
+ one_added = false;
+ one_failed = false;
+ for (unsigned int i=1; i<numIntersEdges; ++i)
+ if (added[i] == false)
+ {
+ Handle(Geom_Curve) curve = intersections[i];
+ Handle(Geom_BoundedCurve) bcurve = Handle(Geom_BoundedCurve)::DownCast(curve);
+ check = convert_bspline.Add(bcurve,tolerance,0,1,0);
+ if (check == false) // If we failed, try again with the reversed curve
+ {
+ curve->Reverse();
+ Handle(Geom_BoundedCurve) bcurve = Handle(Geom_BoundedCurve)::DownCast(curve);
+ check = convert_bspline.Add(bcurve,tolerance,0,1,0);
+ }
+ one_failed = one_failed || (check == false);
+ one_added = one_added || (check == true);
+ added[i] = check;
+ }
+ }
+
+ Assert(one_failed == false,
+ ExcMessage("Joining some of the Edges failed."));
+
+ Handle(Geom_Curve) bspline = convert_bspline.BSplineCurve();
+
+ out_shape = BRepBuilderAPI_MakeEdge(bspline);
+ return out_shape;
+ }
+
+
+ Point<3> axis_intersection(const TopoDS_Shape in_shape,
+ const Point<3> origin,
+ const Point<3> direction,
+ const double tolerance)
{
// translating original Point<dim> to gp point
// shape from the iges model
IntCurvesFace_ShapeIntersector Inters;
Inters.Load(in_shape,tolerance);
-
+
Inters.PerformNearest(line,-RealLast(),+RealLast());
Assert(Inters.IsDone(), ExcMessage("Could not project point."));
return Pnt(Inters.Pnt(1));
}
TopoDS_Edge interpolation_curve(std::vector<Point<3> > &curve_points,
- const Point<3> direction,
- const bool closed,
- const double tolerance)
+ const Point<3> direction,
+ const bool closed,
+ const double tolerance)
{
unsigned int n_vertices = curve_points.size();
- if(direction*direction > 0)
+ if (direction*direction > 0)
{
- std::sort(curve_points.begin(), curve_points.end(),
- boost::bind(&OpenCASCADE::point_compare, _1, _2, direction, tolerance));
+ std::sort(curve_points.begin(), curve_points.end(),
+ boost::bind(&OpenCASCADE::point_compare, _1, _2, direction, tolerance));
}
- // set up array of vertices
+ // set up array of vertices
Handle(TColgp_HArray1OfPnt) vertices = new TColgp_HArray1OfPnt(1,n_vertices);
for (unsigned int vertex=0; vertex<n_vertices; ++vertex)
{
- vertices->SetValue(vertex+1,Pnt(curve_points[vertex]));
+ vertices->SetValue(vertex+1,Pnt(curve_points[vertex]));
}
GeomAPI_Interpolate bspline_generator(vertices, closed, tolerance);
bspline_generator.Perform();
Assert( (bspline_generator.IsDone()), ExcMessage("Interpolated bspline generation failed"));
-
+
Handle(Geom_BSplineCurve) bspline = bspline_generator.Curve();
TopoDS_Edge out_shape = BRepBuilderAPI_MakeEdge(bspline);
return out_shape;
}
- Point<3> closest_point(const TopoDS_Shape in_shape,
- const Point<3> origin,
- TopoDS_Shape &out_shape,
- double &u,
- double &v,
- const double tolerance) {
-
+ Point<3> closest_point(const TopoDS_Shape in_shape,
+ const Point<3> origin,
+ TopoDS_Shape &out_shape,
+ double &u,
+ double &v,
+ const double tolerance)
+ {
+
TopExp_Explorer exp;
gp_Pnt Pproj = Pnt(origin);
double minDistance = 1e7;
gp_Pnt tmp_proj(0.0,0.0,0.0);
-
+
unsigned int counter = 0;
- u=0; v=0;
-
- for(exp.Init(in_shape, TopAbs_FACE); exp.More(); exp.Next()) {
- TopoDS_Face face = TopoDS::Face(exp.Current());
-
- // the projection function needs a surface, so we obtain the
- // surface upon which the face is defined
- Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
-
- ShapeAnalysis_Surface projector(SurfToProj);
- gp_Pnt2d proj_params = projector.ValueOfUV(Pnt(origin), tolerance);
-
- SurfToProj->D0(proj_params.X(),proj_params.Y(),tmp_proj);
-
- double distance = Pnt(tmp_proj).distance(origin);
- if (distance < minDistance)
- {
- minDistance = distance;
- Pproj = tmp_proj;
- out_shape = face;
- u=proj_params.X();
- v=proj_params.Y();
- ++counter;
- }
- }
-
- for(exp.Init(in_shape, TopAbs_EDGE); exp.More(); exp.Next()) {
- TopoDS_Edge edge = TopoDS::Edge(exp.Current());
-
- TopLoc_Location L;
- Standard_Real First;
- Standard_Real Last;
-
- // the projection function needs a Curve, so we obtain the
- // curve upon which the edge is defined
- Handle(Geom_Curve) CurveToProj = BRep_Tool::Curve(edge,L,First,Last);
-
- GeomAPI_ProjectPointOnCurve Proj(Pnt(origin),CurveToProj);
- unsigned int num_proj_points = Proj.NbPoints();
- if ((num_proj_points > 0) && (Proj.LowerDistance() < minDistance))
- {
- minDistance = Proj.LowerDistance();
- Pproj = Proj.NearestPoint();
- out_shape = edge;
- u=Proj.LowerDistanceParameter();
- ++counter;
- }
- }
-
+ u=0;
+ v=0;
+
+ for (exp.Init(in_shape, TopAbs_FACE); exp.More(); exp.Next())
+ {
+ TopoDS_Face face = TopoDS::Face(exp.Current());
+
+ // the projection function needs a surface, so we obtain the
+ // surface upon which the face is defined
+ Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
+
+ ShapeAnalysis_Surface projector(SurfToProj);
+ gp_Pnt2d proj_params = projector.ValueOfUV(Pnt(origin), tolerance);
+
+ SurfToProj->D0(proj_params.X(),proj_params.Y(),tmp_proj);
+
+ double distance = Pnt(tmp_proj).distance(origin);
+ if (distance < minDistance)
+ {
+ minDistance = distance;
+ Pproj = tmp_proj;
+ out_shape = face;
+ u=proj_params.X();
+ v=proj_params.Y();
+ ++counter;
+ }
+ }
+
+ for (exp.Init(in_shape, TopAbs_EDGE); exp.More(); exp.Next())
+ {
+ TopoDS_Edge edge = TopoDS::Edge(exp.Current());
+ if (!BRep_Tool::Degenerated(edge))
+ {
+ TopLoc_Location L;
+ Standard_Real First;
+ Standard_Real Last;
+
+ // the projection function needs a Curve, so we obtain the
+ // curve upon which the edge is defined
+ Handle(Geom_Curve) CurveToProj = BRep_Tool::Curve(edge,L,First,Last);
+
+ GeomAPI_ProjectPointOnCurve Proj(Pnt(origin),CurveToProj);
+ unsigned int num_proj_points = Proj.NbPoints();
+ if ((num_proj_points > 0) && (Proj.LowerDistance() < minDistance))
+ {
+ minDistance = Proj.LowerDistance();
+ Pproj = Proj.NearestPoint();
+ out_shape = edge;
+ u=Proj.LowerDistanceParameter();
+ ++counter;
+ }
+ }
+ }
+
Assert(counter > 0, ExcMessage("Could not find projection points."));
return Pnt(Pproj);
}
+ void create_triangulation(const TopoDS_Face &face,
+ Triangulation<2,3> &tria)
+ {
+ BRepAdaptor_Surface surf(face);
+ const double u0 = surf.FirstUParameter();
+ const double u1 = surf.LastUParameter();
+ const double v0 = surf.FirstVParameter();
+ const double v1 = surf.LastVParameter();
+
+ std::vector<CellData<2> > cells;
+ std::vector<Point<3> > vertices;
+ SubCellData t;
+
+ vertices.push_back(Pnt(surf.Value(u0,v0)));
+ vertices.push_back(Pnt(surf.Value(u1,v0)));
+ vertices.push_back(Pnt(surf.Value(u0,v1)));
+ vertices.push_back(Pnt(surf.Value(u1,v1)));
+
+ CellData<2> cell;
+ for (unsigned int i=0; i<4; ++i)
+ cell.vertices[i] = i;
+
+ cells.push_back(cell);
+ tria.create_triangulation(vertices, cells, t);
+ }
+
} // end namespace
DEAL_II_NAMESPACE_CLOSE
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <deal.II/grid/occ_utilities.h>
+#include <deal.II/grid/occ_boundary_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <TopTools.hxx>
+#include <TopoDS_Shape.hxx>
+#include <Standard_Stream.hxx>
+
+// Create a Triangulation, interpolate its boundary points to a smooth
+// BSpline, and use that as an arlength Boundary Descriptor.
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+
+ // Create a bspline passign 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 edge = interpolation_curve(pts);
+ ArclengthProjectionLineManifold<1,3> manifold(edge);
+
+ Triangulation<1,3> tria;
+ GridGenerator::hyper_cube(tria);
+
+ tria.begin_active()->set_all_manifold_ids(1);
+ tria.set_manifold(1, manifold);
+
+ tria.refine_global(4);
+ GridOut gridout;
+ gridout.write_msh(tria, logfile);
+
+ return 0;
+}
+
--- /dev/null
+$NOD
+17
+1 0 0 0
+2 1 0 0
+3 0.5 1.125 0
+4 -0.213053 0.770812 0
+5 1.21305 0.770812 0
+6 -0.209733 0.361552 0
+7 0.090632 1.04613 0
+8 0.909368 1.04613 0
+9 1.20973 0.361552 0
+10 -0.118931 0.172725 0
+11 -0.251438 0.566341 0
+12 -0.0875601 0.936891 0
+13 0.291313 1.106 0
+14 0.708687 1.106 0
+15 1.08756 0.936891 0
+16 1.25144 0.566341 0
+17 1.11893 0.172725 0
+$ENDNOD
+$ELM
+16
+1 1 0 0 2 1 10
+2 1 0 0 2 10 6
+3 1 0 0 2 6 11
+4 1 0 0 2 11 4
+5 1 0 0 2 4 12
+6 1 0 0 2 12 7
+7 1 0 0 2 7 13
+8 1 0 0 2 13 3
+9 1 0 0 2 3 14
+10 1 0 0 2 14 8
+11 1 0 0 2 8 15
+12 1 0 0 2 15 5
+13 1 0 0 2 5 16
+14 1 0 0 2 16 9
+15 1 0 0 2 9 17
+16 1 0 0 2 17 2
+$ENDELM
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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/grid/occ_utilities.h>
+#include <deal.II/grid/occ_boundary_lib.h>
+
+#include <fstream>
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/grid_out.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_boundary_lib.h>
+
+#include <gp_Pnt.hxx>
+#include <gp_Dir.hxx>
+#include <gp_Ax2.hxx>
+#include <GC_MakeCircle.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+
+using namespace OpenCASCADE;
+
+int main()
+{
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ // A curve passing through the vertices of the unit square.
+ std::vector<Point<3> > vertices;
+ vertices.push_back(Point<3>(0.5,.5-std::sqrt(2.0)/2.0,0));
+ vertices.push_back(Point<3>(1,0,0));
+ vertices.push_back(Point<3>(1,1,0));
+ vertices.push_back(Point<3>(0,1,0));
+ vertices.push_back(Point<3>(0,0,0));
+ TopoDS_Edge edge = interpolation_curve(vertices, Point<3>(), true);
+
+ // Create a boundary projector.
+ ArclengthProjectionLineManifold<2,3> boundary_line(edge);
+
+ // The unit square.
+ Triangulation<2,3> tria;
+ GridGenerator::hyper_cube(tria);
+
+ // Set the exterior boundary
+ tria.set_manifold(1, boundary_line);
+
+ // This is here to ignore the points created in the interior of the
+ // face.
+ tria.begin()->set_all_manifold_ids(1);
+ tria.begin()->set_manifold_id(0);
+
+ // We refine twice, and expect the outer points to end up on a
+ // smooth curve interpolating the square vertices.
+ tria.refine_global(2);
+
+
+ // You can open the generated file with gmsh
+ GridOut gridout;
+ gridout.write_msh (tria, logfile);
+ return 0;
+}
--- /dev/null
+
+$NOD
+25
+1 0.00000 0.00000 0.00000
+2 1.00000 0.00000 0.00000
+3 0.00000 1.00000 0.00000
+4 1.00000 1.00000 0.00000
+5 0.506410 -0.207730 0.00000
+6 -0.166677 0.501960 0.00000
+7 1.16638 0.503363 0.00000
+8 0.499605 1.18595 0.00000
+9 0.500715 0.497943 0.00000
+10 0.231188 -0.153665 0.00000
+11 0.782598 -0.172689 0.00000
+12 -0.124206 0.236965 0.00000
+13 -0.131434 0.767789 0.00000
+14 1.12068 0.239092 0.00000
+15 1.13309 0.769132 0.00000
+16 0.230639 1.14264 0.00000
+17 0.768406 1.14135 0.00000
+18 0.503562 0.145107 0.00000
+19 0.500160 0.841948 0.00000
+20 0.167019 0.499952 0.00000
+21 0.833548 0.500653 0.00000
+22 0.202251 0.190066 0.00000
+23 0.801737 0.188217 0.00000
+24 0.200003 0.804773 0.00000
+25 0.800238 0.805043 0.00000
+$ENDNOD
+$ELM
+16
+1 3 0 0 4 1 10 22 12
+2 3 0 0 4 10 5 18 22
+3 3 0 0 4 12 22 20 6
+4 3 0 0 4 22 18 9 20
+5 3 0 0 4 5 11 23 18
+6 3 0 0 4 11 2 14 23
+7 3 0 0 4 18 23 21 9
+8 3 0 0 4 23 14 7 21
+9 3 0 0 4 6 20 24 13
+10 3 0 0 4 20 9 19 24
+11 3 0 0 4 13 24 16 3
+12 3 0 0 4 24 19 8 16
+13 3 0 0 4 9 21 25 19
+14 3 0 0 4 21 7 15 25
+15 3 0 0 4 19 25 17 8
+16 3 0 0 4 25 15 4 17
+$ENDELM
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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 BSpline surface, and test axis projection.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <deal.II/grid/occ_utilities.h>
+#include <deal.II/grid/occ_boundary_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <TopTools.hxx>
+#include <TopoDS_Shape.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+#include <BRepFill.hxx>
+#include <Standard_Stream.hxx>
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+
+ // Create a bspline passign 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);
+
+
+ DirectionalProjectionBoundary<1,3> manifold(face, Point<3>(0,1,0));
+
+ Triangulation<1,3> tria;
+ GridGenerator::hyper_cube(tria);
+
+ tria.begin_active()->set_all_manifold_ids(1);
+ tria.set_manifold(1, manifold);
+
+ tria.refine_global(4);
+ GridOut gridout;
+ gridout.write_msh(tria, logfile);
+
+ return 0;
+}
--- /dev/null
+$NOD
+17
+1 0 0 0
+2 1 0 0
+3 0.5 1.125 0
+4 0.25 1.09745 0
+5 0.75 1.09745 0
+6 0.125 1.05999 0
+7 0.375 1.11829 0
+8 0.625 1.11829 0
+9 0.875 1.05999 0
+10 0.0625 1.03342 0
+11 0.1875 1.08105 0
+12 0.3125 1.10974 0
+13 0.4375 1.12333 0
+14 0.5625 1.12333 0
+15 0.6875 1.10974 0
+16 0.8125 1.08105 0
+17 0.9375 1.03342 0
+$ENDNOD
+$ELM
+16
+1 1 0 0 2 1 10
+2 1 0 0 2 10 6
+3 1 0 0 2 6 11
+4 1 0 0 2 11 4
+5 1 0 0 2 4 12
+6 1 0 0 2 12 7
+7 1 0 0 2 7 13
+8 1 0 0 2 13 3
+9 1 0 0 2 3 14
+10 1 0 0 2 14 8
+11 1 0 0 2 8 15
+12 1 0 0 2 15 5
+13 1 0 0 2 5 16
+14 1 0 0 2 16 9
+15 1 0 0 2 9 17
+16 1 0 0 2 17 2
+$ENDELM
--- /dev/null
+
+//---------------------------- iges_read.cc ---------------------------
+// $Id: testsuite.html 13373 2006-07-13 13:12:08Z kanschat $
+// Version: $Name$
+//
+// Copyright (C) 2005 by the deal.II authors
+//
+// This file is subject to QPL 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.
+//
+//---------------------------- iges_read.cc ---------------------------
+
+
+// Read goteborg.iges and dump its topological structure to the logfile.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <deal.II/grid/occ_utilities.h>
+#include <deal.II/grid/occ_boundary_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <TopTools.hxx>
+#include <TopoDS_Shape.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+#include <BRepFill.hxx>
+#include <Standard_Stream.hxx>
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+
+ // 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(face, Point<3>(0,0,1));
+
+ Triangulation<2,3> tria;
+ GridGenerator::hyper_cube(tria);
+
+ tria.begin_active()->set_all_manifold_ids(1);
+ tria.set_manifold(1, manifold);
+
+ tria.refine_global(2);
+ GridOut gridout;
+ gridout.write_msh(tria, logfile);
+
+ return 0;
+}
--- /dev/null
+$NOD
+25
+1 0 0 0
+2 1 0 0
+3 0 1 0
+4 1 1 0
+5 0.5 0 0
+6 0 0.5 0
+7 1 0.5 0
+8 0.5 1 1
+9 0.5 0.5 0.5
+10 0.25 0 0
+11 0.75 0 0
+12 0 0.25 0
+13 0 0.75 0
+14 1 0.25 0
+15 1 0.75 0
+16 0.25 1 0.75
+17 0.75 1 0.75
+18 0.5 0.25 0.25
+19 0.5 0.75 0.75
+20 0.25 0.5 0.375
+21 0.75 0.5 0.375
+22 0.25 0.25 0.1875
+23 0.75 0.25 0.1875
+24 0.25 0.75 0.5625
+25 0.75 0.75 0.5625
+$ENDNOD
+$ELM
+16
+1 3 0 0 4 1 10 22 12
+2 3 0 0 4 10 5 18 22
+3 3 0 0 4 12 22 20 6
+4 3 0 0 4 22 18 9 20
+5 3 0 0 4 5 11 23 18
+6 3 0 0 4 11 2 14 23
+7 3 0 0 4 18 23 21 9
+8 3 0 0 4 23 14 7 21
+9 3 0 0 4 6 20 24 13
+10 3 0 0 4 20 9 19 24
+11 3 0 0 4 13 24 16 3
+12 3 0 0 4 24 19 8 16
+13 3 0 0 4 9 21 25 19
+14 3 0 0 4 21 7 15 25
+15 3 0 0 4 19 25 17 8
+16 3 0 0 4 25 15 4 17
+$ENDELM
-
-//---------------------------- occ_circle_projection.cc ---------------------------
-// $Id: testsuite.html 13373 2006-07-13 13:12:08Z kanschat $
-// Version: $Name$
+//-----------------------------------------------------------
//
-// Copyright (C) 2005 by the deal.II authors
+// Copyright (C) 2014 by the deal.II authors
//
-// This file is subject to QPL and may not be distributed
+// 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.
//
-//---------------------------- occ_circle_projection.cc ---------------------------
-
+//-----------------------------------------------------------
// Create a circle, a Triangulation, and try to project normally on
// it.
// Set the exterior boundary
tria.set_boundary(0, boundary_line);
- // Uncomment the following line to
- // see the difference.
-// tria.set_boundary(0, boundary_line_deal);
-
// This is here to ignore the
// points created in the interior
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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.
+//
+//-----------------------------------------------------------
+
// Makes a OpenCASCADE cicrular arc, and project a few points onto it.
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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.
+//
+//-----------------------------------------------------------
+
+// Read the first face of goteborg.iges, attach it to a projector,
+// create a single-cell Triangulation, and refine it with the new
+// projector.
+
+#include "../tests.h"
+
+#include <deal.II/grid/occ_boundary_lib.h>
+
+#include <fstream>
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/grid_out.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_boundary_lib.h>
+
+#include <gp_Pnt.hxx>
+#include <gp_Dir.hxx>
+#include <gp_Ax2.hxx>
+#include <GC_MakeCircle.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+#include <BRepPrimAPI_MakeSphere.hxx>
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+
+ TopoDS_Shape sh = read_IGES(SOURCE_DIR "/iges_files/goteborg.iges");
+ std::vector<TopoDS_Face> faces;
+ std::vector<TopoDS_Edge> edges;
+ std::vector<TopoDS_Vertex> vertices;
+
+ extract_geometrical_shapes(sh, faces, edges, vertices);
+
+ // Create a boundary projector on the first face.
+ NormalProjectionBoundary<2,3> boundary(faces[0], 1e-9);
+
+ // Create a Triangulation with a single cell
+ Triangulation<2,3> tria;
+ create_triangulation(faces[0], tria);
+
+ // Set the boundary
+ tria.set_manifold(1, boundary);
+ tria.begin()->set_all_manifold_ids(1);
+
+ tria.refine_global(2);
+
+ // You can open the generated file with gmsh
+ GridOut gridout;
+ gridout.write_msh (tria, logfile);
+
+ return 0;
+}
+
--- /dev/null
+$NOD
+25
+1 -1.98714 1.02631e-19 -0.149998
+2 -2.63758 -6.93889e-18 -0.249603
+3 -2.00886 0.314666 0.428234
+4 -3.05955 0 0.501418
+5 -2.30472 0 -0.245614
+6 -2.00073 0.167844 0.134255
+7 -2.72534 4.25119e-05 0.204053
+8 -2.54605 0.193924 0.440052
+9 -2.40412 0.0655112 0.139204
+10 -2.14241 8.03383e-19 -0.216122
+11 -2.46936 -3.70788e-20 -0.267918
+12 -1.99609 0.0965804 -0.0146004
+13 -2.00426 0.238868 0.282461
+14 -2.63221 9.34116e-07 -0.0158842
+15 -2.87756 -5.67114e-09 0.370719
+16 -2.27934 0.26212 0.42932
+17 -2.80896 0.111433 0.459585
+18 -2.35424 0.031169 -0.0529932
+19 -2.47235 0.120909 0.295137
+20 -2.20235 0.116409 0.136836
+21 -2.56575 0.0372536 0.170331
+22 -2.17413 0.0607452 -0.0340272
+23 -2.51144 0.0146272 -0.0396781
+24 -2.23928 0.183308 0.286763
+25 -2.68303 0.0677305 0.321656
+$ENDNOD
+$ELM
+16
+1 3 0 0 4 1 10 22 12
+2 3 0 0 4 10 5 18 22
+3 3 0 0 4 12 22 20 6
+4 3 0 0 4 22 18 9 20
+5 3 0 0 4 5 11 23 18
+6 3 0 0 4 11 2 14 23
+7 3 0 0 4 18 23 21 9
+8 3 0 0 4 23 14 7 21
+9 3 0 0 4 6 20 24 13
+10 3 0 0 4 20 9 19 24
+11 3 0 0 4 13 24 16 3
+12 3 0 0 4 24 19 8 16
+13 3 0 0 4 9 21 25 19
+14 3 0 0 4 21 7 15 25
+15 3 0 0 4 19 25 17 8
+16 3 0 0 4 25 15 4 17
+$ENDELM
-
-//---------------------------- iges_read.cc ---------------------------
-// $Id: testsuite.html 13373 2006-07-13 13:12:08Z kanschat $
-// Version: $Name$
+//-----------------------------------------------------------
//
-// Copyright (C) 2005 by the deal.II authors
+// Copyright (C) 2014 by the deal.II authors
//
-// This file is subject to QPL and may not be distributed
+// 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.
//
-//---------------------------- iges_read.cc ---------------------------
-
+//-----------------------------------------------------------
-// Read goteborg.iges and dump its topological structure to the logfile.
+// Create a bspline, and write it as an IGES file.
#include "../tests.h"
#include <fstream>
int main ()
{
- // std::ofstream logfile("output");
- // deallog.attach(logfile);
- // deallog.depth_console(0);
-
// Create a bspline passign through the points
std::vector<Point<3> > pts;
pts.push_back(Point<3>(0,0,0));
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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.
+//
+//-----------------------------------------------------------
+
+// Read goteborg.iges and dump its topological structure to the
+// logfile.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <deal.II/grid/occ_utilities.h>
+#include <TopTools.hxx>
+#include <TopoDS_Shape.hxx>
+#include <Standard_Stream.hxx>
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ TopoDS_Shape sh = read_IGES(SOURCE_DIR "/iges_files/goteborg.iges");
+ std::vector<TopoDS_Face> faces;
+ std::vector<TopoDS_Edge> edges;
+ std::vector<TopoDS_Vertex> vertices;
+
+ extract_geometrical_shapes(sh, faces, edges, vertices);
+
+ unsigned int nf=0, ne=0, nv=0;
+ count_elements(sh, nf, ne, nv);
+ deallog << "Shape contains " << nf << " faces, "
+ << ne << " edges, and "
+ << nv << " vertices." << std::endl;
+
+ if(nf != faces.size())
+ deallog << "Error!" << std::endl;
+ if(ne != edges.size())
+ deallog << "Error!" << std::endl;
+ if(nv != vertices.size())
+ deallog << "Error!" << std::endl;
+
+ return 0;
+}
+
--- /dev/null
+
+DEAL::Shape contains 3 faces, 12 edges, and 24 vertices.
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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.
+//
+//-----------------------------------------------------------
+
+// Read goteborg.iges and dump its topological structure to the
+// logfile.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <deal.II/grid/occ_utilities.h>
+#include <TopTools.hxx>
+#include <TopoDS_Shape.hxx>
+#include <Standard_Stream.hxx>
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ TopoDS_Shape sh = read_IGES(SOURCE_DIR "/iges_files/goteborg.iges");
+ std::vector<TopoDS_Compound> compounds;
+ std::vector<TopoDS_CompSolid> compsolids;
+ std::vector<TopoDS_Solid> solids;
+ std::vector<TopoDS_Shell> shells;
+ std::vector<TopoDS_Wire> wires;
+ extract_compound_shapes(sh, compounds, compsolids, solids, shells, wires);
+
+ deallog << "Shape contains " << compounds.size() << " compounds, "
+ << compsolids.size() << " compound solids, "
+ << solids.size() << " solids, "
+ << shells.size() << " shells, "
+ << wires.size() << " wires." << std::endl;
+ return 0;
+}
+
--- /dev/null
+
+DEAL::Shape contains 1 compounds, 0 compound solids, 0 solids, 0 shells, 3 wires.
--- /dev/null
+ S0000001
+,,31HOpen CASCADE IGES processor 6.7,13HFilename.iges, G0000001
+16HOpen CASCADE 6.7,31HOpen CASCADE IGES processor 6.7,32,308,15,308,15,G0000002
+,1.,1,4HINCH,1,0.000394,15H20140814.114632,3.93701E-09,0.049421, G0000003
+6Hheltai,,11,0,15H20140814.114632,; G0000004
+ 144 1 0 0 0 0 0 000000000D0000001
+ 144 0 0 1 0 0D0000002
+ 122 2 0 0 0 0 0 000010000D0000003
+ 122 0 0 1 0 0D0000004
+ 126 3 0 0 0 0 0 000010000D0000005
+ 126 0 0 6 0 0D0000006
+ 142 9 0 0 0 0 0 000010500D0000007
+ 142 0 0 1 0 0D0000008
+ 102 10 0 0 0 0 0 000010000D0000009
+ 102 0 0 1 0 0D0000010
+ 110 11 0 0 0 0 0 000010000D0000011
+ 110 0 0 1 0 0D0000012
+ 110 12 0 0 0 0 0 000010000D0000013
+ 110 0 0 1 0 0D0000014
+ 110 13 0 0 0 0 0 000010000D0000015
+ 110 0 0 1 0 0D0000016
+ 110 14 0 0 0 0 0 000010000D0000017
+ 110 0 0 1 0 0D0000018
+ 102 15 0 0 0 0 0 000010000D0000019
+ 102 0 0 1 0 0D0000020
+ 110 16 0 0 0 0 0 000010000D0000021
+ 110 0 0 1 0 0D0000022
+ 126 17 0 0 0 0 0 000010000D0000023
+ 126 0 0 6 0 0D0000024
+ 110 23 0 0 0 0 0 000010000D0000025
+ 110 0 0 2 0 0D0000026
+ 126 25 0 0 0 0 0 000010000D0000027
+ 126 0 0 5 0 0D0000028
+144,3,1,0,7; 0000001P0000001
+122,5,0.,0.,0.; 0000003P0000002
+126,5,3,1,0,1,0,0.,0.,0.,0.,1.,2.,3.,3.,3.,3.,1.,1.,1.,1.,1.,1., 0000005P0000003
+0.,0.,3.937007874E-02,-1.531058618E-02,1.968503937E-02, 0000005P0000004
+3.937007874E-02,-6.56167979E-03,4.593175853E-02,3.937007874E-02, 0000005P0000005
+4.593175853E-02,4.593175853E-02,3.937007874E-02,5.468066492E-02, 0000005P0000006
+1.968503937E-02,3.937007874E-02,3.937007874E-02,0., 0000005P0000007
+3.937007874E-02,0.,3.,1.602995811E-17,3.205991622E-17,1.; 0000005P0000008
+142,0,3,9,19,3; 0000007P0000009
+102,4,11,13,15,17; 0000009P0000010
+110,0.,1.,0.,0.,0.,0.; 0000011P0000011
+110,0.,0.,0.,1.,0.,0.; 0000013P0000012
+110,1.,0.,0.,1.,1.,0.; 0000015P0000013
+110,1.,1.,0.,0.,1.,0.; 0000017P0000014
+102,4,21,23,25,27; 0000019P0000015
+110,0.,0.,0.,0.,0.,3.937007874E-02; 0000021P0000016
+126,5,3,1,0,1,0,0.,0.,0.,0.,1.,2.,3.,3.,3.,3.,1.,1.,1.,1.,1.,1., 0000023P0000017
+0.,0.,3.937007874E-02,-1.531058618E-02,1.968503937E-02, 0000023P0000018
+3.937007874E-02,-6.56167979E-03,4.593175853E-02,3.937007874E-02, 0000023P0000019
+4.593175853E-02,4.593175853E-02,3.937007874E-02,5.468066492E-02, 0000023P0000020
+1.968503937E-02,3.937007874E-02,3.937007874E-02,0., 0000023P0000021
+3.937007874E-02,0.,3.,1.602995811E-17,3.205991622E-17,1.; 0000023P0000022
+110,3.937007874E-02,2.185478395E-18,3.937007874E-02, 0000025P0000023
+3.937007874E-02,2.185478395E-18,0.; 0000025P0000024
+126,5,3,1,0,1,0,0.,0.,0.,0.,1.,2.,3.,3.,3.,3.,1.,1.,1.,1.,1.,1., 0000027P0000025
+3.937007874E-02,0.,0.,5.468066492E-02,1.968503937E-02,0., 0000027P0000026
+4.593175853E-02,4.593175853E-02,0.,-6.56167979E-03, 0000027P0000027
+4.593175853E-02,0.,-1.531058618E-02,1.968503937E-02,0.,0.,0.,0., 0000027P0000028
+0.,3.,0.,0.,1.; 0000027P0000029
+S 1G 4D 28P 29 T0000001
-
-//---------------------------- iges_read.cc ---------------------------
-// $Id: testsuite.html 13373 2006-07-13 13:12:08Z kanschat $
-// Version: $Name$
+//-----------------------------------------------------------
//
-// Copyright (C) 2005 by the deal.II authors
+// Copyright (C) 2014 by the deal.II authors
//
-// This file is subject to QPL and may not be distributed
+// 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.
//
-//---------------------------- iges_read.cc ---------------------------
-
+//-----------------------------------------------------------
-// Read goteborg.iges and dump its topological structure to the logfile.
+// Read goteborg.iges and dump its topological structure to the
+// logfile.
#include "../tests.h"
#include <fstream>
<< ne << " edges, and "
<< nv << " vertices." << std::endl;
- // try
- // {
- // TopoDS_Shape sh_fail = read_IGES("iges_files/nonexistent.iges");
- // }
- // catch(...)
- // {
- // deallog << "Missing file correctly threw an exception." << std::endl;
- // }
-
return 0;
}
-
-//---------------------------- iges_read.cc ---------------------------
-// $Id: testsuite.html 13373 2006-07-13 13:12:08Z kanschat $
-// Version: $Name$
+//-----------------------------------------------------------
//
-// Copyright (C) 2005 by the deal.II authors
+// Copyright (C) 2014 by the deal.II authors
//
-// This file is subject to QPL and may not be distributed
+// 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.
//
-//---------------------------- iges_read.cc ---------------------------
-
+//-----------------------------------------------------------
-// Read goteborg.iges and dump its topological structure to the logfile.
+// Read a file in iges format, and write it out again in the same
+// format.
#include "../tests.h"
#include <fstream>
int main ()
{
- // std::ofstream logfile("output");
- // deallog.attach(logfile);
- // deallog.depth_console(0);
-
TopoDS_Shape sh = read_IGES(SOURCE_DIR "/iges_files/wigley.iges");
unsigned int nf=0, ne=0, nv=0;
count_elements(sh, nf, ne, nv);
- // deallog << "Shape contains " << nf << " faces, "
- // << ne << " edges, and "
- // << nv << " vertices." << std::endl;
write_IGES(sh, "output");
-
- // try
- // {
- // TopoDS_Shape sh_fail = read_IGES("iges_files/nonexistent.iges");
- // }
- // catch(...)
- // {
- // deallog << "Missing file correctly threw an exception." << std::endl;
- // }
-
return 0;
}
-
-
-//---------------------------- occ_circle_projection.cc ---------------------------
-// $Id: testsuite.html 13373 2006-07-13 13:12:08Z kanschat $
-// Version: $Name$
+//-----------------------------------------------------------
//
-// Copyright (C) 2005 by the deal.II authors
+// Copyright (C) 2014 by the deal.II authors
//
-// This file is subject to QPL and may not be distributed
+// 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.
//
-//---------------------------- occ_circle_projection.cc ---------------------------
-
+//-----------------------------------------------------------
// Create a Triangulation, interpolate its boundary points to a close
// smooth BSpline, and use that as a Boundary Descriptor.
tria.refine_global(2);
- // You can open the generated file with paraview.
+ // You can open the generated file with gmsh.
GridOut gridout;
gridout.write_msh (tria, logfile);
- std::ofstream of("/tmp/grid.msh");
- gridout.write_msh(tria, of);
-
return 0;
}
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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 circle, a Triangulation, and try to project normally on
+// it.
+
+#include "../tests.h"
+
+#include <deal.II/grid/occ_boundary_lib.h>
+
+#include <fstream>
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/grid_out.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_boundary_lib.h>
+
+#include <gp_Pnt.hxx>
+#include <gp_Dir.hxx>
+#include <gp_Ax2.hxx>
+#include <GC_MakeCircle.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+#include <BRepPrimAPI_MakeSphere.hxx>
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+ std::ofstream logfile("output");
+
+ gp_Pnt center(.5,.5,.5);
+ Standard_Real radius(Point<3>().distance(Pnt(center)));
+
+
+ TopoDS_Face face = BRepPrimAPI_MakeSphere(center, radius);
+
+ // Create a boundary projector.
+ NormalProjectionBoundary<3,3> sphere(face);
+
+
+ // The unit cube.
+ Triangulation<3,3> tria;
+ GridGenerator::hyper_cube(tria);
+
+ // Set the exterior boundary
+ tria.set_manifold(1, sphere);
+
+ // This is here to ignore the
+ // points created in the interior
+ // of the face.
+ tria.begin()->set_all_manifold_ids(1);
+ tria.begin()->set_manifold_id(0);
+
+ tria.refine_global(2);
+
+
+ // You can open the generated file with gmsh
+ GridOut gridout;
+ gridout.write_msh (tria, logfile);
+
+ return 0;
+}
+
--- /dev/null
+$NOD
+125
+1 0 0 0
+2 1 0 0
+3 0 1 0
+4 1 1 0
+5 0 0 1
+6 1 0 1
+7 0 1 1
+8 1 1 1
+9 0.5 -0.112372 -0.112372
+10 -0.112372 0.5 -0.112372
+11 -0.112372 -0.112372 0.5
+12 1.11237 0.5 -0.112372
+13 1.11237 -0.112372 0.5
+14 0.5 1.11237 -0.112372
+15 -0.112372 1.11237 0.5
+16 1.11237 1.11237 0.5
+17 0.5 -0.112372 1.11237
+18 -0.112372 0.5 1.11237
+19 1.11237 0.5 1.11237
+20 0.5 1.11237 1.11237
+21 0.5 -0.366025 0.5
+22 0.5 0.5 -0.366025
+23 -0.366025 0.5 0.5
+24 1.36603 0.5 0.5
+25 0.5 1.36603 0.5
+26 0.5 0.5 1.36603
+27 0.5 0.5 0.5
+28 0.237676 -0.0836036 -0.0836036
+29 0.762324 -0.0836036 -0.0836036
+30 -0.0836036 0.237676 -0.0836036
+31 -0.0836036 0.762324 -0.0836036
+32 -0.0836036 -0.0836036 0.237676
+33 -0.0836036 -0.0836036 0.762324
+34 1.0836 0.237676 -0.0836036
+35 1.0836 0.762324 -0.0836036
+36 1.0836 -0.0836036 0.237676
+37 1.0836 -0.0836036 0.762324
+38 0.237676 1.0836 -0.0836036
+39 0.762324 1.0836 -0.0836036
+40 -0.0836036 1.0836 0.237676
+41 -0.0836036 1.0836 0.762324
+42 1.0836 1.0836 0.237676
+43 1.0836 1.0836 0.762324
+44 0.237676 -0.0836036 1.0836
+45 0.762324 -0.0836036 1.0836
+46 -0.0836036 0.237676 1.0836
+47 -0.0836036 0.762324 1.0836
+48 1.0836 0.237676 1.0836
+49 1.0836 0.762324 1.0836
+50 0.237676 1.0836 1.0836
+51 0.762324 1.0836 1.0836
+52 0.168586 -0.300103 0.5
+53 0.831414 -0.300103 0.5
+54 0.5 -0.300103 0.168586
+55 0.5 -0.300103 0.831414
+56 0.5 0.168586 -0.300103
+57 0.5 0.831414 -0.300103
+58 0.168586 0.5 -0.300103
+59 0.831414 0.5 -0.300103
+60 -0.300103 0.5 0.168586
+61 -0.300103 0.5 0.831414
+62 -0.300103 0.168586 0.5
+63 -0.300103 0.831414 0.5
+64 1.3001 0.5 0.168586
+65 1.3001 0.5 0.831414
+66 1.3001 0.168586 0.5
+67 1.3001 0.831414 0.5
+68 0.168586 1.3001 0.5
+69 0.831414 1.3001 0.5
+70 0.5 1.3001 0.168586
+71 0.5 1.3001 0.831414
+72 0.5 0.168586 1.3001
+73 0.5 0.831414 1.3001
+74 0.168586 0.5 1.3001
+75 0.831414 0.5 1.3001
+76 0.5 0.5 0.933013
+77 0.5 0.5 0.0669873
+78 0.933013 0.5 0.5
+79 0.0669873 0.5 0.5
+80 0.5 0.933013 0.5
+81 0.5 0.0669873 0.5
+82 0.183268 -0.241189 0.183268
+83 0.183268 -0.241189 0.816732
+84 0.816732 -0.241189 0.183268
+85 0.816732 -0.241189 0.816732
+86 0.183268 0.183268 -0.241189
+87 0.816732 0.183268 -0.241189
+88 0.183268 0.816732 -0.241189
+89 0.816732 0.816732 -0.241189
+90 -0.241189 0.183268 0.183268
+91 -0.241189 0.816732 0.183268
+92 -0.241189 0.183268 0.816732
+93 -0.241189 0.816732 0.816732
+94 1.24119 0.183268 0.183268
+95 1.24119 0.816732 0.183268
+96 1.24119 0.183268 0.816732
+97 1.24119 0.816732 0.816732
+98 0.183268 1.24119 0.183268
+99 0.183268 1.24119 0.816732
+100 0.816732 1.24119 0.183268
+101 0.816732 1.24119 0.816732
+102 0.183268 0.183268 1.24119
+103 0.816732 0.183268 1.24119
+104 0.183268 0.816732 1.24119
+105 0.816732 0.816732 1.24119
+106 0.880366 0.880366 0.5
+107 0.119634 0.880366 0.5
+108 0.880366 0.119634 0.5
+109 0.119634 0.119634 0.5
+110 0.880366 0.5 0.880366
+111 0.880366 0.5 0.119634
+112 0.119634 0.5 0.880366
+113 0.119634 0.5 0.119634
+114 0.5 0.880366 0.880366
+115 0.5 0.119634 0.880366
+116 0.5 0.880366 0.119634
+117 0.5 0.119634 0.119634
+118 0.151399 0.151399 0.151399
+119 0.848601 0.151399 0.151399
+120 0.151399 0.848601 0.151399
+121 0.848601 0.848601 0.151399
+122 0.151399 0.151399 0.848601
+123 0.848601 0.151399 0.848601
+124 0.151399 0.848601 0.848601
+125 0.848601 0.848601 0.848601
+$ENDNOD
+$ELM
+64
+1 5 0 0 8 1 28 82 32 30 86 118 90
+2 5 0 0 8 28 9 54 82 86 56 117 118
+3 5 0 0 8 30 86 118 90 10 58 113 60
+4 5 0 0 8 86 56 117 118 58 22 77 113
+5 5 0 0 8 32 82 52 11 90 118 109 62
+6 5 0 0 8 82 54 21 52 118 117 81 109
+7 5 0 0 8 90 118 109 62 60 113 79 23
+8 5 0 0 8 118 117 81 109 113 77 27 79
+9 5 0 0 8 9 29 84 54 56 87 119 117
+10 5 0 0 8 29 2 36 84 87 34 94 119
+11 5 0 0 8 56 87 119 117 22 59 111 77
+12 5 0 0 8 87 34 94 119 59 12 64 111
+13 5 0 0 8 54 84 53 21 117 119 108 81
+14 5 0 0 8 84 36 13 53 119 94 66 108
+15 5 0 0 8 117 119 108 81 77 111 78 27
+16 5 0 0 8 119 94 66 108 111 64 24 78
+17 5 0 0 8 10 58 113 60 31 88 120 91
+18 5 0 0 8 58 22 77 113 88 57 116 120
+19 5 0 0 8 31 88 120 91 3 38 98 40
+20 5 0 0 8 88 57 116 120 38 14 70 98
+21 5 0 0 8 60 113 79 23 91 120 107 63
+22 5 0 0 8 113 77 27 79 120 116 80 107
+23 5 0 0 8 91 120 107 63 40 98 68 15
+24 5 0 0 8 120 116 80 107 98 70 25 68
+25 5 0 0 8 22 59 111 77 57 89 121 116
+26 5 0 0 8 59 12 64 111 89 35 95 121
+27 5 0 0 8 57 89 121 116 14 39 100 70
+28 5 0 0 8 89 35 95 121 39 4 42 100
+29 5 0 0 8 77 111 78 27 116 121 106 80
+30 5 0 0 8 111 64 24 78 121 95 67 106
+31 5 0 0 8 116 121 106 80 70 100 69 25
+32 5 0 0 8 121 95 67 106 100 42 16 69
+33 5 0 0 8 11 52 83 33 62 109 122 92
+34 5 0 0 8 52 21 55 83 109 81 115 122
+35 5 0 0 8 62 109 122 92 23 79 112 61
+36 5 0 0 8 109 81 115 122 79 27 76 112
+37 5 0 0 8 33 83 44 5 92 122 102 46
+38 5 0 0 8 83 55 17 44 122 115 72 102
+39 5 0 0 8 92 122 102 46 61 112 74 18
+40 5 0 0 8 122 115 72 102 112 76 26 74
+41 5 0 0 8 21 53 85 55 81 108 123 115
+42 5 0 0 8 53 13 37 85 108 66 96 123
+43 5 0 0 8 81 108 123 115 27 78 110 76
+44 5 0 0 8 108 66 96 123 78 24 65 110
+45 5 0 0 8 55 85 45 17 115 123 103 72
+46 5 0 0 8 85 37 6 45 123 96 48 103
+47 5 0 0 8 115 123 103 72 76 110 75 26
+48 5 0 0 8 123 96 48 103 110 65 19 75
+49 5 0 0 8 23 79 112 61 63 107 124 93
+50 5 0 0 8 79 27 76 112 107 80 114 124
+51 5 0 0 8 63 107 124 93 15 68 99 41
+52 5 0 0 8 107 80 114 124 68 25 71 99
+53 5 0 0 8 61 112 74 18 93 124 104 47
+54 5 0 0 8 112 76 26 74 124 114 73 104
+55 5 0 0 8 93 124 104 47 41 99 50 7
+56 5 0 0 8 124 114 73 104 99 71 20 50
+57 5 0 0 8 27 78 110 76 80 106 125 114
+58 5 0 0 8 78 24 65 110 106 67 97 125
+59 5 0 0 8 80 106 125 114 25 69 101 71
+60 5 0 0 8 106 67 97 125 69 16 43 101
+61 5 0 0 8 76 110 75 26 114 125 105 73
+62 5 0 0 8 110 65 19 75 125 97 49 105
+63 5 0 0 8 114 125 105 73 71 101 51 20
+64 5 0 0 8 125 97 49 105 101 43 8 51
+$ENDELM
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2014 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 sphere, intersect it with a plane, and use the resulting
+// object as a boundary.
+
+#include "../tests.h"
+
+#include <deal.II/grid/occ_utilities.h>
+#include <deal.II/grid/occ_boundary_lib.h>
+
+#include <fstream>
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/grid_out.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_boundary_lib.h>
+
+#include <gp_Pnt.hxx>
+#include <gp_Dir.hxx>
+#include <gp_Ax2.hxx>
+#include <GC_MakeCircle.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+#include <BRepPrimAPI_MakeSphere.hxx>
+#include <TopExp_Explorer.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
+#include <TopoDS.hxx>
+using namespace OpenCASCADE;
+
+int main()
+{
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ // Create a sphere
+ gp_Pnt center(.5,.5,.5);
+ Standard_Real radius(Point<3>().distance(Pnt(center)));
+
+ TopoDS_Face face = BRepPrimAPI_MakeSphere(center, radius);
+
+ // Intersect the sphere with the plane z = 0
+ TopoDS_Shape intersect = intersect_plane(face, 0, 0, 1, 0);
+ // This shape should be the standard circle. Only one edge should be
+ // in intersect: our edge.
+ TopExp_Explorer exp(intersect, TopAbs_EDGE);
+ TopoDS_Edge edge = TopoDS::Edge(exp.Current());
+
+ ArclengthProjectionLineManifold<2,3> boundary_line(edge);
+
+ // The unit square.
+ Triangulation<2,3> tria;
+ GridGenerator::hyper_cube(tria);
+
+ // Set the exterior boundary
+ tria.set_manifold(1, boundary_line);
+
+ // This is here to ignore the points created in the interior of the
+ // face.
+ tria.begin()->set_all_manifold_ids(1);
+ tria.begin()->set_manifold_id(0);
+
+ // We refine twice, and expect the outer points to end up on a
+ // smooth curve interpolating the square vertices.
+ tria.refine_global(2);
+
+
+ // You can open the generated file with gmsh.
+ GridOut gridout;
+ gridout.write_msh (tria, logfile);
+ return 0;
+}
--- /dev/null
+
+$NOD
+25
+1 0.00000 0.00000 0.00000
+2 1.00000 0.00000 0.00000
+3 0.00000 1.00000 0.00000
+4 1.00000 1.00000 0.00000
+5 0.500000 -0.207107 0.00000
+6 -0.207107 0.500000 0.00000
+7 1.20711 0.500000 0.00000
+8 0.500000 1.20711 0.00000
+9 0.500000 0.500000 0.00000
+10 0.229402 -0.153281 0.00000
+11 0.770598 -0.153281 0.00000
+12 -0.153281 0.229402 0.00000
+13 -0.153281 0.770598 0.00000
+14 1.15328 0.229402 0.00000
+15 1.15328 0.770598 0.00000
+16 0.229402 1.15328 0.00000
+17 0.770598 1.15328 0.00000
+18 0.500000 0.146447 0.00000
+19 0.500000 0.853553 0.00000
+20 0.146447 0.500000 0.00000
+21 0.853553 0.500000 0.00000
+22 0.189433 0.189433 0.00000
+23 0.810567 0.189433 0.00000
+24 0.189433 0.810567 0.00000
+25 0.810567 0.810567 0.00000
+$ENDNOD
+$ELM
+16
+1 3 0 0 4 1 10 22 12
+2 3 0 0 4 10 5 18 22
+3 3 0 0 4 12 22 20 6
+4 3 0 0 4 22 18 9 20
+5 3 0 0 4 5 11 23 18
+6 3 0 0 4 11 2 14 23
+7 3 0 0 4 18 23 21 9
+8 3 0 0 4 23 14 7 21
+9 3 0 0 4 6 20 24 13
+10 3 0 0 4 20 9 19 24
+11 3 0 0 4 13 24 16 3
+12 3 0 0 4 24 19 8 16
+13 3 0 0 4 9 21 25 19
+14 3 0 0 4 21 7 15 25
+15 3 0 0 4 19 25 17 8
+16 3 0 0 4 25 15 4 17
+$ENDELM