Refactoring

src/ifcconvert/IfcConvert.cpp

@@ -38,6 +38,7 @@
 #include "../ifcgeom_schema_agnostic/IfcGeomFilter.h"
 #include "../ifcgeom_schema_agnostic/IfcGeomIterator.h"
 #include "../ifcgeom_schema_agnostic/IfcGeomRenderStyles.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #include "../ifcparse/utils.h"
 
@@ -1632,7 +1633,7 @@ void fix_quantities(IfcParse::IfcFile& f, bool no_progress, bool quiet, bool std
 				auto quantity_count = latebound_access::create(f, "IfcQuantityCount");
 				latebound_access::set(quantity_count, "Name", std::string("Surface Genus"));
 				latebound_access::set(quantity_count, "Description", '#' + boost::lexical_cast<std::string>(part.ItemId()));
-				latebound_access::set(quantity_count, "CountValue", IfcGeom::Kernel::surface_genus(part.Shape()));
+				latebound_access::set(quantity_count, "CountValue", IfcGeom::util::surface_genus(part.Shape()));
 
 				quantities_2->push(quantity_count);				
 			}

src/ifcgeom/IfcAxis2Placement2D.cpp

@@ -24,6 +24,7 @@
 #include <gp_Trsf2d.hxx>
 #include <gp_Ax2d.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -44,7 +45,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcAxis2Placement2D* l, gp_Trsf2d
 
 	gp_Ax2d axis(gp_Pnt2d(P.X(),P.Y()), gp_Dir2d(V.X(),V.Y()));
 	
-	if (!axis_equal(axis, gp_Ax2d(), getValue(GV_PRECISION))) {
+	if (!util::axis_equal(axis, gp_Ax2d(), getValue(GV_PRECISION))) {
 		trsf.SetTransformation(axis, gp_Ax2d());
 	}
 

src/ifcgeom/IfcAxis2Placement3D.cpp

@@ -23,6 +23,7 @@
 #include <gp_Trsf.hxx>
 #include <gp_Ax3.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -66,7 +67,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcAxis2Placement3D* l, gp_Trsf&
 
 	gp_Ax3 ax3(o, axis, refDirection);
 
-	if (!axis_equal(ax3, (gp_Ax3) gp::XOY(), getValue(GV_PRECISION))) {
+	if (!util::axis_equal(ax3, (gp_Ax3) gp::XOY(), getValue(GV_PRECISION))) {
 		trsf.SetTransformation(ax3, gp::XOY());
 	}
 	

src/ifcgeom/IfcBooleanResult.cpp

@@ -20,6 +20,8 @@
 #include <TopoDS_Wire.hxx>
 #include <Standard_Version.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
+#include "../ifcgeom_schema_agnostic/boolean_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -179,6 +181,11 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcBooleanResult* l, TopoDS_Shape
 	
 	bool valid_result;
 
+	util::boolean_settings bst;
+	bst.attempt_2d = getValue(GV_BOOLEAN_ATTEMPT_2D) > 0.;
+	bst.debug = getValue(GV_DEBUG_BOOLEAN) > 0.;
+	bst.precision = getValue(GV_PRECISION);
+
 	if (s1.ShapeType() == TopAbs_COMPOUND && TopoDS_Iterator(s1).More() && util::is_nested_compound_of_solid(s1)) {
 		TopoDS_Compound C;
 		BRep_Builder B;
@@ -187,7 +194,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcBooleanResult* l, TopoDS_Shape
 		valid_result = true;
 		for (; it.More(); it.Next()) {
 			TopoDS_Shape part;
-			if (boolean_operation(it.Value(), second_operand_shapes, occ_op, part)) {
+			if (util::boolean_operation(bst, it.Value(), second_operand_shapes, occ_op, part)) {
 				B.Add(C, part);
 			} else {
 				valid_result = false;
@@ -195,7 +202,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcBooleanResult* l, TopoDS_Shape
 		}
 		shape = C;
 	} else {
-		valid_result = boolean_operation(s1, second_operand_shapes, occ_op, shape);
+		valid_result = util::boolean_operation(bst, s1, second_operand_shapes, occ_op, shape);
 	}
 
 #endif

src/ifcgeom/IfcCartesianTransformationOperator2D.cpp

@@ -24,6 +24,7 @@
 #include <gp_Trsf2d.hxx>
 #include <gp_Ax2d.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -60,7 +61,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcCartesianTransformationOperato
 		trsf.SetScaleFactor(*l->Scale());
 	}
 
-	if (is_identity(trsf, getValue(GV_PRECISION))) {
+	if (util::is_identity(trsf, getValue(GV_PRECISION))) {
 		trsf = gp_Trsf2d();
 	}
 

src/ifcgeom/IfcCartesianTransformationOperator2DnonUniform.cpp

@@ -25,6 +25,7 @@
 #include <gp_Trsf2d.hxx>
 #include <gp_Ax2d.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -61,7 +62,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcCartesianTransformationOperato
 	gtrsf.SetValue(2,2,scale2);
 	gtrsf.Multiply(trsf);
 
-	if (is_identity(gtrsf, getValue(GV_PRECISION))) {
+	if (util::is_identity(gtrsf, getValue(GV_PRECISION))) {
 		gtrsf = gp_GTrsf2d();
 	}
 

src/ifcgeom/IfcCartesianTransformationOperator3D.cpp

@@ -22,6 +22,7 @@
 #include <gp_Trsf.hxx>
 #include <gp_Ax3.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -38,7 +39,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcCartesianTransformationOperato
 	gp_Ax3 ax3 (origin,axis3,axis1);
 	if ( axis2.Dot(ax3.YDirection()) < 0 ) ax3.YReverse();
 	
-	if (!axis_equal(ax3, (gp_Ax3) gp::XOY(), getValue(GV_PRECISION))) {
+	if (!util::axis_equal(ax3, (gp_Ax3) gp::XOY(), getValue(GV_PRECISION))) {
 		trsf.SetTransformation(ax3);
 		trsf.Invert();
 	}

src/ifcgeom/IfcCartesianTransformationOperator3DnonUniform.cpp

@@ -23,6 +23,7 @@
 #include <gp_Trsf.hxx>
 #include <gp_Ax3.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -50,7 +51,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcCartesianTransformationOperato
 	gtrsf.SetValue(3,3,scale3);
 	gtrsf.PreMultiply(trsf);
 
-	if (is_identity(gtrsf, getValue(GV_PRECISION))) {
+	if (util::is_identity(gtrsf, getValue(GV_PRECISION))) {
 		gtrsf = gp_GTrsf();
 	}
 

src/ifcgeom/IfcFace.cpp

@@ -125,14 +125,14 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcFace* l, TopoDS_Shape& result)
 				}
 			}
 
-			for (auto& wire : wires) {
+			for (auto& w : wires) {
 				if (!same_sense) {
-					wire.Reverse();
+					w.Reverse();
 				}
 
-				wire_senses.Bind(wire.Oriented(TopAbs_FORWARD), same_sense ? TopAbs_FORWARD : TopAbs_REVERSED);
+				wire_senses.Bind(w.Oriented(TopAbs_FORWARD), same_sense ? TopAbs_FORWARD : TopAbs_REVERSED);
 
-				fd.wires().emplace_back(TopoDS::Wire(wire));
+				fd.wires().emplace_back(TopoDS::Wire(w));
 			}
 		}
 	}

src/ifcgeom/IfcGeom.h

@@ -167,8 +167,6 @@ private:
 	double boolean_debug_setting;
 	double boolean_attempt_2d;
 
-	size_t operation_counter_ = 0;
-
 	// For stopping PlacementRelTo recursion in convert(const IfcSchema::IfcObjectPlacement* l, gp_Trsf& trsf)
 	const IfcParse::declaration* placement_rel_to_type_;
 	const IfcUtil::IfcBaseEntity* placement_rel_to_instance_;
@@ -273,7 +271,6 @@ public:
 	bool convert_curve(const IfcUtil::IfcBaseInterface* L, Handle(Geom_Curve)& result);
 	bool convert_face(const IfcUtil::IfcBaseInterface* L, TopoDS_Shape& result);
 	bool convert_openings(const IfcSchema::IfcProduct* entity, const IfcSchema::IfcRelVoidsElement::list::ptr& openings, const IfcRepresentationShapeItems& entity_shapes, const gp_Trsf& entity_trsf, IfcRepresentationShapeItems& cut_shapes);
-	bool convert_openings_fast(const IfcSchema::IfcProduct* entity, const IfcSchema::IfcRelVoidsElement::list::ptr& openings, const IfcRepresentationShapeItems& entity_shapes, const gp_Trsf& entity_trsf, IfcRepresentationShapeItems& cut_shapes);
 	void assert_closed_wire(TopoDS_Wire& wire);
 
 	bool convert_layerset(const IfcSchema::IfcProduct*, std::vector<Handle_Geom_Surface>&, std::vector<std::shared_ptr<const SurfaceStyle>>&, std::vector<double>&);
@@ -284,14 +281,6 @@ public:
 	bool split_solid_by_surface(const TopoDS_Shape&, const Handle_Geom_Surface&, TopoDS_Shape&, TopoDS_Shape&);
 	bool split_solid_by_shell(const TopoDS_Shape&, const TopoDS_Shape& s, TopoDS_Shape&, TopoDS_Shape&);
 
-#if OCC_VERSION_HEX < 0x60900
-	bool boolean_operation(const TopoDS_Shape&, const TopTools_ListOfShape&, BOPAlgo_Operation, TopoDS_Shape&);
-	bool boolean_operation(const TopoDS_Shape&, const TopoDS_Shape&, BOPAlgo_Operation, TopoDS_Shape&);
-#else
-	bool boolean_operation(const TopoDS_Shape&, const TopTools_ListOfShape&, BOPAlgo_Operation, TopoDS_Shape&, double fuzziness = -1.);
-	bool boolean_operation(const TopoDS_Shape&, const TopoDS_Shape&, BOPAlgo_Operation, TopoDS_Shape&, double fuzziness = -1.);
-#endif
-
 	bool fit_halfspace(const TopoDS_Shape& a, const TopoDS_Shape& b, TopoDS_Shape& box, double& height);
 
 	const Handle_Geom_Curve intersect(const Handle_Geom_Surface&, const Handle_Geom_Surface&);

src/ifcgeom/IfcGeomIteratorImplementation.h

@@ -303,7 +303,7 @@ namespace IfcGeom {
 				task_result_ptr_initialized = true;
 			}
 
-			progress_ = ++processed_ * 100 / tasks_.size();
+			progress_ = (int) (++processed_ * 100 / tasks_.size());
 		}
 
 		void process_concurrently() {

src/ifcgeom/IfcIndexedPolyCurve.cpp

@@ -50,7 +50,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcIndexedPolyCurve* l, TopoDS_Wi
 			coords.size() < 3 ? 0. : coords[2] * getValue(GV_LENGTH_UNIT)));
 	}
 
