Merge branch 'dev' into lexer

.gitignore

@@ -41,3 +41,6 @@ venv.bak/
 
 # Profiling debris.
 prof/
+
+# MacOS cruft
+.DS_Store

CONTRIBUTING.md

@@ -3,8 +3,8 @@ tests, ensure they build and pass, and ensure that `pylint` and `mypy`
 are happy with your code.
 
 - Install `pip` following their [documentation](https://pip.pypa.io/en/stable/installation/).
-- Install development dependencies: `pip install -e .[development]`
+- Install development dependencies: `pip install -e ".[development]"`
-- Install docs dependencies: `pip install -e .[docs]`
+- Install docs dependencies: `pip install -e ".[docs]"`
 - Install `build123d` in editable mode from current dir:  `pip install -e .`
 - Run tests with: `python -m pytest -n auto`
 - Build docs with: `cd docs && make html`

NOTICE

@@ -0,0 +1,15 @@
+build123d
+Copyright (c) 2022–2025 The build123d Contributors
+
+Licensed under the Apache License, Version 2.0 (the "License");
+you may not use this file except in compliance with the License.
+You may obtain a copy of the License at:
+
+    http://www.apache.org/licenses/LICENSE-2.0
+
+-------------------------------------------------------------------------------
+
+This project was originally derived from portions of the CadQuery codebase
+(https://github.com/CadQuery/cadquery) but has since been extensively
+refactored and restructured into an independent system.
+CadQuery is licensed under the Apache License, Version 2.0.

README.md

@@ -1,5 +1,5 @@
 <p align="center">
-    <img alt="build123d logo" src="docs/assets/build123d_logo/logo-banner.svg">
+    <img alt="build123d logo" src="https://github.com/gumyr/build123d/raw/dev/docs/assets/build123d_logo/logo-banner.svg">
 </p>
 
 [![Documentation Status](https://readthedocs.org/projects/build123d/badge/?version=latest)](https://build123d.readthedocs.io/en/latest/?badge=latest)
@@ -19,9 +19,17 @@
 [![DOI](https://zenodo.org/badge/510925389.svg)](https://doi.org/10.5281/zenodo.14872322)
 
 
-Build123d is a python-based, parametric, [boundary representation (BREP)][BREP] modeling framework for 2D and 3D CAD. It's built on the [Open Cascade] geometric kernel and allows for the creation of complex models using a simple and intuitive python syntax. Build123d can be used to create models for 3D printing, CNC machining, laser cutting, and other manufacturing processes.  Models can be exported to a wide variety of popular CAD tools such as [FreeCAD] and SolidWorks.
+Build123d is a Python-based, parametric [boundary representation (BREP)][BREP] modeling framework for 2D and 3D CAD. Built on the [Open Cascade] geometric kernel, it provides a clean, fully Pythonic interface for creating precise models suitable for 3D printing, CNC machining, laser cutting, and other manufacturing processes. Models can be exported to popular CAD tools such as [FreeCAD] and SolidWorks.
 
-Build123d could be considered as an evolution of [CadQuery] where the somewhat restrictive Fluent API (method chaining) is replaced with stateful context managers - e.g. `with` blocks - thus enabling the full python toolbox: for loops, references to objects, object sorting and filtering, etc.
+Designed for modern, maintainable CAD-as-code, build123d combines clear architecture with expressive, algebraic modeling. It offers:
+- Minimal or no internal state depending on mode,
+- Explicit 1D, 2D, and 3D geometry classes with well-defined operations,
+- Extensibility through subclassing and functional composition—no monkey patching,
+- Standards-compliant code (PEP 8, mypy, pylint) with rich pylance type hints,
+- Deep Python integration—selectors as lists, locations as iterables, and natural conversions (Solid(shell), tuple(Vector)),
+- Operator-driven modeling (obj += sub_obj, Plane.XZ * Pos(X=5) * Rectangle(1, 1)) for algebraic, readable, and composable design logic.
+
+The result is a framework that feels native to Python while providing the full power of OpenCascade geometry underneath.
 
 The documentation for **build123d** can be found at [readthedocs](https://build123d.readthedocs.io/en/latest/index.html).
 
@@ -62,6 +70,10 @@ python3 -m pip install -e .
 
 Further installation instructions are available (e.g. Poetry) see the [installation section on readthedocs](https://build123d.readthedocs.io/en/latest/installation.html).
 
+Attribution:
+
+Build123d was originally derived from portions of the [CadQuery] codebase but has since been extensively refactored and restructured into an independent system.
+
 [BREP]: https://en.wikipedia.org/wiki/Boundary_representation
 [CadQuery]: https://cadquery.readthedocs.io/en/latest/index.html
 [FreeCAD]: https://www.freecad.org/

docs/index.rst

@@ -29,58 +29,35 @@
   :align: center
   :alt: build123d logo
 
-Build123d is a python-based, parametric, boundary representation (BREP) modeling
-framework for 2D and 3D CAD. It's built on the Open Cascade geometric kernel and
-allows for the creation of complex models using a simple and intuitive python
-syntax. Build123d can be used to create models for 3D printing, CNC machining,
-laser cutting, and other manufacturing processes.  Models can be exported to a
-wide variety of popular CAD tools such as FreeCAD and SolidWorks.
-
-Build123d could be considered as an evolution of
-`CadQuery <https://cadquery.readthedocs.io/en/latest/index.html>`_ where the
-somewhat restrictive Fluent API (method chaining) is replaced with stateful
-context managers - i.e. `with` blocks - thus enabling the full python toolbox:
-for loops, references to objects, object sorting and filtering, etc.
-
-Note that this documentation is available in
-`pdf <https://build123d.readthedocs.io/_/downloads/en/latest/pdf/>`_ and
-`epub <https://build123d.readthedocs.io/_/downloads/en/latest/epub/>`_ formats
-for reference while offline.
-
 ########
-Overview
+About
 ########
 
-build123d uses the standard python context manager - i.e. the ``with`` statement often used when
+Build123d is a Python-based, parametric (BREP) modeling framework for 2D and 3D CAD. 
-working with files - as a builder of the object under construction. Once the object is complete
+Built on the Open Cascade geometric kernel, it provides a clean, fully Pythonic interface 
-it can be extracted from the builders and used in other ways: for example exported as a STEP
+for creating precise models suitable for 3D printing, CNC machining, laser cutting, and 
-file or used in an Assembly.  There are three builders available:
+other manufacturing processes. Models can be exported to popular CAD tools such as FreeCAD 
+and SolidWorks.
 
-* **BuildLine**: a builder of one dimensional objects - those with the property
+Designed for modern, maintainable CAD-as-code, build123d combines clear architecture with 
-  of length but not of area or volume - typically used
+expressive, algebraic modeling. It offers:
-  to create complex lines used in sketches or paths.
-* **BuildSketch**: a builder of planar two dimensional objects - those with the property
-  of area but not of volume - typically used to create 2D drawings that are extruded into 3D parts.
-* **BuildPart**: a builder of three dimensional objects - those with the property of volume -
-  used to create individual parts.
 
-The three builders work together in a hierarchy as follows:
+* Minimal or no internal state depending on mode
+* Explicit 1D, 2D, and 3D geometry classes with well-defined operations
+* Extensibility through subclassing and functional composition—no monkey patching
+* Standards-compliant code (PEP 8, mypy, pylint) with rich pylance type hints
+* Deep Python integration—selectors as lists, locations as iterables, and natural 
+  conversions (``Solid(shell)``, ``tuple(Vector)``)
+* Operator-driven modeling (``obj += sub_obj``, ``Plane.XZ * Pos(X=5) * Rectangle(1, 1)``) 
+  for algebraic, readable, and composable design logic
 
 .. code-block:: build123d
 
-    with BuildPart() as my_part:
-        ...
-        with BuildSketch() as my_sketch:
-            ...
-            with BuildLine() as my_line:
-                ...
-            ...
-        ...
 
-where ``my_line`` will be added to ``my_sketch`` once the line is complete and ``my_sketch`` will be
+With build123d, intricate parametric models can be created in just a few lines of readable 
-added to ``my_part`` once the sketch is complete.
+Python code—as demonstrated by the tea cup example below.
 
-As an example, consider the design of a tea cup:
+.. dropdown:: Teacup Example
 
 .. literalinclude:: ../examples/tea_cup.py
     :language: build123d
@@ -92,6 +69,14 @@ As an example, consider the design of a tea cup:
     <script type="module" src="https://unpkg.com/@google/model-viewer/dist/model-viewer.min.js"></script>
     <model-viewer poster="_images/tea_cup.png" src="_static/tea_cup.glb" alt="A tea cup modelled in build123d" auto-rotate camera-controls style="width: 100%; height: 50vh;"></model-viewer>
 
+.. note::
+
+
+  This documentation is available in
+  `pdf <https://build123d.readthedocs.io/_/downloads/en/latest/pdf/>`_ and
+  `epub <https://build123d.readthedocs.io/_/downloads/en/latest/epub/>`_ formats
+  for reference while offline.
+
 .. note::
 
   There is a `Discord <https://discord.com/invite/Bj9AQPsCfx>`_ server (shared with CadQuery) where

examples/extrude.py

@@ -68,8 +68,7 @@ with BuildPart() as ex26:
     with BuildSketch() as ex26_sk:
         with Locations((0, rev)):
             Circle(rad)
-    revolve(axis=Axis.X, revolution_arc=90)
+    revolve(axis=Axis.X, revolution_arc=180)
-    mirror(about=Plane.XZ)
     with BuildSketch() as ex26_sk2:
         Rectangle(rad, rev)
     ex26_target = ex26.part

examples/extrude_algebra.py

@@ -26,8 +26,8 @@ rad, rev = 3, 25
 
 # Extrude last
 circle = Pos(0, rev) * Circle(rad)
-ex26_target = revolve(circle, Axis.X, revolution_arc=90)
+ex26_target = revolve(circle, Axis.X, revolution_arc=180)
-ex26_target = ex26_target + mirror(ex26_target, Plane.XZ)
+ex26_target = ex26_target
 
 rect = Rectangle(rad, rev)
 

pyproject.toml

@@ -68,7 +68,7 @@ development = [
     "black",
     "mypy",
     "pylint",
-    "pytest",
+    "pytest==8.4.2", # TODO: unpin on resolution of pytest-dev/pytest-xdist/issues/1273
     "pytest-benchmark",
     "pytest-cov",
     "pytest-xdist",

src/build123d/build_common.py

@@ -466,7 +466,7 @@ class Builder(ABC, Generic[ShapeT]):
                 elif mode == Mode.INTERSECT:
                     if self._obj is None:
                         raise RuntimeError("Nothing to intersect with")
-                    combined = self._obj.intersect(*typed[self._shape])
+                    combined = self._obj.intersect(Compound(typed[self._shape]))
                 elif mode == Mode.REPLACE:
                     combined = self._sub_class(list(typed[self._shape]))
 

src/build123d/drafting.py

@@ -453,7 +453,7 @@ class DimensionLine(BaseSketchObject):
             if self_intersection is None:
                 self_intersection_area = 0.0
             else:
-                self_intersection_area = self_intersection.area
+                self_intersection_area = sum(f.area for f in self_intersection.faces())
             d_line += placed_label
             bbox_size = d_line.bounding_box().diagonal
 
@@ -467,7 +467,7 @@ class DimensionLine(BaseSketchObject):
                 if line_intersection is None:
                     common_area = 0.0
                 else:
-                    common_area = line_intersection.area
+                    common_area = sum(f.area for f in line_intersection.faces())
             common_area += self_intersection_area
             score = (d_line.area - 10 * common_area) / bbox_size
             d_lines[d_line] = score

src/build123d/exporters.py

@@ -758,7 +758,7 @@ class ExportDXF(Export2D):
         )
 
         # need to apply the transform on the geometry level
-        if edge.wrapped is None or edge.location is None:
+        if not edge or edge.location is None:
             raise ValueError(f"Edge is empty {edge}.")
         t = edge.location.wrapped.Transformation()
         spline.Transform(t)
@@ -1345,7 +1345,7 @@ class ExportSVG(Export2D):
         u2 = adaptor.LastParameter()
 
         # Apply the shape location to the geometry.
-        if edge.wrapped is None or edge.location is None:
+        if not edge or edge.location is None:
             raise ValueError(f"Edge is empty {edge}.")
         t = edge.location.wrapped.Transformation()
         spline.Transform(t)
@@ -1411,7 +1411,7 @@ class ExportSVG(Export2D):
     }
 
     def _edge_segments(self, edge: Edge, reverse: bool) -> list[PathSegment]:
-        if edge.wrapped is None:
+        if not edge:
             raise ValueError(f"Edge is empty {edge}.")
         edge_reversed = edge.wrapped.Orientation() == TopAbs_Orientation.TopAbs_REVERSED
         geom_type = edge.geom_type

src/build123d/mesher.py

@@ -295,7 +295,7 @@ class Mesher:
                 ocp_mesh_vertices.append(pnt)
 
             # Store the triangles from the triangulated faces
-            if facet.wrapped is None:
+            if not facet:
                 continue
             facet_reversed = facet.wrapped.Orientation() == ta.TopAbs_REVERSED
             order = [1, 3, 2] if facet_reversed else [1, 2, 3]

src/build123d/objects_curve.py

@@ -29,6 +29,7 @@ license:
 from __future__ import annotations
 
 import copy as copy_module
+import warnings
 import numpy as np
 import sympy  # type: ignore
 from collections.abc import Iterable
@@ -792,7 +793,7 @@ class FilletPolyline(BaseLineObject):
 
     Args:
         pts (VectorLike | Iterable[VectorLike]): sequence of two or more points
-        radius (float): fillet radius
+        radius (float | Iterable[float]): radius to fillet at each vertex or a single value for all vertices
         close (bool, optional): close end points with extra Edge and corner fillets.
             Defaults to False
         mode (Mode, optional): combination mode. Defaults to Mode.ADD
@@ -807,7 +808,7 @@ class FilletPolyline(BaseLineObject):
     def __init__(
         self,
         *pts: VectorLike | Iterable[VectorLike],
-        radius: float,
+        radius: float | Iterable[float],
         close: bool = False,
         mode: Mode = Mode.ADD,
     ):
@@ -818,8 +819,18 @@ class FilletPolyline(BaseLineObject):
 
         if len(points) < 2:
             raise ValueError("FilletPolyline requires two or more pts")
-        if radius <= 0:
+
-            raise ValueError("radius must be positive")
+        if isinstance(radius, (int, float)):
+            radius_list = [radius] * len(points)  # Single radius for all points
+        else:
+            radius_list = list(radius)
+            if len(radius_list) != len(points) - int(not close) * 2:
+                raise ValueError(
+                    f"radius list length ({len(radius_list)}) must match angle count ({ len(points) - int(not close) * 2})"
+                )
+        for r in radius_list:
+            if r <= 0:
+                raise ValueError(f"radius {r} must be positive")
 
         lines_pts = WorkplaneList.localize(*points)
 
@@ -851,12 +862,14 @@ class FilletPolyline(BaseLineObject):
 
         # For each corner vertex create a new fillet Edge
         fillets = []
-        for vertex, edges in vertex_to_edges.items():
+        for i, (vertex, edges) in enumerate(vertex_to_edges.items()):
             if len(edges) != 2:
                 continue
             other_vertices = {ve for e in edges for ve in e.vertices() if ve != vertex}
             third_edge = Edge.make_line(*[v for v in other_vertices])
-            fillet_face = Face(Wire(edges + [third_edge])).fillet_2d(radius, [vertex])
+            fillet_face = Face(Wire(edges + [third_edge])).fillet_2d(
+                radius_list[i - int(not close)], [vertex]
+            )
             fillets.append(fillet_face.edges().filter_by(GeomType.CIRCLE)[0])
 
         # Create the Edges that join the fillets
@@ -1362,6 +1375,12 @@ class PointArcTangentLine(BaseEdgeObject):
         mode (Mode, optional): combination mode. Defaults to Mode.ADD
     """
 
+    warnings.warn(
+        "The 'PointArcTangentLine' object is deprecated and will be removed in a future version.",
+        DeprecationWarning,
+        stacklevel=2,
+    )
+
     _applies_to = [BuildLine._tag]
 
     def __init__(
@@ -1441,6 +1460,12 @@ class PointArcTangentArc(BaseEdgeObject):
         RuntimeError: No tangent arc found
     """
 
+    warnings.warn(
+        "The 'PointArcTangentArc' object is deprecated and will be removed in a future version.",
+        DeprecationWarning,
+        stacklevel=2,
+    )
+
     _applies_to = [BuildLine._tag]
 
     def __init__(
@@ -1585,6 +1610,12 @@ class ArcArcTangentLine(BaseEdgeObject):
         mode (Mode, optional): combination mode. Defaults to Mode.ADD
     """
 
+    warnings.warn(
+        "The 'ArcArcTangentLine' object is deprecated and will be removed in a future version.",
+        DeprecationWarning,
+        stacklevel=2,
+    )
+
     _applies_to = [BuildLine._tag]
 
     def __init__(
@@ -1685,6 +1716,12 @@ class ArcArcTangentArc(BaseEdgeObject):
         mode (Mode, optional): combination mode. Defaults to Mode.ADD
     """
 
+    warnings.warn(
+        "The 'ArcArcTangentArc' object is deprecated and will be removed in a future version.",
+        DeprecationWarning,
+        stacklevel=2,
+    )
+
     _applies_to = [BuildLine._tag]
 
     def __init__(

src/build123d/operations_generic.py

@@ -365,7 +365,7 @@ def chamfer(
 
     if target._dim == 1:
         if isinstance(target, BaseLineObject):
-            if target.wrapped is None:
+            if not target:
                 target = Wire([])  # empty wire
             else:
                 target = Wire(target.wrapped)
@@ -465,7 +465,7 @@ def fillet(
 
     if target._dim == 1:
         if isinstance(target, BaseLineObject):
-            if target.wrapped is None:
+            if not target:
                 target = Wire([])  # empty wire
             else:
                 target = Wire(target.wrapped)

src/build123d/topology/composite.py

@@ -58,13 +58,12 @@ import copy
 import os
 import sys
 import warnings
-from itertools import combinations
-from typing import Type, Union
-
 from collections.abc import Iterable, Iterator, Sequence
+from itertools import combinations
+from typing_extensions import Self
 
 import OCP.TopAbs as ta
-from OCP.BRepAlgoAPI import BRepAlgoAPI_Fuse
+from OCP.BRepAlgoAPI import BRepAlgoAPI_Common, BRepAlgoAPI_Fuse, BRepAlgoAPI_Section
 from OCP.Font import (
     Font_FA_Bold,
     Font_FA_BoldItalic,
@@ -107,7 +106,6 @@ from build123d.geometry import (
     VectorLike,
     logger,
 )
-from typing_extensions import Self
 
 from .one_d import Edge, Wire, Mixin1D
 from .shape_core import (
@@ -130,7 +128,7 @@ from .utils import (
 from .zero_d import Vertex
 
 
-class Compound(Mixin3D, Shape[TopoDS_Compound]):
+class Compound(Mixin3D[TopoDS_Compound]):
     """A Compound in build123d is a topological entity representing a collection of
     geometric shapes grouped together within a single structure. It serves as a
     container for organizing diverse shapes like edges, faces, or solids. This
@@ -455,7 +453,7 @@ class Compound(Mixin3D, Shape[TopoDS_Compound]):
         will be a Wire, otherwise a Shape.
         """
         if self._dim == 1:
-            curve = Curve() if self.wrapped is None else Curve(self.wrapped)
+            curve = Curve() if self._wrapped is None else Curve(self.wrapped)
             sum1d: Edge | Wire | ShapeList[Edge] = curve + other
             if isinstance(sum1d, ShapeList):
                 result1d: Curve | Wire = Curve(sum1d)
@@ -517,7 +515,7 @@ class Compound(Mixin3D, Shape[TopoDS_Compound]):
         Check if empty.
         """
 
-        return TopoDS_Iterator(self.wrapped).More()
+        return self._wrapped is not None and TopoDS_Iterator(self.wrapped).More()
 
     def __iter__(self) -> Iterator[Shape]:
         """
@@ -534,7 +532,7 @@ class Compound(Mixin3D, Shape[TopoDS_Compound]):
     def __len__(self) -> int:
         """Return the number of subshapes"""
         count = 0
-        if self.wrapped is not None:
+        if self._wrapped is not None:
             for _ in self:
                 count += 1
         return count
@@ -602,7 +600,7 @@ class Compound(Mixin3D, Shape[TopoDS_Compound]):
 
     def compounds(self) -> ShapeList[Compound]:
         """compounds - all the compounds in this Shape"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             return ShapeList()
         if isinstance(self.wrapped, TopoDS_Compound):
             # pylint: disable=not-an-iterable
@@ -651,11 +649,7 @@ class Compound(Mixin3D, Shape[TopoDS_Compound]):
                     children[child_index_pair[1]]
                 )
                 if obj_intersection is not None:
-                    common_volume = (
+                    common_volume = sum(s.volume for s in obj_intersection.solids())
-                        0.0
-                        if isinstance(obj_intersection, list)
-                        else obj_intersection.volume
-                    )
                     if common_volume > tolerance:
                         return (
                             True,
@@ -711,6 +705,148 @@ class Compound(Mixin3D, Shape[TopoDS_Compound]):
 
         return results
 
+    def intersect(
+        self, *to_intersect: Shape | Vector | Location | Axis | Plane
+    ) -> None | ShapeList[Vertex | Edge | Face | Solid]:
+        """Intersect Compound with Shape or geometry object
+
+        Args:
+            to_intersect (Shape | Vector | Location | Axis | Plane): objects to intersect
+
+        Returns:
+            ShapeList[Vertex | Edge | Face | Solid] | None: ShapeList of vertices, edges,
+                faces, and/or solids.
+        """
+
+        def to_vector(objs: Iterable) -> ShapeList:
+            return ShapeList([Vector(v) if isinstance(v, Vertex) else v for v in objs])
+
+        def to_vertex(objs: Iterable) -> ShapeList:
+            return ShapeList([Vertex(v) if isinstance(v, Vector) else v for v in objs])
+
+        def bool_op(
+            args: Sequence,
+            tools: Sequence,
+            operation: BRepAlgoAPI_Common | BRepAlgoAPI_Section,
+        ) -> ShapeList:
+            # Wrap Shape._bool_op for corrected output
+            intersections: Shape | ShapeList = Shape()._bool_op(args, tools, operation)
+            if isinstance(intersections, ShapeList):
+                return intersections
+            if isinstance(intersections, Shape) and not intersections.is_null:
+                return ShapeList([intersections])
+            return ShapeList()
+
+        def expand_compound(compound: Compound) -> ShapeList:
+            shapes = ShapeList(compound.children)
+            for shape_type in [Vertex, Edge, Wire, Face, Shell, Solid]:
+                shapes.extend(compound.get_type(shape_type))
+            return shapes
+
+        def filter_shapes_by_order(shapes: ShapeList, orders: list) -> ShapeList:
+            # Remove lower order shapes from list which *appear* to be part of
+            # a higher order shape using a lazy distance check
+            # (sufficient for vertices, may be an issue for higher orders)
+            order_groups = []
+            for order in orders:
+                order_groups.append(
+                    ShapeList([s for s in shapes if isinstance(s, order)])
+                )
+
+            filtered_shapes = order_groups[-1]
+            for i in range(len(order_groups) - 1):
+                los = order_groups[i]
+                his: list = sum(order_groups[i + 1 :], [])
+                filtered_shapes.extend(
+                    ShapeList(
+                        lo
+                        for lo in los
+                        if all(lo.distance_to(hi) > TOLERANCE for hi in his)
+                    )
+                )
+
+            return filtered_shapes
+
+        common_set: ShapeList[Shape] = expand_compound(self)
+        target: ShapeList | Shape
+        for other in to_intersect:
+            # Conform target type
+            match other:
+                case Axis():
+                    # BRepAlgoAPI_Section seems happier if Edge isnt infinite
+                    bbox = self.bounding_box()
+                    dist = self.distance_to(other.position)
+                    dist = dist if dist >= 1 else 1
+                    target = Edge.make_line(
+                        other.position - other.direction * bbox.diagonal * dist,
+                        other.position + other.direction * bbox.diagonal * dist,
+                    )
+                case Plane():
+                    target = Face(other)
+                case Vector():
+                    target = Vertex(other)
+                case Location():
+                    target = Vertex(other.position)
+                case Compound():
+                    target = expand_compound(other)
+                case _ if issubclass(type(other), Shape):
+                    target = other
+                case _:
+                    raise ValueError(f"Unsupported type to_intersect: {type(other)}")
+
+            # Find common matches
+            common: list[Vertex | Edge | Wire | Face | Shell | Solid] = []
+            result: ShapeList
+            for obj in common_set:
+                if isinstance(target, Shape):
+                    target = ShapeList([target])
+                result = ShapeList()
+                for t in target:
+                    operation = BRepAlgoAPI_Section()
+                    result.extend(bool_op((obj,), (t,), operation))
+                    if (
+                        not isinstance(obj, Edge | Wire)
+                        and not isinstance(t, Edge | Wire)
+                    ) or (
+                        isinstance(obj, Solid | Compound)
+                        or isinstance(t, Solid | Compound)
+                    ):
+                        # Face + Edge combinations may produce an intersection
+                        # with Common but always with Section.
+                        # No easy way to deduplicate
+                        # Many Solid + Edge combinations need Common
+                        operation = BRepAlgoAPI_Common()
+                        result.extend(bool_op((obj,), (t,), operation))
+
+                if result:
+                    common.extend(result)
+
+            expanded: ShapeList = ShapeList()
+            if common:
+                for shape in common:
+                    if isinstance(shape, Compound):
+                        expanded.extend(expand_compound(shape))
+                    else:
+                        expanded.append(shape)
+
+            if expanded:
+                common_set = ShapeList()
+                for shape in expanded:
+                    if isinstance(shape, Wire):
+                        common_set.extend(shape.edges())
+                    elif isinstance(shape, Shell):
+                        common_set.extend(shape.faces())
+                    else:
+                        common_set.append(shape)
+                common_set = to_vertex(set(to_vector(common_set)))
+                common_set = filter_shapes_by_order(
+                    common_set, [Vertex, Edge, Face, Solid]
+                )
+            else:
+                return None
+
+        return ShapeList(common_set)
+
     def unwrap(self, fully: bool = True) -> Self | Shape:
         """Strip unnecessary Compound wrappers
 

src/build123d/topology/constrained_lines.py

@@ -174,7 +174,7 @@ def _as_gcc_arg(obj: Edge | Vector, constaint: Tangency) -> tuple[
     - Edge -> (QualifiedCurve, h2d, first, last, True)
     - Vector -> (CartesianPoint, None, None, None, False)
     """
-    if obj.wrapped is None:
+    if not obj:
         raise TypeError("Can't create a qualified curve from empty edge")
 
     if isinstance(obj.wrapped, TopoDS_Edge):

src/build123d/topology/one_d.py

@@ -52,12 +52,11 @@ license:
 from __future__ import annotations
 
 import copy
-import numpy as np
 import warnings
-from collections.abc import Iterable
+from collections.abc import Iterable, Sequence
 from itertools import combinations
 from math import atan2, ceil, copysign, cos, floor, inf, isclose, pi, radians
-from typing import TYPE_CHECKING, Literal, TypeAlias, overload
+from typing import TYPE_CHECKING, Literal, overload
 from typing import cast as tcast
 
 import numpy as np
@@ -217,6 +216,7 @@ from build123d.geometry import (
 )
 
 from .shape_core import (
+    TOPODS,
     Shape,
     ShapeList,
     SkipClean,
@@ -250,7 +250,7 @@ if TYPE_CHECKING:  # pragma: no cover
     from .two_d import Face, Shell  # pylint: disable=R0801
 
 
-class Mixin1D(Shape):
+class Mixin1D(Shape[TOPODS]):
     """Methods to add to the Edge and Wire classes"""
 
     # ---- Properties ----
@@ -263,14 +263,14 @@ class Mixin1D(Shape):
     @property
     def is_closed(self) -> bool:
         """Are the start and end points equal?"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't determine if empty Edge or Wire is closed")
         return BRep_Tool.IsClosed_s(self.wrapped)
 
     @property
     def is_forward(self) -> bool:
         """Does the Edge/Wire loop forward or reverse"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't determine direction of empty Edge or Wire")
         return self.wrapped.Orientation() == TopAbs_Orientation.TopAbs_FORWARD
 
@@ -388,8 +388,7 @@ class Mixin1D(Shape):
                 shape
                 # for o in (other if isinstance(other, (list, tuple)) else [other])
                 for o in ([other] if isinstance(other, Shape) else other)
-                if o is not None
+                for shape in get_top_level_topods_shapes(o.wrapped if o else None)
-                for shape in get_top_level_topods_shapes(o.wrapped)
             ]
         # If there is nothing to add return the original object
         if not topods_summands:
@@ -404,7 +403,7 @@ class Mixin1D(Shape):
         )
         summand_edges = [e for summand in summands for e in summand.edges()]
 
-        if self.wrapped is None:  # an empty object
+        if self._wrapped is None:  # an empty object
             if len(summands) == 1:
                 sum_shape: Edge | Wire | ShapeList[Edge] = summands[0]
             else:
@@ -452,7 +451,7 @@ class Mixin1D(Shape):
         Returns:
             Vector: center
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find center of empty edge/wire")
 
         if center_of == CenterOf.GEOMETRY:
@@ -578,7 +577,7 @@ class Mixin1D(Shape):
             >>> show(my_wire, Curve(comb))
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't create curvature_comb for empty curve")
         pln = self.common_plane()
         if pln is None or not isclose(abs(pln.z_dir.Z), 1.0, abs_tol=TOLERANCE):
@@ -729,122 +728,103 @@ class Mixin1D(Shape):
         def to_vertex(objs: Iterable) -> ShapeList:
             return ShapeList([Vertex(v) if isinstance(v, Vector) else v for v in objs])
 
-        common_set: ShapeList[Vertex | Edge] = ShapeList(self.edges())
+        def bool_op(
-        target: ShapeList | Shape | Plane
+            args: Sequence,
+            tools: Sequence,
+            operation: BRepAlgoAPI_Section | BRepAlgoAPI_Common,
+        ) -> ShapeList:
+            # Wrap Shape._bool_op for corrected output
+            intersections: Shape | ShapeList = Shape()._bool_op(args, tools, operation)
+            if isinstance(intersections, ShapeList):
+                return intersections or ShapeList()
+            if isinstance(intersections, Shape) and not intersections.is_null:
+                return ShapeList([intersections])
+            return ShapeList()
+
+        def filter_shapes_by_order(shapes: ShapeList, orders: list) -> ShapeList:
+            # Remove lower order shapes from list which *appear* to be part of
+            # a higher order shape using a lazy distance check
+            # (sufficient for vertices, may be an issue for higher orders)
+            order_groups = []
+            for order in orders:
+                order_groups.append(
+                    ShapeList([s for s in shapes if isinstance(s, order)])
+                )
+
+            filtered_shapes = order_groups[-1]
+            for i in range(len(order_groups) - 1):
+                los = order_groups[i]
+                his: list = sum(order_groups[i + 1 :], [])
+                filtered_shapes.extend(
+                    ShapeList(
+                        lo
+                        for lo in los
+                        if all(lo.distance_to(hi) > TOLERANCE for hi in his)
+                    )
+                )
+
+            return filtered_shapes
+
+        common_set: ShapeList[Vertex | Edge | Wire] = ShapeList([self])
+        target: Shape | Plane
         for other in to_intersect:
             # Conform target type
-            # Vertices need to be Vector for set()
             match other:
                 case Axis():
-                    target = ShapeList([Edge(other)])
+                    # BRepAlgoAPI_Section seems happier if Edge isnt infinite
+                    bbox = self.bounding_box()
+                    dist = self.distance_to(other.position)
+                    dist = dist if dist >= 1 else 1
+                    target = Edge.make_line(
+                        other.position - other.direction * bbox.diagonal * dist,
+                        other.position + other.direction * bbox.diagonal * dist,
+                    )
                 case Plane():
                     target = other
                 case Vector():
                     target = Vertex(other)
                 case Location():
                     target = Vertex(other.position)
-                case Edge():
-                    target = ShapeList([other])
-                case Wire():
-                    target = ShapeList(other.edges())
                 case _ if issubclass(type(other), Shape):
                     target = other
                 case _:
                     raise ValueError(f"Unsupported type to_intersect: {type(other)}")
 
             # Find common matches
-            common: list[Vector | Edge] = []
+            common: list[Vertex | Edge | Wire] = []
-            result: ShapeList | Shape | None
+            result: ShapeList | None
             for obj in common_set:
                 match (obj, target):
-                    case obj, Shape() as target:
+                    case (_, Plane()):
-                        # Find Shape with Edge/Wire
+                        target = Shape(BRepBuilderAPI_MakeFace(other.wrapped).Face())
-                        if isinstance(target, Vertex):
+                        operation = BRepAlgoAPI_Section()
-                            result = Shape.intersect(obj, target)
+                        result = bool_op((obj,), (target,), operation)
-                        else:
+                        operation = BRepAlgoAPI_Common()
-                            result = target.intersect(obj)
+                        result.extend(bool_op((obj,), (target,), operation))
 
-                        if result:
+                    case (_, Vertex() | Edge() | Wire()):
-                            if not isinstance(result, list):
+                        operation = BRepAlgoAPI_Section()
-                                result = ShapeList([result])
+                        section = bool_op((obj,), (target,), operation)
-                            common.extend(to_vector(result))
+                        result = section
+                        if not section:
+                            operation = BRepAlgoAPI_Common()
+                            result.extend(bool_op((obj,), (target,), operation))
 
-                    case Vertex() as obj, target:
+                    case _ if issubclass(type(target), Shape):
-                        if not isinstance(target, ShapeList):
+                        result = target.intersect(obj)
-                            target = ShapeList([target])
 
-                        for tar in target:
+                if result:
-                            if isinstance(tar, Edge):
+                    common.extend(result)
-                                result = Shape.intersect(obj, tar)
-                            else:
-                                result = obj.intersect(tar)
-
-                            if result:
-                                if not isinstance(result, list):
-                                    result = ShapeList([result])
-                                common.extend(to_vector(result))
-
-                    case Edge() as obj, ShapeList() as targets:
-                        # Find any edge / edge intersection points
-                        for tar in targets:
-                            # Find crossing points
-                            try:
-                                intersection_points = obj.find_intersection_points(tar)
-                                common.extend(intersection_points)
-                            except ValueError:
-                                pass
-
-                            # Find common end points
-                            obj_end_points = set(Vector(v) for v in obj.vertices())
-                            tar_end_points = set(Vector(v) for v in tar.vertices())
-                            points = set.intersection(obj_end_points, tar_end_points)
-                            common.extend(points)
-
-                        # Find Edge/Edge overlaps
-                        result = obj._bool_op(
-                            (obj,), targets, BRepAlgoAPI_Common()
-                        ).edges()
-                        common.extend(result if isinstance(result, list) else [result])
-
-                    case Edge() as obj, Plane() as plane:
-                        # Find any edge / plane intersection points & edges
-                        # Find point intersections
-                        if obj.wrapped is None:
-                            continue
-                        geom_line = BRep_Tool.Curve_s(
-                            obj.wrapped, obj.param_at(0), obj.param_at(1)
-                        )
-                        geom_plane = Geom_Plane(plane.local_coord_system)
-                        intersection_calculator = GeomAPI_IntCS(geom_line, geom_plane)
-                        plane_intersection_points: list[Vector] = []
-                        if intersection_calculator.IsDone():
-                            plane_intersection_points = [
-                                Vector(intersection_calculator.Point(i + 1))
-                                for i in range(intersection_calculator.NbPoints())
-                            ]
-                        common.extend(plane_intersection_points)
-
-                        # Find edge intersections
-                        if all(
-                            plane.contains(v)
-                            for v in obj.positions(i / 7 for i in range(8))
-                        ):  # is a 2D edge
-                            common.append(obj)
 
             if common:
-                common_set = to_vertex(set(common))
+                common_set = ShapeList()
-                # Remove Vertex intersections coincident to Edge intersections
+                for shape in common:
-                vts = common_set.vertices()
+                    if isinstance(shape, Wire):
-                eds = common_set.edges()
+                        common_set.extend(shape.edges())
-                if vts and eds:
+                    else:
-                    filtered_vts = ShapeList(
+                        common_set.append(shape)
-                        [
+                common_set = to_vertex(set(to_vector(common_set)))
-                            v
+                common_set = filter_shapes_by_order(common_set, [Vertex, Edge])
-                            for v in vts
-                            if all(v.distance_to(e) > TOLERANCE for e in eds)
-                        ]
-                    )
-                    common_set = filtered_vts + eds
             else:
                 return None
 
@@ -991,7 +971,7 @@ class Mixin1D(Shape):
         Returns:
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find normal of empty edge/wire")
 
         curve = self.geom_adaptor()
@@ -1225,7 +1205,7 @@ class Mixin1D(Shape):
         Returns:
 
         """
-        if self.wrapped is None or face.wrapped is None:
+        if self._wrapped is None or not face:
             raise ValueError("Can't project an empty Edge or Wire onto empty Face")
 
         bldr = BRepProj_Projection(
@@ -1297,7 +1277,7 @@ class Mixin1D(Shape):
 
             return edges
 
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't project empty edge/wire")
 
         # Setup the projector
@@ -1400,7 +1380,7 @@ class Mixin1D(Shape):
             - **Keep.BOTH**: Returns a tuple `(inside, outside)` where each element is
               either a `Self` or `list[Self]`, or `None` if no corresponding part is found.
         """
-        if self.wrapped is None or tool.wrapped is None:
+        if self._wrapped is None or not tool:
             raise ValueError("Can't split an empty edge/wire/tool")
 
         shape_list = TopTools_ListOfShape()
@@ -1566,7 +1546,7 @@ class Mixin1D(Shape):
         return Shape.get_shape_list(self, "Wire")
 
 
-class Edge(Mixin1D, Shape[TopoDS_Edge]):
+class Edge(Mixin1D[TopoDS_Edge]):
     """An Edge in build123d is a fundamental element in the topological data structure
     representing a one-dimensional geometric entity within a 3D model. It encapsulates
     information about a curve, which could be a line, arc, or other parametrically
@@ -1647,7 +1627,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         Returns:
             Edge: extruded shape
         """
-        if obj.wrapped is None:
+        if not obj:
             raise ValueError("Can't extrude empty vertex")
         return Edge(TopoDS.Edge_s(_extrude_topods_shape(obj.wrapped, direction)))
 
@@ -2538,7 +2518,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         extension_factor: float = 0.1,
     ):
         """Helper method to slightly extend an edge that is bound to a surface"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't extend empty spline")
         if self.geom_type != GeomType.BSPLINE:
             raise TypeError("_extend_spline only works with splines")
@@ -2595,7 +2575,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         Returns:
             ShapeList[Vector]: list of intersection points
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find intersections of empty edge")
 
         # Convert an Axis into an edge at least as large as self and Axis start point
@@ -2723,7 +2703,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
 
     def geom_adaptor(self) -> BRepAdaptor_Curve:
         """Return the Geom Curve from this Edge"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find adaptor for empty edge")
         return BRepAdaptor_Curve(self.wrapped)
 
@@ -2811,7 +2791,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
             float: Normalized parameter in [0.0, 1.0] corresponding to the point's
             closest location on the edge.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find param on empty edge")
 
         pnt = Vector(point)
@@ -2945,7 +2925,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         Returns:
             Edge: reversed
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("An empty edge can't be reversed")
 
         assert isinstance(self.wrapped, TopoDS_Edge)
@@ -3025,7 +3005,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         # if start_u >= end_u:
         #     raise ValueError(f"start ({start_u}) must be less than end ({end_u})")
 
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't trim empty edge")
 
         self_copy = copy.deepcopy(self)
@@ -3060,7 +3040,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         Returns:
             Edge: trimmed edge
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't trim empty edge")
 
         start_u = Mixin1D._to_param(self, start, "start")
@@ -3089,7 +3069,7 @@ class Edge(Mixin1D, Shape[TopoDS_Edge]):
         return Edge(new_edge)
 
 
-class Wire(Mixin1D, Shape[TopoDS_Wire]):
+class Wire(Mixin1D[TopoDS_Wire]):
     """A Wire in build123d is a topological entity representing a connected sequence
     of edges forming a continuous curve or path in 3D space. Wires are essential
     components in modeling complex objects, defining boundaries for surfaces or
@@ -3623,7 +3603,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
         Returns:
             Wire: chamfered wire
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't chamfer empty wire")
 
         reference_edge = edge
@@ -3638,7 +3618,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
         )
 
         for v in vertices:
-            if v.wrapped is None:
+            if not v:
                 continue
             edge_list = vertex_edge_map.FindFromKey(v.wrapped)
 
@@ -3695,7 +3675,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
         Returns:
             Wire: filleted wire
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't fillet an empty wire")
 
         # Create a face to fillet
@@ -3723,7 +3703,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
         Returns:
             Wire: fixed wire
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't fix an empty edge")
 
         sf_w = ShapeFix_Wireframe(self.wrapped)
@@ -3735,7 +3715,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
 
     def geom_adaptor(self) -> BRepAdaptor_CompCurve:
         """Return the Geom Comp Curve for this Wire"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't get geom adaptor of empty wire")
 
         return BRepAdaptor_CompCurve(self.wrapped)
@@ -3779,7 +3759,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
             float: Normalized parameter in [0.0, 1.0] representing the relative
             position of the projected point along the wire.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find point on empty wire")
 
         point_on_curve = Vector(point)
@@ -3932,7 +3912,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
 
         """
         # pylint: disable=too-many-branches
-        if self.wrapped is None or target_object.wrapped is None:
+        if self._wrapped is None or not target_object:
             raise ValueError("Can't project empty Wires or to empty Shapes")
 
         if direction is not None and center is None:
@@ -4021,7 +4001,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
         Returns:
             Wire: stitched wires
         """
-        if self.wrapped is None or other.wrapped is None:
+        if self._wrapped is None or not other:
             raise ValueError("Can't stitch empty wires")
 
         wire_builder = BRepBuilderAPI_MakeWire()
@@ -4065,7 +4045,7 @@ class Wire(Mixin1D, Shape[TopoDS_Wire]):
         """
         # Build a single Geom_BSplineCurve from the wire, in *topological order*
         builder = GeomConvert_CompCurveToBSplineCurve()
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't convert an empty wire")
         wire_explorer = BRepTools_WireExplorer(self.wrapped)
 
@@ -4217,9 +4197,9 @@ def topo_explore_connected_edges(
     parent = parent if parent is not None else edge.topo_parent
     if parent is None:
         raise ValueError("edge has no valid parent")
-    given_topods_edge = edge.wrapped
+    if not edge:
-    if given_topods_edge is None:
         raise ValueError("edge is empty")
+    given_topods_edge = edge.wrapped
     connected_edges = set()
 
     # Find all the TopoDS_Edges for this Shape
@@ -4262,11 +4242,11 @@ def topo_explore_connected_faces(
 ) -> list[TopoDS_Face]:
     """Given an edge extracted from a Shape, return the topods_faces connected to it"""
 
-    if edge.wrapped is None:
+    if not edge:
         raise ValueError("Can't explore from an empty edge")
 
     parent = parent if parent is not None else edge.topo_parent
-    if parent is None or parent.wrapped is None:
+    if not parent:
         raise ValueError("edge has no valid parent")
 
     # make a edge --> faces mapping

src/build123d/topology/shape_core.py

@@ -287,7 +287,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         color: ColorLike | None = None,
         parent: Compound | None = None,
     ):
-        self.wrapped: TOPODS | None = (
+        self._wrapped: TOPODS | None = (
             tcast(Optional[TOPODS], downcast(obj)) if obj is not None else None
         )
         self.for_construction = False
@@ -304,6 +304,18 @@ class Shape(NodeMixin, Generic[TOPODS]):
 
     # pylint: disable=too-many-instance-attributes, too-many-public-methods
 
+    @property
+    def wrapped(self):
+        assert self._wrapped
+        return self._wrapped
+
+    @wrapped.setter
+    def wrapped(self, shape: TOPODS):
+        self._wrapped = shape
+
+    def __bool__(self):
+        return self._wrapped is not None
+
     @property
     @abstractmethod
     def _dim(self) -> int | None:
@@ -312,7 +324,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
     @property
     def area(self) -> float:
         """area -the surface area of all faces in this Shape"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             return 0.0
         properties = GProp_GProps()
         BRepGProp.SurfaceProperties_s(self.wrapped, properties)
@@ -351,7 +363,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             GeomType: The geometry type of the shape
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot determine geometry type of an empty shape")
 
         shape: TopAbs_ShapeEnum = shapetype(self.wrapped)
@@ -380,7 +392,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             bool: is the shape manifold or water tight
         """
         # Extract one or more (if a Compound) shape from self
-        if self.wrapped is None:
+        if self._wrapped is None:
             return False
         shape_stack = get_top_level_topods_shapes(self.wrapped)
 
@@ -431,12 +443,12 @@ class Shape(NodeMixin, Generic[TOPODS]):
         underlying shape with the potential to be given a location and an
         orientation.
         """
-        return self.wrapped is None or self.wrapped.IsNull()
+        return self._wrapped is None or self.wrapped.IsNull()
 
     @property
     def is_planar_face(self) -> bool:
         """Is the shape a planar face even though its geom_type may not be PLANE"""
-        if self.wrapped is None or not isinstance(self.wrapped, TopoDS_Face):
+        if self._wrapped is None or not isinstance(self.wrapped, TopoDS_Face):
             return False
         surface = BRep_Tool.Surface_s(self.wrapped)
         is_face_planar = GeomLib_IsPlanarSurface(surface, TOLERANCE)
@@ -448,7 +460,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         subshapes. See the OCCT docs on BRepCheck_Analyzer::IsValid for a full
         description of what is checked.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return True
         chk = BRepCheck_Analyzer(self.wrapped)
         chk.SetParallel(True)
@@ -474,7 +486,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
     @property
     def location(self) -> Location:
         """Get this Shape's Location"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't find the location of an empty shape")
         return Location(self.wrapped.Location())
 
@@ -518,7 +530,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             - It is commonly used in structural analysis, mechanical simulations,
               and physics-based motion calculations.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't calculate matrix for empty shape")
         properties = GProp_GProps()
         BRepGProp.VolumeProperties_s(self.wrapped, properties)
@@ -546,7 +558,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
     @property
     def position(self) -> Vector:
         """Get the position component of this Shape's Location"""
-        if self.wrapped is None or self.location is None:
+        if self._wrapped is None or self.location is None:
             raise ValueError("Can't find the position of an empty shape")
         return self.location.position
 
@@ -575,7 +587,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             (Vector(0, 1, 0), 1000.0),
             (Vector(0, 0, 1), 300.0)]
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't calculate properties for empty shape")
 
         properties = GProp_GProps()
@@ -615,7 +627,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             (150.0, 200.0, 50.0)
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't calculate moments for empty shape")
 
         properties = GProp_GProps()
@@ -785,7 +797,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if obj.wrapped is None:
+        if not obj:
             return 0.0
 
         properties = GProp_GProps()
@@ -805,7 +817,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         ],
     ) -> ShapeList:
         """Helper to extract entities of a specific type from a shape."""
-        if shape.wrapped is None:
+        if not shape:
             return ShapeList()
         shape_list = ShapeList(
             [shape.__class__.cast(i) for i in shape.entities(entity_type)]
@@ -859,7 +871,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         if not all(summand._dim == addend_dim for summand in summands):
             raise ValueError("Only shapes with the same dimension can be added")
 
-        if self.wrapped is None:  # an empty object
+        if self._wrapped is None:  # an empty object
             if len(summands) == 1:
                 sum_shape = summands[0]
             else:
@@ -876,7 +888,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         """intersect shape with self operator &"""
         others = other if isinstance(other, (list, tuple)) else [other]
 
-        if self.wrapped is None or (isinstance(other, Shape) and other.wrapped is None):
+        if not self or (isinstance(other, Shape) and not other):
             raise ValueError("Cannot intersect shape with empty compound")
         new_shape = self.intersect(*others)
 
@@ -948,7 +960,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
 
     def __hash__(self) -> int:
         """Return hash code"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             return 0
         return hash(self.wrapped)
 
@@ -966,7 +978,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
     def __sub__(self, other: None | Shape | Iterable[Shape]) -> Self | ShapeList[Self]:
         """cut shape from self operator -"""
 
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot subtract shape from empty compound")
 
         # Convert `other` to list of base objects and filter out None values
@@ -1014,7 +1026,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             BoundBox: A box sized to contain this Shape
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return BoundBox(Bnd_Box())
         tolerance = TOLERANCE if tolerance is None else tolerance
         return BoundBox.from_topo_ds(self.wrapped, tolerance=tolerance, optimal=optimal)
@@ -1033,7 +1045,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Shape: Original object with extraneous internal edges removed
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return self
         upgrader = ShapeUpgrade_UnifySameDomain(self.wrapped, True, True, True)
         upgrader.AllowInternalEdges(False)
@@ -1112,7 +1124,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None or other.wrapped is None:
+        if self._wrapped is None or not other:
             raise ValueError("Cannot calculate distance to or from an empty shape")
 
         return BRepExtrema_DistShapeShape(self.wrapped, other.wrapped).Value()
@@ -1125,7 +1137,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         self, other: Shape | VectorLike
     ) -> tuple[float, Vector, Vector]:
         """Minimal distance between two shapes and the points on each shape"""
-        if self.wrapped is None or (isinstance(other, Shape) and other.wrapped is None):
+        if self._wrapped is None or (isinstance(other, Shape) and not other):
             raise ValueError("Cannot calculate distance to or from an empty shape")
 
         if isinstance(other, Shape):
@@ -1155,14 +1167,14 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot calculate distance to or from an empty shape")
 
         dist_calc = BRepExtrema_DistShapeShape()
         dist_calc.LoadS1(self.wrapped)
 
         for other_shape in others:
-            if other_shape.wrapped is None:
+            if not other_shape:
                 raise ValueError("Cannot calculate distance to or from an empty shape")
             dist_calc.LoadS2(other_shape.wrapped)
             dist_calc.Perform()
@@ -1181,7 +1193,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
 
     def entities(self, topo_type: Shapes) -> list[TopoDS_Shape]:
         """Return all of the TopoDS sub entities of the given type"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             return []
         return _topods_entities(self.wrapped, topo_type)
 
@@ -1209,7 +1221,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             list[Face]: A list of intersected faces sorted by distance from axis.position
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return ShapeList()
 
         line = gce_MakeLin(axis.wrapped).Value()
@@ -1239,7 +1251,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
 
     def fix(self) -> Self:
         """fix - try to fix shape if not valid"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             return self
         if not self.is_valid:
             shape_copy: Shape = copy.deepcopy(self, None)
@@ -1281,7 +1293,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
     #     self, child_type: Shapes, parent_type: Shapes
     # ) -> Dict[Shape, list[Shape]]:
     #     """This function is very slow on M1 macs and is currently unused"""
-    #     if self.wrapped is None:
+    #     if self._wrapped is None:
     #         return {}
 
     #     res = TopTools_IndexedDataMapOfShapeListOfShape()
@@ -1319,7 +1331,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             (e.g., edges, vertices) and other compounds, the method returns a list
             of only the simple shapes directly contained at the top level.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return ShapeList()
         return ShapeList(
             self.__class__.cast(s) for s in get_top_level_topods_shapes(self.wrapped)
@@ -1327,7 +1339,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
 
     def intersect(
         self, *to_intersect: Shape | Vector | Location | Axis | Plane
-    ) -> None | Self | ShapeList[Self]:
+    ) -> None | ShapeList[Self]:
         """Intersection of the arguments and this shape
 
         Args:
@@ -1335,8 +1347,8 @@ class Shape(NodeMixin, Generic[TOPODS]):
                 intersect with
 
         Returns:
-            Self | ShapeList[Self]: Resulting object may be of a different class than self
+            None | ShapeList[Self]: Resulting ShapeList may contain different class
-                or a ShapeList if multiple non-Compound object created
+                than self
         """
 
         def _to_vertex(vec: Vector) -> Vertex:
@@ -1380,15 +1392,12 @@ class Shape(NodeMixin, Generic[TOPODS]):
 
         # Find the shape intersections
         intersect_op = BRepAlgoAPI_Common()
-        shape_intersections = self._bool_op((self,), objs, intersect_op)
+        intersections = self._bool_op((self,), objs, intersect_op)
-        if isinstance(shape_intersections, ShapeList) and not shape_intersections:
+        if isinstance(intersections, ShapeList):
-            return None
+            return intersections or None
-        if (
+        if isinstance(intersections, Shape) and not intersections.is_null:
-            not isinstance(shape_intersections, ShapeList)
+            return ShapeList([intersections])
-            and shape_intersections.is_null
+        return None
-        ):
-            return None
-        return shape_intersections
 
     def is_equal(self, other: Shape) -> bool:
         """Returns True if two shapes are equal, i.e. if they share the same
@@ -1401,7 +1410,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None or other.wrapped is None:
+        if self._wrapped is None or not other:
             return False
         return self.wrapped.IsEqual(other.wrapped)
 
@@ -1416,7 +1425,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None or other.wrapped is None:
+        if self._wrapped is None or not other:
             return False
         return self.wrapped.IsSame(other.wrapped)
 
@@ -1429,7 +1438,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot locate an empty shape")
         if loc.wrapped is None:
             raise ValueError("Cannot locate a shape at an empty location")
@@ -1448,7 +1457,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Shape: copy of Shape at location
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot locate an empty shape")
         if loc.wrapped is None:
             raise ValueError("Cannot locate a shape at an empty location")
@@ -1466,7 +1475,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot mesh an empty shape")
 
         if not BRepTools.Triangulation_s(self.wrapped, tolerance):
@@ -1487,7 +1496,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         if not mirror_plane:
             mirror_plane = Plane.XY
 
-        if self.wrapped is None:
+        if self._wrapped is None:
             return self
         transformation = gp_Trsf()
         transformation.SetMirror(
@@ -1505,7 +1514,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
 
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot move an empty shape")
         if loc.wrapped is None:
             raise ValueError("Cannot move a shape at an empty location")
@@ -1525,7 +1534,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Shape: copy of Shape moved to relative location
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot move an empty shape")
         if loc.wrapped is None:
             raise ValueError("Cannot move a shape at an empty location")
@@ -1539,7 +1548,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             OrientedBoundBox: A box oriented and sized to contain this Shape
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return OrientedBoundBox(Bnd_OBB())
         return OrientedBoundBox(self)
 
@@ -1641,7 +1650,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             - The radius of gyration is computed based on the shape’s mass properties.
             - It is useful for evaluating structural stability and rotational behavior.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't calculate radius of gyration for empty shape")
 
         properties = GProp_GProps()
@@ -1660,7 +1669,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             DeprecationWarning,
             stacklevel=2,
         )
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot relocate an empty shape")
         if loc.wrapped is None:
             raise ValueError("Cannot relocate a shape at an empty location")
@@ -1855,7 +1864,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
                 "keep must be one of Keep.INSIDE, Keep.OUTSIDE, or Keep.BOTH"
             )
 
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot split an empty shape")
 
         # Process the perimeter
@@ -1863,7 +1872,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
             raise ValueError("perimeter must be a closed Wire or Edge")
         perimeter_edges = TopTools_SequenceOfShape()
         for perimeter_edge in perimeter.edges():
-            if perimeter_edge.wrapped is None:
+            if not perimeter_edge:
                 continue
             perimeter_edges.Append(perimeter_edge.wrapped)
 
@@ -1871,7 +1880,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         lefts: list[Shell] = []
         rights: list[Shell] = []
         for target_shell in self.shells():
-            if target_shell.wrapped is None:
+            if not target_shell:
                 continue
             constructor = BRepFeat_SplitShape(target_shell.wrapped)
             constructor.Add(perimeter_edges)
@@ -1900,7 +1909,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         self, tolerance: float, angular_tolerance: float = 0.1
     ) -> tuple[list[Vector], list[tuple[int, int, int]]]:
         """General triangulated approximation"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot tessellate an empty shape")
 
         self.mesh(tolerance, angular_tolerance)
@@ -1962,7 +1971,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Self: Approximated shape
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot approximate an empty shape")
 
         params = ShapeCustom_RestrictionParameters()
@@ -1999,7 +2008,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Shape: a copy of the object, but with geometry transformed
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return self
         new_shape = copy.deepcopy(self, None)
         transformed = downcast(
@@ -2022,7 +2031,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Shape: copy of transformed shape with all objects keeping their type
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return self
         new_shape = copy.deepcopy(self, None)
         transformed = downcast(
@@ -2095,7 +2104,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             Shape: copy of transformed Shape
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return self
         shape_copy: Shape = copy.deepcopy(self, None)
         transformed_shape = BRepBuilderAPI_Transform(
@@ -2126,7 +2135,11 @@ class Shape(NodeMixin, Generic[TOPODS]):
         args = list(args)
         tools = list(tools)
         # Find the highest order class from all the inputs Solid > Vertex
-        order_dict = {type(s): type(s).order for s in [self] + args + tools}
+        order_dict = {
+            type(s): type(s).order
+            for s in [self] + args + tools
+            if hasattr(type(s), "order")
+        }
         highest_order = sorted(order_dict.items(), key=lambda item: item[1])[-1]
 
         # The base of the operation
@@ -2200,7 +2213,7 @@ class Shape(NodeMixin, Generic[TOPODS]):
         Returns:
             tuple[ShapeList[Vertex], ShapeList[Edge]]: section results
         """
-        if self.wrapped is None or other.wrapped is None:
+        if self._wrapped is None or not other:
             return (ShapeList(), ShapeList())
 
         section = BRepAlgoAPI_Section(self.wrapped, other.wrapped)
@@ -2701,15 +2714,16 @@ class ShapeList(list[T]):
                     tol_digits,
                 )
 
-        elif hasattr(group_by, "wrapped"):
+        elif not group_by:
-            if group_by.wrapped is None:
+            raise ValueError("Cannot group by an empty object")
-                raise ValueError("Cannot group by an empty object")
 
-            if isinstance(group_by.wrapped, (TopoDS_Edge, TopoDS_Wire)):
+        elif hasattr(group_by, "wrapped") and isinstance(
+            group_by.wrapped, (TopoDS_Edge, TopoDS_Wire)
+        ):
 
-                def key_f(obj):
+            def key_f(obj):
-                    pnt1, _pnt2 = group_by.closest_points(obj.center())
+                pnt1, _pnt2 = group_by.closest_points(obj.center())
-                    return round(group_by.param_at_point(pnt1), tol_digits)
+                return round(group_by.param_at_point(pnt1), tol_digits)
 
         elif isinstance(group_by, SortBy):
             if group_by == SortBy.LENGTH:
@@ -2815,22 +2829,22 @@ class ShapeList(list[T]):
                 ).position.Z,
                 reverse=reverse,
             )
-        elif hasattr(sort_by, "wrapped"):
+        elif not sort_by:
-            if sort_by.wrapped is None:
+            raise ValueError("Cannot sort by an empty object")
-                raise ValueError("Cannot sort by an empty object")
+        elif hasattr(sort_by, "wrapped") and isinstance(
+            sort_by.wrapped, (TopoDS_Edge, TopoDS_Wire)
+        ):
 
-            if isinstance(sort_by.wrapped, (TopoDS_Edge, TopoDS_Wire)):
+            def u_of_closest_center(obj) -> float:
+                """u-value of closest point between object center and sort_by"""
+                assert not isinstance(sort_by, SortBy)
+                pnt1, _pnt2 = sort_by.closest_points(obj.center())
+                return sort_by.param_at_point(pnt1)
 
-                def u_of_closest_center(obj) -> float:
+            # pylint: disable=unnecessary-lambda
-                    """u-value of closest point between object center and sort_by"""
+            objects = sorted(
-                    assert not isinstance(sort_by, SortBy)
+                self, key=lambda o: u_of_closest_center(o), reverse=reverse
-                    pnt1, _pnt2 = sort_by.closest_points(obj.center())
+            )
-                    return sort_by.param_at_point(pnt1)
-
-                # pylint: disable=unnecessary-lambda
-                objects = sorted(
-                    self, key=lambda o: u_of_closest_center(o), reverse=reverse
-                )
 
         elif isinstance(sort_by, SortBy):
             if sort_by == SortBy.LENGTH:

src/build123d/topology/three_d.py

@@ -56,13 +56,13 @@ from __future__ import annotations
 
 import platform
 import warnings
+from collections.abc import Iterable, Sequence
 from math import radians, cos, tan
-from typing import Union, TYPE_CHECKING
+from typing import TYPE_CHECKING
-
+from typing_extensions import Self
-from collections.abc import Iterable
 
 import OCP.TopAbs as ta
-from OCP.BRepAlgoAPI import BRepAlgoAPI_Cut
+from OCP.BRepAlgoAPI import BRepAlgoAPI_Common, BRepAlgoAPI_Cut, BRepAlgoAPI_Section
 from OCP.BRepBuilderAPI import BRepBuilderAPI_MakeSolid
 from OCP.BRepClass3d import BRepClass3d_SolidClassifier
 from OCP.BRepFeat import BRepFeat_MakeDPrism
@@ -95,6 +95,7 @@ from OCP.gp import gp_Ax2, gp_Pnt
 from build123d.build_enums import CenterOf, GeomType, Kind, Transition, Until
 from build123d.geometry import (
     DEG2RAD,
+    TOLERANCE,
     Axis,
     BoundBox,
     Color,
@@ -104,10 +105,9 @@ from build123d.geometry import (
     Vector,
     VectorLike,
 )
-from typing_extensions import Self
 
 from .one_d import Edge, Wire, Mixin1D
-from .shape_core import Shape, ShapeList, Joint, downcast, shapetype
+from .shape_core import TOPODS, Shape, ShapeList, Joint, downcast, shapetype
 from .two_d import sort_wires_by_build_order, Mixin2D, Face, Shell
 from .utils import (
     _extrude_topods_shape,
@@ -122,7 +122,7 @@ if TYPE_CHECKING:  # pragma: no cover
     from .composite import Compound, Curve, Sketch, Part  # pylint: disable=R0801
 
 
-class Mixin3D(Shape):
+class Mixin3D(Shape[TOPODS]):
     """Additional methods to add to 3D Shape classes"""
 
     project_to_viewport = Mixin1D.project_to_viewport
@@ -420,6 +420,130 @@ class Mixin3D(Shape):
 
         return return_value
 
+    def intersect(
+        self, *to_intersect: Shape | Vector | Location | Axis | Plane
+    ) -> None | ShapeList[Vertex | Edge | Face | Solid]:
+        """Intersect Solid with Shape or geometry object
+
+        Args:
+            to_intersect (Shape | Vector | Location | Axis | Plane): objects to intersect
+
+        Returns:
+            ShapeList[Vertex | Edge | Face | Solid] | None: ShapeList of vertices, edges,
+                faces, and/or solids.
+        """
+
+        def to_vector(objs: Iterable) -> ShapeList:
+            return ShapeList([Vector(v) if isinstance(v, Vertex) else v for v in objs])
+
+        def to_vertex(objs: Iterable) -> ShapeList:
+            return ShapeList([Vertex(v) if isinstance(v, Vector) else v for v in objs])
+
+        def bool_op(
+            args: Sequence,
+            tools: Sequence,
+            operation: BRepAlgoAPI_Common | BRepAlgoAPI_Section,
+        ) -> ShapeList:
+            # Wrap Shape._bool_op for corrected output
+            intersections: Shape | ShapeList = Shape()._bool_op(args, tools, operation)
+            if isinstance(intersections, ShapeList):
+                return intersections or ShapeList()
+            if isinstance(intersections, Shape) and not intersections.is_null:
+                return ShapeList([intersections])
+            return ShapeList()
+
+        def filter_shapes_by_order(shapes: ShapeList, orders: list) -> ShapeList:
+            # Remove lower order shapes from list which *appear* to be part of
+            # a higher order shape using a lazy distance check
+            # (sufficient for vertices, may be an issue for higher orders)
+            order_groups = []
+            for order in orders:
+                order_groups.append(
+                    ShapeList([s for s in shapes if isinstance(s, order)])
+                )
+
+            filtered_shapes = order_groups[-1]
+            for i in range(len(order_groups) - 1):
+                los = order_groups[i]
+                his: list = sum(order_groups[i + 1 :], [])
+                filtered_shapes.extend(
+                    ShapeList(
+                        lo
+                        for lo in los
+                        if all(lo.distance_to(hi) > TOLERANCE for hi in his)
+                    )
+                )
+
+            return filtered_shapes
+
+        common_set: ShapeList[Vertex | Edge | Face | Solid] = ShapeList([self])
+        target: Shape
+        for other in to_intersect:
+            # Conform target type
+            match other:
+                case Axis():
+                    # BRepAlgoAPI_Section seems happier if Edge isnt infinite
+                    bbox = self.bounding_box()
+                    dist = self.distance_to(other.position)
+                    dist = dist if dist >= 1 else 1
+                    target = Edge.make_line(
+                        other.position - other.direction * bbox.diagonal * dist,
+                        other.position + other.direction * bbox.diagonal * dist,
+                    )
+                case Plane():
+                    target = Face(other)
+                case Vector():
+                    target = Vertex(other)
+                case Location():
+                    target = Vertex(other.position)
+                case _ if issubclass(type(other), Shape):
+                    target = other
+                case _:
+                    raise ValueError(f"Unsupported type to_intersect: {type(other)}")
+
+            # Find common matches
+            common: list[Vertex | Edge | Wire | Face | Shell | Solid] = []
+            result: ShapeList | None
+            for obj in common_set:
+                match (obj, target):
+                    case (_, Vertex() | Edge() | Wire() | Face() | Shell() | Solid()):
+                        operation = BRepAlgoAPI_Section()
+                        result = bool_op((obj,), (target,), operation)
+                        if (
+                            not isinstance(obj, Edge | Wire)
+                            and not isinstance(target, (Edge | Wire))
+                        ) or (isinstance(obj, Solid) or isinstance(target, Solid)):
+                            # Face + Edge combinations may produce an intersection
+                            # with Common but always with Section.
+                            # No easy way to deduplicate
+                            # Many Solid + Edge combinations need Common
+                            operation = BRepAlgoAPI_Common()
+                            result.extend(bool_op((obj,), (target,), operation))
+
+                    case _ if issubclass(type(target), Shape):
+                        result = target.intersect(obj)
+
+                if result:
+                    common.extend(result)
+
+            if common:
+                common_set = ShapeList()
+                for shape in common:
+                    if isinstance(shape, Wire):
+                        common_set.extend(shape.edges())
+                    elif isinstance(shape, Shell):
+                        common_set.extend(shape.faces())
+                    else:
+                        common_set.append(shape)
+                common_set = to_vertex(set(to_vector(common_set)))
+                common_set = filter_shapes_by_order(
+                    common_set, [Vertex, Edge, Face, Solid]
+                )
+            else:
+                return None
+
+        return ShapeList(common_set)
+
     def is_inside(self, point: VectorLike, tolerance: float = 1.0e-6) -> bool:
         """Returns whether or not the point is inside a solid or compound
         object within the specified tolerance.
@@ -590,7 +714,7 @@ class Mixin3D(Shape):
         return Shape.get_shape_list(self, "Solid")
 
 
-class Solid(Mixin3D, Shape[TopoDS_Solid]):
+class Solid(Mixin3D[TopoDS_Solid]):
     """A Solid in build123d represents a three-dimensional solid geometry
     in a topological structure. A solid is a closed and bounded volume, enclosing
     a region in 3D space. It comprises faces, edges, and vertices connected in a
@@ -1269,7 +1393,7 @@ class Solid(Mixin3D, Shape[TopoDS_Solid]):
             outer_wire = section
             inner_wires = inner_wires if inner_wires else []
 
-        shapes = []
+        shapes: list[Mixin3D[TopoDS_Shape]] = []
         for wire in [outer_wire] + inner_wires:
             builder = BRepOffsetAPI_MakePipeShell(Wire(path).wrapped)
 

src/build123d/topology/two_d.py

@@ -66,7 +66,7 @@ import OCP.TopAbs as ta
 from OCP.BRep import BRep_Builder, BRep_Tool
 from OCP.BRepAdaptor import BRepAdaptor_Curve, BRepAdaptor_Surface
 from OCP.BRepAlgo import BRepAlgo
-from OCP.BRepAlgoAPI import BRepAlgoAPI_Common
+from OCP.BRepAlgoAPI import BRepAlgoAPI_Common, BRepAlgoAPI_Section
 from OCP.BRepBuilderAPI import (
     BRepBuilderAPI_MakeEdge,
     BRepBuilderAPI_MakeFace,
@@ -139,6 +139,7 @@ from build123d.geometry import (
 
 from .one_d import Edge, Mixin1D, Wire
 from .shape_core import (
+    TOPODS,
     Shape,
     ShapeList,
     SkipClean,
@@ -165,7 +166,7 @@ if TYPE_CHECKING:  # pragma: no cover
 T = TypeVar("T", Edge, Wire, "Face")
 
 
-class Mixin2D(ABC, Shape):
+class Mixin2D(ABC, Shape[TOPODS]):
     """Additional methods to add to Face and Shell class"""
 
     project_to_viewport = Mixin1D.project_to_viewport
@@ -213,7 +214,7 @@ class Mixin2D(ABC, Shape):
 
     def __neg__(self) -> Self:
         """Reverse normal operator -"""
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Invalid Shape")
         new_surface = copy.deepcopy(self)
         new_surface.wrapped = downcast(self.wrapped.Complemented())
@@ -244,7 +245,7 @@ class Mixin2D(ABC, Shape):
         Returns:
             list[tuple[Vector, Vector]]: Point and normal of intersection
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return []
 
         intersection_line = gce_MakeLin(other.wrapped).Value()
@@ -278,6 +279,126 @@ class Mixin2D(ABC, Shape):
 
         return result
 
+    def intersect(
+        self, *to_intersect: Shape | Vector | Location | Axis | Plane
+    ) -> None | ShapeList[Vertex | Edge | Face]:
+        """Intersect Face with Shape or geometry object
+
+        Args:
+            to_intersect (Shape | Vector | Location | Axis | Plane): objects to intersect
+
+        Returns:
+            ShapeList[Vertex | Edge | Face] | None: ShapeList of vertices, edges, and/or
+                faces.
+        """
+
+        def to_vector(objs: Iterable) -> ShapeList:
+            return ShapeList([Vector(v) if isinstance(v, Vertex) else v for v in objs])
+
+        def to_vertex(objs: Iterable) -> ShapeList:
+            return ShapeList([Vertex(v) if isinstance(v, Vector) else v for v in objs])
+
+        def bool_op(
+            args: Sequence,
+            tools: Sequence,
+            operation: BRepAlgoAPI_Section | BRepAlgoAPI_Common,
+        ) -> ShapeList:
+            # Wrap Shape._bool_op for corrected output
+            intersections: Shape | ShapeList = Shape()._bool_op(args, tools, operation)
+            if isinstance(intersections, ShapeList):
+                return intersections or ShapeList()
+            if isinstance(intersections, Shape) and not intersections.is_null:
+                return ShapeList([intersections])
+            return ShapeList()
+
+        def filter_shapes_by_order(shapes: ShapeList, orders: list) -> ShapeList:
+            # Remove lower order shapes from list which *appear* to be part of
+            # a higher order shape using a lazy distance check
+            # (sufficient for vertices, may be an issue for higher orders)
+            order_groups = []
+            for order in orders:
+                order_groups.append(
+                    ShapeList([s for s in shapes if isinstance(s, order)])
+                )
+
+            filtered_shapes = order_groups[-1]
+            for i in range(len(order_groups) - 1):
+                los = order_groups[i]
+                his: list = sum(order_groups[i + 1 :], [])
+                filtered_shapes.extend(
+                    ShapeList(
+                        lo
+                        for lo in los
+                        if all(lo.distance_to(hi) > TOLERANCE for hi in his)
+                    )
+                )
+
+            return filtered_shapes
+
+        common_set: ShapeList[Vertex | Edge | Face | Shell] = ShapeList([self])
+        target: Shape
+        for other in to_intersect:
+            # Conform target type
+            match other:
+                case Axis():
+                    # BRepAlgoAPI_Section seems happier if Edge isnt infinite
+                    bbox = self.bounding_box()
+                    dist = self.distance_to(other.position)
+                    dist = dist if dist >= 1 else 1
+                    target = Edge.make_line(
+                        other.position - other.direction * bbox.diagonal * dist,
+                        other.position + other.direction * bbox.diagonal * dist,
+                    )
+                case Plane():
+                    target = Face(other)
+                case Vector():
+                    target = Vertex(other)
+                case Location():
+                    target = Vertex(other.position)
+                case _ if issubclass(type(other), Shape):
+                    target = other
+                case _:
+                    raise ValueError(f"Unsupported type to_intersect: {type(other)}")
+
+            # Find common matches
+            common: list[Vertex | Edge | Wire | Face | Shell] = []
+            result: ShapeList | None
+            for obj in common_set:
+                match (obj, target):
+                    case (_, Vertex() | Edge() | Wire() | Face() | Shell()):
+                        operation = BRepAlgoAPI_Section()
+                        result = bool_op((obj,), (target,), operation)
+                        if not isinstance(obj, Edge | Wire) and not isinstance(
+                            target, (Edge | Wire)
+                        ):
+                            # Face + Edge combinations may produce an intersection
+                            # with Common but always with Section.
+                            # No easy way to deduplicate
+                            operation = BRepAlgoAPI_Common()
+                            result.extend(bool_op((obj,), (target,), operation))
+
+                    case _ if issubclass(type(target), Shape):
+                        result = target.intersect(obj)
+
+                if result:
+                    common.extend(result)
+
+            if common:
+                common_set = ShapeList()
+                for shape in common:
+                    if isinstance(shape, Wire):
+                        common_set.extend(shape.edges())
+                    elif isinstance(shape, Shell):
+                        common_set.extend(shape.faces())
+                    else:
+                        common_set.append(shape)
+                common_set = to_vertex(set(to_vector(common_set)))
+                common_set = filter_shapes_by_order(common_set, [Vertex, Edge, Face])
+            else:
+                return None
+
+        return ShapeList(common_set)
+
     @abstractmethod
     def location_at(self, *args: Any, **kwargs: Any) -> Location:
         """A location from a face or shell"""
@@ -350,7 +471,7 @@ class Mixin2D(ABC, Shape):
                 world_point, world_point - target_object_center
             )
 
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't wrap around an empty face")
 
         # Initial setup
@@ -411,7 +532,7 @@ class Mixin2D(ABC, Shape):
             raise RuntimeError(
                 f"Length error of {length_error:.6f} exceeds tolerance {tolerance}"
             )
-        if wrapped_edge.wrapped is None or not wrapped_edge.is_valid:
+        if not wrapped_edge or not wrapped_edge.is_valid:
             raise RuntimeError("Wrapped edge is invalid")
 
         if not snap_to_face:
@@ -434,7 +555,7 @@ class Mixin2D(ABC, Shape):
         return projected_edge
 
 
-class Face(Mixin2D, Shape[TopoDS_Face]):
+class Face(Mixin2D[TopoDS_Face]):
     """A Face in build123d represents a 3D bounded surface within the topological data
     structure. It encapsulates geometric information, defining a face of a 3D shape.
     These faces are integral components of complex structures, such as solids and
@@ -449,7 +570,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
     @overload
     def __init__(
         self,
-        obj: TopoDS_Face,
+        obj: TopoDS_Face | Plane,
         label: str = "",
         color: Color | None = None,
         parent: Compound | None = None,
@@ -457,7 +578,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         """Build a Face from an OCCT TopoDS_Shape/TopoDS_Face
 
         Args:
-            obj (TopoDS_Shape, optional): OCCT Face.
+            obj (TopoDS_Shape | Plane, optional): OCCT Face or Plane.
             label (str, optional): Defaults to ''.
             color (Color, optional): Defaults to None.
             parent (Compound, optional): assembly parent. Defaults to None.
@@ -487,7 +608,9 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
 
         if args:
             l_a = len(args)
-            if isinstance(args[0], TopoDS_Shape):
+            if isinstance(args[0], Plane):
+                obj = args[0]
+            elif isinstance(args[0], TopoDS_Shape):
                 obj, label, color, parent = args[:4] + (None,) * (4 - l_a)
             elif isinstance(args[0], Wire):
                 outer_wire, inner_wires, label, color, parent = args[:5] + (None,) * (
@@ -516,6 +639,9 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         color = kwargs.get("color", color)
         parent = kwargs.get("parent", parent)
 
+        if isinstance(obj, Plane):
+            obj = BRepBuilderAPI_MakeFace(obj.wrapped).Face()
+
         if outer_wire is not None:
             inner_topods_wires = (
                 [w.wrapped for w in inner_wires] if inner_wires is not None else []
@@ -545,7 +671,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
             float: The total surface area, including the area of holes. Returns 0.0 if
             the face is empty.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             return 0.0
 
         return self.without_holes().area
@@ -605,7 +731,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
             ValueError: If the face or its underlying representation is empty.
             ValueError: If the face is not planar.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Can't determine axes_of_symmetry of empty face")
 
         if not self.is_planar_face:
@@ -671,15 +797,13 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
             ).sort_by(Axis(cog, cross_dir))
 
             bottom_area = sum(f.area for f in bottom_list)
-            intersect_area = 0.0
             for flipped_face, bottom_face in zip(top_flipped_list, bottom_list):
                 intersection = flipped_face.intersect(bottom_face)
-                if intersection is None or isinstance(intersection, list):
+                if intersection is None:
                     intersect_area = -1.0
                     break
                 else:
-                    assert isinstance(intersection, Face)
+                    intersect_area = sum(f.area for f in intersection.faces())
-                    intersect_area += intersection.area
 
             if intersect_area == -1.0:
                 continue
@@ -871,7 +995,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         Returns:
             Face: extruded shape
         """
-        if obj.wrapped is None:
+        if not obj:
             raise ValueError("Can't extrude empty object")
         return Face(TopoDS.Face_s(_extrude_topods_shape(obj.wrapped, direction)))
 
@@ -981,7 +1105,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
                 )
                 return single_point_curve
 
-            if shape.wrapped is None:
+            if not shape:
                 raise ValueError("input Edge cannot be empty")
 
             adaptor = BRepAdaptor_Curve(shape.wrapped)
@@ -1015,6 +1139,12 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         plane: Plane = Plane.XY,
     ) -> Face:
         """Create a unlimited size Face aligned with plane"""
+        warnings.warn(
+            "The 'make_plane' method is deprecated and will be removed in a future version.",
+            DeprecationWarning,
+            stacklevel=2,
+        )
+
         pln_shape = BRepBuilderAPI_MakeFace(plane.wrapped).Face()
         return cls(pln_shape)
 
@@ -1104,7 +1234,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
             raise ValueError("exterior must be a Wire or list of Edges")
 
         for edge in outside_edges:
-            if edge.wrapped is None:
+            if not edge:
                 raise ValueError("exterior contains empty edges")
             surface.Add(edge.wrapped, GeomAbs_C0)
 
@@ -1135,7 +1265,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         if interior_wires:
             makeface_object = BRepBuilderAPI_MakeFace(surface_face.wrapped)
             for wire in interior_wires:
-                if wire.wrapped is None:
+                if not wire:
                     raise ValueError("interior_wires contain an empty wire")
                 makeface_object.Add(wire.wrapped)
             try:
@@ -1329,7 +1459,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
             ) from err
 
         result = result.fix()
-        if not result.is_valid or result.wrapped is None:
+        if not result.is_valid or not result:
             raise RuntimeError("Non planar face is invalid")
 
         return result
@@ -1940,7 +2070,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
             DeprecationWarning,
             stacklevel=2,
         )
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot approximate an empty shape")
 
         return self.__class__.cast(BRepAlgo.ConvertFace_s(self.wrapped, tolerance))
@@ -1953,7 +2083,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         Returns:
             Face: A new Face instance identical to the original but without any holes.
         """
-        if self.wrapped is None:
+        if self._wrapped is None:
             raise ValueError("Cannot remove holes from an empty face")
 
         if not (inner_wires := self.inner_wires()):
@@ -2327,7 +2457,7 @@ class Face(Mixin2D, Shape[TopoDS_Face]):
         return wrapped_wire
 
 
-class Shell(Mixin2D, Shape[TopoDS_Shell]):
+class Shell(Mixin2D[TopoDS_Shell]):
     """A Shell is a fundamental component in build123d's topological data structure
     representing a connected set of faces forming a closed surface in 3D space. As
     part of a geometric model, it defines a watertight enclosure, commonly encountered
@@ -2359,7 +2489,7 @@ class Shell(Mixin2D, Shape[TopoDS_Shell]):
             obj = obj_list[0]
 
         if isinstance(obj, Face):
-            if obj.wrapped is None:
+            if not obj:
                 raise ValueError(f"Can't create a Shell from empty Face")
             builder = BRep_Builder()
             shell = TopoDS_Shell()

src/build123d/topology/utils.py

@@ -263,7 +263,10 @@ def _make_topods_face_from_wires(
     for inner_wire in inner_wires:
         if not BRep_Tool.IsClosed_s(inner_wire):
             raise ValueError("Cannot build face(s): inner wire is not closed")
-        face_builder.Add(inner_wire)
+        sf_s = ShapeFix_Shape(inner_wire)
+        sf_s.Perform()
+        fixed_inner_wire = TopoDS.Wire_s(sf_s.Shape())
+        face_builder.Add(fixed_inner_wire)
 
     face_builder.Build()
 

src/build123d/vtk_tools.py

@@ -80,7 +80,7 @@ def to_vtk_poly_data(
     if not HAS_VTK:
         warnings.warn("VTK not supported", stacklevel=2)
 
-    if obj.wrapped is None:
+    if not obj:
         raise ValueError("Cannot convert an empty shape")
 
     vtk_shape = IVtkOCC_Shape(obj.wrapped)

tests/test_build_line.py

@@ -183,6 +183,59 @@ class BuildLineTests(unittest.TestCase):
         self.assertEqual(len(p.edges().filter_by(GeomType.CIRCLE)), 2)
         self.assertEqual(len(p.edges().filter_by(GeomType.LINE)), 3)
 
+        with self.assertRaises(ValueError):
+            p = FilletPolyline(
+                (0, 0),
+                (10, 0),
+                (10, 10),
+                (0, 10),
+                radius=(1, 2, 3, 0),
+                close=True,
+            )
+
+        with self.assertRaises(ValueError):
+            p = FilletPolyline(
+                (0, 0),
+                (10, 0),
+                (10, 10),
+                (0, 10),
+                radius=(1, 2, 3, 4),
+                close=False,
+            )
+
+        with self.assertRaises(ValueError):
+            p = FilletPolyline(
+                (0, 0),
+                (10, 0),
+                (10, 10),
+                (0, 10),
+                radius=-1,
+                close=True,
+            )
+
+        with self.assertRaises(ValueError):
+            p = FilletPolyline(
+                (0, 0),
+                (10, 0),
+                (10, 10),
+                (0, 10),
+                radius=(1, 2),
+                close=True,
+            )
+
+        with BuildLine(Plane.YZ):
+            p = FilletPolyline(
+                (0, 0),
+                (10, 0),
+                (10, 10),
+                (0, 10),
+                radius=(1, 2, 3, 4),
+                close=True,
+            )
+        self.assertEqual(len(p.edges()), 8)
+        self.assertEqual(len(p.edges().filter_by(GeomType.CIRCLE)), 4)
+        self.assertEqual(len(p.edges().filter_by(GeomType.LINE)), 4)
+
         with BuildLine(Plane.YZ):
             p = FilletPolyline(
                 (0, 0, 0), (0, 0, 10), (10, 2, 10), (10, 0, 0), radius=2, close=True

tests/test_direct_api/test_face.py

@@ -130,8 +130,8 @@ class TestFace(unittest.TestCase):
                 distance=1, distance2=2, vertices=[vertex], edge=other_edge
             )
 
-    def test_make_rect(self):
+    def test_plane_as_face(self):
-        test_face = Face.make_plane()
+        test_face = Face(Plane.XY)
         self.assertAlmostEqual(test_face.normal_at(), (0, 0, 1), 5)
 
     def test_length_width(self):

tests/test_direct_api/test_intersection.py

@@ -152,6 +152,7 @@ def test_shape_0d(obj, target, expected):
     run_test(obj, target, expected)
 
 
+# 1d Shapes
 ed1 = Line((0, 0), (5, 0)).edge()
 ed2 = Line((0, -1), (5, 1)).edge()
 ed3 = Line((0, 0, 5), (5, 0, 5)).edge()
@@ -220,6 +221,195 @@ def test_shape_1d(obj, target, expected):
     run_test(obj, target, expected)
 
 
+# 2d Shapes
+fc1 = Rectangle(5, 5).face()
+fc2 = Pos(Z=5) * Rectangle(5, 5).face()
+fc3 = Rot(Y=90) * Rectangle(5, 5).face()
+fc4 = Rot(Z=45) * Rectangle(5, 5).face()
+fc5 = Pos(2.5, 2.5, 2.5) * Rot(0, 90) * Rectangle(5, 5).face()
+fc6 = Pos(2.5, 2.5) * Rot(0, 90, 45, Extrinsic.XYZ) * Rectangle(5, 5).face()
+fc7 = (Rot(90) * Cylinder(2, 4)).faces().filter_by(GeomType.CYLINDER)[0]
+
+fc11 = Rectangle(4, 4).face()
+fc22 = sweep(Rot(90) * CenterArc((0, 0), 2, 0, 180), Line((0, 2), (0, -2)))
+sh1 = Shell([Pos(-4) * fc11, fc22])
+sh2 = Pos(Z=1) * sh1
+sh3 = Shell([Pos(-4) * fc11, fc22, Pos(2, 0, -2) * Rot(0, 90) * fc11])
+sh4 = Shell([Pos(-4) * fc11, fc22, Pos(4) * fc11])
+sh5 = Pos(Z=1) * Shell([Pos(-2, 0, -2) * Rot(0, -90) * fc11, fc22, Pos(2, 0, -2) * Rot(0, 90) * fc11])
+
+shape_2d_matrix = [
+    Case(fc1, vl2, None, "non-coincident", None),
+    Case(fc1, vl1, [Vertex], "coincident", None),
+
+    Case(fc1, lc2, None, "non-coincident", None),
+    Case(fc1, lc1, [Vertex], "coincident", None),
+
+    Case(fc2, ax1, None, "parallel/skew", None),
+    Case(fc3, ax1, [Vertex], "intersecting", None),
+    Case(fc1, ax1, [Edge], "collinear", None),
+
+    Case(fc1, pl3, None, "parallel/skew", None),
+    Case(fc1, pl1, [Edge], "intersecting", None),
+    Case(fc1, pl2, [Face], "collinear", None),
+    Case(fc7, pl1, [Edge, Edge], "multi intersect", None),
+
+    Case(fc1, vt2, None, "non-coincident", None),
+    Case(fc1, vt1, [Vertex], "coincident", None),
+
+    Case(fc1, ed3, None, "parallel/skew", None),
+    Case(Pos(1) * fc3, ed1, [Vertex], "intersecting", None),
+    Case(fc1, ed1, [Edge], "collinear", None),
+    Case(Pos(1.1) * fc3, ed4, [Vertex, Vertex], "multi intersect", None),
+
+    Case(fc1, wi6, None, "parallel/skew", None),
+    Case(Pos(1) * fc3, wi4, [Vertex], "intersecting", None),
+    Case(fc1, wi1, [Edge, Edge], "2 collinear", None),
+    Case(Rot(90) * fc4, wi5, [Vertex, Vertex], "multi intersect", None),
+    Case(Rot(90) * fc4, wi2, [Vertex, Edge], "intersect + collinear", None),
+
+    Case(fc1, fc2, None, "parallel/skew", None),
+    Case(fc1, fc3, [Edge], "intersecting", None),
+    Case(fc1, fc4, [Face], "coplanar", None),
+    Case(fc1, fc5, [Edge], "intersecting edge", None),
+    Case(fc1, fc6, [Vertex], "intersecting vertex", None),
+    Case(fc1, fc7, [Edge, Edge], "multi-intersecting", None),
+    Case(fc7, Pos(Y=2) * fc7, [Face], "cyl intersecting", None),
+
+    Case(sh2, fc1, None, "parallel/skew", None),
+    Case(Pos(Z=1) * sh3, fc1, [Edge], "intersecting", None),
+    Case(sh1, fc1, [Face, Edge], "coplanar + intersecting", None),
+    Case(sh4, fc1, [Face, Face], "2 coplanar", None),
+    Case(sh5, fc1, [Edge, Edge], "2 intersecting", None),
+
+    Case(fc1, [fc4, Pos(2, 2) * fc1], [Face], "multi to_intersect, intersecting", None),
+    Case(fc1, [ed1, Pos(2.5, 2.5) * fc1], [Edge], "multi to_intersect, intersecting", None),
+    Case(fc7, [wi5, fc1], [Vertex], "multi to_intersect, intersecting", None),
+]
+
+@pytest.mark.parametrize("obj, target, expected", make_params(shape_2d_matrix))
+def test_shape_2d(obj, target, expected):
+    run_test(obj, target, expected)
+
+# 3d Shapes
+sl1 = Box(2, 2, 2).solid()
+sl2 = Pos(Z=5) * Box(2, 2, 2).solid()
+sl3 = Cylinder(2, 1).solid() - Cylinder(1.5, 1).solid()
+
+wi7 = Wire([l1 := sl3.faces().sort_by(Axis.Z)[-1].edges()[0].trim(.3, .4),
+          l2 := l1.trim(2, 3),
+          RadiusArc(l1 @ 1, l2 @ 0, 1, short_sagitta=False)
+          ])
+
+shape_3d_matrix = [
+    Case(sl2, vl1, None, "non-coincident", None),
+    Case(Pos(2) * sl1, vl1, [Vertex], "contained", None),
+    Case(Pos(1, 1, -1) * sl1, vl1, [Vertex], "coincident", None),
+
+    Case(sl2, lc1, None, "non-coincident", None),
+    Case(Pos(2) * sl1, lc1, [Vertex], "contained", None),
+    Case(Pos(1, 1, -1) * sl1, lc1, [Vertex], "coincident", None),
+
+    Case(sl2, ax1, None, "non-coincident", None),
+    Case(sl1, ax1, [Edge], "intersecting", None),
+    Case(Pos(1, 1, 1) * sl1, ax2, [Edge], "coincident", None),
+
+    Case(sl1, pl3, None, "non-coincident", None),
+    Case(sl1, pl2, [Face], "intersecting", None),
+
+    Case(sl2, vt1, None, "non-coincident", None),
+    Case(Pos(2) * sl1, vt1, [Vertex], "contained", None),
+    Case(Pos(1, 1, -1) * sl1, vt1, [Vertex], "coincident", None),
+
+    Case(sl1, ed3, None, "non-coincident", None),
+    Case(sl1, ed1, [Edge], "intersecting", None),
+    Case(sl1, Pos(0, 1, 1) * ed1, [Edge], "edge collinear", "duplicate edges, BRepAlgoAPI_Common and _Section both return edge"),
+    Case(sl1, Pos(1, 1, 1) * ed1, [Vertex], "corner coincident", None),
+    Case(Pos(2.1, 1) * sl1, ed4, [Edge, Edge], "multi-intersect", None),
+
+    Case(Pos(2, .5, -1) * sl1, wi6, None, "non-coincident", None),
+    Case(Pos(2, .5, 1) * sl1, wi6, [Edge, Edge], "multi-intersecting", None),
+    Case(sl3, wi7, [Edge, Edge], "multi-coincident, is_equal check", None),
+
+    Case(sl2, fc1, None, "non-coincident", None),
+    Case(sl1, fc1, [Face], "intersecting", None),
+    Case(Pos(3.5, 0, 1) * sl1, fc1, [Edge], "edge collinear", None),
+    Case(Pos(3.5, 3.5) * sl1, fc1, [Vertex], "corner coincident", None),
+    Case(Pos(.9) * sl1, fc7, [Face, Face], "multi-intersecting", None),
+
+    Case(sl2, sh1, None, "non-coincident", None),
+    Case(Pos(-2) * sl1, sh1, [Face, Face], "multi-intersecting", None),
+
+    Case(sl1, sl2, None, "non-coincident", None),
+    Case(sl1, Pos(1, 1, 1) * sl1, [Solid], "intersecting", None),
+    Case(sl1, Pos(2, 2, 1) * sl1, [Edge], "edge collinear", None),
+    Case(sl1, Pos(2, 2, 2) * sl1, [Vertex], "corner coincident", None),
+    Case(sl1, Pos(.45) * sl3, [Solid, Solid], "multi-intersect", None),
+
+    Case(Pos(1.5, 1.5) * sl1, [sl3, Pos(.5, .5) * sl1], [Solid], "multi to_intersect, intersecting", None),
+    Case(Pos(1.5, 1.5) * sl1, [sl3, Pos(Z=.5) * fc1], [Face], "multi to_intersect, intersecting", None),
+]
+
+@pytest.mark.parametrize("obj, target, expected", make_params(shape_3d_matrix))
+def test_shape_3d(obj, target, expected):
+    run_test(obj, target, expected)
+
+# Compound Shapes
+cp1 = Compound(GridLocations(5, 0, 2, 1) * Vertex())
+cp2 = Compound(GridLocations(5, 0, 2, 1) * Line((0, -1), (0, 1)))
+cp3 = Compound(GridLocations(5, 0, 2, 1) * Rectangle(2, 2))
+cp4 = Compound(GridLocations(5, 0, 2, 1) * Box(2, 2, 2))
+
+cv1 = Curve() + [ed1, ed2, ed3]
+sk1 = Sketch() + [fc1, fc2, fc3]
+pt1 = Part() + [sl1, sl2, sl3]
+
+
+shape_compound_matrix = [
+    Case(cp1, vl1, None, "non-coincident", None),
+    Case(Pos(-.5) * cp1, vl1, [Vertex], "intersecting", None),
+
+    Case(cp2, lc1, None, "non-coincident", None),
+    Case(Pos(-.5) * cp2, lc1, [Vertex], "intersecting", None),
+
+    Case(Pos(Z=1) * cp3, ax1, None, "non-coincident", None),
+    Case(cp3, ax1, [Edge, Edge], "intersecting", None),
+
+    Case(Pos(Z=3) * cp4, pl2, None, "non-coincident", None),
+    Case(cp4, pl2, [Face, Face], "intersecting", None),
+
+    Case(cp1, vt1, None, "non-coincident", None),
+    Case(Pos(-.5) * cp1, vt1, [Vertex], "intersecting", None),
+
+    Case(Pos(Z=1) * cp2, ed1, None, "non-coincident", None),
+    Case(cp2, ed1, [Vertex], "intersecting", None),
+
+    Case(Pos(Z=1) * cp3, fc1, None, "non-coincident", None),
+    Case(cp3, fc1, [Face, Face], "intersecting", None),
+
+    Case(Pos(Z=5) * cp4, sl1, None, "non-coincident", None),
+    Case(Pos(2) * cp4, sl1, [Solid], "intersecting", None),
+
+    Case(cp1, Pos(Z=1) * cp1, None, "non-coincident", None),
+    Case(cp1, cp2, [Vertex, Vertex], "intersecting", None),
+    Case(cp2, cp3, [Edge, Edge], "intersecting", None),
+    Case(cp3, cp4, [Face, Face], "intersecting", None),
+
+    Case(cp1, Compound(children=cp1.get_type(Vertex)), [Vertex, Vertex], "mixed child type", None),
+    Case(cp4, Compound(children=cp3.get_type(Face)), [Face, Face], "mixed child type", None),
+
+    Case(cp2, [cp3, cp4], [Edge, Edge], "multi to_intersect, intersecting", None),
+
+    Case(cv1, cp3, [Edge, Edge], "intersecting", "duplicate edges, BRepAlgoAPI_Common and _Section both return edge"),
+    Case(sk1, cp3, [Face, Face], "intersecting", None),
+    Case(pt1, cp3, [Face, Face], "intersecting", None),
+
+]
+
+@pytest.mark.parametrize("obj, target, expected", make_params(shape_compound_matrix))
+def test_shape_compound(obj, target, expected):
+    run_test(obj, target, expected)
+
 # FreeCAD issue example
 c1 = CenterArc((0, 0), 10, 0, 360).edge()
 c2 = CenterArc((19, 0), 10, 0, 360).edge()
@@ -240,7 +430,7 @@ freecad_matrix = [
     Case(vert, e1, [Vertex], "vertical, ellipse, tangent", None),
     Case(horz, e1, [Vertex], "horizontal, ellipse, tangent", None),
 
-    Case(c1, c2, [Vertex, Vertex], "circle, skew, intersect", "Should return 2 Vertices"),
+    Case(c1, c2, [Vertex, Vertex], "circle, skew, intersect", None),
     Case(c1, horz, [Vertex], "circle, horiz, tangent", None),
     Case(c2, horz, [Vertex], "circle, horiz, tangent", None),
     Case(c1, vert, [Vertex], "circle, vert, tangent", None),
@@ -263,11 +453,11 @@ w1 = Wire.make_circle(0.5)
 f1 = Face(Wire.make_circle(0.5))
 
 issues_matrix = [
-    Case(t, t, [Face, Face], "issue #1015", "Returns Compound"),
+    Case(t, t, [Face, Face], "issue #1015", None),
-    Case(l, s, [Edge], "issue #945", "Edge.intersect only takes 1D"),
+    Case(l, s, [Edge], "issue #945", None),
-    Case(a, b, [Edge], "issue #918", "Returns empty Compound"),
+    Case(a, b, [Edge], "issue #918", None),
-    Case(e1, w1, [Vertex, Vertex], "issue #697"),
+    Case(e1, w1, [Vertex, Vertex], "issue #697", None),
-    Case(e1, f1, [Edge], "issue #697", "Edge.intersect only takes 1D"),
+    Case(e1, f1, [Edge], "issue #697", None),
 ]
 
 @pytest.mark.parametrize("obj, target, expected", make_params(issues_matrix))
@@ -279,6 +469,9 @@ def test_issues(obj, target, expected):
 exception_matrix = [
     Case(vt1, Color(), None, "Unsupported type", None),
     Case(ed1, Color(), None, "Unsupported type", None),
+    Case(fc1, Color(), None, "Unsupported type", None),
+    Case(sl1, Color(), None, "Unsupported type", None),
+    Case(cp1, Color(), None, "Unsupported type", None),
 ]
 
 @pytest.mark.skip

tests/test_direct_api/test_oriented_bound_box.py

@@ -229,13 +229,15 @@ class TestOrientedBoundBox(unittest.TestCase):
             obb = OrientedBoundBox(rect)
             corners = obb.corners
             poly = Polygon(*corners, align=None)
-            self.assertAlmostEqual(rect.intersect(poly).area, rect.area, 5)
+            area = sum(f.area for f in rect.intersect(poly).faces())
+            self.assertAlmostEqual(area, rect.area, 5)
 
         for face in Box(1, 2, 3).faces():
             obb = OrientedBoundBox(face)
             corners = obb.corners
             poly = Polygon(*corners, align=None)
-            self.assertAlmostEqual(face.intersect(poly).area, face.area, 5)
+            area = sum(f.area for f in face.intersect(poly).faces())
+            self.assertAlmostEqual(area, face.area, 5)
 
     def test_line_corners(self):
         """

tests/test_direct_api/test_shape.py

@@ -299,7 +299,8 @@ class TestShape(unittest.TestCase):
         predicted_location = Location(offset) * Rotation(*rotation)
         located_shape = Solid.make_box(1, 1, 1).locate(predicted_location)
         intersect = shape.intersect(located_shape)
-        self.assertAlmostEqual(intersect.volume, 1, 5)
+        volume = sum(s.volume for s in intersect.solids())
+        self.assertAlmostEqual(volume, 1, 5)
 
     def test_position_and_orientation(self):
         box = Solid.make_box(1, 1, 1).locate(Location((1, 2, 3), (10, 20, 30)))
@@ -475,7 +476,7 @@ class TestShape(unittest.TestCase):
         self.assertAlmostEqual(Vector(verts[0]), (1, 2, 0), 5)
         self.assertListEqual(edges, [])
 
-        verts, edges = Vertex(1, 2, 0)._ocp_section(Face.make_plane(Plane.XY))
+        verts, edges = Vertex(1, 2, 0)._ocp_section(Face(Plane.XY))
         self.assertAlmostEqual(Vector(verts[0]), (1, 2, 0), 5)
         self.assertListEqual(edges, [])
 
@@ -493,7 +494,7 @@ class TestShape(unittest.TestCase):
         self.assertEqual(len(edges1), 1)
         self.assertAlmostEqual(edges1[0].length, 20, 5)
 
-        vertices2, edges2 = cylinder._ocp_section(Face.make_plane(pln))
+        vertices2, edges2 = cylinder._ocp_section(Face(pln))
         self.assertEqual(len(vertices2), 1)
         self.assertEqual(len(edges2), 1)
         self.assertAlmostEqual(Vector(vertices2[0]), (5, 0, 0), 5)
@@ -588,7 +589,7 @@ class TestShape(unittest.TestCase):
             empty.distance_to_with_closest_points(Vector(1, 1, 1))
         with self.assertRaises(ValueError):
             empty.distance_to(Vector(1, 1, 1))
-        with self.assertRaises(ValueError):
+        with self.assertRaises(AttributeError):
             box.intersect(empty_loc)
         self.assertEqual(empty._ocp_section(Vertex(1, 1, 1)), ([], []))
         self.assertEqual(empty.faces_intersected_by_axis(Axis.Z), ShapeList())

tests/test_direct_api/test_solid.py

@@ -153,7 +153,9 @@ class TestSolid(unittest.TestCase):
         self.assertAlmostEqual(twist.volume, 1, 5)
         top = twist.faces().sort_by(Axis.Z)[-1].rotate(Axis.Z, 45)
         bottom = twist.faces().sort_by(Axis.Z)[0]
-        self.assertAlmostEqual(top.translate((0, 0, -1)).intersect(bottom).area, 1, 5)
+        intersect = top.translate((0, 0, -1)).intersect(bottom)
+        area = sum(f.area for f in intersect.faces())
+        self.assertAlmostEqual(area, 1, 5)
         # Wire
         base = Wire.make_rect(1, 1)
         twist = Solid.extrude_linear_with_rotation(
@@ -162,7 +164,9 @@ class TestSolid(unittest.TestCase):
         self.assertAlmostEqual(twist.volume, 1, 5)
         top = twist.faces().sort_by(Axis.Z)[-1].rotate(Axis.Z, 45)
         bottom = twist.faces().sort_by(Axis.Z)[0]
-        self.assertAlmostEqual(top.translate((0, 0, -1)).intersect(bottom).area, 1, 5)
+        intersect = top.translate((0, 0, -1)).intersect(bottom)
+        area = sum(f.area for f in intersect.faces())
+        self.assertAlmostEqual(area, 1, 5)
 
     def test_make_loft(self):
         loft = Solid.make_loft(

tests/test_drafting.py

@@ -231,7 +231,8 @@ class DimensionLineTestCase(unittest.TestCase):
                 ],
                 draft=metric,
             )
-        self.assertGreater(hole.intersect(d_line).area, 0)
+        area = sum(f.area for f in hole.intersect(d_line).faces())
+        self.assertGreater(area, 0)
 
     def test_outside_arrows(self):
         d_line = DimensionLine([(0, 0, 0), (15, 0, 0)], draft=metric)