area(::Spherical, geom) infers Any: boxed multi-method local closure (_polygon_area) in NaiveTriangulatedSphericalArea

[asinghvi17]

Summary

area(::Spherical, geom) is type-unstable: the internal applyreduce in the NaiveTriangulatedSphericalArea method returns Any, because its local helper _polygon_area is a multi-method local function captured by a wrapping lambda, which Julia lowers into a Core.Box{Any}. This forces a heap allocation + dynamic dispatch per call, and leaks out as area(...)::Any whenever the trailing Float64(...) conversion can't re-pin the type.

I hit this downstream in ConservativeRegridding.jl, where it boxed ~4000 Float64s per call inside a quadrature loop (a ::Float64 call-site barrier worked around it — ~36× fewer allocations in that assembly).

MWE

using GeometryOps, GeoInterface, StaticArrays
const GO = GeometryOps; const GI = GeoInterface
const USP = GO.UnitSpherical.UnitSphericalPoint

poly = GI.Polygon([GI.LinearRing([
    USP(0.0, 0.0, 1.0), USP(1.0, 0.0, 0.0), USP(0.0, 1.0, 0.0), USP(0.0, 0.0, 1.0),
])])

Base.return_types(GO.area, (typeof(GO.Spherical()), typeof(poly)))
# 1-element Vector{Any}:
#  Any          # <- expected Float64

GO.area(GO.Spherical(), poly)   # value is correct: 6.375823512160898e13

Independent of the ring's point container — Vector- and SVector{4}-backed LinearRings both infer Any.

Root cause

src/methods/area.jl, the NaiveTriangulatedSphericalArea method (≈ area.jl:259):

function area(alg::NaiveTriangulatedSphericalArea, geom, ::Type{T} = Float64; ...) where T <: AbstractFloat
    function _polygon_area(trait::GI.PolygonTrait, alg::SphericalTriangleAreaMethod, poly)   # method 1
        ...
    end
    _polygon_area(::GI.PointTrait, alg, geom) = zero(T)                                       # method 2

    unit_area = applyreduce(
        WithTrait((trait, g) -> _polygon_area(trait, alg.method, g)),   # <- captures multi-method local
        +, TraitTarget{Union{GI.PolygonTrait, GI.PointTrait}}(), geom; init=zero(T), ...,
    )
    return T(unit_area * manifold(alg).radius^2)
end

_polygon_area has two methods and is captured by the (trait, g) -> … lambda. Julia lowers a multi-method local function that is also captured into a Core.Box (:contents::Any); @code_warntype shows _polygon_area::Core.Box, so the applyreduce(WithTrait(...), …) call returns Any and unit_area::Any.

Minimal demonstration that the trigger is specifically multi-method local + capture — not T, alg, the radius, or the ring math (_naive_triangulated_spherical_ring_area infers Float64, and Spherical().radius::Float64):

function probe(geom, ::Type{T}) where T<:AbstractFloat
    f(::Val{1}, x) = one(T)
    f(::Val{2}, x) = zero(T)            # 2nd method makes `f` a *multi-method* local
    applyreduce(g -> f(Val(1), g), +, GI.PolygonTrait(), geom; init=zero(T))
end
Base.return_types(probe, (typeof(poly), Type{Float64}))   # Any[Any]
# Delete the Val(2) method (single-method local), OR pass `f` directly with no
# wrapping lambda, and it infers Float64.

Two caveats on reproducing:

• The trailing T(unit_area * radius²) at ≈ area.jl:282 is meant to pin the result but cannot recover once unit_area::Any. In environments where another loaded package contributes Float64 constructor methods (in my case a ClimaCore.Utilities.AutoBroadcaster method shows up in Base.return_types(Float64, (Any,))), Float64(::Any) is itself non-inferrable, so the conversion is no safety net at all and the Any leaks all the way out (as in the MWE above).
• Even in environments where that conversion does mask the final return type as Float64, the internal Core.Box allocation and dynamic dispatch through the boxed _polygon_area remain — so it's a per-call performance bug regardless.

Suggested fix

De-box the closure — either is verified to restore Float64 with identical numeric output:

1. Lift _polygon_area to a top-level @inline helper taking method / ::Type{T} as arguments (no closure); or
2. Collapse _polygon_area to a single local that branches on the trait internally (trait isa GI.PolygonTrait ? … : zero(T)) and pass it directly to WithTrait (drop the (trait, g) -> … wrapper).

Versions

• GeometryOps v0.1.40, GeometryOpsCore v0.1.10
• Julia 1.12.6