diff --git a/cadquery/assembly.py b/cadquery/assembly.py
index 6974e5873..bed1fdc8a 100644
--- a/cadquery/assembly.py
+++ b/cadquery/assembly.py
@@ -33,10 +33,10 @@
     exportVRML,
     exportGLTF,
     STEPExportModeLiterals,
-    ExportModes,
 )
 
 from .selectors import _expression_grammar as _selector_grammar
+from .utils import deprecate
 
 # type definitions
 AssemblyObjects = Union[Shape, Workplane, None]
@@ -451,6 +451,7 @@ def solve(self, verbosity: int = 0) -> "Assembly":
 
         return self
 
+    @deprecate()
     def save(
         self,
         path: str,
@@ -507,6 +508,62 @@ def save(
 
         return self
 
+    def export(
+        self,
+        path: str,
+        exportType: Optional[ExportLiterals] = None,
+        mode: STEPExportModeLiterals = "default",
+        tolerance: float = 0.1,
+        angularTolerance: float = 0.1,
+        **kwargs,
+    ) -> "Assembly":
+        """
+        Save assembly to a file.
+
+        :param path: Path and filename for writing.
+        :param exportType: export format (default: None, results in format being inferred form the path)
+        :param mode: STEP only - See :meth:`~cadquery.occ_impl.exporters.assembly.exportAssembly`.
+        :param tolerance: the deflection tolerance, in model units. Only used for glTF, VRML. Default 0.1.
+        :param angularTolerance: the angular tolerance, in radians. Only used for glTF, VRML. Default 0.1.
+        :param \\**kwargs: Additional keyword arguments.  Only used for STEP, glTF and STL.
+            See :meth:`~cadquery.occ_impl.exporters.assembly.exportAssembly`.
+        :param ascii: STL only - Sets whether or not STL export should be text or binary
+        :type ascii: bool
+        """
+
+        # Make sure the export mode setting is correct
+        if mode not in get_args(STEPExportModeLiterals):
+            raise ValueError(f"Unknown assembly export mode {mode} for STEP")
+
+        if exportType is None:
+            t = path.split(".")[-1].upper()
+            if t in ("STEP", "XML", "VRML", "VTKJS", "GLTF", "GLB", "STL"):
+                exportType = cast(ExportLiterals, t)
+            else:
+                raise ValueError("Unknown extension, specify export type explicitly")
+
+        if exportType == "STEP":
+            exportAssembly(self, path, mode, **kwargs)
+        elif exportType == "XML":
+            exportCAF(self, path)
+        elif exportType == "VRML":
+            exportVRML(self, path, tolerance, angularTolerance)
+        elif exportType == "GLTF" or exportType == "GLB":
+            exportGLTF(self, path, None, tolerance, angularTolerance)
+        elif exportType == "VTKJS":
+            exportVTKJS(self, path)
+        elif exportType == "STL":
+            # Handle the ascii setting for STL export
+            export_ascii = False
+            if "ascii" in kwargs:
+                export_ascii = bool(kwargs.get("ascii"))
+
+            self.toCompound().exportStl(path, tolerance, angularTolerance, export_ascii)
+        else:
+            raise ValueError(f"Unknown format: {exportType}")
+
+        return self
+
     @classmethod
     def load(cls, path: str) -> "Assembly":
 
diff --git a/cadquery/cq.py b/cadquery/cq.py
index 76a0f98ad..d3d25ba0b 100644
--- a/cadquery/cq.py
+++ b/cadquery/cq.py
@@ -33,9 +33,10 @@
     List,
     cast,
     Dict,
+    Iterator,
 )
 from typing_extensions import Literal
-from inspect import Parameter, Signature, isbuiltin
+from inspect import Parameter, Signature
 
 
 from .occ_impl.geom import Vector, Plane, Location
@@ -51,8 +52,9 @@
 )
 
 from .occ_impl.exporters.svg import getSVG, exportSVG
+from .occ_impl.exporters import export
 
