# Arc length for a spline

**URL:** <https://discourse.julialang.org/t/arc-length-for-a-spline/36161>\
**Category:** General Usage\
**Tags:** question, quadgk, splines\
**Created:** [March 18, 2020, 9:31pm UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161 "2020-03-18T21:31:08Z")\
**Posts on this page:** 8\
**Page:** 2

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [February 5, 2025, 12:13am UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/21 "2025-02-05T00:13:03Z")

</div>

> [@danielwe](#):
>
> Can’t you just call `quadgk(f, 0, knots..., pi; kws...)`?

You can, but splatting a large array is problematic too.

~~I’ve been meaning to add an API for this so that you at least don’t need to call the undocumented `QuadGK.Segment` constructor, e.g. a new method for `alloc_segbuf` that lets you pass an array of points.~~ There’s already a documented API `quadgk(f, knots; kws...)` (no splatting), see below.

---

<div class="post-metadata">

**Author:** ![danielwe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielwe/32/35657_2.png) [@danielwe](https://discourse.julialang.org/u/danielwe)\
**Post date:** [February 5, 2025, 12:14am UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/22 "2025-02-05T00:14:01Z")

</div>

I think you already did. Looks like `quadgk(f, knots; kws...)` (no splatting) should work too: [QuadGK.jl/src/api.jl at master · JuliaMath/QuadGK.jl · GitHub](https://github.com/JuliaMath/QuadGK.jl/blob/master/src/api.jl#L240-L268)

You called the helper `to_segbuf` rather than adding a method to `alloc_segbuf`, but since users can pass the array to `quadgk` directly that doesn’t seem important.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [February 5, 2025, 2:01am UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/23 "2025-02-05T02:01:02Z")

</div>

> [@danielwe](#):
>
> I think you already did.

Doh, I forgot! It’s [even documented](https://github.com/JuliaMath/QuadGK.jl/blob/89a3b8fbeb8f6861617e39e510647bd8111dc584/src/api.jl#L27-L30).

---

<div class="post-metadata">

**Author:** ![yakir12](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/yakir12/32/297_2.png) [@yakir12](https://discourse.julialang.org/u/yakir12)\
**Post date:** [February 5, 2025, 7:40am UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/24 "2025-02-05T07:40:07Z")

</div>

I guess even more improvements might come once that PR gets cleared and merged (thank you for that!), but here are the corrected benchmarks:

```julia
julia> @btime method1($spl, extrema(t)...)
  5.940 ms (60 allocations: 9.16 MiB)
70.44577023545651

julia> @btime method2($spl, extrema(t)...)
  1.424 ms (51046 allocations: 1.51 MiB)
70.44577026541424

julia> @btime method2a($spl)
  134.311 μs (4983 allocations: 151.69 KiB)
70.44577024751408

julia> @btime method3($spl, $(gauss(15))...)
  83.455 μs (3354 allocations: 101.50 KiB)
70.44577024750694

julia> @btime method4($spl, order=4, rtol=1e-5)
  69.240 μs (2467 allocations: 75.44 KiB)
70.44577023146906

```

These have been compared so the relative accuracy is ~1e-10 (and I tried to fix `method4` following your discussion, hope I didn’t mess it up):

```julia
function method1(spl, t1, t2)
    tl = LinRange(t1, t2, 100_000)
    xy = spl(tl)
    Δxy = xy[:,2:end] - xy[:,1:end-1]
    sum(norm.(eachcol(Δxy)))
end
function method2(spl, t1, t2)
    res, _ = quadgk(t -> norm(derivative(spl, t)), t1, t2, rtol = 1e-8)
    res
end
function method2a(spl)
    knots = Dierckx.get_knots(spl)
    s = 0.0
    for (k1, k2) in zip(knots[1:end-1], knots[2:end])
        res, _ = quadgk(t -> norm(derivative(spl, t)), k1, k2)
        s += res
    end
    return s
end
function method3(spl, x_gauss, w_gauss)
    knots = Dierckx.get_knots(spl)
    integral = sum(eachindex(knots)[2:end]) do i
        t1, t2 = knots[i-1], knots[i]
        scale = (t2 - t1) / 2
        sum(zip(x_gauss, w_gauss)) do ((xg, wg))
            t = (xg + 1) * scale + t1 # map from (-1,1) to (t1,t2)
            norm(derivative(spl, t)) * wg
        end * scale
    end
end
function method4(spl; kws...)
    knots = Dierckx.get_knots(spl)
    s, _ = quadgk(t -> norm(derivative(spl, t)), knots; kws...)
    return s
end

row = df[5, :]
spl = row.spl
t = row.t
s, _ = quadgk(t -> norm(derivative(spl, t)), extrema(t)..., rtol = 1e-9)
rel(x) = (x - s)/s

rel(method1(spl, extrema(t)...))
rel(method2(spl, extrema(t)...))
rel(method2a(spl))
rel(method3(spl, gauss(15)...))
rel(method4(spl, order=4, rtol=1e-5))

```

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [February 5, 2025, 10:15am UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/25 "2025-02-05T10:15:22Z")

</div>

For future reference, if you need a higher-level interface for integrals over geometries, including `BezierCurve` and `ParametrizedCurve`, consider the MeshIntegrals.jl package developed and maintained by @mike.ingold , @JoshuaLampert and @kylebeggs

> **[About · MeshIntegrals.jl](https://juliageometry.github.io/MeshIntegrals.jl/stable/)**
>
> Documentation for MeshIntegrals.jl.

It relies on QuadGK.jl, FastGaussQuadrature.jl and HCubature.jl as backends, and can facilitate the lives of developers who need to perform integration over curves, surfaces or volumes.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [February 5, 2025, 12:35pm UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/26 "2025-02-05T12:35:03Z")

</div>

> [@yakir12](#):
>
> `@btime method3($spl, $gauss(15)...)`

I think you need more parens `@btime method3($spl, $(gauss(15))...)`, otherwise you are only interpolating the function name.

---

<div class="post-metadata">

**Author:** ![yakir12](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/yakir12/32/297_2.png) [@yakir12](https://discourse.julialang.org/u/yakir12)\
**Post date:** [February 5, 2025, 12:43pm UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/27 "2025-02-05T12:43:57Z")

</div>

Oops, you’re right, I edited it now.

---

<div class="post-metadata">

**Author:** ![yakir12](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/yakir12/32/297_2.png) [@yakir12](https://discourse.julialang.org/u/yakir12)\
**Post date:** [February 7, 2025, 8:04am UTC](https://discourse.julialang.org/t/arc-length-for-a-spline/36161/28 "2025-02-07T08:04:11Z")

</div>

Just to address your excellent point: I need a _smoothing_ spline, hence the use of `Dierckx.jl` (and not `Interpolations.jl` for instance).

In other news, I’ve updated the winning implementation, “`method4`”, to allow for bounded arc-lengths:

```julia
function arclength(spl, t1, t2; kws...)
    knots = get_knots(spl)
    filter!(t -> t1 < t < t2, knots)
    pushfirst!(knots, t1)
    push!(knots, t2)
    s, _ = quadgk(t -> norm(derivative(spl, t)), knots; kws...)
    return s
end

```

[Previous page](https://discourse.julialang.org/t/arc-length-for-a-spline/36161.md?page=1)