-	int max_index = points.size();
+	int max_index = (int) points.size();
 
 	BRepBuilderAPI_MakeWire w;
 

src/ifcgeom/IfcPolygonalFaceSet.cpp

@@ -41,12 +41,12 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcPolygonalFaceSet* pfs, TopoDS_
 	
 	for (auto& f : *polygonal_faces) {
 		loop_grouping.emplace_back();
-		loop_grouping.back().push_back(indices.size());
+		loop_grouping.back().push_back((int) indices.size());
 		indices.push_back(f->CoordIndex());
 		if (f->as<IfcSchema::IfcIndexedPolygonalFaceWithVoids>()) {
 			auto inner_coordinates = f->as<IfcSchema::IfcIndexedPolygonalFaceWithVoids>()->InnerCoordIndices();
 			for (auto& x : inner_coordinates) {
-				loop_grouping.back().push_back(indices.size());
+				loop_grouping.back().push_back((int) indices.size());
 				indices.push_back(x);
 			}
 		}

src/ifcgeom/IfcRepresentation.cpp

@@ -33,9 +33,9 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcRepresentation* l, IfcRepresen
 				TopoDS_Shape s;
 				if (convert_shape(representation_item, s)) {
 					if (s.ShapeType() == TopAbs_COMPOUND && TopoDS_Iterator(s).More() && TopoDS_Iterator(s).Value().ShapeType() == TopAbs_SOLID) {
-						TopoDS_Iterator it(s);
-						for (; it.More(); it.Next()) {
-							shapes.push_back(IfcRepresentationShapeItem(representation_item->data().id(), it.Value(), get_style(representation_item)));
+						TopoDS_Iterator topo_it(s);
+						for (; topo_it.More(); topo_it.Next()) {
+							shapes.push_back(IfcRepresentationShapeItem(representation_item->data().id(), topo_it.Value(), get_style(representation_item)));
 						}
 					} else {
 						shapes.push_back(IfcRepresentationShapeItem(representation_item->data().id(), s, get_style(representation_item)));

src/ifcgeom/IfcSurfaceCurveSweptAreaSolid.cpp

@@ -34,6 +34,7 @@
 #include <ShapeFix_Edge.hxx>
 #include <ShapeAnalysis_Surface.hxx>
 #include "../ifcgeom/IfcGeom.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -51,7 +52,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcSurfaceCurveSweptAreaSolid* l,
 			Logger::Error("Failed to convert reference surface", l);
 			return false;
 		}
-		if (count(surface_shell, TopAbs_FACE) != 1) {
+		if (util::count(surface_shell, TopAbs_FACE) != 1) {
 			Logger::Error("Non-continuous reference surface", l);
 			return false;
 		}

src/ifcgeom/IfcSweptDiskSolid.cpp

@@ -37,6 +37,7 @@
 #include "../ifcgeom_schema_agnostic/sweep_utils.h"
 #include "../ifcgeom_schema_agnostic/wire_utils.h"
 #include "../ifcgeom_schema_agnostic/face_definition.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #define Kernel MAKE_TYPE_NAME(Kernel)
 
@@ -81,7 +82,7 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcSweptDiskSolid* l, TopoDS_Shap
 	ep = l->EndParam();
 #endif
 
-	if (count(wire, TopAbs_EDGE) == 1 && sp && ep) {
+	if (util::count(wire, TopAbs_EDGE) == 1 && sp && ep) {
 		TopoDS_Vertex v0, v1;
 		TopExp::Vertices(wire, v0, v1);
 		if (v0.IsSame(v1)) {
@@ -153,11 +154,13 @@ bool IfcGeom::Kernel::convert(const IfcSchema::IfcSweptDiskSolid* l, TopoDS_Shap
 						gp_Pln plane(c1->Position().Location(), n);
 						auto face = BRepBuilderAPI_MakeFace(plane).Face();
 
-						TopoDS_Wire wire;
-						BB.MakeWire(wire);
-						BB.Add(wire, a);
-						BB.Add(wire, b);
-						BB.Add(face, wire);
+						{
+							TopoDS_Wire w;
+							BB.MakeWire(w);
+							BB.Add(w, a);
+							BB.Add(w, b);
+							BB.Add(face, w);
+						}
 						TopExp::CommonVertex(a, b, V);
 
 						BRepFilletAPI_MakeFillet2d mf2d(face);

src/ifcgeom_schema_agnostic/IfcGeomRepresentation.cpp

@@ -29,6 +29,7 @@
 
 #include "../ifcparse/IfcLogger.h"
 #include "../ifcgeom_schema_agnostic/Kernel.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 IfcGeom::Representation::Serialization::Serialization(const BRep& brep)
 	: Representation(brep.settings())
@@ -214,7 +215,7 @@ bool IfcGeom::Representation::BRep::calculate_volume(double& volume) const {
 		volume = 0.;
 
 		for (IfcGeom::IfcRepresentationShapeItems::const_iterator it = begin(); it != end(); ++it) {
-			if (Kernel::is_manifold(it->Shape())) {
+			if (util::is_manifold(it->Shape())) {
 				GProp_GProps prop;
 				BRepGProp::VolumeProperties(it->Shape(), prop);
 				volume += prop.Mass();
@@ -238,7 +239,7 @@ bool IfcGeom::Representation::BRep::calculate_projected_surface_area(const gp_Ax
 			double x, y, z;
 			surface_area_along_direction(settings().deflection_tolerance(), it->Shape(), ax, x, y, z);
 
-			if (Kernel::is_manifold(it->Shape())) {
+			if (util::is_manifold(it->Shape())) {
 				x /= 2.;
 				y /= 2.;
 				z /= 2.;

src/ifcgeom_schema_agnostic/IfcGeomTree.h

@@ -410,7 +410,7 @@ namespace IfcGeom {
 					
 					// Assumption is that the number of styles is small, so the linear lookup time is not significant.
 					auto sit = std::find(styles_.begin(), styles_.end(), *adaptor);
-					int index;
+					size_t index;
 					if (sit == styles_.end()) {
 						index = styles_.size();
 						styles_.push_back(*adaptor);
@@ -420,7 +420,7 @@ namespace IfcGeom {
 
 					TopExp_Explorer exp(it.Value(), TopAbs_FACE);
 					for (; exp.More(); exp.Next()) {
-						face_styles_.Bind(exp.Current(), index);
+						face_styles_.Bind(exp.Current(), (int) index);
 					}
 				}
 			}

src/ifcgeom_schema_agnostic/Kernel.cpp

@@ -23,33 +23,6 @@ IfcGeom::Kernel::Kernel(IfcParse::IfcFile* file) {
 	}
 }
 
-int IfcGeom::Kernel::count(const TopoDS_Shape& s, TopAbs_ShapeEnum t, bool unique) {
-	if (unique) {
-		TopTools_IndexedMapOfShape map;
-		TopExp::MapShapes(s, t, map);
-		return map.Extent();
-	} else {
-		int i = 0;
-		TopExp_Explorer exp(s, t);
-		for (; exp.More(); exp.Next()) {
-			++i;
-		}
-		return i;
-	}
-}
-
-
-int IfcGeom::Kernel::surface_genus(const TopoDS_Shape& s) {
-	int nv = count(s, TopAbs_VERTEX, true);
-	int ne = count(s, TopAbs_EDGE, true);
-	int nf = count(s, TopAbs_FACE, true);
-
-	const int euler = nv - ne + nf;
-	const int genus = (2 - euler) / 2;
-
-	return genus;
-}
-
 IfcGeom::impl::KernelFactoryImplementation& IfcGeom::impl::kernel_implementations() {
 	static KernelFactoryImplementation impl;
 	return impl;
@@ -247,102 +220,3 @@ std::map<std::string, IfcUtil::IfcBaseEntity*> IfcGeom::Kernel::get_layers(IfcUt
 
 	throw IfcParse::IfcException("Unexpected entity " + inst->declaration().name());
 }
-
-bool IfcGeom::Kernel::is_manifold(const TopoDS_Shape& a) {
-	if (a.ShapeType() == TopAbs_COMPOUND || a.ShapeType() == TopAbs_SOLID) {
-		TopoDS_Iterator it(a);
-		for (; it.More(); it.Next()) {
-			if (!is_manifold(it.Value())) {
-				return false;
-			}
-		}
-		return true;
-	} else {
-		TopTools_IndexedDataMapOfShapeListOfShape map;
-		TopExp::MapShapesAndAncestors(a, TopAbs_EDGE, TopAbs_FACE, map);
-
-		for (int i = 1; i <= map.Extent(); ++i) {
-			const TopoDS_Edge& e = TopoDS::Edge(map.FindKey(i));
-
-			TopoDS_Vertex v0, v1;
-			TopExp::Vertices(e, v0, v1);
-			const bool degenerate = !v0.IsNull() && !v1.IsNull() && v0.IsSame(v1);
-
-			if (degenerate) {
-				continue;
-			}
-
-			if (map.FindFromIndex(i).Extent() != 2) {
-				return false;
-			}
-		}
-
-		return true;
-	}
-}
-
-bool IfcGeom::util::is_nested_compound_of_solid(const TopoDS_Shape& s, int depth) {
-	if (s.ShapeType() == TopAbs_COMPOUND) {
-		TopoDS_Iterator it(s);
-		for (; it.More(); it.Next()) {
-			if (!is_nested_compound_of_solid(it.Value(), depth + 1)) {
-				return false;
-			}
-		}
-		return true;
-	} else if (s.ShapeType() == TopAbs_SOLID) {
-		return depth > 0;
-	} else {
-		return false;
-	}
-}
-
-namespace {
-	template <typename T> struct dimension_count {};
-	template <>           struct dimension_count <gp_Trsf2d > { static const int n = 2; };
-	template <>           struct dimension_count <gp_GTrsf2d> { static const int n = 2; };
-	template <>           struct dimension_count < gp_Trsf  > { static const int n = 3; };
-	template <>           struct dimension_count < gp_GTrsf > { static const int n = 3; };
-
-	template <typename T>
-	bool is_identity_helper(const T& t, double tolerance) {
-		// Note the {1, n+1} range due to Open Cascade's 1-based indexing
-		// Note the {1, n+2} range due to the translation part of the matrix
-		for (int i = 1; i < dimension_count<T>::n + 2; ++i) {
-			for (int j = 1; j < dimension_count<T>::n + 1; ++j) {
-				const double iden_value = i == j ? 1. : 0.;
-				const double trsf_value = t.Value(j, i);
-				if (fabs(trsf_value - iden_value) > tolerance) {
-					return false;
-				}
-			}
-		}
-		return true;
-	}
-}
-
-bool IfcGeom::Kernel::is_identity(const gp_Trsf2d& t, double tolerance) {
-	return is_identity_helper(t, tolerance);
-}
-
-bool IfcGeom::Kernel::is_identity(const gp_GTrsf2d& t, double tolerance) {
-	return is_identity_helper(t, tolerance);
-}
-
-bool IfcGeom::Kernel::is_identity(const gp_Trsf& t, double tolerance) {
-	return is_identity_helper(t, tolerance);
-}
-
-bool IfcGeom::Kernel::is_identity(const gp_GTrsf& t, double tolerance) {
-	return is_identity_helper(t, tolerance);
-}
-
-gp_Trsf IfcGeom::Kernel::combine_offset_and_rotation(const gp_Vec & offset, const gp_Quaternion & rotation) {
-	auto offset_transform = gp_Trsf{};
-	offset_transform.SetTranslation(offset);
-
-	auto rotation_transform = gp_Trsf{};
-	rotation_transform.SetRotation(rotation);
-
-	return rotation_transform * offset_transform;
-}

src/ifcgeom_schema_agnostic/Kernel.h

@@ -124,38 +124,10 @@ namespace IfcGeom {
 			return implementation_->convert_placement(item, trsf);
 		}
 
-		IFC_PARSE_API static int count(const TopoDS_Shape&, TopAbs_ShapeEnum, bool unique = false);
-		IFC_PARSE_API static int surface_genus(const TopoDS_Shape&);
-
-		IFC_PARSE_API static bool is_manifold(const TopoDS_Shape& a);
+		
 		IFC_PARSE_API static IfcUtil::IfcBaseEntity* get_decomposing_entity(IfcUtil::IfcBaseEntity*, bool include_openings = true);
 		IFC_PARSE_API static std::map<std::string, IfcUtil::IfcBaseEntity*> get_layers(IfcUtil::IfcBaseEntity*);
 
-		// For axis placements detect equality early in order for the
-		// relatively computionaly expensive gp_Trsf calculation to be skipped
-		IFC_PARSE_API static bool axis_equal(const gp_Ax3& a, const gp_Ax3& b, double tolerance) {
-			if (!a.Location().IsEqual(b.Location(), tolerance)) return false;
-			// Note that the tolerance below is angular, above is linear. Since architectural
-			// objects are about 1m'ish in scale, it should be somewhat equivalent. Besides,
-			// this is mostly a filter for NULL or default values in the placements.
-			if (!a.Direction().IsEqual(b.Direction(), tolerance)) return false;
-			if (!a.XDirection().IsEqual(b.XDirection(), tolerance)) return false;
-			if (!a.YDirection().IsEqual(b.YDirection(), tolerance)) return false;
-			return true;
-		}
-
-		IFC_PARSE_API static bool axis_equal(const gp_Ax2d& a, const gp_Ax2d& b, double tolerance) {
-			if (!a.Location().IsEqual(b.Location(), tolerance)) return false;
-			if (!a.Direction().IsEqual(b.Direction(), tolerance)) return false;
-			return true;
-		}
-
-		IFC_PARSE_API static bool is_identity(const gp_Trsf2d& t, double tolerance);
-		IFC_PARSE_API static bool is_identity(const gp_GTrsf2d& t, double tolerance);
-		IFC_PARSE_API static bool is_identity(const gp_Trsf& t, double tolerance);
-		IFC_PARSE_API static bool is_identity(const gp_GTrsf& t, double tolerance);
-
-		IFC_PARSE_API static gp_Trsf combine_offset_and_rotation(const gp_Vec &offset, const gp_Quaternion& rotation);
 	};
 
 	namespace impl {
@@ -171,10 +143,6 @@ namespace IfcGeom {
 		KernelFactoryImplementation& kernel_implementations();
 	}
 
-	namespace util {
-		bool is_nested_compound_of_solid(const TopoDS_Shape& s, int depth = 0);
-	}
-
 	class IFC_GEOM_API geometry_exception : public std::exception {
 	protected:
 		std::string message;

src/ifcgeom_schema_agnostic/base_utils.cpp

@@ -0,0 +1,155 @@
+#include "base_utils.h"
+
+#include <TopTools_IndexedMapOfShape.hxx>
+#include <TopExp.hxx>
+#include <TopExp_Explorer.hxx>
+#include <TopoDS.hxx>
+#include <TopoDS_Vertex.hxx>
+#include <gp_GTrsf.hxx>
+#include <gp_GTrsf2d.hxx>
+
+// For axis placements detect equality early in order for the
+// relatively computionaly expensive gp_Trsf calculation to be skipped
+bool IfcGeom::util::axis_equal(const gp_Ax3 & a, const gp_Ax3 & b, double tolerance) {
+	if (!a.Location().IsEqual(b.Location(), tolerance)) return false;
+	// Note that the tolerance below is angular, above is linear. Since architectural
+	// objects are about 1m'ish in scale, it should be somewhat equivalent. Besides,
+	// this is mostly a filter for NULL or default values in the placements.
+	if (!a.Direction().IsEqual(b.Direction(), tolerance)) return false;
+	if (!a.XDirection().IsEqual(b.XDirection(), tolerance)) return false;
+	if (!a.YDirection().IsEqual(b.YDirection(), tolerance)) return false;
+	return true;
+}
+
+bool IfcGeom::util::axis_equal(const gp_Ax2d & a, const gp_Ax2d & b, double tolerance) {
+	if (!a.Location().IsEqual(b.Location(), tolerance)) return false;
+	if (!a.Direction().IsEqual(b.Direction(), tolerance)) return false;
+	return true;
+}
+
+int IfcGeom::util::count(const TopoDS_Shape& s, TopAbs_ShapeEnum t, bool unique) {
+	if (unique) {
+		TopTools_IndexedMapOfShape map;
+		TopExp::MapShapes(s, t, map);
+		return map.Extent();
+	} else {
+		int i = 0;
+		TopExp_Explorer exp(s, t);
+		for (; exp.More(); exp.Next()) {
+			++i;
+		}
+		return i;
+	}
+}
+
+
+int IfcGeom::util::surface_genus(const TopoDS_Shape& s) {
+	int nv = count(s, TopAbs_VERTEX, true);
+	int ne = count(s, TopAbs_EDGE, true);
+	int nf = count(s, TopAbs_FACE, true);
+
+	const int euler = nv - ne + nf;
+	const int genus = (2 - euler) / 2;
+
+	return genus;
+}
+
+bool IfcGeom::util::is_manifold(const TopoDS_Shape& a) {
+	if (a.ShapeType() == TopAbs_COMPOUND || a.ShapeType() == TopAbs_SOLID) {
+		TopoDS_Iterator it(a);
+		for (; it.More(); it.Next()) {
+			if (!is_manifold(it.Value())) {
+				return false;
+			}
+		}
+		return true;
+	} else {
+		TopTools_IndexedDataMapOfShapeListOfShape map;
+		TopExp::MapShapesAndAncestors(a, TopAbs_EDGE, TopAbs_FACE, map);
+
+		for (int i = 1; i <= map.Extent(); ++i) {
+			const TopoDS_Edge& e = TopoDS::Edge(map.FindKey(i));
+
+			TopoDS_Vertex v0, v1;
+			TopExp::Vertices(e, v0, v1);
+			const bool degenerate = !v0.IsNull() && !v1.IsNull() && v0.IsSame(v1);
+
+			if (degenerate) {
+				continue;
+			}
+
+			if (map.FindFromIndex(i).Extent() != 2) {
+				return false;
+			}
+		}
+
+		return true;
+	}
+}
+
+
+bool IfcGeom::util::is_nested_compound_of_solid(const TopoDS_Shape& s, int depth) {
+	if (s.ShapeType() == TopAbs_COMPOUND) {
+		TopoDS_Iterator it(s);
+		for (; it.More(); it.Next()) {
+			if (!is_nested_compound_of_solid(it.Value(), depth + 1)) {
+				return false;
+			}
+		}
+		return true;
+	} else if (s.ShapeType() == TopAbs_SOLID) {
+		return depth > 0;
+	} else {
+		return false;
+	}
+}
+
+namespace {
+	template <typename T> struct dimension_count {};
+	template <>           struct dimension_count <gp_Trsf2d > { static const int n = 2; };
+	template <>           struct dimension_count <gp_GTrsf2d> { static const int n = 2; };
+	template <>           struct dimension_count < gp_Trsf  > { static const int n = 3; };
+	template <>           struct dimension_count < gp_GTrsf > { static const int n = 3; };
+
+	template <typename T>
+	bool is_identity_helper(const T& t, double tolerance) {
+		// Note the {1, n+1} range due to Open Cascade's 1-based indexing
+		// Note the {1, n+2} range due to the translation part of the matrix
+		for (int i = 1; i < dimension_count<T>::n + 2; ++i) {
+			for (int j = 1; j < dimension_count<T>::n + 1; ++j) {
+				const double iden_value = i == j ? 1. : 0.;
+				const double trsf_value = t.Value(j, i);
+				if (fabs(trsf_value - iden_value) > tolerance) {
+					return false;
+				}
+			}
+		}
+		return true;
+	}
+}
+
+bool IfcGeom::util::is_identity(const gp_Trsf2d& t, double tolerance) {
+	return is_identity_helper(t, tolerance);
+}
+
+bool IfcGeom::util::is_identity(const gp_GTrsf2d& t, double tolerance) {
+	return is_identity_helper(t, tolerance);
+}
+
+bool IfcGeom::util::is_identity(const gp_Trsf& t, double tolerance) {
+	return is_identity_helper(t, tolerance);
+}
+
+bool IfcGeom::util::is_identity(const gp_GTrsf& t, double tolerance) {
+	return is_identity_helper(t, tolerance);
+}
+
+gp_Trsf IfcGeom::util::combine_offset_and_rotation(const gp_Vec & offset, const gp_Quaternion & rotation) {
+	auto offset_transform = gp_Trsf{};
+	offset_transform.SetTranslation(offset);
+
+	auto rotation_transform = gp_Trsf{};
+	rotation_transform.SetRotation(rotation);
+
+	return rotation_transform * offset_transform;
+}

src/ifcgeom_schema_agnostic/base_utils.h

@@ -0,0 +1,32 @@
+#ifndef BASE_UTILS_H
+#define BASE_UTILS_H
+
+#include <TopoDS_Shape.hxx>
+#include <gp_Ax3.hxx>
+
+namespace IfcGeom {
+	namespace util {
+
+		int count(const TopoDS_Shape&, TopAbs_ShapeEnum, bool unique = false);
+		int surface_genus(const TopoDS_Shape&);
+
+		bool is_manifold(const TopoDS_Shape& a);
+
+		// For axis placements detect equality early in order for the
+		// relatively computionaly expensive gp_Trsf calculation to be skipped
+		bool axis_equal(const gp_Ax3& a, const gp_Ax3& b, double tolerance);
+
+		bool axis_equal(const gp_Ax2d& a, const gp_Ax2d& b, double tolerance);
+
+		bool is_identity(const gp_Trsf2d& t, double tolerance);
+		bool is_identity(const gp_GTrsf2d& t, double tolerance);
+		bool is_identity(const gp_Trsf& t, double tolerance);
+		bool is_identity(const gp_GTrsf& t, double tolerance);
+
+		gp_Trsf combine_offset_and_rotation(const gp_Vec &offset, const gp_Quaternion& rotation);
+
+		bool is_nested_compound_of_solid(const TopoDS_Shape& s, int depth = 0);
+	}
+}
+
+#endif

src/ifcgeom_schema_agnostic/boolean_utils.cpp

@@ -1,6 +1,7 @@
 #include "boolean_utils.h"
 
 #include "../ifcgeom_schema_agnostic/IfcGeomTree.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #include <BRepBuilderAPI_Copy.hxx>
 #include <TopExp_Explorer.hxx>
@@ -17,8 +18,18 @@
 #include <ShapeAnalysis_Surface.hxx>
 #include <BRepBuilderAPI_MakeFace.hxx>
 #include <Standard_Version.hxx>
+#include <BRepAlgoAPI_Fuse.hxx>
+#include <BRepPrimAPI_MakePrism.hxx>
+#include <BOPAlgo_PaveFiller.hxx>
+#include <BOPAlgo_Alerts.hxx>
+#include <ShapeFix_Shape.hxx>
+#include <BRepCheck_Analyzer.hxx>
+#include <BRepCheck_ListIteratorOfListOfStatus.hxx>
+#include <BRepCheck.hxx>
+#include <ShapeAnalysis_Edge.hxx>
 
 #include <vector>
+#include <thread>
 
 void IfcGeom::util::copy_operand(const TopTools_ListOfShape & l, TopTools_ListOfShape & r) {
 #if OCC_VERSION_HEX < 0x70000
@@ -609,6 +620,18 @@ bool IfcGeom::util::boolean_subtraction_2d_using_builder(const TopoDS_Shape & a_
 		}
 	}
 
+	{
+		TopoDS_Compound C;
+		BRep_Builder BB;
+		BB.MakeCompound(C);
+		
+		for (auto& w : wires) {
+			BB.Add(C, w);
+		}
+
+		BRepTools::Write(C, "debug.brep");
+	}
+
 	shape_index = 0;
 	edge_index = 0;
 
@@ -786,4 +809,551 @@ bool IfcGeom::util::points_on_planar_face_generator::operator()(gp_Pnt& p) {
 	}
 
 	return false;
-}
+}
+
+
+bool IfcGeom::util::boolean_operation(const boolean_settings& settings, const TopoDS_Shape& a_input, const TopTools_ListOfShape& b_input, BOPAlgo_Operation op, TopoDS_Shape& result, double fuzziness) {
+	using namespace std::string_literals;
+
+	const bool do_unify = true;
+	const bool do_subtraction_eliminate_disjoint_bbox = true;
+	const bool do_subtraction_eliminate_touching = true;
+	const bool do_attempt_2d_boolean = settings.attempt_2d;
+	const bool debug = settings.debug;
+
+	std::string debug_identifier;
+	if (debug) {
+		static my_thread_local size_t operation_counter_ = 0;
+		std::stringstream ss;
+		ss << "bool-" << std::this_thread::get_id() << "-" << (operation_counter_++);
+		debug_identifier = ss.str();
+		Logger::Notice("Boolean debug identifier: " + debug_identifier);
+	}
+
+	if (fuzziness < 0.) {
+		fuzziness = settings.precision / 10.;
+	}
+
+	// @todo, it does seem a bit odd, we first triangulate non-planar faces
+	// to later unify them again. Can we make this a bit more intelligent?
+	TopoDS_Shape a;
+	TopTools_ListOfShape b;
+
+	if (do_unify) {
+		PERF("boolean operation: unifying operands");
+
+		a = unify(a_input, fuzziness * 1000.);
+
+		Logger::Message(
+			Logger::LOG_DEBUG,
+			"Simplified operand A from "s +
+			std::to_string(count(a_input, TopAbs_FACE)) +
+			" to "s +
+			std::to_string(count(a, TopAbs_FACE))
+		);
+
+		{
+			TopTools_ListIteratorOfListOfShape it(b_input);
+			for (; it.More(); it.Next()) {
+				b.Append(unify(it.Value(), fuzziness));
+				Logger::Message(
+					Logger::LOG_DEBUG,
+					"Simplified operand B from "s +
+					std::to_string(count(it.Value(), TopAbs_FACE)) +
+					" to "s +
+					std::to_string(count(b.Last(), TopAbs_FACE))
+				);
+			}
+		}
+	} else {
+		a = a_input;
+		b = b_input;
+	}
+
+	bool success = false;
+	BRepAlgoAPI_BooleanOperation* builder;
+	TopTools_ListOfShape b_tmp;
+
+	if (op == BOPAlgo_CUT) {
+		builder = new BRepAlgoAPI_Cut();
+
+		if (do_subtraction_eliminate_disjoint_bbox) {
+			PERF("boolean subtraction: eliminate disjoint bbox");
+
+			auto N = bounding_box_overlap(fuzziness, a, b, b_tmp);
+			if (N) {
+				Logger::Notice("Eliminated " + std::to_string(N) + " disjoint operands");
+				std::swap(b, b_tmp);
+			}
+		}
+
+		if (do_subtraction_eliminate_touching) {
+			PERF("boolean subtraction: eliminate touching");
+
+			b_tmp.Clear();
+			auto N = eliminate_touching_operands(fuzziness, a, b, b_tmp);
+			if (N) {
+				Logger::Notice("Eliminated " + std::to_string(N) + " touching operands");
+				std::swap(b, b_tmp);
+			}
+		}
+
+	} else if (op == BOPAlgo_COMMON) {
+		builder = new BRepAlgoAPI_Common();
+	} else if (op == BOPAlgo_FUSE) {
+		builder = new BRepAlgoAPI_Fuse();
+	} else {
+		return false;
+	}
+
+	if (b.Extent() == 0) {
+		result = a;
+		return true;
+	}
+
+	if (Logger::LOG_NOTICE >= Logger::Verbosity()) {
+		PERF("preliminary manifoldness check");
+
+		Logger::Notice("Operand A is " + (is_manifold(a) ? ""s : "non-"s) + "manifold");
+
+		TopTools_ListIteratorOfListOfShape it(b);
+		for (int i = 0; it.More(); it.Next(), ++i) {
+			Logger::Notice("Operand B " + std::to_string(i) + " is " + (is_manifold(it.Value()) ? ""s : "non-"s) + "manifold");
+		}
+	}
+
+	// Find a sensible value for the fuzziness, based on precision
+	// and limited by edge lengths and vertex-edge distances.
+	double min_length_orig;
+
+	{
+		PERF("boolean operation: min edge length");
+
+		min_length_orig = min_edge_length(a);
+		TopTools_ListIteratorOfListOfShape it(b);
+		for (; it.More(); it.Next()) {
+			double d = min_edge_length(it.Value());
+			if (d < min_length_orig) {
+				min_length_orig = d;
+			}
+		}
+	}
+
+	{
+		PERF("boolean operation: min vertex-edge dist");
+
+		double d = min_vertex_edge_distance(a, settings.precision, min_length_orig);
+		if (d < min_length_orig) {
+			min_length_orig = d;
+		}
+
+		TopTools_ListIteratorOfListOfShape it(b);
+		for (; it.More(); it.Next()) {
+			d = min_vertex_edge_distance(it.Value(), settings.precision, min_length_orig);
+			if (d < min_length_orig) {
+				min_length_orig = d;
+			}
+		}
+	}
+
+	const double fuzz = (std::min)(min_length_orig / 3., fuzziness);
+
+	Logger::Notice("Used fuzziness: " + std::to_string(fuzz));
+
+	TopTools_ListOfShape s1s;
+	s1s.Append(copy_operand(a));
+
+	if (debug) {
+		TopTools_ListOfShape* lists[2] = { &s1s, &b };
+		static std::string operand_names[2] = { "a", "b" };
+		for (int i = 0; i < 2; ++i) {
+			TopTools_ListIteratorOfListOfShape it(*lists[i]);
+			for (int j = 0; it.More(); it.Next(), ++j) {
+				std::string fn = debug_identifier + "-" + operand_names[i] + "-" + std::to_string(j) + ".brep";
+				BRepTools::Write(it.Value(), fn.c_str());
+			}
+		}
+	}
+
+	if (op == BOPAlgo_CUT) {
+		TopoDS_Face a_face;
+		std::pair<double, double> a_interval;
+
+		TopTools_ListOfShape b_faces, b_remainder_3d;
+
+		bool is_extrusion_a = false;
+		if (do_attempt_2d_boolean) {
+			PERF("boolean subtraction: extrusion check");
+
+			is_extrusion_a = is_extrusion(gp::DY(), a, a_face, a_interval);
+		}
+
+		if (is_extrusion_a) {
+			Logger::Notice("Operand A 1/1 is an extrusion");
+
+			TopTools_ListIteratorOfListOfShape it(b);
+			for (int nb = 1; it.More(); it.Next(), ++nb) {
+				bool process_2d = false;
+				TopoDS_Face b_face;
+				std::pair<double, double> b_interval;
+
+				bool is_extrusion_b;
+				{
+					PERF("boolean subtraction: extrusion check");
+
+					is_extrusion_b = is_extrusion(gp::DY(), it.Value(), b_face, b_interval);
+				}
+
+				if (is_extrusion_b) {
+					Logger::Notice("Operand B " + std::to_string(nb) + "/" + std::to_string(b.Extent()) + " is an extrusion");
+
+					if (b_interval.first < a_interval.first + fuzz && b_interval.second > a_interval.second - fuzz) {
+						Logger::Notice("Operand B creates a through hole");
+
+						// Align b with a operand
+						gp_Trsf trsf;
+						trsf.SetTranslation(gp_Vec(gp::DY()) * (a_interval.first - b_interval.first));
+
+						b_faces.Append(b_face.Moved(trsf));
+						process_2d = true;
+					}
+				}
+
+				if (!process_2d) {
+					b_remainder_3d.Append(it.Value());
+				}
+			}
+
+			if (b_faces.Extent()) {
+				TopoDS_Shape face_result;
+
+				bool boolean_op_2d_success;
+				{
+					PERF("boolean operation: 2d builder");
+					// First try using face builder
+
+					boolean_op_2d_success = boolean_subtraction_2d_using_builder(a_face, b_faces, face_result, fuzziness);
+				}
+
+				if (!boolean_op_2d_success) {
+					PERF("boolean operation: 2d");
+					// Retry using generic 2d using boolean algo on faces
+
+					boolean_op_2d_success = boolean_operation(settings, a_face, b_faces, op, face_result, fuzziness);
+				}
+
+				if (boolean_op_2d_success) {
+					PERF("boolean operation: 2d to 3d");
+
+					BRepPrimAPI_MakePrism mp(face_result, gp_Vec(gp::DY()) * (a_interval.second - a_interval.first));
+					if (mp.IsDone()) {
+						if (b_remainder_3d.Extent()) {
+							Logger::Notice(std::to_string(b_remainder_3d.Extent()) + " operands remaining to process in 3D");
+							b = b_remainder_3d;
+							s1s.Clear();
+							s1s.Append(mp.Shape());
+						} else {
+							Logger::Notice("Processed fully in 2D");
+							result = mp.Shape();
+							return true;
+						}
+					} else {
+						Logger::Notice("Failed to extrude 2D boolean result. Retrying in 3D.");
+					}
+				} else {
+					Logger::Notice("Failed to perform 2D boolean operation. Retrying in 3D.");
+				}
+			} else {
+				Logger::Notice("No second operands can be processed as 2D inner bounds. Retrying in 3D.");
+			}
+		}
+	}
+
+#if OCC_VERSION_HEX >= 0x70000
+	builder->SetNonDestructive(true);
+#endif
+	builder->SetFuzzyValue(fuzz);
+	builder->SetArguments(s1s);
+	copy_operand(b, b_tmp);
+	std::swap(b, b_tmp);
+	builder->SetTools(b);
+	{
+		PERF("boolean operation: build");
+
+		builder->Build();
+	}
+	if (builder->IsDone()) {
+		if (builder->DSFiller()->HasWarning(STANDARD_TYPE(BOPAlgo_AlertAcquiredSelfIntersection))) {
+			Logger::Notice("Builder reports self-intersection in output");
+			success = false;
+		} else {
+			TopoDS_Shape r = *builder;
+
+			{
+				PERF("boolean operation: shape healing");
+
+				ShapeFix_Shape fix(r);
+				try {
+					fix.SetMaxTolerance(fuzz);
+					fix.Perform();
+					r = fix.Shape();
+				} catch (...) {
+					Logger::Error("Shape healing failed on boolean result");
+				}
+			}
+
+			{
+				PERF("boolean operation: shape analysis");
+
+				BRepCheck_Analyzer ana(r);
+				success = ana.IsValid() != 0;
+
+				if (!success) {
+					Logger::Notice("Boolean operation yields invalid result");
+
+					std::stringstream str;
+					bool any_emitted = false;
+
+					std::function<void(const TopoDS_Shape&)> dump;
+					dump = [&ana, &str, &dump, &any_emitted](const TopoDS_Shape& s) {
+						if (!ana.Result(s).IsNull()) {
+							BRepCheck_ListIteratorOfListOfStatus itl;
+							itl.Initialize(ana.Result(s)->Status());
+							for (; itl.More(); itl.Next()) {
+								if (itl.Value() != BRepCheck_NoError) {
+									if (any_emitted) {
+										str << ", ";
+									}
+									BRepCheck::Print(itl.Value(), str);
+									str.seekp(str.tellp() - (std::streamoff)1);
+									str << " on ";
+									TopAbs::Print(s.ShapeType(), str);
+									any_emitted = true;
+								}
+							}
+						}
+						for (TopoDS_Iterator it(s); it.More(); it.Next()) {
+							dump(it.Value());
+						}
+					};
+
+					dump(r);
+
+					Logger::Notice(str.str());
+				}
+			}
+
+			if (success) {
+
+				{
+					PERF("boolean operation: manifoldness check");
+
+					success = !is_manifold(a) || is_manifold(r);
+				}
+
+				if (!success) {
+					PERF("boolean operation: manifoldness check excemption");
+
+					// An excemption for the requirement to be manifold: When the cut operands have overlapping edge belonging to faces that do not overlap.
+					bool operands_nonmanifold = false;
+					if (op == BOPAlgo_CUT) {
+						TopTools_IndexedMapOfShape edges;
+						TopTools_IndexedDataMapOfShapeListOfShape map;
+						for (TopTools_ListIteratorOfListOfShape it2(b); it2.More(); it2.Next()) {
+							auto& bb = it2.Value();
+							TopExp::MapShapes(bb, TopAbs_EDGE, edges);
+							TopExp::MapShapesAndAncestors(bb, TopAbs_EDGE, TopAbs_FACE, map);
+						}
+						IfcGeom::impl::tree<int> tree;
+						for (int i = 1; i <= edges.Extent(); ++i) {
+							tree.add(i, edges.FindKey(i));
+						}
+						for (int i = 1; i <= edges.Extent(); ++i) {
+							const TopoDS_Edge& ei = TopoDS::Edge(edges.FindKey(i));
+							Bnd_Box bb;
+							BRepBndLib::Add(ei, bb);
+							bb.Enlarge(fuzziness);
+							auto ii = tree.select_box(bb, false);
+							for (int j : ii) {
+								if (j != i) {
+									const TopoDS_Edge& ej = TopoDS::Edge(edges.FindKey(j));
+									ShapeAnalysis_Edge sae;
+									double f = fuzziness;
+									bool edges_overlapping = sae.CheckOverlapping(ei, ej, f, 0.) ||
+										sae.CheckOverlapping(ej, ei, f, 0.);
+
+									if (edges_overlapping) {
+										auto faces_i = map.FindFromKey(edges.FindKey(i));
+										auto faces_j = map.FindFromKey(edges.FindKey(j));
+										bool overlap = false;
+										for (TopTools_ListIteratorOfListOfShape it4(faces_i); it4.More(); it4.Next()) {
+											auto& fi = it4.Value();
+											for (TopTools_ListIteratorOfListOfShape it2(faces_j); it2.More(); it2.Next()) {
+												auto& fj = it2.Value();
+												if (faces_overlap(TopoDS::Face(fi), TopoDS::Face(fj))) {
+													overlap = true;
+												}
+											}
+											if (overlap) {
+												break;
+											}
+										}
+										operands_nonmanifold = !overlap;
+										break;
+									}
+								}
+							}
+							if (operands_nonmanifold) {
+								break;
+							}
+						}
+					}
+					success = operands_nonmanifold;
+				}
+
+				if (success) {
+
+					bool all_faces_included_in_result = true;
+					bool has_open_shells = false;
+
+					if (op == BOPAlgo_CUT) {
+						PERF("boolean operation: open shell face adition check");
+
+						for (TopExp_Explorer exp(a, TopAbs_SHELL); exp.More(); exp.Next()) {
+							if (!exp.Current().Closed()) {
+								// This 'face addition check' is only done when the first operand
+								// contains open shells (which was initially the aim of this check
+								// see #1472).
+								// Later in #1914 we found that the logic to apply openings in groups
+								// of similar edge lengths can create a situation of inner voids, which
+								// trigger a false positive in this check. This could have also been
+								// solved below by checking whether the opening(s) are included as a
+								// unmodified (interior) shell within a solid of multiple shells.
+								// Checking for open shells in first operand was quicker and more
+								// straightforward. The question still is whether in cases like #1472
+								// we need to first try the boolean union as solid/solid interference
+								// to trigger this case or whether we can immediately proceed to a face/
+								// solid operation.
+								has_open_shells = true;
+								break;
+							}
+						}
+
+						if (has_open_shells) {
+							TopTools_IndexedMapOfShape faces;
+							TopExp::MapShapes(r, TopAbs_FACE, faces);
+							for (TopExp_Explorer exp(a, TopAbs_FACE); exp.More(); exp.Next()) {
+								auto& f = TopoDS::Face(exp.Current());
+								if (!faces.Contains(f)) {
+									all_faces_included_in_result = false;
+									break;
+								}
+							}
+						} else {
+							all_faces_included_in_result = false;
+						}
+					}
+
+					int result_n_faces = count(r, TopAbs_FACE);
+					int first_op_n_faces = count(a, TopAbs_FACE);
+
+					if (op == BOPAlgo_CUT && has_open_shells && all_faces_included_in_result && result_n_faces > first_op_n_faces) {
+						success = false;
+						Logger::Notice("Boolean result discarded because subtractions results in only the addition of faces");
+					} else {
+						// when there are edges or vertex-edge distances close to the used fuzziness, the
+						// output is not trusted and the operation is attempted with a higher fuzziness.
+						int reason = 0;
+						double v;
+
+						{
+							PERF("boolean operation: result min edge length check");
+
+							if ((v = min_edge_length(r)) < fuzziness * 3.) {
+								reason = 0;
+								success = false;
+
+								goto skip_further_checks;
+							}
+						}
+
+						{
+							PERF("boolean operation: result min vertex-edge dist check");
+
+							if ((v = min_vertex_edge_distance(r, settings.precision, fuzziness * 3.)) < fuzziness * 3.) {
+								reason = 1;
+								success = false;
+
+								goto skip_further_checks;
+							}
+						}
+
+						{
+							PERF("boolean operation: result min face-face dist check");
+
+							if ((v = min_face_face_distance(r, 1.e-4)) < 1.e-4) {
+								// #2095 Check if this distance wasn't already realized in the input first operand.
+								if (v < min_face_face_distance(a, 1.e-4)) {
+									reason = 2;
+									success = false;
+								}
+							}
+						}
+
+					skip_further_checks:
+						if (!success) {
+							static const char* const reason_strings[] = { "edge length", "vertex-edge", "face-face" };
+							std::stringstream str;
+							str << "Boolean operation result failing " << reason_strings[reason] << " interference check, with fuzziness " << fuzziness << " with length " << v;
+							Logger::Notice(str.str());
+						}
+					}
+
+					if (success) {
+						result = r;
+					}
+
+				} else {
+					Logger::Notice("Boolean operation yields non-manifold result");
+				}
+			}
+		}
+	} else {
+		std::stringstream str;
+
+#if OCC_VERSION_HEX >= 0x70200
+
+		if (builder->HasError(STANDARD_TYPE(BOPAlgo_AlertBOPNotAllowed))) {
+			Logger::Error("Invalid operands. Using first operand");
+			result = a;
+			success = true;
+		}
+#endif
+
+#if OCC_VERSION_HEX >= 0x70000
+		builder->DumpErrors(str);
+#else
+		str << "Error code: " << builder->ErrorStatus();
+#endif
+		std::string str_str = str.str();
+		if (str_str.size()) {
+			Logger::Notice(str_str);
+		}
+	}
+	delete builder;
+	if (!success) {
+		const double new_fuzziness = fuzziness * 10.;
+		if (new_fuzziness - 1e-15 <= settings.precision * 10000. && new_fuzziness < min_length_orig) {
+			return boolean_operation(settings, a, b, op, result, new_fuzziness);
+		} else {
+			Logger::Notice("No longer attempting boolean operation with higher fuzziness");
+		}
+	}
+	return success && !result.IsNull();
+}
+
+bool IfcGeom::util::boolean_operation(const boolean_settings& settings, const TopoDS_Shape& a, const TopoDS_Shape& b, BOPAlgo_Operation op, TopoDS_Shape& result, double fuzziness) {
+	TopTools_ListOfShape bs;
+	bs.Append(b);
+	return boolean_operation(settings, a, bs, op, result, fuzziness);
+}

src/ifcgeom_schema_agnostic/boolean_utils.h

@@ -28,6 +28,7 @@
 #include <BRep_Tool.hxx>
 #include <BRepTools.hxx>
 #include <TopTools_IndexedMapOfShape.hxx>
+#include <BOPAlgo_Operation.hxx>
 
 namespace IfcGeom {
 	namespace util {
@@ -84,6 +85,14 @@ namespace IfcGeom {
 
 		bool boolean_subtraction_2d_using_builder(const TopoDS_Shape& a_input, const TopTools_ListOfShape& b_input, TopoDS_Shape& result, double eps);
 
+		struct boolean_settings {
+			bool debug, attempt_2d;
+			double precision;
+		};
+
+		bool boolean_operation(const boolean_settings& settings, const TopoDS_Shape&, const TopTools_ListOfShape&, BOPAlgo_Operation, TopoDS_Shape&, double fuzziness = -1.);
+
+		bool boolean_operation(const boolean_settings& settings, const TopoDS_Shape&, const TopoDS_Shape&, BOPAlgo_Operation, TopoDS_Shape&, double fuzziness = -1.);
 	}
 }
 

src/ifcgeom_schema_agnostic/sweep_utils.cpp

@@ -1,8 +1,5 @@
 #include "sweep_utils.h"
 
-#include "../ifcparse/IfcLogger.h"
-#include "../ifcgeom_schema_agnostic/Kernel.h"
-
 #include <gp_Ax2.hxx>
 #include <Geom_Line.hxx>
 #include <Geom_Circle.hxx>
@@ -24,6 +21,10 @@
 #include <BRepBuilderAPI_MakeFace.hxx>
 #include <BRepOffsetAPI_MakePipeShell.hxx>
 
+#include "../ifcparse/IfcLogger.h"
+#include "../ifcgeom_schema_agnostic/Kernel.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
+
 bool IfcGeom::util::wire_is_c1_continuous(const TopoDS_Wire & w, double tol) {
 	// NB Note that c0 continuity is NOT checked!
 
@@ -193,7 +194,7 @@ void IfcGeom::util::sort_edges(const TopoDS_Wire & wire, std::vector<TopoDS_Edge
 
 	std::set<TopoDS_TShape*> seen;
 
-	auto num_edges = IfcGeom::Kernel::count(wire, TopAbs_EDGE);
+	auto num_edges = count(wire, TopAbs_EDGE);
 
 	TopoDS_Vertex v0, v1;
 	// @todo this creates the ancestor map twice

src/ifcgeom_schema_agnostic/wire_utils.cpp

@@ -2,6 +2,7 @@
 
 #include "../ifcparse/IfcLogger.h"
 #include "../ifcgeom_schema_agnostic/Kernel.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 #include "../ifcgeom_schema_agnostic/IfcGeomTree.h"
 
 #include <TopExp.hxx>
@@ -84,8 +85,8 @@ bool IfcGeom::util::approximate_plane_through_wire(const TopoDS_Wire& wire, gp_P
 
 	exp.Init(wire);
 	for (; exp.More(); exp.Next()) {
-		const TopoDS_Vertex& v = exp.CurrentVertex();
-		current = BRep_Tool::Pnt(v);
+		const TopoDS_Vertex& vrt = exp.CurrentVertex();
+		current = BRep_Tool::Pnt(vrt);
 		if (plane.SquareDistance(current) > eps2) {
 			return false;
 		}
@@ -355,7 +356,7 @@ bool IfcGeom::util::wire_intersections(const TopoDS_Wire& wire, TopTools_ListOfS
 		return false;
 	}
 
-	int n = IfcGeom::Kernel::count(wire, TopAbs_EDGE);
+	int n = util::count(wire, TopAbs_EDGE);
 	if (n < 3) {
 		wires.Append(wire);
 		return false;

src/serializers/SvgSerializer.cpp

@@ -78,6 +78,7 @@
 #include <Extrema_ExtPElS.hxx>
 
 #include "../ifcparse/IfcGlobalId.h"
+#include "../ifcgeom_schema_agnostic/base_utils.h"
 
 #include <boost/format.hpp>
 #include <boost/tokenizer.hpp>
@@ -439,7 +440,7 @@ namespace {
 			return boost::none;
 		}
 
-		if (IfcGeom::Kernel::count(shell, TopAbs_FACE) != 6) {
+		if (IfcGeom::util::count(shell, TopAbs_FACE) != 6) {
 			return boost::none;
 		}
 
@@ -1442,7 +1443,7 @@ void SvgSerializer::write(const geometry_data& data) {
 
 				}
 				
-				if (file && data.product->declaration().is("IfcBuildingStorey") && storey_height_display_ != SH_NONE && wires->Length() == 1 && IfcGeom::Kernel::count(wire, TopAbs_EDGE) == 1) {
+				if (file && data.product->declaration().is("IfcBuildingStorey") && storey_height_display_ != SH_NONE && wires->Length() == 1 && IfcGeom::util::count(wire, TopAbs_EDGE) == 1) {
 					
 					std::string elev_str;