DFBDF performs poorly when integrating sin(t)

[topolarity]

Describe the bug 🐞

DFBDF appears to solve the system dx ~ y, y ~ sin(t) very poorly.

Minimal Reproducible Example 👇

julia> using OrdinaryDiffEq, Plots
julia> dae_f(du, u, _, t) = [du[1] - u[2], u[2] - sin(t)]
julia> prob = DAEProblem(DAEFunction(dae_f), [0.0, 0.0], [0.0, 0.0], (0., 2π), nothing, differential_vars=[true,false], abstol=1e-6, reltol=1e-6)
julia> sol = solve(prob, DFBDF(); maxiters=10_000_000)
retcode: Success
Interpolation: 1st order linear
t: 1359394-element Vector{Float64}:
 0.0
 6.283185307179586e-6
 7.924885500379142e-6
 8.703751558664616e-6
 ⋮
u: 1359394-element Vector{Vector{Float64}}:
 [0.0, 0.0]
 [0.0, 0.0]
 [1.0315106532637465e-11, 6.283185307138244e-6]
 [1.6487530864616867e-11, 7.92488550029619e-6]
  ⋮

julia> t′ = 0.0:1e-3:2π
julia> plot(t′, -cos.(t′) .+ 1.0 .- [u[1] for u in sol(t′)]; label="u[1] error")
julia> plot!(t′, sin.(t′) - [u[2] for u in sol(t′)]; label="u[2] error")

The solution with abstol = 1e-3, reltol = 1e-3 is also noticeably bad:

A few things jump out to me here:

• 1.3 million points is... a lot of time steps
• u[2] seems to be solved very precisely but for the previous timestep (for the first few time points at least)
• the error in u is about [2e-5, 5e-5] which is ~one and a half orders of magnitude worse than the tolerances

Environment:

• Output of using Pkg; Pkg.status()

(tmp) pkg> status
Status `~/repos/tmp/tmp/Project.toml`
  [459566f4] DiffEqCallbacks v2.36.1
  [1dea7af3] OrdinaryDiffEq v6.70.0
  [c3572dad] Sundials v4.23.2

[topolarity]

For comparison, IDA solves this in just 115 timesteps with much better results: