diff --git a/benchmark/benchmarks.jl b/benchmark/benchmarks.jl
index 14896a03..6a750efa 100644
--- a/benchmark/benchmarks.jl
+++ b/benchmark/benchmarks.jl
@@ -16,20 +16,22 @@ integrands = (
     (name = "Vector", f = p -> fill(norm(to(p)), 3))
 )
 rules = (
-    (name = "GaussLegendre", rule = GaussLegendre(100)),
-    (name = "GaussKronrod", rule = GaussKronrod()),
-    (name = "HAdaptiveCubature", rule = HAdaptiveCubature())
+    GaussLegendre(100),
+    GaussKronrod(),
+    HAdaptiveCubature()
 )
 geometries = (
-    (name = "Segment", item = Segment(Point(0, 0, 0), Point(1, 1, 1))),
-    (name = "Sphere", item = Sphere(Point(0, 0, 0), 1.0))
+    let ẑ = Vec(0, 0, 1)
+        CylinderSurface(Plane(Point(0, 0, 0), ẑ), Plane(Point(0, 0, 3), ẑ), 2.5)
+    end,
+    Segment(Point(0, 0, 0), Point(1, 1, 1)),
+    Sphere(Point(0, 0, 0), 1.0)
 )
 
 SUITE["Integrals"] = let s = BenchmarkGroup()
     for (int, rule, geometry) in Iterators.product(integrands, rules, geometries)
-        n1 = geometry.name
-        n2 = "$(int.name) $(rule.name)"
-        s[n1][n2] = @benchmarkable integral($int.f, $geometry.item, $rule.rule)
+        name = "$(nameof(typeof(geometry))), $(int.name), $(nameof(typeof(rule)))"
+        s[name] = @benchmarkable integral($int.f, $geometry, $rule)
     end
     s
 end
