18#ifdef DEAL_II_WITH_OPENCASCADE
24# include <IGESControl_Controller.hxx>
25# include <IGESControl_Reader.hxx>
26# include <IGESControl_Writer.hxx>
27# include <STEPControl_Controller.hxx>
28# include <STEPControl_Reader.hxx>
29# include <STEPControl_Writer.hxx>
30# include <Standard_Version.hxx>
31# include <TopExp_Explorer.hxx>
33# include <TopoDS_Edge.hxx>
34# include <TopoDS_Face.hxx>
35# include <TopoDS_Shape.hxx>
40# if (OCC_VERSION_MAJOR < 7)
41# include <Handle_Standard_Transient.hxx>
43# include <Standard_Transient.hxx>
46# include <BRepAdaptor_Curve.hxx>
47# include <BRepAdaptor_HCompCurve.hxx>
48# include <BRepAdaptor_HCurve.hxx>
49# include <BRepAdaptor_Surface.hxx>
50# include <BRepAlgo_Section.hxx>
51# include <BRepBndLib.hxx>
52# include <BRepBuilderAPI_MakeEdge.hxx>
53# include <BRepBuilderAPI_Sewing.hxx>
54# include <BRepBuilderAPI_Transform.hxx>
55# include <BRepMesh_IncrementalMesh.hxx>
56# include <BRepTools.hxx>
57# include <BRep_Builder.hxx>
58# include <GCPnts_AbscissaPoint.hxx>
59# include <GeomAPI_Interpolate.hxx>
60# include <GeomAPI_ProjectPointOnCurve.hxx>
61# include <GeomAPI_ProjectPointOnSurf.hxx>
62# include <GeomConvert_CompCurveToBSplineCurve.hxx>
63# include <GeomLProp_SLProps.hxx>
64# include <Geom_BoundedCurve.hxx>
65# include <Geom_Plane.hxx>
66# include <IntCurvesFace_ShapeIntersector.hxx>
67# include <Poly_Triangulation.hxx>
68# include <ShapeAnalysis_Surface.hxx>
69# include <StlAPI_Reader.hxx>
70# include <StlAPI_Writer.hxx>
71# include <TColStd_HSequenceOfTransient.hxx>
72# include <TColStd_SequenceOfTransient.hxx>
73# include <TColgp_HArray1OfPnt.hxx>
74# include <TopLoc_Location.hxx>
87 std::tuple<unsigned int, unsigned int, unsigned int>
91 unsigned int n_faces = 0, n_edges = 0, n_vertices = 0;
92 for (
exp.Init(shape, TopAbs_FACE);
exp.More();
exp.Next(), ++n_faces)
95 for (
exp.Init(shape, TopAbs_EDGE);
exp.More();
exp.Next(), ++n_edges)
98 for (
exp.Init(shape, TopAbs_VERTEX);
exp.More();
exp.Next(), ++n_vertices)
101 return std::tuple<unsigned int, unsigned int, unsigned int>(n_faces,
108 std::vector<TopoDS_Face> & faces,
109 std::vector<TopoDS_Edge> & edges,
110 std::vector<TopoDS_Vertex> &
vertices)
117 for (
exp.Init(shape, TopAbs_FACE);
exp.More();
exp.Next())
119 faces.push_back(TopoDS::Face(
exp.Current()));
121 for (
exp.Init(shape, TopAbs_EDGE);
exp.More();
exp.Next())
123 edges.push_back(TopoDS::Edge(
exp.Current()));
125 for (
exp.Init(shape, TopAbs_VERTEX);
exp.More();
exp.Next())
134 std::vector<TopoDS_Compound> & compounds,
135 std::vector<TopoDS_CompSolid> &compsolids,
136 std::vector<TopoDS_Solid> & solids,
137 std::vector<TopoDS_Shell> & shells,
138 std::vector<TopoDS_Wire> & wires)
141 compsolids.resize(0);
147 for (
exp.Init(shape, TopAbs_COMPOUND);
exp.More();
exp.Next())
149 compounds.push_back(TopoDS::Compound(
exp.Current()));
151 for (
exp.Init(shape, TopAbs_COMPSOLID);
exp.More();
exp.Next())
153 compsolids.push_back(TopoDS::CompSolid(
exp.Current()));
155 for (
exp.Init(shape, TopAbs_SOLID);
exp.More();
exp.Next())
157 solids.push_back(TopoDS::Solid(
exp.Current()));
159 for (
exp.Init(shape, TopAbs_SHELL);
exp.More();
exp.Next())
161 shells.push_back(TopoDS::Shell(
exp.Current()));
163 for (
exp.Init(shape, TopAbs_WIRE);
exp.More();
exp.Next())
165 wires.push_back(TopoDS::Wire(
exp.Current()));
169 template <
int spacedim>
176 return gp_Pnt(p[0], 0, 0);
178 return gp_Pnt(p[0], p[1], 0);
180 return gp_Pnt(p[0], p[1], p[2]);
186 template <
int spacedim>
188 point(
const gp_Pnt &p,
const double tolerance)
196 "Cannot convert OpenCASCADE point to 1d if p.Y() != 0."));
199 "Cannot convert OpenCASCADE point to 1d if p.Z() != 0."));
204 "Cannot convert OpenCASCADE point to 2d if p.Z() != 0."));
218 const double tolerance)
220 const double rel_tol =
222 if (direction.
norm() > 0.0)
223 return (p1 * direction < p2 * direction - rel_tol);
225 for (
int d = dim;
d >= 0; --
d)
226 if (p1[
d] < p2[
d] - rel_tol)
228 else if (p2[
d] < p1[
d] - rel_tol)
238 read_IGES(
const std::string &filename,
const double scale_factor)
240 IGESControl_Reader reader;
241 IFSelect_ReturnStatus stat;
242 stat = reader.ReadFile(filename.c_str());
245 Standard_Boolean failsonly = Standard_False;
246 IFSelect_PrintCount mode = IFSelect_ItemsByEntity;
247 reader.PrintCheckLoad(failsonly, mode);
249 Standard_Integer nRoots = reader.TransferRoots();
258 scale.SetScale(Origin, scale_factor);
260 TopoDS_Shape sh = reader.OneShape();
261 BRepBuilderAPI_Transform trans(sh,
scale);
263 return trans.Shape();
267 write_IGES(
const TopoDS_Shape &shape,
const std::string &filename)
269 IGESControl_Controller::Init();
270 IGESControl_Writer ICW(
"MM", 0);
271 Standard_Boolean ok = ICW.AddShape(shape);
274 Standard_Boolean OK = ICW.Write(filename.c_str());
282 StlAPI_Reader reader;
284 reader.Read(shape, filename.c_str());
291 const std::string & filename,
292 const double deflection,
293 const bool sew_different_faces,
294 const double sewer_tolerance,
295 const bool is_relative,
296 const double angular_deflection,
297 const bool in_parallel)
300 std::vector<TopoDS_Vertex>
vertices;
301 std::vector<TopoDS_Edge> edges;
302 std::vector<TopoDS_Face> faces;
304 const bool mesh_is_present =
305 std::none_of(faces.begin(), faces.end(), [&Loc](
const TopoDS_Face &face) {
306 Handle(Poly_Triangulation) theTriangulation =
307 BRep_Tool::Triangulation(face, Loc);
308 return theTriangulation.IsNull();
310 TopoDS_Shape shape_to_be_written = shape;
311 if (!mesh_is_present)
313 if (sew_different_faces)
315 BRepBuilderAPI_Sewing sewer(sewer_tolerance);
316 sewer.Add(shape_to_be_written);
318 shape_to_be_written = sewer.SewedShape();
321 shape_to_be_written = shape;
325 BRepMesh_IncrementalMesh mesh_im(shape_to_be_written,
332 StlAPI_Writer writer;
334# if ((OCC_VERSION_MAJOR * 100 + OCC_VERSION_MINOR * 10) >= 690)
336 const auto error = writer.Write(shape_to_be_written, filename.c_str());
339# if ((OCC_VERSION_MAJOR * 100 + OCC_VERSION_MINOR * 10) < 720)
350 writer.Write(shape_to_be_written, filename.c_str());
355 read_STEP(
const std::string &filename,
const double scale_factor)
357 STEPControl_Reader reader;
358 IFSelect_ReturnStatus stat;
359 stat = reader.ReadFile(filename.c_str());
362 Standard_Boolean failsonly = Standard_False;
363 IFSelect_PrintCount mode = IFSelect_ItemsByEntity;
364 reader.PrintCheckLoad(failsonly, mode);
366 Standard_Integer nRoots = reader.TransferRoots();
375 scale.SetScale(Origin, scale_factor);
377 TopoDS_Shape sh = reader.OneShape();
378 BRepBuilderAPI_Transform trans(sh,
scale);
380 return trans.Shape();
384 write_STEP(
const TopoDS_Shape &shape,
const std::string &filename)
386 STEPControl_Controller::Init();
387 STEPControl_Writer SCW;
388 IFSelect_ReturnStatus status;
389 status = SCW.Transfer(shape, STEPControl_AsIs);
391 ExcMessage(
"Failed to add shape to STEP controller."));
393 status = SCW.Write(filename.c_str());
396 ExcMessage(
"Failed to write translated shape to STEP file."));
402 double tolerance = 0.0;
404 std::vector<TopoDS_Face> faces;
405 std::vector<TopoDS_Edge> edges;
406 std::vector<TopoDS_Vertex>
vertices;
411 tolerance = std::fmax(tolerance, BRep_Tool::Tolerance(vertex));
413 for (
const auto &edge : edges)
414 tolerance = std::fmax(tolerance, BRep_Tool::Tolerance(edge));
416 for (
const auto &face : faces)
417 tolerance = std::fmax(tolerance, BRep_Tool::Tolerance(face));
431 Handle(Geom_Plane) plane =
new Geom_Plane(c_x, c_y, c_z, c);
432 BRepAlgo_Section section(in_shape, plane);
433 TopoDS_Shape edges = section.Shape();
438 join_edges(
const TopoDS_Shape &in_shape,
const double tolerance)
440 TopoDS_Edge out_shape;
441 const TopoDS_Shape & edges = in_shape;
442 std::vector<Handle_Geom_BoundedCurve> intersections;
446 gp_Pnt PIn(0.0, 0.0, 0.0);
447 gp_Pnt PFin(0.0, 0.0, 0.0);
448 gp_Pnt PMid(0.0, 0.0, 0.0);
449 TopExp_Explorer edgeExplorer(edges, TopAbs_EDGE);
451 while (edgeExplorer.More())
453 edge = TopoDS::Edge(edgeExplorer.Current());
454 Handle(Geom_Curve) curve = BRep_Tool::Curve(edge,
L, First, Last);
455 intersections.push_back(Handle(Geom_BoundedCurve)::DownCast(curve));
461 unsigned int numIntersEdges = intersections.size();
464 GeomConvert_CompCurveToBSplineCurve convert_bspline(intersections[0]);
466 bool check =
false, one_added =
true, one_failed =
true;
467 std::vector<bool> added(numIntersEdges,
false);
469 while (one_added ==
true)
473 for (
unsigned int i = 1; i < numIntersEdges; ++i)
474 if (added[i] ==
false)
476 Handle(Geom_Curve) curve = intersections[i];
477 Handle(Geom_BoundedCurve) bcurve =
478 Handle(Geom_BoundedCurve)::DownCast(curve);
479 check = convert_bspline.Add(bcurve, tolerance,
false,
true, 0);
484 Handle(Geom_BoundedCurve) bcurve =
485 Handle(Geom_BoundedCurve)::DownCast(curve);
487 convert_bspline.Add(bcurve, tolerance,
false,
true, 0);
489 one_failed = one_failed || (check ==
false);
490 one_added = one_added || (check ==
true);
495 Assert(one_failed ==
false,
496 ExcMessage(
"Joining some of the Edges failed."));
498 Handle(Geom_Curve) bspline = convert_bspline.BSplineCurve();
500 out_shape = BRepBuilderAPI_MakeEdge(bspline);
509 const double tolerance)
513 gp_Pnt P0 =
point(origin);
516 dim > 1 ? direction[1] : 0,
517 dim > 2 ? direction[2] : 0));
521 gp_Pnt Pproj(0.0, 0.0, 0.0);
525 IntCurvesFace_ShapeIntersector Inters;
526 Inters.Load(in_shape, tolerance);
530 Inters.Perform(line, -RealLast(), +RealLast());
533 double minDistance = 1e7;
535 for (
int i = 0; i < Inters.NbPnt(); ++i)
537 const double distance =
point(origin).Distance(Inters.Pnt(i + 1));
540 if (distance < minDistance)
542 minDistance = distance;
543 result = point<dim>(Inters.Pnt(i + 1));
555 const double tolerance)
557 unsigned int n_vertices = curve_points.size();
559 if (direction * direction > 0)
561 std::sort(curve_points.begin(),
564 return OpenCASCADE::point_compare(p1,
572 Handle(TColgp_HArray1OfPnt)
vertices =
573 new TColgp_HArray1OfPnt(1, n_vertices);
574 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
580 GeomAPI_Interpolate bspline_generator(
vertices, closed, tolerance);
581 bspline_generator.Perform();
582 Assert((bspline_generator.IsDone()),
583 ExcMessage(
"Interpolated bspline generation failed"));
585 Handle(Geom_BSplineCurve) bspline = bspline_generator.Curve();
586 TopoDS_Edge out_shape = BRepBuilderAPI_MakeEdge(bspline);
587 out_shape.Closed(closed);
593 template <
int spacedim>
594 std::vector<TopoDS_Edge>
602 std::map<unsigned int, std::pair<unsigned int, unsigned int>> vert_to_faces;
603 std::map<unsigned int, std::pair<unsigned int, unsigned int>> face_to_verts;
604 std::map<unsigned int, bool> visited_faces;
605 std::map<unsigned int, Point<spacedim>> vert_to_point;
607 unsigned int face_index;
609 for (
const auto &cell :
triangulation.active_cell_iterators())
611 if (cell->face(f)->at_boundary())
614 face_index = cell->face(f)->index();
615 const unsigned int v0 = cell->face(f)->vertex_index(0);
616 const unsigned int v1 = cell->face(f)->vertex_index(1);
617 face_to_verts[face_index].first =
v0;
618 face_to_verts[face_index].second =
v1;
619 visited_faces[face_index] =
false;
624 f, 0,
true,
false,
false)];
626 f, 1,
true,
false,
false)];
629 if (vert_to_faces.find(
v0) == vert_to_faces.end())
631 vert_to_faces[
v0].first = face_index;
635 vert_to_faces[
v0].second = face_index;
637 if (vert_to_faces.find(
v1) == vert_to_faces.end())
639 vert_to_faces[
v1].first = face_index;
643 vert_to_faces[
v1].second = face_index;
649 std::vector<TopoDS_Edge> interpolation_curves;
650 bool finished = (face_to_verts.size() == 0);
651 face_index = finished ? 0 : face_to_verts.begin()->first;
653 while (finished ==
false)
655 const unsigned int start_point_index = face_to_verts[face_index].first;
656 unsigned int point_index = start_point_index;
659 std::vector<Point<spacedim>> pointlist;
662 visited_faces[face_index] =
true;
663 auto current_point = vert_to_point[point_index];
664 pointlist.push_back(current_point);
667 if (face_to_verts[face_index].
first != point_index)
668 point_index = face_to_verts[face_index].first;
670 point_index = face_to_verts[face_index].second;
673 if (vert_to_faces[point_index].
first != face_index)
674 face_index = vert_to_faces[point_index].first;
676 face_index = vert_to_faces[point_index].second;
678 while (point_index != start_point_index);
680 interpolation_curves.push_back(
684 for (
const auto &f : visited_faces)
685 if (f.second ==
false)
687 face_index = f.first;
692 return interpolation_curves;
697 std::tuple<Point<dim>, TopoDS_Shape, double,
double>
700 const double tolerance)
703 gp_Pnt Pproj =
point(origin);
705 double minDistance = 1e7;
706 gp_Pnt tmp_proj(0.0, 0.0, 0.0);
708 unsigned int counter = 0;
709 unsigned int face_counter = 0;
711 TopoDS_Shape out_shape;
715 for (
exp.Init(in_shape, TopAbs_FACE);
exp.More();
exp.Next())
717 TopoDS_Face face = TopoDS::Face(
exp.Current());
721 Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
723 ShapeAnalysis_Surface projector(SurfToProj);
724 gp_Pnt2d proj_params = projector.ValueOfUV(
point(origin), tolerance);
726 SurfToProj->D0(proj_params.X(), proj_params.Y(), tmp_proj);
728 double distance = point<dim>(tmp_proj).distance(origin);
729 if (distance < minDistance)
731 minDistance = distance;
745 if (face_counter == 0)
746 for (
exp.Init(in_shape, TopAbs_EDGE);
exp.More();
exp.Next())
748 TopoDS_Edge edge = TopoDS::Edge(
exp.Current());
749 if (!BRep_Tool::Degenerated(edge))
757 Handle(Geom_Curve) CurveToProj =
758 BRep_Tool::Curve(edge,
L, First, Last);
760 GeomAPI_ProjectPointOnCurve Proj(
point(origin), CurveToProj);
761 unsigned int num_proj_points = Proj.NbPoints();
762 if ((num_proj_points > 0) && (Proj.LowerDistance() < minDistance))
764 minDistance = Proj.LowerDistance();
765 Pproj = Proj.NearestPoint();
767 u = Proj.LowerDistanceParameter();
774 return std::tuple<Point<dim>, TopoDS_Shape, double,
double>(
775 point<dim>(Pproj), out_shape, u, v);
783 const double tolerance)
785 std::tuple<Point<dim>, TopoDS_Shape, double,
double> ref =
787 return std::get<0>(ref);
793 const double tolerance)
796 std::tuple<Point<3>, TopoDS_Shape, double,
double> shape_and_params =
799 TopoDS_Shape &out_shape = std::get<1>(shape_and_params);
800 double & u = std::get<2>(shape_and_params);
801 double & v = std::get<3>(shape_and_params);
805 std::tuple<unsigned int, unsigned int, unsigned int>
numbers =
812 "Could not find normal: the shape containing the closest point has 0 faces."));
816 "Could not find normal: the shape containing the closest point has more than 1 face."));
820 exp.Init(out_shape, TopAbs_FACE);
821 TopoDS_Face face = TopoDS::Face(
exp.Current());
827 push_forward(
const TopoDS_Shape &in_shape,
const double u,
const double v)
829 switch (in_shape.ShapeType())
833 BRepAdaptor_Surface surf(TopoDS::Face(in_shape));
834 return point<dim>(surf.Value(u, v));
838 BRepAdaptor_Curve curve(TopoDS::Edge(in_shape));
839 return point<dim>(curve.Value(u));
853 Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
854 GeomLProp_SLProps props(SurfToProj, u, v, 1, 1
e-7);
855 gp_Pnt Value = props.Value();
857 gp_Dir Normal = props.Normal();
858 Assert(props.IsCurvatureDefined(),
859 ExcMessage(
"Curvature is not well defined!"));
860 Standard_Real Min_Curvature = props.MinCurvature();
861 Standard_Real Max_Curvature = props.MaxCurvature();
868 if (face.Orientation() == TopAbs_REVERSED)
875 return std::tuple<Point<3>,
Tensor<1, 3>, double,
double>(point<3>(Value),
883 template <
int spacedim>
888 BRepAdaptor_Surface surf(face);
889 const double u0 = surf.FirstUParameter();
890 const double u1 = surf.LastUParameter();
891 const double v0 = surf.FirstVParameter();
892 const double v1 = surf.LastVParameter();
894 std::vector<CellData<2>> cells;
895 std::vector<Point<spacedim>>
vertices;
898 vertices.push_back(point<spacedim>(surf.Value(u0,
v0)));
899 vertices.push_back(point<spacedim>(surf.Value(u1,
v0)));
900 vertices.push_back(point<spacedim>(surf.Value(u0,
v1)));
901 vertices.push_back(point<spacedim>(surf.Value(u1,
v1)));
904 for (
unsigned int i = 0; i < 4; ++i)
907 cells.push_back(cell);
911# include "utilities.inst"
Abstract base class for mapping classes.
virtual boost::container::small_vector< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
numbers::NumberTraits< Number >::real_type norm() const
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcUnsupportedShape()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
graph_traits< Graph >::vertex_descriptor Vertex
TopoDS_Edge interpolation_curve(std::vector< Point< dim > > &curve_points, const Tensor< 1, dim > &direction=Tensor< 1, dim >(), const bool closed=false, const double tolerance=1e-7)
std::tuple< Point< 3 >, Tensor< 1, 3 >, double, double > push_forward_and_differential_forms(const TopoDS_Face &face, const double u, const double v, const double tolerance=1e-7)
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)
Point< dim > push_forward(const TopoDS_Shape &in_shape, const double u, const double v)
std::tuple< unsigned int, unsigned int, unsigned int > count_elements(const TopoDS_Shape &shape)
void write_STEP(const TopoDS_Shape &shape, const std::string &filename)
std::tuple< Point< 3 >, Tensor< 1, 3 >, double, double > closest_point_and_differential_forms(const TopoDS_Shape &in_shape, const Point< 3 > &origin, const double tolerance=1e-7)
void extract_geometrical_shapes(const TopoDS_Shape &shape, std::vector< TopoDS_Face > &faces, std::vector< TopoDS_Edge > &edges, std::vector< TopoDS_Vertex > &vertices)
TopoDS_Edge join_edges(const TopoDS_Shape &in_shape, const double tolerance=1e-7)
TopoDS_Shape read_STEP(const std::string &filename, const double scale_factor=1e-3)
bool point_compare(const Point< dim > &p1, const Point< dim > &p2, const Tensor< 1, dim > &direction=Tensor< 1, dim >(), const double tolerance=1e-10)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Point< dim > closest_point(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
std::tuple< Point< dim >, TopoDS_Shape, double, double > project_point_and_pull_back(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
Point< dim > line_intersection(const TopoDS_Shape &in_shape, const Point< dim > &origin, const Tensor< 1, dim > &direction, const double tolerance=1e-7)
void write_STL(const TopoDS_Shape &shape, const std::string &filename, const double deflection, const bool sew_different_faces=false, const double sewer_tolerance=1e-6, const bool is_relative=false, const double angular_deflection=0.5, const bool in_parallel=false)
void create_triangulation(const TopoDS_Face &face, Triangulation< 2, spacedim > &tria)
TopoDS_Shape read_STL(const std::string &filename)
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)
std::vector< TopoDS_Edge > create_curves_from_triangulation_boundary(const Triangulation< 2, spacedim > &triangulation, const Mapping< 2, spacedim > &mapping=StaticMappingQ1< 2, spacedim >::mapping)
double get_shape_tolerance(const TopoDS_Shape &shape)
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)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
::VectorizedArray< Number, width > exp(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
std::vector< unsigned int > vertices
static unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)