diff --git a/AC-Vector.cabal b/AC-Vector.cabal
--- a/AC-Vector.cabal
+++ b/AC-Vector.cabal
@@ -1,6 +1,6 @@
 Cabal-Version: >= 1.6
 Name:          AC-Vector
-Version:       2.3.0
+Version:       2.3.1
 Stability:     Experimental
 Synopsis:      Efficient geometric vectors and transformations.
 
@@ -14,9 +14,8 @@
   just contains data structures that are useful for graphics
   work.)
 
-  The @Vector@ class now implies @Num@ and @Fractional@.
-  Hopefully this will be the last breaking change for a while.
-  Also added @vpromote@, @(\/|)@ and @(|\/)@.
+  Added the @unions@ function. (More efficient than
+  @foldr1 union@.)
 
 Category:      Data, Math, Numerical, Graphics
 License:       BSD3
diff --git a/Data/BoundingBox/B1.hs b/Data/BoundingBox/B1.hs
--- a/Data/BoundingBox/B1.hs
+++ b/Data/BoundingBox/B1.hs
@@ -6,10 +6,10 @@
 
 import Data.Vector.Class
 import Data.Vector.V1
-import Data.BoundingBox.Range as R
+import qualified Data.BoundingBox.Range as R
 
 -- | The 'BBox1' type is basically a 'Range', but all the operations over it work with 'Vector1' (which is really 'Scalar'). While it's called a bounding /box/, a 1-dimensional box is in truth a simple line interval, just like 'Range'.
-newtype BBox1 = BBox1 {range :: Range} deriving (Eq, Show)
+newtype BBox1 = BBox1 {range :: R.Range} deriving (Eq, Show)
 
 -- | Given two vectors, construct a bounding box (swapping the endpoints if necessary).
 bound_corners :: Vector1 -> Vector1 -> BBox1
@@ -40,3 +40,7 @@
 isect (BBox1 r0) (BBox1 r1) = do
   r <- (r0 `R.isect` r1)
   return (BBox1 r)
+
+-- | Efficiently compute the union of a list of bounding boxes.
+unions :: [BBox1] -> BBox1
+unions = BBox1 . R.unions . map range
diff --git a/Data/BoundingBox/B2.hs b/Data/BoundingBox/B2.hs
--- a/Data/BoundingBox/B2.hs
+++ b/Data/BoundingBox/B2.hs
@@ -6,22 +6,22 @@
 
 import Data.Vector.Class
 import Data.Vector.V2