@@ -41,11 +43,8 @@ end
 spec = (
     f = p -> norm(to(p)),
     f_exp = p::Point -> exp(-norm(to(p))^2 / u"m^2"),
-    g = (
+    geometries = (
         bezier = BezierCurve([Point(t, sin(t), 0) for t in range(-pi, pi, length = 361)]),
-        line = Line(Point(0, 0, 0), Point(1, 1, 1)),
-        plane = Plane(Point(0, 0, 0), Vec(0, 0, 1)),
-        ray = Ray(Point(0, 0, 0), Vec(0, 0, 1)),
         rope = Rope([Point(t, t, t) for t in 1:32]...),
         triangle = Triangle(Point(1, 0, 0), Point(0, 1, 0), Point(0, 0, 1)),
         tetrahedron = let
@@ -56,19 +55,29 @@ spec = (
             Tetrahedron(a, b, c, a + ẑ)
         end
     ),
-    rule_gl = GaussLegendre(100),
-    rule_gk = GaussKronrod(),
-    rule_hc = HAdaptiveCubature()
+    geometries_exp = (
+        line = Line(Point(0, 0, 0), Point(1, 1, 1)),
+        plane = Plane(Point(0, 0, 0), Vec(0, 0, 1)),
+        ray = Ray(Point(0, 0, 0), Vec(0, 0, 1))
+    ),
+    rules = (
+        GaussLegendre(100),
+        HAdaptiveCubature()
+    )
 )
 
-SUITE["Specializations/Scalar GaussLegendre"] = let s = BenchmarkGroup()
-    s["BezierCurve"] = @benchmarkable integral($spec.f, $spec.g.bezier, $spec.rule_gl)
-    s["Line"] = @benchmarkable integral($spec.f_exp, $spec.g.line, $spec.rule_gl)
-    s["Plane"] = @benchmarkable integral($spec.f_exp, $spec.g.plane, $spec.rule_gl)
-    s["Ray"] = @benchmarkable integral($spec.f_exp, $spec.g.ray, $spec.rule_gl)
-    s["Rope"] = @benchmarkable integral($spec.f, $spec.g.rope, $spec.rule_gl)
-    s["Triangle"] = @benchmarkable integral($spec.f, $spec.g.triangle, $spec.rule_gl)
-    s["Tetrahedron"] = @benchmarkable integral($spec.f, $spec.g.tetrahedron, $spec.rule_gl)
+SUITE["Specializations"] = let s = BenchmarkGroup()
+    # Benchmark most specialization geometries
+    for rule in spec.rules, geometry in spec.geometries
+        name = "$(nameof(typeof(geometry))), Scalar, $(nameof(typeof(rule)))"
+        s[name] = @benchmarkable integral($spec.f, $geometry, $rule)
+    end
+
+    # Geometries that span an infinite domain use exp function instead
+    for rule in spec.rules, geometry in spec.geometries_exp
+        name = "$(nameof(typeof(geometry))), Scalar, $(nameof(typeof(rule)))"
+        s[name] = @benchmarkable integral($spec.f_exp, $geometry, $rule)
+    end
     s
 end
 
diff --git a/src/specializations/CylinderSurface.jl b/src/specializations/CylinderSurface.jl
index c359e3af..f6375edd 100644
--- a/src/specializations/CylinderSurface.jl
+++ b/src/specializations/CylinderSurface.jl
@@ -31,3 +31,33 @@ function integral(
 
     return sides + top + bottom
 end
+
+function integral(
+    f,
+    cyl::Meshes.CylinderSurface,
+    rule::HAdaptiveCubature;
+    FP::Type{T} = Float64,
+    diff_method::DM = _default_diff_method(cyl, FP),
+    kwargs...
+) where {T <: AbstractFloat, DM <: DifferentiationMethod}
+    _check_diff_method_support(cyl, diff_method)
+
+    # Use a sample integrand to develop and append a buffer to the given rule
+    f_sample(ts) = Unitful.ustrip.(f(cyl(ts...)) * differential(cyl, ts))
+    N = Meshes.paramdim(cyl)
+    buffer = HCubature.hcubature_buffer(f_sample, _zeros(FP, N), _ones(FP, N))
+    rule = HAdaptiveCubature(rule.kwargs..., buffer = buffer)
+
+    # The generic method only parametrizes the sides
+    sides = _integral(f, cyl, rule; diff_method = diff_method, FP = FP, kwargs...)
+
+    # Integrate the Disk at the top
+    disk_top = Meshes.Disk(cyl.top, cyl.radius)
+    top = _integral(f, disk_top, rule; diff_method = diff_method, FP = FP, kwargs...)
+
+    # Integrate the Disk at the bottom
+    disk_bottom = Meshes.Disk(cyl.bot, cyl.radius)
+    bottom = _integral(f, disk_bottom, rule; diff_method = diff_method, FP = FP, kwargs...)
+
+    return sides + top + bottom
+end
diff --git a/src/specializations/Rope.jl b/src/specializations/Rope.jl
index e917708b..abd63936 100644
--- a/src/specializations/Rope.jl
+++ b/src/specializations/Rope.jl
@@ -35,3 +35,27 @@ function integral(
     # Convert the Rope into Segments, sum the integrals of those
     return sum(segment -> integral(f, segment, rule; kwargs...), Meshes.segments(rope))
 end
+
+# Use HCubature.hcubature_buffer to reduce allocations
+function integral(
+    f,
+    rope::Meshes.Rope,
+    rule::HAdaptiveCubature;
+    FP::Type{T} = Float64,
+    kwargs...
+) where {T <: AbstractFloat}
+    # Geometry information
+    N = Meshes.paramdim(rope)
+    segments = Meshes.segments(rope)
+
+    # Use a sample integrand to develop and append a buffer to the given rule
+    sample = first(segments)
+    integrand(ts) = f(sample(ts...)) * differential(sample, ts)
+    uintegrand(ts) = Unitful.ustrip.(integrand(ts))
+    buffer = HCubature.hcubature_buffer(uintegrand, _zeros(FP, N), _ones(FP, N))
+    rule = HAdaptiveCubature(rule.kwargs..., buffer = buffer)
+
+    # Convert the Rope into Segments, sum the integrals of those
+    _subintegral(seg) = _integral(f, seg, rule; FP = FP, kwargs...)
+    return sum(_subintegral, segments)
+end