Obsolete use of nested GaussKronrod rules

CHANGELOG.md

@@ -7,6 +7,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Removed
+
+- Removed previously-deprecated support for use of `GaussKronrod` rules on geometries with more than one parametric dimension.
+
 
 ## [0.16.4] - 2025-08-10
 

Project.toml

@@ -1,6 +1,6 @@
 name = "MeshIntegrals"
 uuid = "dadec2fd-bbe0-4da4-9dbe-476c782c8e47"
-version = "0.16.4"
+version = "0.17.0-DEV"
 
 [deps]
 CliffordNumbers = "3998ac73-6bd4-4031-8035-f167dd3ed523"

benchmark/benchmarks.jl

@@ -21,12 +21,19 @@ rules = (
     (name = "HAdaptiveCubature", rule = HAdaptiveCubature())
 )
 geometries = (
-    (name = "Segment", item = Segment(Point(0, 0, 0), Point(1, 1, 1))),
-    (name = "Sphere", item = Sphere(Point(0, 0, 0), 1.0))
+    (name = "Segment", item = Segment(Point(0, 0, 0), Point(1, 1, 1))), # 1D
+    (name = "Sphere", item = Sphere(Point(0, 0), 1.0)), # 2D
+    (name = "Sphere", item = Sphere(Point(0, 0, 0), 1.0)) # 3D
 )
 
 SUITE["Integrals"] = let s = BenchmarkGroup()
     for (int, rule, geometry) in Iterators.product(integrands, rules, geometries)
+        # Skip unsupported applications of GaussKronrod
+        if (Meshes.paramdim(geometry.item) > 1) && (rule.rule == GaussKronrod())
+            continue
+        end
+
+        # Construct benchmark and add it to test suite
         n1 = geometry.name
         n2 = "$(int.name) $(rule.name)"
         s[n1][n2] = @benchmarkable integral($int.f, $geometry.item, $rule.rule)

docs/src/support.md