-import Data.BoundingBox.Range as R
+import qualified Data.BoundingBox.Range as R
 
 -- | A 'BBox2' is a 2D bounding box (aligned to the coordinate axies).
 data BBox2 = BBox2 {minX, minY, maxX, maxY :: {-# UNPACK #-} !Scalar} deriving (Eq, Show)
 
 -- | Return the X-range that this bounding box covers.
-rangeX :: BBox2 -> Range
-rangeX b = Range (minX b) (maxX b)
+rangeX :: BBox2 -> R.Range
+rangeX b = R.Range (minX b) (maxX b)
 
 -- | Return the Y-range that this bounding box covers.
-rangeY :: BBox2 -> Range
-rangeY b = Range (minY b) (maxY b)
+rangeY :: BBox2 -> R.Range
+rangeY b = R.Range (minY b) (maxY b)
 
 -- | Given ranges for each coordinate axis, construct a bounding box.
-rangeXY :: Range -> Range -> BBox2
-rangeXY (Range x0 x1) (Range y0 y1) = BBox2 x0 y0 x1 y1
+rangeXY :: R.Range -> R.Range -> BBox2
+rangeXY (R.Range x0 x1) (R.Range y0 y1) = BBox2 x0 y0 x1 y1
 
 -- | Given a pair of corner points, construct a bounding box. (The points must be from opposite corners, but it doesn't matter /which/ corners nor which order they are given in.)
 bound_corners :: Vector2 -> Vector2 -> BBox2
@@ -63,3 +63,14 @@
   rx <- (rangeX b0) `R.isect` (rangeX b1)
   ry <- (rangeY b0) `R.isect` (rangeY b1)
   return (rangeXY rx ry)
+
+-- | Efficiently compute the union of a list of bounding boxes.
+unions :: [BBox2] -> BBox2
+unions bs =
+  let
+    minP = map min_point bs
+    maxP = map max_point bs
+  in
+    BBox2
+      (minimum $ map v2x minP) (minimum $ map v2y minP)
+      (maximum $ map v2x maxP) (maximum $ map v2y maxP)
diff --git a/Data/BoundingBox/B3.hs b/Data/BoundingBox/B3.hs
--- a/Data/BoundingBox/B3.hs
+++ b/Data/BoundingBox/B3.hs
@@ -6,26 +6,26 @@
 
 import Data.Vector.Class
 import Data.Vector.V3
-import Data.BoundingBox.Range as R
+import qualified Data.BoundingBox.Range as R
 
 -- | A 'BBox3' is a 3D bounding box (aligned to the coordinate axies).
 data BBox3 = BBox3 {minX, minY, minZ, maxX, maxY, maxZ :: {-# UNPACK #-} !Scalar} deriving (Eq, Show)
 
 -- | Return the X-range that this bounding box covers.
-rangeX :: BBox3 -> Range
-rangeX b = Range (minX b) (maxX b)
+rangeX :: BBox3 -> R.Range
+rangeX b = R.Range (minX b) (maxX b)
 
 -- | Return the Y-range that this bounding box covers.
-rangeY :: BBox3 -> Range
-rangeY b = Range (minY b) (maxY b)
+rangeY :: BBox3 -> R.Range
+rangeY b = R.Range (minY b) (maxY b)
 
 -- | Return the Z-range that this bounding box covers.
-rangeZ :: BBox3 -> Range
-rangeZ b = Range (minZ b) (maxZ b)
+rangeZ :: BBox3 -> R.Range
+rangeZ b = R.Range (minZ b) (maxZ b)
 
 -- | Given ranges for each coordinate axis, construct a bounding box.
-rangeXYZ :: Range -> Range -> Range -> BBox3
-rangeXYZ (Range x0 x1) (Range y0 y1) (Range z0 z1) = BBox3 x0 y0 z0 x1 y1 z1
+rangeXYZ :: R.Range -> R.Range -> R.Range -> BBox3
+rangeXYZ (R.Range x0 x1) (R.Range y0 y1) (R.Range z0 z1) = BBox3 x0 y0 z0 x1 y1 z1
 
 -- | Given a pair of corner points, construct a bounding box. (The points must be from opposite corners, but it doesn't matter /which/ corners nor which order they are given in.)
 bound_corners :: Vector3 -> Vector3 -> BBox3
@@ -71,3 +71,14 @@
   ry <- (rangeY b0) `R.isect` (rangeY b1)
   rz <- (rangeZ b0) `R.isect` (rangeZ b1)
   return (rangeXYZ rx ry rz)
+
+-- | Efficiently compute the union of a list of bounding boxes.
+unions :: [BBox3] -> BBox3
+unions bs =
+  let
+    minP = map min_point bs
+    maxP = map max_point bs
+  in
+    BBox3
+      (minimum $ map v3x minP) (minimum $ map v3y minP) (minimum $ map v3z minP)
+      (maximum $ map v3x maxP) (maximum $ map v3y maxP) (maximum $ map v3z maxP)
diff --git a/Data/BoundingBox/B4.hs b/Data/BoundingBox/B4.hs
--- a/Data/BoundingBox/B4.hs
+++ b/Data/BoundingBox/B4.hs
@@ -6,30 +6,30 @@
 
 import Data.Vector.Class
 import Data.Vector.V4
-import Data.BoundingBox.Range as R
+import qualified Data.BoundingBox.Range as R
 
 -- | A 'BBox4' is a 4D bounding box (aligned to the coordinate axies).
 data BBox4 = BBox4 {minX, minY, minZ, minW, maxX, maxY, maxZ, maxW :: {-# UNPACK #-} !Scalar} deriving (Eq, Show)
 
 -- | Return the X-range that this bounding box covers.
-rangeX :: BBox4 -> Range
-rangeX b = Range (minX b) (maxX b)
+rangeX :: BBox4 -> R.Range
+rangeX b = R.Range (minX b) (maxX b)
 
 -- | Return the Y-range that this bounding box covers.
-rangeY :: BBox4 -> Range
-rangeY b = Range (minY b) (maxY b)
+rangeY :: BBox4 -> R.Range
+rangeY b = R.Range (minY b) (maxY b)
 
 -- | Return the Z-range that this bounding box covers.
-rangeZ :: BBox4 -> Range
-rangeZ b = Range (minZ b) (maxZ b)
+rangeZ :: BBox4 -> R.Range
+rangeZ b = R.Range (minZ b) (maxZ b)
 
 -- | Return the W-range (4th coordinate) that this bounding box covers.
-rangeW :: BBox4 -> Range
-rangeW b = Range (minW b) (maxW b)
+rangeW :: BBox4 -> R.Range
+rangeW b = R.Range (minW b) (maxW b)
 
 -- | Given ranges for each coordinate axis, construct a bounding box.
-rangeXYZW :: Range -> Range -> Range -> Range -> BBox4
-rangeXYZW (Range x0 x1) (Range y0 y1) (Range z0 z1) (Range w0 w1) = BBox4 x0 y0 z0 w0 x1 y1 z1 w1
+rangeXYZW :: R.Range -> R.Range -> R.Range -> R.Range -> BBox4
+rangeXYZW (R.Range x0 x1) (R.Range y0 y1) (R.Range z0 z1) (R.Range w0 w1) = BBox4 x0 y0 z0 w0 x1 y1 z1 w1
 
 -- | Given a pair of corner points, construct a bounding box. (The points must be from opposite corners, but it doesn't matter /which/ corners nor which order they are given in.)
 bound_corners :: Vector4 -> Vector4 -> BBox4
@@ -80,3 +80,14 @@
   rz <- (rangeZ b0) `R.isect` (rangeZ b1)
   rw <- (rangeW b0) `R.isect` (rangeW b1)
   return (rangeXYZW rx ry rz rw)
+
+-- | Efficiently compute the union of a list of bounding boxes.
+unions :: [BBox4] -> BBox4
+unions bs =
+  let
+    minP = map min_point bs
+    maxP = map max_point bs
+  in
+    BBox4
+      (minimum $ map v4x minP) (minimum $ map v4y minP) (minimum $ map v4z minP) (minimum $ map v4w minP)
+      (maximum $ map v4x maxP) (maximum $ map v4y maxP) (maximum $ map v4z maxP) (maximum $ map v4w maxP)
diff --git a/Data/BoundingBox/Range.hs b/Data/BoundingBox/Range.hs
--- a/Data/BoundingBox/Range.hs
+++ b/Data/BoundingBox/Range.hs
@@ -33,3 +33,7 @@
   if x00 < x10
     then if x10 < x01 then Just (Range x10 (min x01 x11)) else Nothing
     else if x00 < x11 then Just (Range x00 (min x01 x11)) else Nothing
+
+-- | Efficiently compute the union of a list of ranges.
+unions :: [Range] -> Range
+unions rs = Range (minimum $ map min_point rs) (maximum $ map max_point rs)