diff --git a/src/Approximations/box_approximation.jl b/src/Approximations/box_approximation.jl
index 63eb291f58..15baf2bbe4 100644
--- a/src/Approximations/box_approximation.jl
+++ b/src/Approximations/box_approximation.jl
@@ -27,7 +27,7 @@ function load_staticarrays_directions()
 end  # quote / load_staticarrays_directions
 
 """
-    box_approximation(S::LazySet{N}) where {N}
+    box_approximation(S::LazySet)
 
 Overapproximate a set by a tight hyperrectangle.
 
@@ -48,7 +48,11 @@ An alias for this function is `interval_hull`.
 The center and radius of the hyperrectangle are obtained by averaging the low
 and high coordinates of `S` computed with the `extrema` function.
 """
-function box_approximation(S::LazySet{N}) where {N}
+function box_approximation(S::LazySet)
+    return _box_approximation_extrema(S)
+end
+
+function _box_approximation_extrema(S::LazySet{N}) where {N}
     n = dim(S)
     c = Vector{N}(undef, n)
     r = Vector{N}(undef, n)
@@ -73,6 +77,19 @@ function box_approximation(S::LazySet{N}) where {N}
     return Hyperrectangle(c, r)
 end
 
+function box_approximation(cap::Intersection{N,T1,T2}) where {N,T1,T2}
+    if ispolyhedral(cap.X) && ispolyhedral(cap.Y)
+        if isboundedtype(T1) || isboundedtype(T2)
+            S = convert(HPolytope, cap)
+        else
+            S = convert(HPolyhedron, cap)
+        end
+    else
+        S = cap
+    end
+    return _box_approximation_extrema(S)
+end
+
 """
     interval_hull
 
diff --git a/test/Approximations/box_approximation.jl b/test/Approximations/box_approximation.jl
index 843980f2cb..4a92e03928 100644
--- a/test/Approximations/box_approximation.jl
+++ b/test/Approximations/box_approximation.jl
@@ -41,6 +41,18 @@ for N in [Float64, Rational{Int}, Float32]
     # alias constructors
     @test interval_hull(b) == □(b) == box_approximation(b)
 
+    # box_approximation of lazy intersection
+    # both arguments are bounded
+    X = Ball1(zeros(N, 2), N(1))
+    Y = Ball1(N[1//2, 0], N(1))
+    Z = box_approximation(X ∩ Y)
+    @test isequivalent(Z, Hyperrectangle(N[1//4, 0], N[3//4, 3//4]))
+    # one argument is unbounded
+    X = Hyperrectangle(N[0], N[1]) × Universe{N}(1)
+    Y = Hyperrectangle(N[0, 0], N[1, 1])
+    Z = box_approximation(X ∩ Y)
+    @test isequivalent(Z, Y)
+
     # ===================================================================
     # Testing box_approximation_symmetric (= symmetric interval hull)
     # ===================================================================