-from .utils import deprecate, deprecate_kwarg_name
+from .utils import deprecate, deprecate_kwarg_name, get_arity
 
 from .selectors import (
     Selector,
@@ -270,7 +272,7 @@ def split(self: T, keepTop: bool = False, keepBottom: bool = False) -> T:
         ...
 
     @overload
-    def split(self: T, splitter: Union[T, Shape]) -> T:
+    def split(self: T, splitter: Union["Workplane", Shape]) -> T:
         ...
 
     def split(self: T, *args, **kwargs) -> T:
@@ -383,9 +385,7 @@ def combineSolids(
             raise ValueError("Cannot Combine: at least one solid required!")
 
         # get context solid and we don't want to find our own objects
-        ctxSolid = self._findType(
-            (Solid, Compound), searchStack=False, searchParents=True
-        )
+        ctxSolid = self._findType((Solid,), searchStack=False, searchParents=True)
         if ctxSolid is None:
             ctxSolid = toCombine.pop(0)
 
@@ -748,8 +748,20 @@ def end(self, n: int = 1) -> "Workplane":
     def _findType(self, types, searchStack=True, searchParents=True):
 
         if searchStack:
-            rv = [s for s in self.objects if isinstance(s, types)]
-            if rv and types == (Solid, Compound):
+            rv = []
+            for obj in self.objects:
+                if isinstance(obj, types):
+                    rv.append(obj)
+                # unpack compounds in a special way when looking for Solids
+                elif isinstance(obj, Compound) and types == (Solid,):
+                    for T in types:
+                        # _entities(...) needed due to weird behavior with shelled object unpacking
+                        rv.extend(T(el) for el in obj._entities(T.__name__))
+                # otherwise unpack compounds normally
+                elif isinstance(obj, Compound):
+                    rv.extend(el for el in obj if isinstance(el, type))
+
+            if rv and types == (Solid,):
                 return Compound.makeCompound(rv)
             elif rv:
                 return rv[0]
@@ -781,7 +793,7 @@ def findSolid(
         results with an object already on the stack.
         """
 
-        found = self._findType((Solid, Compound), searchStack, searchParents)
+        found = self._findType((Solid,), searchStack, searchParents)
 
         if found is None:
             message = "on the stack or " if searchStack else ""
@@ -802,7 +814,7 @@ def findFace(self, searchStack: bool = True, searchParents: bool = True) -> Face
         :returns: A face or None if no face is found.
         """
 
-        found = self._findType(Face, searchStack, searchParents)
+        found = self._findType((Face,), searchStack, searchParents)
 
         if found is None:
             message = "on the stack or " if searchStack else ""
@@ -1460,8 +1472,8 @@ def rarray(
         If you want to position the array at another point, create another workplane
         that is shifted to the position you would like to use as a reference
 
-        :param xSpacing: spacing between points in the x direction ( must be > 0)
-        :param ySpacing: spacing between points in the y direction ( must be > 0)
+        :param xSpacing: spacing between points in the x direction ( must be >= 0)
+        :param ySpacing: spacing between points in the y direction ( must be >= 0)
         :param xCount: number of points ( > 0 )
         :param yCount: number of points ( > 0 )
         :param center: If True, the array will be centered around the workplane center.
@@ -1470,8 +1482,8 @@ def rarray(
           centering along each axis.
         """
 
-        if xSpacing <= 0 or ySpacing <= 0 or xCount < 1 or yCount < 1:
-            raise ValueError("Spacing and count must be > 0 ")
+        if (xSpacing <= 0 and ySpacing <= 0) or xCount < 1 or yCount < 1:
+            raise ValueError("Spacing and count must be > 0 in at least one direction")
 
         if isinstance(center, bool):
             center = (center, center)
@@ -3310,9 +3322,7 @@ def _fuseWithBase(self: T, obj: Shape) -> T:
         :return: a new object that represents the result of combining the base object with obj,
            or obj if one could not be found
         """
-        baseSolid = self._findType(
-            (Solid, Compound), searchStack=True, searchParents=True
-        )
+        baseSolid = self._findType((Solid,), searchStack=True, searchParents=True)
         r = obj
         if baseSolid is not None:
             r = baseSolid.fuse(obj)
@@ -3328,7 +3338,7 @@ def _cutFromBase(self: T, obj: Shape) -> T:
         :return: a new object that represents the result of combining the base object with obj,
            or obj if one could not be found
         """
-        baseSolid = self._findType((Solid, Compound), True, True)
+        baseSolid = self._findType((Solid,), True, True)
 
         r = obj
         if baseSolid is not None:
@@ -3399,9 +3409,7 @@ def union(
 
         # now combine with existing solid, if there is one
         # look for parents to cut from
-        solidRef = self._findType(
-            (Solid, Compound), searchStack=True, searchParents=True
-        )
+        solidRef = self._findType((Solid,), searchStack=True, searchParents=True)
         if solidRef is not None:
             r = solidRef.fuse(*newS, glue=glue, tol=tol)
         elif len(newS) > 1:
@@ -3414,7 +3422,8 @@ def union(
 
         return self.newObject([r])
 
-    def __or__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
+    @deprecate()
+    def __or__(self: T, other: Union["Workplane", Solid, Compound]) -> T:
         """
         Syntactic sugar for union.
 
@@ -3426,15 +3435,15 @@ def __or__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
             Sphere = Workplane("XY").sphere(1)
             result = Box | Sphere
         """
-        return self.union(toUnion)
+        return self.union(other)
 
-    def __add__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
+    def __add__(self: T, other: Union["Workplane", Solid, Compound]) -> T:
         """
         Syntactic sugar for union.
 
         Notice that :code:`r = a + b` is equivalent to :code:`r = a.union(b)` and :code:`r = a | b`.
         """
-        return self.union(toUnion)
+        return self.union(other)
 
     def cut(
         self: T,
@@ -3472,7 +3481,7 @@ def cut(
 
         return self.newObject([newS])
 
-    def __sub__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
+    def __sub__(self: T, other: Union["Workplane", Solid, Compound]) -> T:
         """
         Syntactic sugar for cut.
 
@@ -3484,7 +3493,7 @@ def __sub__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
             Sphere = Workplane("XY").sphere(1)
             result = Box - Sphere
         """
-        return self.cut(toUnion)
+        return self.cut(other)
 
     def intersect(
         self: T,
@@ -3522,7 +3531,8 @@ def intersect(
 
         return self.newObject([newS])
 
-    def __and__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
+    @deprecate()
+    def __and__(self: T, other: Union["Workplane", Solid, Compound]) -> T:
         """
         Syntactic sugar for intersect.
 
@@ -3534,7 +3544,38 @@ def __and__(self: T, toUnion: Union["Workplane", Solid, Compound]) -> T:
             Sphere = Workplane("XY").sphere(1)
             result = Box & Sphere
         """
-        return self.intersect(toUnion)
+
+        return self.intersect(other)
+
+    def __mul__(self: T, other: Union["Workplane", Solid, Compound]) -> T:
+        """
+        Syntactic sugar for intersect.
+
+        Notice that :code:`r = a * b` is equivalent to :code:`r = a.intersect(b)`.
+
+        Example::
+
+            Box = Workplane("XY").box(1, 1, 1, centered=(False, False, False))
+            Sphere = Workplane("XY").sphere(1)
+            result = Box * Sphere
+        """
+
+        return self.intersect(other)
+
+    def __truediv__(self: T, other: Union["Workplane", Solid, Compound]) -> T:
+        """
+        Syntactic sugar for intersect.
+
+        Notice that :code:`r = a / b` is equivalent to :code:`r = a.split(b)`.
+
+        Example::
+
+            Box = Workplane("XY").box(1, 1, 1, centered=(False, False, False))
+            Sphere = Workplane("XY").sphere(1)
+            result = Box / Sphere
+        """
+
+        return self.split(other)
 
     def cutBlind(
         self: T,
@@ -3680,6 +3721,10 @@ def _getFaces(self) -> List[Face]:
         for el in self.objects:
             if isinstance(el, Sketch):
                 rv.extend(el)
+            elif isinstance(el, Face):
+                rv.append(el)
+            elif isinstance(el, Compound):
+                rv.extend(subel for subel in el if isinstance(subel, Face))
 
         if not rv:
             rv.extend(wiresToFaces(self.ctx.popPendingWires()))
@@ -4441,6 +4486,19 @@ def __getitem__(self: T, item: Union[int, Sequence[int], slice]) -> T:
 
         return rv
 
+    def __iter__(self: T) -> Iterator[Shape]:
+        """
+        Special method for iterating over Shapes in objects
+        """
+
+        for el in self.objects:
+            if isinstance(el, Compound):
+                yield from el
+            elif isinstance(el, Shape):
+                yield el
+            elif isinstance(el, Sketch):
+                yield from el
+
     def filter(self: T, f: Callable[[CQObject], bool]) -> T:
         """
         Filter items using a boolean predicate.
@@ -4488,11 +4546,7 @@ def invoke(
         :return: Workplane object.
         """
 
-        if isbuiltin(f):
-            arity = 0  # assume 0 arity for builtins; they cannot be introspected
-        else:
-            arity = f.__code__.co_argcount  # NB: this is not understood by mypy
-
+        arity = get_arity(f)
         rv = self
 
         if arity == 0:
@@ -4506,6 +4560,29 @@ def invoke(
 
         return rv
 
+    def export(
+        self: T,
+        fname: str,
+        tolerance: float = 0.1,
+        angularTolerance: float = 0.1,
+        opt: Optional[Dict[str, Any]] = None,
+    ) -> T:
+        """
+        Export Workplane to file.
+        
+        :param path: Filename.
+        :param tolerance: the deflection tolerance, in model units. Default 0.1.
+        :param angularTolerance: the angular tolerance, in radians. Default 0.1.
+        :param opt: additional options passed to the specific exporter. Default None.
+        :return: Self.
+        """
+
+        export(
+            self, fname, tolerance=tolerance, angularTolerance=angularTolerance, opt=opt
+        )
+
+        return self
+
 
 # alias for backward compatibility
 CQ = Workplane
diff --git a/cadquery/occ_impl/exporters/__init__.py b/cadquery/occ_impl/exporters/__init__.py
index e79da1a9f..985bba629 100644
--- a/cadquery/occ_impl/exporters/__init__.py
+++ b/cadquery/occ_impl/exporters/__init__.py
@@ -2,14 +2,13 @@
 import os
 import io as StringIO
 
-from typing import IO, Optional, Union, cast, Dict, Any
+from typing import IO, Optional, Union, cast, Dict, Any, Iterable
 from typing_extensions import Literal
 
 from OCP.VrmlAPI import VrmlAPI
 
-from ...cq import Workplane
 from ...utils import deprecate
-from ..shapes import Shape
+from ..shapes import Shape, compound
 
 from .svg import getSVG
 from .json import JsonMesh
@@ -17,7 +16,6 @@
 from .threemf import ThreeMFWriter
 from .dxf import exportDXF, DxfDocument
 from .vtk import exportVTP
-from .utils import toCompound
 
 
 class ExportTypes:
@@ -30,15 +28,16 @@ class ExportTypes:
     VRML = "VRML"
     VTP = "VTP"
     THREEMF = "3MF"
+    BREP = "BREP"
 
 
 ExportLiterals = Literal[
-    "STL", "STEP", "AMF", "SVG", "TJS", "DXF", "VRML", "VTP", "3MF"
+    "STL", "STEP", "AMF", "SVG", "TJS", "DXF", "VRML", "VTP", "3MF", "BREP"
 ]
 
 
 def export(
-    w: Union[Shape, Workplane],
+    w: Union[Shape, Iterable[Shape]],
     fname: str,
     exportType: Optional[ExportLiterals] = None,
     tolerance: float = 0.1,
@@ -49,7 +48,7 @@ def export(
     """
     Export Workplane or Shape to file. Multiple entities are converted to compound.
 
-    :param w:  Shape or Workplane to be exported.
+    :param w:  Shape or Iterable[Shape] (e.g. Workplane) to be exported.
     :param fname: output filename.
     :param exportType: the exportFormat to use. If None will be inferred from the extension. Default: None.
     :param tolerance: the deflection tolerance, in model units. Default 0.1.
@@ -63,10 +62,10 @@ def export(
     if not opt:
         opt = {}
 
-    if isinstance(w, Workplane):
-        shape = toCompound(w)
-    else:
+    if isinstance(w, Shape):
         shape = w
+    else:
+        shape = compound(*w)
 
     if exportType is None:
         t = fname.split(".")[-1].upper()
@@ -106,10 +105,7 @@ def export(
             tmfw.write3mf(f)
 
     elif exportType == ExportTypes.DXF:
-        if isinstance(w, Workplane):
-            exportDXF(w, fname, **opt)
-        else:
-            raise ValueError("Only Workplanes can be exported as DXF")
+        exportDXF(w, fname, **opt)
 
     elif exportType == ExportTypes.STEP:
         shape.exportStep(fname, **opt)
@@ -129,6 +125,9 @@ def export(
     elif exportType == ExportTypes.VTP:
         exportVTP(shape, fname, tolerance, angularTolerance)
 
+    elif exportType == ExportTypes.BREP:
+        shape.exportBrep(fname)
+
     else:
         raise ValueError("Unknown export type")
 
@@ -142,7 +141,7 @@ def toString(shape, exportType, tolerance=0.1, angularTolerance=0.05):
 
 @deprecate()
 def exportShape(
-    w: Union[Shape, Workplane],
+    w: Union[Shape, Iterable[Shape]],
     exportType: ExportLiterals,
     fileLike: IO,
     tolerance: float = 0.1,
@@ -164,10 +163,10 @@ def tessellate(shape, angularTolerance):
         return shape.tessellate(tolerance, angularTolerance)
 
     shape: Shape
-    if isinstance(w, Workplane):
-        shape = toCompound(w)
-    else:
+    if isinstance(w, Shape):
         shape = w
+    else:
+        shape = compound(*w)
 
     if exportType == ExportTypes.TJS:
         tess = tessellate(shape, angularTolerance)
diff --git a/cadquery/occ_impl/exporters/dxf.py b/cadquery/occ_impl/exporters/dxf.py
index 8300d3fc2..7a711731b 100644
--- a/cadquery/occ_impl/exporters/dxf.py
+++ b/cadquery/occ_impl/exporters/dxf.py
@@ -1,6 +1,17 @@
 """DXF export utilities."""
 
-from typing import Any, Dict, List, Literal, Optional, Tuple, Union
+from typing import (
+    Any,
+    Dict,
+    List,
+    Literal,
+    Optional,
+    Tuple,
+    Union,
+    Iterable,
+    Protocol,
+    runtime_checkable,
+)
 
 import ezdxf
 from ezdxf import units, zoom
@@ -10,10 +21,10 @@
 from OCP.GC import GC_MakeArcOfEllipse
 from typing_extensions import Self
 
-from ...cq import Face, Plane, Workplane
 from ...units import RAD2DEG
-from ..shapes import Edge, Shape, Compound
-from .utils import toCompound
+from ..shapes import Face, Edge, Shape, Compound, compound
+from ..geom import Plane
+
 
 ApproxOptions = Literal["spline", "arc"]
 DxfEntityAttributes = Tuple[
@@ -21,6 +32,16 @@
 ]
 
 
+@runtime_checkable
+class WorkplaneLike(Protocol):
+    @property
+    def plane(self) -> Plane:
+        ...
+
+    def __iter__(self) -> Iterable[Shape]:
+        ...
+
+
 class DxfDocument:
     """Create DXF document from CadQuery objects.
 
@@ -125,14 +146,18 @@ def add_layer(
 
         return self
 
-    def add_shape(self, workplane: Workplane, layer: str = "") -> Self:
+    def add_shape(self, shape: Union[WorkplaneLike, Shape], layer: str = "") -> Self:
         """Add CadQuery shape to a DXF layer.
 
-        :param workplane: CadQuery Workplane
+        :param s: CadQuery Workplane or Shape
         :param layer: layer definition name
         """
-        plane = workplane.plane
-        shape = toCompound(workplane).transformShape(plane.fG)
+
+        if isinstance(shape, WorkplaneLike):
+            plane = shape.plane
+            shape_ = compound(*shape).transformShape(plane.fG)
+        else:
+            shape_ = shape
 
         general_attributes = {}
         if layer:
@@ -141,20 +166,20 @@ def add_shape(self, workplane: Workplane, layer: str = "") -> Self:
         if self.approx == "spline":
             edges = [
                 e.toSplines() if e.geomType() == "BSPLINE" else e
-                for e in self._ordered_edges(shape)
+                for e in self._ordered_edges(shape_)
             ]
 
         elif self.approx == "arc":
             edges = []
 
             # this is needed to handle free wires
-            for el in shape.Wires():
+            for el in shape_.Wires():
                 edges.extend(
                     self._ordered_edges(Face.makeFromWires(el).toArcs(self.tolerance))
                 )
 
         else:
-            edges = self._ordered_edges(shape)
+            edges = self._ordered_edges(shape_)
 
         for edge in edges:
             converter = self._DISPATCH_MAP.get(edge.geomType(), None)
@@ -342,7 +367,7 @@ def _dxf_spline(cls, edge: Edge, plane: Plane) -> DxfEntityAttributes:
 
 
 def exportDXF(
-    w: Workplane,
+    w: Union[WorkplaneLike, Shape, Iterable[Shape]],
     fname: str,
     approx: Optional[ApproxOptions] = None,
     tolerance: float = 1e-3,
@@ -362,5 +387,11 @@ def exportDXF(
     """
 
     dxf = DxfDocument(approx=approx, tolerance=tolerance, doc_units=doc_units)
-    dxf.add_shape(w)
+
+    if isinstance(w, (WorkplaneLike, Shape)):
+        dxf.add_shape(w)
+    else:
+        for s in w:
+            dxf.add_shape(s)
+
     dxf.document.saveas(fname)
diff --git a/cadquery/occ_impl/exporters/threemf.py b/cadquery/occ_impl/exporters/threemf.py
index cc81ebc7b..e33ed5d64 100644
--- a/cadquery/occ_impl/exporters/threemf.py
+++ b/cadquery/occ_impl/exporters/threemf.py
@@ -4,7 +4,7 @@
 from typing import IO, List, Literal, Tuple, Union
 from zipfile import ZipFile, ZIP_DEFLATED, ZIP_STORED
 
-from ...cq import Compound, Shape, Vector
+from ..shapes import Compound, Shape, Vector
 
 
 class CONTENT_TYPES(object):
@@ -48,8 +48,8 @@ def __init__(
     def write3mf(
         self, outfile: Union[PathLike, str, IO[bytes]],
     ):
-        """ 
-        Write to the given file. 
+        """
+        Write to the given file.
         """
 
         try:
diff --git a/cadquery/occ_impl/exporters/utils.py b/cadquery/occ_impl/exporters/utils.py
deleted file mode 100644
index 32129941b..000000000
--- a/cadquery/occ_impl/exporters/utils.py
+++ /dev/null
@@ -1,7 +0,0 @@
-from ...cq import Workplane
-from ..shapes import Compound, Shape
-
-
-def toCompound(shape: Workplane) -> Compound:
-
-    return Compound.makeCompound(val for val in shape.vals() if isinstance(val, Shape))
diff --git a/cadquery/occ_impl/shapes.py b/cadquery/occ_impl/shapes.py
index ca1d2d68e..3ca248e96 100644
--- a/cadquery/occ_impl/shapes.py
+++ b/cadquery/occ_impl/shapes.py
@@ -1576,6 +1576,23 @@ def __truediv__(self, other: "Shape") -> "Shape":
 
         return split(self, other)
 
+    def export(
+        self: T,
+        fname: str,
+        tolerance: float = 0.1,
+        angularTolerance: float = 0.1,
+        opt: Optional[Dict[str, Any]] = None,
+    ):
+        """
+        Export Shape to file.
+        """
+
+        from .exporters import export  # imported here to prevent circular imports
+
+        export(
+            self, fname, tolerance=tolerance, angularTolerance=angularTolerance, opt=opt
+        )
+
 
 class ShapeProtocol(Protocol):
     @property
@@ -4968,7 +4985,7 @@ def sweep(s: Shape, path: Shape, cap: bool = False) -> Shape:
 @sweep.register
 def sweep(s: Sequence[Shape], path: Shape, cap: bool = False) -> Shape:
     """
-    Sweep edges, wires or faces along a path, multiple sections are supported. 
+    Sweep edges, wires or faces along a path, multiple sections are supported.
     For faces cap has no effect. Do not mix faces with other types.
     """
 
diff --git a/cadquery/sketch.py b/cadquery/sketch.py
index 5ffddf4c0..f76c90dc3 100644
--- a/cadquery/sketch.py
+++ b/cadquery/sketch.py
@@ -12,7 +12,9 @@
     TypeVar,
     cast as tcast,
     Literal,
+    overload,
 )
+
 from math import tan, sin, cos, pi, radians, remainder
 from itertools import product, chain
 from multimethod import multimethod
@@ -21,9 +23,21 @@
 from .hull import find_hull
 from .selectors import StringSyntaxSelector, Selector
 from .types import Real
-
-from .occ_impl.shapes import Shape, Face, Edge, Wire, Compound, Vertex, edgesToWires
+from .utils import get_arity
+
+from .occ_impl.shapes import (
+    Shape,
+    Face,
+    Edge,
+    Wire,
+    Compound,
+    Vertex,
+    edgesToWires,
+    compound,
+    VectorLike,
+)
 from .occ_impl.geom import Location, Vector
+from .occ_impl.exporters import export
 from .occ_impl.importers.dxf import _importDXF
 from .occ_impl.sketch_solver import (
     SketchConstraintSolver,
@@ -35,7 +49,9 @@
     arc_point,
 )
 
-Modes = Literal["a", "s", "i", "c"]  # add, subtract, intersect, construct
+#%% types
+
+Modes = Literal["a", "s", "i", "c", "r"]  # add, subtract, intersect, construct, replace
 Point = Union[Vector, Tuple[Real, Real]]
 TOL = 1e-6
 
@@ -43,6 +59,31 @@
 SketchVal = Union[Shape, Location]
 
 
+# %% utilities
+
+
+def _sanitize_for_bool(obj: SketchVal) -> Shape:
+    """
+    Make sure that a Shape is selected
+    """
+
+    if isinstance(obj, Location):
+        raise ValueError("Location was provided, Shape is required.")
+
+    return obj
+
+
+def _to_compound(obj: Shape) -> Compound:
+
+    if isinstance(obj, Compound):
+        return obj
+    else:
+        return Compound.makeCompound((obj,))
+
+
+# %% Constraint
+
+
 class Constraint(object):
 
     tags: Tuple[str, ...]
@@ -86,6 +127,9 @@ def __init__(
         self.param = tcast(Any, converter)(param) if converter else param
 
 
+# %% Sketch
+
+
 class Sketch(object):
     """
     2D sketch. Supports faces, edges and edges with constraints based construction.
@@ -95,17 +139,21 @@ class Sketch(object):
     locs: List[Location]
 
     _faces: Compound
-    _wires: List[Wire]
     _edges: List[Edge]
 
-    _selection: List[SketchVal]
+    _selection: Optional[List[SketchVal]]
     _constraints: List[Constraint]
 
     _tags: Dict[str, Sequence[SketchVal]]
 
     _solve_status: Optional[Dict[str, Any]]
 
-    def __init__(self: T, parent: Any = None, locs: Iterable[Location] = (Location(),)):
+    def __init__(
+        self: T,
+        parent: Any = None,
+        locs: Iterable[Location] = (Location(),),
+        obj: Optional[Compound] = None,
+    ):
         """
         Construct an empty sketch.
         """
@@ -113,11 +161,10 @@ def __init__(self: T, parent: Any = None, locs: Iterable[Location] = (Location()
         self.parent = parent
         self.locs = list(locs)
 
-        self._faces = Compound.makeCompound(())
-        self._wires = []
+        self._faces = obj if obj else Compound.makeCompound(())
         self._edges = []
 
-        self._selection = []
+        self._selection = None
         self._constraints = []
 
         self._tags = {}
@@ -126,19 +173,23 @@ def __init__(self: T, parent: Any = None, locs: Iterable[Location] = (Location()
 
     def __iter__(self) -> Iterator[Face]:
         """
-        Iterate over faces-locations combinations.
+        Iterate over faces-locations combinations. If not faces are present
+        iterate over edges:
         """
 
-        return iter(f for l in self.locs for f in self._faces.moved(l).Faces())
+        if self._faces:
+            return iter(f for l in self.locs for f in self._faces.moved(l).Faces())
+        else:
+            return iter(e.moved(l) for l in self.locs for e in self._edges)
 
-    def _tag(self: T, val: Sequence[Union[Shape, Location]], tag: str):
+    def _tag(self: T, val: Sequence[SketchVal], tag: str):
 
         self._tags[tag] = val
 
     # face construction
     def face(
         self: T,
-        b: Union[Wire, Iterable[Edge], Compound, T],
+        b: Union[Wire, Iterable[Edge], Shape, T],
         angle: Real = 0,
         mode: Modes = "a",
         tag: Optional[str] = None,
@@ -152,9 +203,11 @@ def face(
 
         if isinstance(b, Wire):
             res = Face.makeFromWires(b)
-        elif isinstance(b, (Sketch, Compound)):
+        elif isinstance(b, Sketch):
             res = b
-        elif isinstance(b, Iterable) and not isinstance(b, Shape):
+        elif isinstance(b, Shape):
+            res = compound(b.Faces())
+        elif isinstance(b, Iterable):
             wires = edgesToWires(tcast(Iterable[Edge], b))
             res = Face.makeFromWires(*(wires[0], wires[1:]))
         else:
@@ -488,6 +541,8 @@ def each(
             self._faces = self._faces.cut(*res)
         elif mode == "i":
             self._faces = self._faces.intersect(*res)
+        elif mode == "r":
+            self._faces = compound(res)
         elif mode == "c":
             if not tag:
                 raise ValueError("No tag specified - the geometry will be unreachable")
@@ -506,10 +561,8 @@ def hull(self: T, mode: Modes = "a", tag: Optional[str] = None) -> T:
             rv = find_hull(el for el in self._selection if isinstance(el, Edge))
         elif self._faces:
             rv = find_hull(el for el in self._faces.Edges())
-        elif self._edges or self._wires:
-            rv = find_hull(
-                chain(self._edges, chain.from_iterable(w.Edges() for w in self._wires))
-            )
+        elif self._edges:
+            rv = find_hull(self._edges)
         else:
             raise ValueError("No objects available for hull construction")
 
@@ -522,10 +575,15 @@ def offset(self: T, d: Real, mode: Modes = "a", tag: Optional[str] = None) -> T:
         Offset selected wires or edges.
         """
 
-        rv = (el.offset2D(d) for el in self._selection if isinstance(el, Wire))
+        if self._selection:
+            rv = (el.offset2D(d) for el in self._selection if isinstance(el, Wire))
 
-        for el in chain.from_iterable(rv):
-            self.face(el, mode=mode, tag=tag, ignore_selection=bool(self._selection))
+            for el in chain.from_iterable(rv):
+                self.face(
+                    el, mode=mode, tag=tag, ignore_selection=bool(self._selection)
+                )
+        else:
+            raise ValueError("Selection is needed to offset")
 
         return self
 
@@ -533,12 +591,17 @@ def _matchFacesToVertices(self) -> Dict[Face, List[Vertex]]:
 
         rv = {}
 
-        for f in self._faces.Faces():
-
-            f_vertices = f.Vertices()
-            rv[f] = [
-                v for v in self._selection if isinstance(v, Vertex) and v in f_vertices
-            ]
+        if self._selection:
+            for f in self._faces.Faces():
+
+                f_vertices = f.Vertices()
+                rv[f] = [
+                    v
+                    for v in self._selection
+                    if isinstance(v, Vertex) and v in f_vertices
+                ]
+        else:
+            raise ValueError("Selection is needed to match vertices to faces")
 
         return rv
 
@@ -622,7 +685,10 @@ def tag(self: T, tag: str) -> T:
         Tag current selection.
         """
 
-        self._tags[tag] = list(self._selection)
+        if self._selection:
+            self._tags[tag] = list(self._selection)
+        else:
+            raise ValueError("Selection is needed to tag")
 
         return self
 
@@ -679,7 +745,7 @@ def reset(self: T) -> T:
         Reset current selection.
         """
 
-        self._selection = []
+        self._selection = None
         return self
 
     def delete(self: T) -> T:
@@ -687,15 +753,18 @@ def delete(self: T) -> T:
         Delete selected object.
         """
 
-        for obj in self._selection:
-            if isinstance(obj, Face):
-                self._faces.remove(obj)
-            elif isinstance(obj, Wire):
-                self._wires.remove(obj)
-            elif isinstance(obj, Edge):
-                self._edges.remove(obj)
+        if self._selection:
+            for obj in self._selection:
+                if isinstance(obj, Face):
+                    self._faces.remove(obj)
+                elif isinstance(obj, Edge):
+                    self._edges.remove(obj)
+                else:
+                    raise ValueError(f"Deletion of {obj} not supported")
+        else:
+            raise ValueError("Selection is needed to delete")
 
-        self._selection = []
+        self.reset()
 
         return self
 
@@ -871,7 +940,7 @@ def bezier(
         The edge will pass through the last points, and the inner points
         are bezier control points.
         """
-        p1 = self._endPoint()
+
         val = Edge.makeBezier([Vector(*p) for p in pts])
 
         return self.edge(val, tag, forConstruction)
@@ -1011,13 +1080,49 @@ def copy(self: T) -> T:
 
         return rv
 
+    @overload
     def moved(self: T, loc: Location) -> T:
+        ...
+
+    @overload
+    def moved(self: T, loc1: Location, loc2: Location, *locs: Location) -> T:
+        ...
+
+    @overload
+    def moved(self: T, locs: Sequence[Location]) -> T:
+        ...
+
+    @overload
+    def moved(
+        self: T,
+        x: Real = 0,
+        y: Real = 0,
+        z: Real = 0,
+        rx: Real = 0,
+        ry: Real = 0,
+        rz: Real = 0,
+    ) -> T:
+        ...
+
+    @overload
+    def moved(self: T, loc: VectorLike) -> T:
+        ...
+
+    @overload
+    def moved(self: T, loc1: VectorLike, loc2: VectorLike, *locs: VectorLike) -> T:
+        ...
+
+    @overload
+    def moved(self: T, loc: Sequence[VectorLike]) -> T:
+        ...
+
+    def moved(self: T, *args, **kwargs) -> T:
         """
         Create a partial copy of the sketch with moved _faces.
         """
 
         rv = self.__class__()
-        rv._faces = self._faces.moved(loc)
+        rv._faces = self._faces.moved(*args, **kwargs)
 
         return rv
 
@@ -1040,14 +1145,203 @@ def finalize(self) -> Any:
 
     def val(self: T) -> SketchVal:
         """
-        Return the first selected item or Location().
+        Return the first selected item, underlying compound or first edge.
         """
 
-        return self._selection[0] if self._selection else Location()
+        if self._selection is not None:
+            rv = self._selection[0]
+        elif not self._faces and self._edges:
+            rv = self._edges[0]
+        else:
+            rv = self._faces
+
+        return rv
 
     def vals(self: T) -> List[SketchVal]:
         """
-        Return the list of selected items.
+        Return all selected items, underlying compound or all edges.
+        """
+
+        rv: List[SketchVal]
+
+        if self._selection is not None:
+            rv = list(self._selection)
+        elif not self._faces and self._edges:
+            rv = list(self._edges)
+        else:
+            rv = list(self._faces)
+
+        return rv
+
+    def add(self: T) -> T:
+        """
+        Add selection to the underlying faces.
+        """
+
+        self._faces += compound(self._selection).faces()
+
+        return self
+
+    def subtract(self: T) -> T:
         """
+        Subtract selection from the underlying faces.
+        """
+
+        self._faces -= compound(self._selection).faces()
 
-        return self._selection
+        return self
+
+    def replace(self: T) -> T:
+        """
+        Replace the underlying faces with the selection.
+        """
+
+        self._faces = compound(self._selection).faces()
+
+        return self
+
+    def __add__(self: T, other: "Sketch") -> T:
+        """
+        Fuse self and other.
+        """
+
+        res = _sanitize_for_bool(self.val()) + _sanitize_for_bool(other.val())
+
+        return self.__class__(obj=_to_compound(res))
+
+    def __sub__(self: T, other: "Sketch") -> T:
+        """
+        Subtract other from self.
+        """
+
+        res = _sanitize_for_bool(self.val()) - _sanitize_for_bool(other.val())
+
+        return self.__class__(obj=_to_compound(res))
+
+    def __mul__(self: T, other: "Sketch") -> T:
+        """
+        Intersect self and other.
+        """
+
+        res = _sanitize_for_bool(self.val()) * _sanitize_for_bool(other.val())
+
+        return self.__class__(obj=_to_compound(res))
+
+    def __truediv__(self: T, other: "Sketch") -> T:
+        """
+        Split self with other.
+        """
+
+        res = _sanitize_for_bool(self.val()) / _sanitize_for_bool(other.val())
+
+        return self.__class__(obj=_to_compound(res))
+
+    def __getitem__(self: T, item: Union[int, Sequence[int], slice]) -> T:
+
+        vals = self.vals()
+
+        if isinstance(item, Iterable):
+            self._selection = [vals[i] for i in item]
+        elif isinstance(item, slice):
+            self._selection = vals[item]
+        else:
+            self._selection = [vals[item]]
+
+        return self
+
+    def filter(self: T, f: Callable[[SketchVal], bool]) -> T:
+        """
+        Filter items using a boolean predicate.
+        
+        :param f: Callable to be used for filtering.
+        :return: Sketch object with filtered items.
+        """
+
+        self._selection = list(filter(f, self.vals()))
+
+        return self
+
+    def map(self: T, f: Callable[[SketchVal], SketchVal]):
+        """
+        Apply a callable to every item separately.
+        
+        :param f: Callable to be applied to every item separately.
+        :return: Sketch object with f applied to all items.
+        """
+
+        self._selection = list(map(f, self.vals()))
+
+        return self
+
+    def apply(self: T, f: Callable[[Iterable[SketchVal]], Iterable[SketchVal]]):
+        """
+        Apply a callable to all items at once.
+        
+        :param f: Callable to be applied.
+        :return: Sketch object with f applied to all items.
+        """
+
+        self._selection = list(f(self.vals()))
+
+        return self
+
+    def sort(self: T, key: Callable[[SketchVal], Any]) -> T:
+        """
+        Sort items using a callable.
+        
+        :param key: Callable to be used for sorting.
+        :return: Sketch object with items sorted.
+        """
+
+        self._selection = list(sorted(self.vals(), key=key))
+
+        return self
+
+    def invoke(
+        self: T, f: Union[Callable[[T], T], Callable[[T], None], Callable[[], None]]
+    ):
+        """
+        Invoke a callable mapping Sketch to Sketch or None. Supports also
+        callables that take no arguments such as breakpoint. Returns self if callable
+        returns None.
+        
+        :param f: Callable to be invoked.
+        :return: Sketch object.
+        """
+
+        arity = get_arity(f)
+        rv = self
+
+        if arity == 0:
+            f()  # type: ignore
+        elif arity == 1:
+            res = f(self)  # type: ignore
+            if res is not None:
+                rv = res
+        else:
+            raise ValueError("Provided function {f} accepts too many arguments")
+
+        return rv
+
+    def export(
+        self: T,
+        fname: str,
+        tolerance: float = 0.1,
+        angularTolerance: float = 0.1,
+        opt: Optional[Dict[str, Any]] = None,
+    ) -> T:
+        """
+        Export Sketch to file.
+        
+        :param path: Filename.
+        :param tolerance: the deflection tolerance, in model units. Default 0.1.
+        :param angularTolerance: the angular tolerance, in radians. Default 0.1.
+        :param opt: additional options passed to the specific exporter. Default None.
+        :return: Self.
+        """
+
+        export(
+            self, fname, tolerance=tolerance, angularTolerance=angularTolerance, opt=opt
+        )
+
+        return self
diff --git a/cadquery/utils.py b/cadquery/utils.py
index b72fe2ed7..efa9e494f 100644
--- a/cadquery/utils.py
+++ b/cadquery/utils.py
@@ -1,5 +1,5 @@
 from functools import wraps
-from inspect import signature
+from inspect import signature, isbuiltin
 from typing import TypeVar, Callable, cast
 from warnings import warn
 
@@ -71,3 +71,15 @@ def wrapped(*args, **kwargs):
             return f(*args, **kwargs)
 
         return wrapped
+
+
+def get_arity(f: TCallable) -> int:
+
+    if isbuiltin(f):
+        rv = 0  # assume 0 arity for builtins; they cannot be introspected
+    else:
+        # NB: this is not understood by mypy
+        n_defaults = len(f.__defaults__) if f.__defaults__ else 0
+        rv = f.__code__.co_argcount - n_defaults
+
+    return rv
diff --git a/cadquery/vis.py b/cadquery/vis.py
index 5e5019f86..352b36e24 100644
--- a/cadquery/vis.py
+++ b/cadquery/vis.py
@@ -21,9 +21,7 @@ def _to_assy(*objs: Union[Shape, Workplane, Assembly, Sketch]) -> Assembly:
         if isinstance(obj, (Shape, Workplane, Assembly)):
             assy.add(obj)
         elif isinstance(obj, Sketch):
-            assy.add(obj._faces)
-            assy.add(Compound.makeCompound(obj._edges))
-            assy.add(Compound.makeCompound(obj._wires))
+            assy.add(Compound.makeCompound(obj))
         elif isinstance(obj, TopoDS_Shape):
             assy.add(Shape(obj))
         else:
diff --git a/doc/extending.rst b/doc/extending.rst
index a0fb0cd18..b43a69609 100644
--- a/doc/extending.rst
+++ b/doc/extending.rst
@@ -200,12 +200,15 @@ Extending CadQuery: Special Methods
 -----------------------------------
 
 The above-mentioned approach has one drawback, it requires monkey-patching or subclassing. To avoid this
-one can also use the following special methods of :py:class:`cadquery.Workplane`
+one can also use the following special methods of :py:class:`cadquery.Workplane` and :py:class:`cadquery.Sketch`
 and write plugins in a more functional style.
 
     * :py:meth:`cadquery.Workplane.map`
     * :py:meth:`cadquery.Workplane.apply`
     * :py:meth:`cadquery.Workplane.invoke`
+    * :py:meth:`cadquery.Sketch.map`
+    * :py:meth:`cadquery.Sketch.apply`
+    * :py:meth:`cadquery.Sketch.invoke`
 
 Here is the same plugin rewritten using one of those methods.
 
diff --git a/doc/importexport.rst b/doc/importexport.rst
index b035e22bf..66133de65 100644
--- a/doc/importexport.rst
+++ b/doc/importexport.rst
@@ -8,10 +8,10 @@ Importing and Exporting Files
 Introduction
 #############
 
-The purpose of this section is to explain how to import external file formats into CadQuery, and export files from 
-it as well. While the external file formats can be used to interchange CAD model data with other software, CadQuery 
-does not support any formats that carry parametric data with them at this time. The only format that is fully 
-parametric is CadQuery's own Python format. Below are lists of the import and export file formats that CadQuery 
+The purpose of this section is to explain how to import external file formats into CadQuery, and export files from
+it as well. While the external file formats can be used to interchange CAD model data with other software, CadQuery
+does not support any formats that carry parametric data with them at this time. The only format that is fully
+parametric is CadQuery's own Python format. Below are lists of the import and export file formats that CadQuery
 supports.
 
 Import Formats
@@ -64,7 +64,7 @@ Importing a DXF profile with default settings and using it within a CadQuery scr
     )
 
 Note the use of the :meth:`Workplane.wires` and :meth:`Workplane.toPending` methods to make the DXF profile
-ready for use during subsequent operations. Calling ``toPending()`` tells CadQuery to make the edges/wires available 
+ready for use during subsequent operations. Calling ``toPending()`` tells CadQuery to make the edges/wires available
 to the next modelling operation that is called in the chain.
 
 Importing STEP
@@ -96,7 +96,7 @@ since it will be determined from the file extension. Below is an example.
     box = cq.Workplane().box(10, 10, 10)
 
     # Export the box
-    cq.exporters.export(box, "/path/to/step/box.step")
+    box.export("/path/to/step/box.step")
 
 Non-Default File Extensions
 ----------------------------
@@ -110,10 +110,10 @@ not recognize the file extension. In that case the export type has to be specifi
     box = cq.Workplane().box(10, 10, 10)
 
     # Export the box
-    cq.exporters.export(box, "/path/to/step/box.stp", cq.exporters.ExportTypes.STEP)
+    box.export("/path/to/step/box.stp", cq.exporters.ExportTypes.STEP)
 
     # The export type may also be specified as a literal
-    cq.exporters.export(box, "/path/to/step/box2.stp", "STEP")
+    box.export("/path/to/step/box2.stp", "STEP")
 
 Setting Extra Options
 ----------------------
@@ -128,10 +128,10 @@ or the :meth:`Assembly.exportAssembly`` method.
     box = cq.Workplane().box(10, 10, 10)
 
     # Export the box, provide additional options with the opt dict
-    cq.exporters.export(box, "/path/to/step/box.step", opt={"write_pcurves": False})
+    box.export("/path/to/step/box.step", opt={"write_pcurves": False})
 
     # or equivalently when exporting a lower level Shape object
-    box.val().exportStep("/path/to/step/box2.step", write_pcurves=False)
+    box.val().export("/path/to/step/box2.step", opt={"write_pcurves": False})
 
 
 Exporting Assemblies to STEP
@@ -159,7 +159,7 @@ export with all defaults is shown below.
     assy.add(pin, color=cq.Color(0, 1, 0), name="pin")
 
     # Save the assembly to STEP
-    assy.save("out.step")
+    assy.export("out.step")
 
 This will produce a STEP file that is nested with auto-generated object names. The colors of each assembly object will be
 preserved, but the names that were set for each will not.
@@ -183,10 +183,10 @@ fused solids.
     assy.add(pin, color=cq.Color(0, 1, 0), name="pin")
 
     # Save the assembly to STEP
-    assy.save("out.stp", "STEP", mode="fused")
+    assy.export("out.stp", "STEP", mode="fused")
 
     # Specify additional options such as glue as keyword arguments
-    assy.save("out_glue.step", mode="fused", glue=True, write_pcurves=False)
+    assy.export("out_glue.step", mode="fused", glue=True, write_pcurves=False)
 
 Naming
 -------
@@ -197,7 +197,7 @@ This is done by setting the name property of the assembly before calling :meth:`
 .. code-block:: python
 
     assy = Assembly(name="my_assembly")
-    assy.save(
+    assy.export(
         "out.stp",
         cq.exporters.ExportTypes.STEP,
         mode=cq.exporters.assembly.ExportModes.FUSED,
@@ -224,13 +224,13 @@ export with all defaults is shown below. To export to a binary glTF file, change
     pin = cq.Workplane().center(2, 2).cylinder(radius=2, height=20)
     assy.add(pin, color=cq.Color(0, 1, 0), name="pin")
 
-    # Save the assembly to STEP
-    assy.save("out.gltf")
+    # Save the assembly to GLTF
+    assy.export("out.gltf")
 
 Exporting SVG
 ###############
 
-The SVG exporter has several options which can be useful for achieving the desired final output. Those 
+The SVG exporter has several options which can be useful for achieving the desired final output. Those
 options are as follows.
 
 * *width* - Width of the resulting image (None to fit based on height).
@@ -245,7 +245,7 @@ options are as follows.
 * *showHidden* - Whether or not to show hidden lines.
 * *focus* - If specified, creates a perspective SVG with the projector at the distance specified.
 
-The options are passed to the exporter in a dictionary, and can be left out to force the SVG to be created with default options. 
+The options are passed to the exporter in a dictionary, and can be left out to force the SVG to be created with default options.
 Below are examples with and without options set.
 
 Without options:
@@ -257,13 +257,13 @@ Without options:
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(result, "/path/to/file/box.svg")
+    result.export("/path/to/file/box.svg")
 
 Results in:
 
 ..  image:: _static/importexport/box_default_options.svg
 
-Note that the exporters API figured out the format type from the file extension. The format 
+Note that the exporters API figured out the format type from the file extension. The format
 type can be set explicitly by using :py:class:`exporters.ExportTypes`.
 
 The following is an example of using options to alter the resulting SVG output by passing in the ``opt`` parameter.
@@ -275,8 +275,7 @@ The following is an example of using options to alter the resulting SVG output b
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(
-        result,
+    result.export(
         "/path/to/file/box_custom_options.svg",
         opt={
             "width": 300,
@@ -308,7 +307,7 @@ The STL exporter is capable of adjusting the quality of the resulting mesh, and
 .. automethod::
     cadquery.occ_impl.shapes.Shape.exportStl
 
-For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the 
+For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the
 optimum values that will produce an acceptable mesh.
 
 .. code-block:: python
@@ -318,7 +317,7 @@ optimum values that will produce an acceptable mesh.
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(result, "/path/to/file/mesh.stl")
+    result.export("/path/to/file/mesh.stl")
 
 Exporting AMF and 3MF
 ######################
@@ -329,7 +328,7 @@ The AMF and 3MF exporters are capable of adjusting the quality of the resulting
 * ``tolerance`` - A linear deflection setting which limits the distance between a curve and its tessellation. Setting this value too low will result in large meshes that can consume computing resources. Setting the value too high can result in meshes with a level of detail that is too low. Default is 0.1, which is good starting point for a range of cases.
 * ``angularTolerance`` - Angular deflection setting which limits the angle between subsequent segments in a polyline. Default is 0.1.
 
-For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the 
+For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the
 optimum values that will produce an acceptable mesh. Note that parameters for color and material are absent.
 
 .. code-block:: python
@@ -339,7 +338,7 @@ optimum values that will produce an acceptable mesh. Note that parameters for co
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(result, "/path/to/file/mesh.amf", tolerance=0.01, angularTolerance=0.1)
+    result.export("/path/to/file/mesh.amf", tolerance=0.01, angularTolerance=0.1)
 
 
 Exporting TJS
@@ -351,7 +350,7 @@ The TJS (ThreeJS) exporter produces a file in JSON format that describes a scene
 * ``tolerance`` - A linear deflection setting which limits the distance between a curve and its tessellation. Setting this value too low will result in large meshes that can consume computing resources. Setting the value too high can result in meshes with a level of detail that is too low. Default is 0.1, which is good starting point for a range of cases.
 * ``angularTolerance`` - Angular deflection setting which limits the angle between subsequent segments in a polyline. Default is 0.1.
 
-For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the 
+For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the
 optimum values that will produce an acceptable mesh.
 
 .. code-block:: python
@@ -361,15 +360,14 @@ optimum values that will produce an acceptable mesh.
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(
-        result,
+    result.export(
         "/path/to/file/mesh.json",
         tolerance=0.01,
         angularTolerance=0.1,
         exportType=exporters.ExportTypes.TJS,
     )
 
-Note that the export type was explicitly specified as ``TJS`` because the extension that was used for the file name was ``.json``. If the extension ``.tjs`` 
+Note that the export type was explicitly specified as ``TJS`` because the extension that was used for the file name was ``.json``. If the extension ``.tjs``
 had been used, CadQuery would have understood to use the ``TJS`` export format.
 
 Exporting VRML
@@ -381,7 +379,7 @@ The VRML exporter is capable of adjusting the quality of the resulting mesh, and
 * ``tolerance`` - A linear deflection setting which limits the distance between a curve and its tessellation. Setting this value too low will result in large meshes that can consume computing resources. Setting the value too high can result in meshes with a level of detail that is too low. Default is 0.1, which is good starting point for a range of cases.
 * ``angularTolerance`` - Angular deflection setting which limits the angle between subsequent segments in a polyline. Default is 0.1.
 
-For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the 
+For more complex objects, some experimentation with ``tolerance`` and ``angularTolerance`` may be required to find the
 optimum values that will produce an acceptable mesh.
 
 .. code-block:: python
@@ -391,8 +389,8 @@ optimum values that will produce an acceptable mesh.
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(
-        result, "/path/to/file/mesh.vrml", tolerance=0.01, angularTolerance=0.1
+    result.export(
+        "/path/to/file/mesh.vrml", tolerance=0.01, angularTolerance=0.1
     )
 
 Exporting DXF
@@ -423,16 +421,26 @@ Options
     See `Units`_.
 
 .. code-block:: python
-    :caption: DXF document without options.
+    :caption: DXF of workplanes.
 
     import cadquery as cq
-    from cadquery import exporters
 
-    result = cq.Workplane().box(10, 10, 10)
+    result = cq.Workplane().box(10, 10, 10).section()
 
     exporters.exportDXF(result, "/path/to/file/object.dxf")
     # or
-    exporters.export(result, "/path/to/file/object.dxf")
+    result.export("/path/to/file/object.dxf")
+
+Sketches can also be directly exported to DXF.
+
+.. code-block:: python
+    :caption: DXF export of sketches.
+
+    import cadquery as cq
+
+    result = cq.Sketch().rect(1,1)
+
+    result.export("/path/to/file/object.dxf")
 
 
 Units
@@ -460,7 +468,7 @@ Document units can be set to any :doc:`unit supported by ezdxf <ezdxf-stable:con
     import cadquery as cq
     from cadquery import exporters
 
-    result = cq.Workplane().box(10, 10, 10)
+    result = cq.Workplane().box(10, 10, 10).section()
 
     exporters.exportDXF(
         result,
@@ -470,8 +478,7 @@ Document units can be set to any :doc:`unit supported by ezdxf <ezdxf-stable:con
 
     # or
 
-    exporters.export(
-        result,
+    result.export(
         "/path/to/file/object.dxf",
         opt={"doc_units": 6},  # set DXF document units to meters
     )
@@ -496,7 +503,7 @@ By default, the DXF exporter will output splines exactly as they are represented
 Exporting Other Formats
 ########################
 
-The remaining export formats do not accept any additional parameters other than file name, and can be exported 
+The remaining export formats do not accept any additional parameters other than file name, and can be exported
 using the following structure.
 
 .. code-block:: python
@@ -506,9 +513,9 @@ using the following structure.
 
     result = cq.Workplane().box(10, 10, 10)
 
-    exporters.export(result, "/path/to/file/object.[file_extension]")
+    result.export("/path/to/file/object.[file_extension]")
 
-Be sure to use the correct file extension so that CadQuery can determine the export format. If in doubt, fall 
+Be sure to use the correct file extension so that CadQuery can determine the export format. If in doubt, fall
 back to setting the type explicitly by using :py:class:`exporters.ExportTypes`.
 
 For example:
@@ -518,6 +525,6 @@ For example:
     import cadquery as cq
     from cadquery import exporters
 
-    result = cq.Workplane().box(10, 10, 10)
+    result = cq.Workplane().box(10, 10, 10).section()
 
-    exporters.export(result, "/path/to/file/object.dxf", exporters.ExportTypes.DXF)
+    result.export("/path/to/file/object.dxf", exporters.ExportTypes.DXF)
diff --git a/doc/selectors.rst b/doc/selectors.rst
index 9670526bd..323826cd0 100644
--- a/doc/selectors.rst
+++ b/doc/selectors.rst
@@ -191,13 +191,16 @@ objects. This includes chaining and combining.
 Additional special methods
 --------------------------
 
-:py:class:`cadquery.Workplane` provides the following special methods that can be used
+:py:class:`cadquery.Workplane` and :py:class:`cadquery.Sketch` provide the following special methods that can be used
 for quick prototyping of selectors when implementing a complete selector via subclassing of
 :py:class:`cadquery.Selector` is not desirable.
 
     * :py:meth:`cadquery.Workplane.filter`
     * :py:meth:`cadquery.Workplane.sort`
     * :py:meth:`cadquery.Workplane.__getitem__`
+    * :py:meth:`cadquery.Sketch.filter`
+    * :py:meth:`cadquery.Sketch.sort`
+    * :py:meth:`cadquery.Sketch.__getitem__`
 
 For example, one could use those methods for selecting objects within a certain range of volumes.
 
diff --git a/doc/sketch.rst b/doc/sketch.rst
index 988ee4a24..6dcc84b7f 100644
--- a/doc/sketch.rst
+++ b/doc/sketch.rst
@@ -13,7 +13,7 @@ approaches.
 Face-based API
 ==============
 
-The main approach for constructing sketches is based on constructing faces and 
+The main approach for constructing sketches is based on constructing faces and
 combining them using boolean operations.
 
 .. cadquery::
@@ -40,7 +40,14 @@ class does not implement history and all modifications happen in-place.
 Modes
 ^^^^^
 
-Every operation from the face API accepts a mode parameter to define how to combine the created object with existing ones. It can be fused (``mode='a'``), cut (``mode='s'``), intersected (``mode='i'``) or just stored for construction (``mode='c'``). In the last case, it is mandatory to specify a ``tag`` in order to be able to refer to the object later on. By default faces are fused together. Note the usage of the subtractive and additive modes in the example above. The additional two are demonstrated below.
+Every operation from the face API accepts a mode parameter to define
+how to combine the created object with existing ones. It can be fused (``mode='a'``),
+cut (``mode='s'``), intersected (``mode='i'``), replaced (``mode='r'``)
+or just stored for construction (``mode='c'``).
+In the last case, it is mandatory to specify a ``tag`` in order to be able to
+refer to the object later on. By default faces are fused together.
+Note the usage of the subtractive and additive modes in the example above.
+The additional two are demonstrated below.
 
 .. cadquery::
     :height: 600px
@@ -285,10 +292,42 @@ Two sketches on different workplanes are needed when using :meth:`~cadquery.Work
 
     s2 = Sketch().rect(2, 1).vertices().fillet(0.2)
 
-    result = Workplane().placeSketch(s1, s2.moved(Location(Vector(0, 0, 3)))).loft()
+    result = Workplane().placeSketch(s1, s2.moved(z=3)).loft()
 
 When lofting only outer wires are taken into account and inner wires are silently ignored. Note that only sketches on the top of stack are considered for the current operation (i.e. there are no pending sketches), so when lofting or sweeping all relevant sketches have to be added in one `placeSketch` call.
 
+
+Combining sketches
+==================
+
+Sketches can be combined using :meth:`~cadquery.Sketch.face`.
+
+.. cadquery::
+   :height: 600px
+
+   import cadquery as cq
+
+   s1 = cq.Sketch().rect(2, 2)
+   s2 = cq.Sketch().circle(0.5)
+
+   result = s1.face(s2, mode='s')
+
+
+It is also possible to use boolean operations to achieve the same effect.
+
+.. cadquery::
+   :height: 600px
+
+   import cadquery as cq
+
+   s1 = cq.Sketch().rect(2, 2).vertices().fillet(0.25).reset()
+   s2 = cq.Sketch().rect(1, 1, angle=45).vertices().chamfer(0.1).reset()
+
+   result = s1 - s2
+
+Boolean operations are selection sensitive, so in this example
+:meth:`~cadquery.Sketch.reset` call is needed.
+
 Offsets made easy
 =================
 
@@ -309,3 +348,31 @@ Conveniently, it is possible to reuse a sketch to create an :meth:`~cadquery.Ske
 
    result = cq.Workplane("front").placeSketch(sketch_offset).extrude(1.0)
    result = result.faces(">Z").workplane().placeSketch(sketch).cutBlind(-0.50)
+
+
+It is obviously possible to use negative offsets, but it requires being more careful with the mode
+of the offset operation. Usually one wants to replace the original face, hence ``mode='r'``.
+
+.. cadquery::
+   :height: 600px
+
+   import cadquery as cq
+
+   sketch  = (cq.Sketch()
+   .rect(1.0, 4.0)
+   .circle(1.0)
+   .clean()
+   )
+
+   sketch_offset = sketch.copy().wires().offset(-0.25, mode='r')
+
+   result = cq.Workplane("front").placeSketch(sketch).extrude(1.0)
+   result = result.faces(">Z").workplane().placeSketch(sketch_offset).cutBlind(-0.50)
+
+
+Exporting and importing
+=======================
+
+It is possible to export sketches using :meth:`~cadquery.Sketch.export`.
+See :ref:`importexport` for more details.
+Importing of DXF files is supported as well using :meth:`~cadquery.Sketch.importDXF`.
diff --git a/tests/test_assembly.py b/tests/test_assembly.py
index 4f4829f2d..e1754d861 100644
--- a/tests/test_assembly.py
+++ b/tests/test_assembly.py
@@ -37,7 +37,7 @@
 from OCP.TopAbs import TopAbs_ShapeEnum
 
 
-@pytest.fixture(scope="module")
+@pytest.fixture(scope="function")
 def tmpdir(tmp_path_factory):
     return tmp_path_factory.mktemp("assembly")
 
@@ -689,6 +689,50 @@ def test_save(extension, args, nested_assy, nested_assy_sphere):
     assert os.path.exists(filename)
 
 
+@pytest.mark.parametrize(
+    "extension, args, kwargs",
+    [
+        ("step", (), {}),
+        ("xml", (), {}),
+        ("vrml", (), {}),
+        ("gltf", (), {}),
+        ("glb", (), {}),
+        ("stl", (), {"ascii": False}),
+        ("stl", (), {"ascii": True}),
+        ("stp", ("STEP",), {}),
+        ("caf", ("XML",), {}),
+        ("wrl", ("VRML",), {}),
+        ("stl", ("STL",), {}),
+    ],
+)
+def test_export(extension, args, kwargs, tmpdir, nested_assy):
+
+    filename = "nested." + extension
+
+    with tmpdir:
+        nested_assy.export(filename, *args, **kwargs)
+        assert os.path.exists(filename)
+
+
+def test_export_vtkjs(tmpdir, nested_assy):
+
+    with tmpdir:
+        nested_assy.export("nested.vtkjs")
+        assert os.path.exists("nested.vtkjs.zip")
+
+
+def test_export_errors(nested_assy):
+
+    with pytest.raises(ValueError):
+        nested_assy.export("nested.1234")
+
+    with pytest.raises(ValueError):
+        nested_assy.export("nested.stl", "1234")
+
+    with pytest.raises(ValueError):
+        nested_assy.export("nested.step", mode="1234")
+
+
 def test_save_stl_formats(nested_assy_sphere):
 
     # Binary export
@@ -1368,7 +1412,6 @@ def resulting_plane(shape0):
     nonplanar_spline = cq.Edge.makeSpline(points1, periodic=True)
     fail_this(nonplanar_spline)
 
-    # planar wire should succeed
     # make a triangle in the XZ plane
     points2 = [cq.Vector(x) for x in [(-1, 0, -1), (0, 0, 1), (1, 0, -1)]]
     points2.append(points2[0])
diff --git a/tests/test_cadquery.py b/tests/test_cadquery.py
index 0437b752c..57726aa3e 100644
--- a/tests/test_cadquery.py
+++ b/tests/test_cadquery.py
@@ -15,6 +15,7 @@
 
 from cadquery import *
 from cadquery import occ_impl
+from cadquery.occ_impl.shapes import *
 from tests import (
     BaseTest,
     writeStringToFile,
@@ -5756,7 +5757,7 @@ def test_map_apply_filter_sort(self):
 
     def test_getitem(self):
 
-        w = Workplane().rarray(2, 1, 5, 1).box(1, 1, 1, combine=False)
+        w = Workplane().rarray(2, 0, 5, 1).box(1, 1, 1, combine=False)
 
         assert w[0].solids().size() == 1
         assert w[-2:].solids().size() == 2
@@ -5764,7 +5765,7 @@ def test_getitem(self):
 
     def test_invoke(self):
 
-        w = Workplane().rarray(2, 1, 5, 1).box(1, 1, 1, combine=False)
+        w = Workplane().rarray(2, 0, 5, 1).box(1, 1, 1, combine=False)
 
         # builtin
         assert w.invoke(print).size() == 5
@@ -5790,3 +5791,41 @@ def test_tessellate(self):
         verts, _ = Face.makePlane(1e-9, 1e-9).tessellate(1e-3)
 
         assert len(verts) == 0
+
+    def test_export(self):
+
+        w = Workplane().box(1, 1, 1).export("box.brep")
+
+        assert (w - Shape.importBrep("box.brep")).val().Volume() == approx(0)
+
+    def test_bool_operators(self):
+
+        w1 = Workplane().box(1, 1, 2)
+        w2 = Workplane().box(2, 2, 1)
+
+        assert (w1 + w2).val().Volume() == approx(5)
+        assert (w1 - w2).val().Volume() == approx(1)
+        assert (w1 * w2).val().Volume() == approx(1)
+        assert (w1 / w2).solids().size() == 3
+
+    def test_extrude_face(self):
+
+        f = face(rect(1, 1))
+        c = compound(f)
+
+        # face
+        assert Workplane().add(f).extrude(1).val().Volume() == approx(1)
+        # compound with face
+        assert Workplane().add(c).extrude(1).val().Volume() == approx(1)
+
+    def test_workplane_iter(self):
+
+        s = Workplane().sketch().rarray(5, 0, 5, 1).rect(1, 1).finalize()
+        w1 = Workplane().pushPoints([(-10, 0), (10, 0)])
+        w2 = w1.box(1, 1, 1)  # NB this results in Compound of two Solids
+        w3 = w1.box(1, 1, 1, combine=False)
+
+        assert len(list(s)) == 5
+        assert len(list(w1)) == 0
+        assert len(list(w2)) == 2  # 2 beacuase __iter__ unpacks Compounds
+        assert len(list(w3)) == 2
diff --git a/tests/test_exporters.py b/tests/test_exporters.py
index f75d0e202..40fc54235 100644
--- a/tests/test_exporters.py
+++ b/tests/test_exporters.py
@@ -11,6 +11,8 @@
 import pytest
 import ezdxf
 
+from uuid import uuid1
+
 from pytest import approx
 
 # my modules
@@ -27,8 +29,9 @@
     Vector,
     Color,
 )
+
+from cadquery.occ_impl.shapes import rect, face, compound
 from cadquery.occ_impl.exporters.dxf import DxfDocument
-from cadquery.occ_impl.exporters.utils import toCompound
 from tests import BaseTest
 from OCP.GeomConvert import GeomConvert
 from OCP.BRepBuilderAPI import BRepBuilderAPI_MakeEdge
@@ -49,6 +52,12 @@ def box123():
     return Workplane().box(1, 2, 3)
 
 
+@pytest.fixture()
+def fname(tmpdir):
+
+    return tmpdir / str(uuid1())
+
+
 def test_step_options(tmpdir):
     """
     Exports a box using the options to decrease STEP file size and
@@ -298,7 +307,7 @@ def test_line(self):
         workplane = Workplane().line(1, 1)
 
         plane = workplane.plane
-        shape = toCompound(workplane).transformShape(plane.fG)
+        shape = compound(*workplane).transformShape(plane.fG)
         edges = shape.Edges()
 
         result = DxfDocument._dxf_line(edges[0])
@@ -311,7 +320,7 @@ def test_circle(self):
         workplane = Workplane().circle(1)
 
         plane = workplane.plane
-        shape = toCompound(workplane).transformShape(plane.fG)
+        shape = compound(*workplane).transformShape(plane.fG)
         edges = shape.Edges()
 
         result = DxfDocument._dxf_circle(edges[0])
@@ -324,7 +333,7 @@ def test_arc(self):
         workplane = Workplane().radiusArc((1, 1), 1)
 
         plane = workplane.plane
-        shape = toCompound(workplane).transformShape(plane.fG)
+        shape = compound(*workplane).transformShape(plane.fG)
         edges = shape.Edges()
 
         result_type, result_attributes = DxfDocument._dxf_circle(edges[0])
@@ -352,7 +361,7 @@ def test_ellipse(self):
         workplane = Workplane().ellipse(2, 1, 0)
 
         plane = workplane.plane
-        shape = toCompound(workplane).transformShape(plane.fG)
+        shape = compound(*workplane).transformShape(plane.fG)
         edges = shape.Edges()
 
         result_type, result_attributes = DxfDocument._dxf_ellipse(edges[0])
@@ -388,7 +397,7 @@ def test_spline(self):
         )
 
         plane = workplane.plane
-        shape = toCompound(workplane).transformShape(plane.fG)
+        shape = compound(*workplane).transformShape(plane.fG)
         edges = shape.Edges()
 
         result_type, result_attributes = DxfDocument._dxf_spline(edges[0], plane)
@@ -644,9 +653,6 @@ def testDXF(self):
 
         exporters.export(self._box().section(), "out.dxf")
 
-        with self.assertRaises(ValueError):
-            exporters.export(self._box().val(), "out.dxf")
-
         s1 = (
             Workplane("XZ")
             .polygon(10, 10)
@@ -945,3 +951,23 @@ def test_dxf_ellipse_arc(tmpdir):
 
     assert w2.val().isValid()
     assert w2.val().Volume() == approx(math.pi * r ** 2 / 4)
+
+
+def test_dxf_sketch(fname):
+
+    s = Sketch().rect(1, 2)
+    exporters.exportDXF(s, fname)
+
+    s_imported = Sketch().importDXF(fname)
+
+    assert (s._faces - s_imported._faces).Volume() == 0
+
+
+def test_dxf_shape(fname):
+
+    s = face(rect(10, 0.5))
+    exporters.exportDXF(s, fname)
+
+    s_imported = Sketch().importDXF(fname).val()
+
+    assert (s - s_imported).Volume() == 0
diff --git a/tests/test_free_functions.py b/tests/test_free_functions.py
index bbf89b221..2092ce24b 100644
--- a/tests/test_free_functions.py
+++ b/tests/test_free_functions.py
@@ -547,3 +547,14 @@ def test_loft():
     assert r4.Volume() > 0
     assert r5.Area() == approx(1)
     assert len(r6.Faces()) == 16
+
+
+# %% export
+def test_export():
+
+    b1 = box(1, 1, 1)
+    b1.export("box.brep")
+
+    b2 = Shape.importBrep("box.brep")
+
+    assert (b1 - b2).Volume() == approx(0)
diff --git a/tests/test_sketch.py b/tests/test_sketch.py
index dc2a2788e..2e21ae482 100644
--- a/tests/test_sketch.py
+++ b/tests/test_sketch.py
@@ -2,8 +2,9 @@
 
 from cadquery.sketch import Sketch, Vector, Location
 from cadquery.selectors import LengthNthSelector
+from cadquery import Edge, Vertex
 
-from pytest import approx, raises
+from pytest import approx, raises, fixture
 from math import pi, sqrt
 
 testdataDir = os.path.join(os.path.dirname(__file__), "testdata")
@@ -47,39 +48,54 @@ def test_face_interface():
     assert len(s7.vertices()._selection) == 3
     assert s7._faces.Area() == approx(0.5)
 
+    s8 = Sketch().face(Sketch().rect(1, 1).val())
+    assert s8._faces.Area() == approx(1)
+
+    s9 = Sketch().face(Sketch().rect(1, 1))
+    assert s9._faces.Area() == approx(1)
+
     with raises(ValueError):
-        Sketch().face(Sketch().rect(1, 1)._faces)
+        Sketch().face(1234)
 
 
 def test_modes():
 
+    # additive mode
     s1 = Sketch().rect(2, 2).rect(1, 1, mode="a")
 
     assert s1._faces.Area() == approx(4)
     assert len(s1._faces.Faces()) == 2
 
+    # subtraction mode
     s2 = Sketch().rect(2, 2).rect(1, 1, mode="s")
 
     assert s2._faces.Area() == approx(4 - 1)
     assert len(s2._faces.Faces()) == 1
 
+    # intersection mode
     s3 = Sketch().rect(2, 2).rect(1, 1, mode="i")
 
     assert s3._faces.Area() == approx(1)
     assert len(s3._faces.Faces()) == 1
 
+    # construction mode
     s4 = Sketch().rect(2, 2).rect(1, 1, mode="c", tag="t")
 
     assert s4._faces.Area() == approx(4)
     assert len(s4._faces.Faces()) == 1
     assert s4._tags["t"][0].Area() == approx(1)
 
+    # construction mode requires tagging
     with raises(ValueError):
         Sketch().rect(2, 2).rect(1, 1, mode="c")
 
     with raises(ValueError):
         Sketch().rect(2, 2).rect(1, 1, mode="dummy")
 
+    # replace mode
+    s5 = Sketch().rect(1, 1).wires().offset(-0.1, mode="r").reset()
+    assert s5.val().Area() == approx(0.8 ** 2)
+
 
 def test_distribute():
 
@@ -372,12 +388,15 @@ def test_delete():
 
     assert len(s2._edges) == 2
 
+    with raises(ValueError):
+        s2.vertices().delete()
+
 
 def test_selectors():
 
     s = Sketch().push([(-2, 0), (2, 0)]).rect(1, 1).rect(0.5, 0.5, mode="s").reset()
 
-    assert len(s._selection) == 0
+    assert s._selection is None
 
     s.vertices()
 
@@ -385,7 +404,7 @@ def test_selectors():
 
     s.reset()
 
-    assert len(s._selection) == 0
+    assert s._selection is None
 
     s.edges()
 
@@ -409,7 +428,7 @@ def test_selectors():
 
     s.tag("test").reset()
 
-    assert len(s._selection) == 0
+    assert s._selection is None
 
     s.select("test")
 
@@ -427,8 +446,9 @@ def test_selectors():
     s.reset().vertices("<X and <Y").val()
     assert s.val().toTuple() == approx((-2.5, -0.5, 0.0))
 
-    s.reset().vertices(">>X[1] and <Y").val()
-    assert s.val().toTuple()[0] == approx((0, 0, 0))
+    with raises(IndexError):
+        s.reset().vertices(">>X[1] and <Y").val()
+    assert len(s._selection) == 0
 
 
 def test_edge_interface():
@@ -745,3 +765,195 @@ def test_dxf_import():
     s5 = Sketch().importDXF(filename, tol=1e-3, exclude=["1"])
 
     assert s5._faces.isValid()
+
+
+def test_val():
+
+    s1 = Sketch().segment((0, 0), (0, 1))
+
+    assert isinstance(s1.val(), Edge)
+
+    s1.vertices()
+
+    assert isinstance(s1.val(), Vertex)
+
+    s2 = Sketch().circle(1)
+
+    assert len(s2.val().Faces()) == 1
+
+
+def test_vals():
+
+    s1 = Sketch().segment((0, 0), (0, 1))
+
+    assert len(s1.vals()) == 1
+
+    s1.vertices()
+
+    assert len(s1.vals()) == 2
+
+    s2 = Sketch().circle(1)
+
+    assert len(s2.vals()) == 1
+
+
+def test_bool_ops():
+
+    s1 = Sketch().rect(0.5, 2)
+    s2 = Sketch().rect(1, 1)
+    s3 = Sketch().segment((-1, 0), (1, 0))
+
+    # sum
+    assert (s1 + s2).val().Area() == approx(1.5)
+    # diff
+    assert (s1 - s2).val().Area() == approx(0.5)
+    assert len((s1 - s2).val().Wires()) == 2
+    # common
+    assert (s1 * s2).val().Area() == approx(0.5)
+    assert len((s1 * s2).val().Wires()) == 1
+    # split
+    assert len((s1 / s3).val().Faces()) == 2
+    assert len((s1 / s2).val().Faces()) == 3
+    assert (s1 / s2).val().Area() == approx(1)
+
+
+def test_export():
+
+    s1 = Sketch().rect(1, 1).export("sketch.dxf")
+    s2 = Sketch().importDXF("sketch.dxf")
+
+    assert (s1 - s2).val().Area() == approx(0)
+
+
+@fixture
+def s1():
+    """
+    Sketch with 2 faces
+
+    """
+
+    return Sketch().rect(1, 1).push([(5, 0)]).rect(1, 4).reset()
+
+
+@fixture
+def s2():
+    """
+    Sketch with 5 faces
+
+    """
+
+    return Sketch().rarray(2, 1, 5, 1).rect(1, 1).reset()
+
+
+@fixture
+def f():
+    """
+    Single face
+    """
+
+    return Sketch().rect(0.1, 1).val()
+
+
+def test_filter(s1, f):
+
+    assert len(s1.filter(lambda x: x.Area() > 2).vals()) == 1
+    assert len(s1.reset().filter(lambda x: x.Area() >= 1).vals()) == 2
+    assert len(s1.filter(lambda x: x.Area() < 1).vals()) == 0
+
+
+def test_sort(s1):
+
+    assert s1.sort(lambda x: -x.Area())[-1].val().Area() == approx(1)
+
+
+def test_apply(s1, f):
+
+    assert s1.apply(lambda x: [f]).val().Area() == approx(0.1)
+
+
+def test_map(s1, f):
+    assert s1.map(
+        lambda x: f.moved(x.Center())
+    ).replace().reset().val().Area() == approx(0.2)
+
+
+def test_getitem(s2):
+
+    assert len(s2[0].vals()) == 1
+    assert len(s2.reset()[-2:].vals()) == 2
+    assert len(s2.reset()[[0, 1]].vals()) == 2
+
+
+def test_invoke(s1, s2):
+
+    # builtin
+    assert len(s2.invoke(print).vals()) == 5
+    # arity 0
+    assert len(s2.invoke(lambda: 1).vals()) == 5
+    # arity 1 and no return
+    assert len(s2.invoke(lambda x: None).vals()) == 5
+    # arity 1
+    assert len(s2.invoke(lambda x: s1).vals()) == 2
+    # test exception with wrong arity
+    with raises(ValueError):
+        s2.invoke(lambda x, y: 1)
+
+
+@fixture
+def s3():
+    """
+    Simple sketch with one face.
+    """
+
+    return Sketch().rect(1, 1)
+
+
+def test_replace(s3, f):
+
+    assert s3.map(lambda x: f).replace().reset().val().Area() == approx(0.1)
+
+
+def test_add(s3, f):
+
+    # we do not clean, so adding will increase the number of edges
+    assert len(s3.map(lambda x: f).add().reset().val().Edges()) == 10
+
+
+def test_subtract(s3, f):
+
+    assert s3.map(lambda x: f).subtract().reset().val().Area() == approx(0.9)
+
+
+def test_iter(s1):
+
+    # __iter__ ofer face
+    assert len(list(s1)) == 2
+
+    # __iter__ over edges
+    assert len(list(Sketch().segment((0, 0), (1, 0)))) == 1
+
+
+def test_sanitize():
+
+    # it does not make sense to fuse faces and Locations
+    with raises(ValueError):
+        Sketch().rect(1, 1) + Sketch().rarray(1, 1, 1, 1)
+
+
+def test_missing_selection(s1):
+
+    # offset requires selected wires
+    with raises(ValueError):
+        s1.offset(0.1)
+
+    # fillet requires selected vertices
+    with raises(ValueError):
+        s1.fillet(0.1)
+
+    # cannot delete without selection
+    with raises(ValueError):
+        s1.delete()
+
+    # cannot tag without selection
+    with raises(ValueError):
+        s1.tag("name")