@@ -2,24 +2,13 @@
 
 This library aims to enable users to calculate the value of integrals over all
 [**Meshes.jl**](https://github.com/JuliaGeometry/Meshes.jl) geometry types using
-a number of numerical integration rules and techniques. However, some combinations
-of geometry types and integration rules are ill-suited (and a few are simply not
-yet implemented).
-
-## General Recommendations
+a number of numerical integration rules and techniques.
 
 In general, `GaussKronrod` integration rules are recommended (and the default) for
 geometries with one parametric dimension. For geometries with more than one
 parametric dimension, e.g. surfaces and volumes, `HAdaptiveCubature` rules are
 recommended (and the default).
 
-While it is currently possible to apply nested `GaussKronrod` rules to numerically
-integrate surfaces, this produces results that are strictly inferior to using an
-equivalent `HAdaptiveCubature` rule, so support for this usage has been deprecated.
-In version 16.x of MeshIntegrals.jl, using a `GaussKronrod` rule for a surface
-will work but will yield a deprecation warning. Beginning with a future version
-17.0, this combination will simply be unsupported and throw an error.
-
 ## The Support Matrix
 
 The following Support Matrix captures the current state of support for all geometry/rule
@@ -28,46 +17,45 @@ designed to check for accuracy.
 
 | `Meshes.Geometry/Domain` | `GaussKronrod` | `GaussLegendre` | `HAdaptiveCubature` |
 |----------|----------------|---------------|---------------------|
-| `Ball` in `𝔼{2}` | ⚠️ | ✅ | ✅ |
+| `Ball` in `𝔼{2}` | 🛑 | ✅ | ✅ |
 | `Ball` in `𝔼{3}` | 🛑 | ✅ | ✅ |
 | `BezierCurve` | ✅ | ✅ | ✅ |
 | `Box` in `𝔼{1}` | ✅ | ✅ | ✅ |
-| `Box` in `𝔼{2}` | ⚠️ | ✅ | ✅ |
+| `Box` in `𝔼{2}` | 🛑 | ✅ | ✅ |
 | `Box` in `𝔼{≥3}` | 🛑 | ✅ | ✅ |
 | `CartesianGrid` | ✅ | ✅ | ✅ |
 | `Circle` | ✅ | ✅ | ✅ |
 | `Cone` | 🛑 | ✅ | ✅ |
-| `ConeSurface` | ⚠️ | ✅ | ✅ |
+| `ConeSurface` | 🛑 | ✅ | ✅ |
 | `Cylinder` | 🛑 | ✅ | ✅ |
-| `CylinderSurface` | ⚠️ | ✅ | ✅ |
-| `Disk` | ⚠️ | ✅ | ✅ |
+| `CylinderSurface` | 🛑 | ✅ | ✅ |
+| `Disk` | 🛑 | ✅ | ✅ |
 | `Ellipsoid` | ✅ | ✅ | ✅ |
-| `Frustum` | ⚠️ | ✅ | ✅ |
-| `FrustumSurface` | ⚠️ | ✅ | ✅ |
+| `Frustum` | 🛑 | ✅ | ✅ |
+| `FrustumSurface` | 🛑 | ✅ | ✅ |
 | `Hexahedron` | ✅ | ✅ | ✅ |
 | `Line` | ✅ | ✅ | ✅ |
-| `ParaboloidSurface` | ⚠️ | ✅ | ✅ |
+| `ParaboloidSurface` | 🛑 | ✅ | ✅ |
 | `ParametrizedCurve` | ✅ | ✅ | ✅ |
 | `Plane` | ✅ | ✅ | ✅ |
-| `PolyArea` | ⚠️ | ✅ | ✅ |
-| `Pyramid` | ⚠️ | ✅ | ✅ |
-| `Quadrangle` | ⚠️ | ✅ | ✅ |
+| `PolyArea` | 🛑 | ✅ | ✅ |
+| `Pyramid` | 🛑 | ✅ | ✅ |
+| `Quadrangle` | 🛑 | ✅ | ✅ |
 | `Ray` | ✅ | ✅ | ✅ |
 | `RegularGrid` | ✅ | ✅ | ✅ |
 | `Ring` | ✅ | ✅ | ✅ |
 | `Rope` | ✅ | ✅ | ✅ |
 | `Segment` | ✅ | ✅ | ✅ |
-| `SimpleMesh` | ⚠️ | ✅ | ✅ |
+| `SimpleMesh` | 🛑 | ✅ | ✅ |
 | `Sphere` in `𝔼{2}` | ✅ | ✅ | ✅ |
-| `Sphere` in `𝔼{3}` | ⚠️ | ✅ | ✅ |
+| `Sphere` in `𝔼{3}` | 🛑 | ✅ | ✅ |
 | `StructuredGrid` | ✅ | ✅ | ✅ |
-| `Tetrahedron` | ⚠️ | ✅ | ✅ |
+| `Tetrahedron` | 🛑 | ✅ | ✅ |
 | `Triangle` | ✅ | ✅ | ✅ |
-| `Torus` | ⚠️ | ✅ | ✅ |
-| `Wedge` | ⚠️ | ✅ | ✅ |
+| `Torus` | 🛑 | ✅ | ✅ |
+| `Wedge` | 🛑 | ✅ | ✅ |
 
 | Symbol | Support Level |
 |--------|---------|
 | ✅ | Supported |
-| ⚠️ | Deprecated |
 | 🛑 | Not supported |

src/integral.jl

@@ -63,24 +63,14 @@ function _integral(
 ) where {DM <: DifferentiationMethod, T <: AbstractFloat}
     _check_diff_method_support(geometry, diff_method)
 
-    # Implementation depends on number of parametric dimensions over which to integrate
-    N = Meshes.paramdim(geometry)
-    if N == 1
-        integrand(t) = f(geometry(t)) * differential(geometry, (t,), diff_method)
-        return QuadGK.quadgk(integrand, zero(FP), one(FP); rule.kwargs...)[1]
-    elseif N == 2
-        # Issue deprecation warning
-        Base.depwarn("Use `HAdaptiveCubature` instead of \
-                     `GaussKronrod` for surfaces.", :integral)
-
-        # Nested integration
-        integrand2D(u, v) = f(geometry(u, v)) * differential(geometry, (u, v), diff_method)
-        ∫₁(v) = QuadGK.quadgk(u -> integrand2D(u, v), zero(FP), one(FP); rule.kwargs...)[1]
-        return QuadGK.quadgk(v -> ∫₁(v), zero(FP), one(FP); rule.kwargs...)[1]
-    else
-        _error_unsupported_combination("geometry with more than two parametric dimensions",
-            "GaussKronrod")
+    # Only supported for 1D geometries
+    if Meshes.paramdim(geometry) != 1
+        msg = "GaussKronrod rules only supported in one parametric dimension."
+        throw(ArgumentError(msg))
     end
+
+    integrand(t) = f(geometry(t)) * differential(geometry, (t,), diff_method)
+    return QuadGK.quadgk(integrand, zero(FP), one(FP); rule.kwargs...)[1]
 end
 
 # GaussLegendre

src/integration_rules.jl

@@ -12,10 +12,9 @@ abstract type IntegrationRule end
     GaussKronrod(kwargs...)
 
 The h-adaptive Gauss-Kronrod quadrature rule implemented by
-[QuadGK.jl](https://github.com/JuliaMath/QuadGK.jl). All standard
-`QuadGK.quadgk` keyword arguments are supported. This rule works natively for one
-dimensional geometries; some two- and three-dimensional geometries are additionally
-supported using nested integral solvers with the specified `kwarg` settings.
+[QuadGK.jl](https://github.com/JuliaMath/QuadGK.jl) which can be used for any
+single-dimensional geometry. All standard `QuadGK.quadgk` keyword arguments are
+supported.
 """
 struct GaussKronrod <: IntegrationRule
     kwargs::Base.Pairs

src/utils.jl

@@ -1,13 +1,3 @@
-################################################################################
-#                           Misc. Internal Tools
-################################################################################
-
-# Common error message structure
-function _error_unsupported_combination(geometry, rule)
-    msg = "Integrating a $geometry using a $rule rule not supported."
-    throw(ArgumentError(msg))
-end
-
 ################################################################################
 #                           DifferentiationMethod
 ################################################################################

test/combinations.jl

@@ -63,7 +63,7 @@ This file includes tests for:
         elseif N == 2
             # surface
             aliases = Bool.((0, 1, 0))
-            rules = Bool.((1, 1, 1))
+            rules = Bool.((0, 1, 1))
             return SupportStatus(aliases..., rules..., autoenzyme)
         elseif N == 3
             # volume
@@ -92,7 +92,7 @@ This file includes tests for:
             if getfield(supports, nameof(alias))
                 @test alias(testable.integrand, testable.geometry)≈testable.solution rtol=rtol
             else
-                @test_throws "not supported" alias(testable.integrand, testable.geometry)
+                @test_throws "supported" alias(testable.integrand, testable.geometry)
             end
         end # for
 
@@ -119,7 +119,7 @@ This file includes tests for:
             else
                 f = testable.integrand
                 geometry = testable.geometry
-                @test_throws "not supported" integral(f, geometry, rule)
+                @test_throws "supported" integral(f, geometry, rule)
             end
         end # for
 
@@ -137,7 +137,7 @@ This file includes tests for:
             if supported
                 @test integral(f, geometry; diff_method = method)≈sol rtol=rtol
             else
-                @test_throws "not supported" integral(f, geometry; diff_method = method)
+                @test_throws "supported" integral(f, geometry; diff_method = method)
             end
         end # for
     end # function