DiscreteCallback is not called after VectorContinuousCallback is called

[levyung]

Describe the bug 🐞

I have a VectorContinuousCallback and a DiscreteCallback. The DiscreteCallback condition is always true; I set its save_positions = (false, false) so that I do not get all time steps doubled. The DiscreteCallback affect is to increase a counter by one, which initially is set to zero. By the end of the evolution, the size of sol.t and counter differ by the number of VectorContinuousCallback events + 1.

Expected behavior

The docs say that DiscreteCallbacks are called after continuous callbacks are called. Thus, I expected the counter to be equal to the size of sol.t array minus one (as I understand the first time step is ignored by discrete callbacks). In the concrete example below I expect the counter to be equal to 25, as length(sol.t) is 26.

Minimal Reproducible Example 👇

using DifferentialEquations

# Define the ODE system (harmonic oscillator)
function harmonic!(du, u, p, t)
    du[1] = u[2]     # velocity
    du[2] = -u[1]    # acceleration
end

# Define multiple conditions to check
function conditions(out, u, t, integrator)
    out[1] = u[1] - 0.5  # Position crosses 0.5
end

# Define what happens at condition crossings
function affect!(integrator, event_index)
    println("Event at t = $(integrator.t)")
end

# Discrete callback 
function discrete_condition(u, t, integrator)
    true  
end

# Counter for calling discrete callback
counter = 0 

# Increasing it every discrete callback
function discrete_affect!(integrator)
    global counter += 1
end

# Setup callbacks
vcb = VectorContinuousCallback(conditions, affect!, 2)
dcb = DiscreteCallback(discrete_condition, discrete_affect!, 
    save_positions = (false, false))

# Combine callbacks
cb = CallbackSet(vcb, dcb)

# Initial conditions and solve
u0 = [1.0, 0.0]
tspan = (0.0, 10.0)
prob = ODEProblem(harmonic!, u0, tspan)
sol = solve(prob, Tsit5(), callback=cb)

sol.t |> length |> println #26
counter |> println #22
# number of continuous events = 3

Error & Stacktrace ⚠️

No error.

Environment (please complete the following information):

• Output of using Pkg; Pkg.status()

Status `C:\Users\levyu\.julia\environments\frg\Project.toml`
  [459566f4] DiffEqCallbacks v4.3.0
  [0c46a032] DifferentialEquations v7.16.0
  [7c1d4256] DynamicPolynomials v0.6.1
  [f6369f11] ForwardDiff v0.10.38
  [f213a82b] HomotopyContinuation v2.12.0
  [de52edbc] Integrals v4.5.0
⌃ [d1acc4aa] IntervalArithmetic v0.22.22
  [d2bf35a9] IntervalRootFinding v0.6.0
⌃ [961ee093] ModelingToolkit v9.61.0
⌃ [8913a72c] NonlinearSolve v4.3.0
  [91a5bcdd] Plots v1.40.9
⌃ [f27b6e38] Polynomials v4.0.13
  [438e738f] PyCall v1.96.4
  [1a520dc8] SimpleRationalFunctions v0.0.2
⌃ [24249f21] SymPy v2.2.1
⌃ [0c5d862f] Symbolics v6.24.0
Info Packages marked with ⌃ have new versions available and may be upgradable.

• Output of using Pkg; Pkg.status(; mode = PKGMODE_MANIFEST)

Status `C:\Users\levyu\.julia\environments\frg\Manifest.toml`
  [47edcb42] ADTypes v1.13.0
  [398f06c4] AbstractLattices v0.3.1
  [1520ce14] AbstractTrees v0.4.5
  [7d9f7c33] Accessors v0.1.41
  [79e6a3ab] Adapt v4.1.1
  [66dad0bd] AliasTables v1.1.3
  [a95523ee] AlmostBlockDiagonals v0.1.10
  [fb37089c] Arblib v1.2.1
  [ec485272] ArnoldiMethod v0.4.0
  [4fba245c] ArrayInterface v7.18.0
  [4c555306] ArrayLayouts v1.11.0
  [aae01518] BandedMatrices v1.9.1
  [e2ed5e7c] Bijections v0.1.9
  [d1d4a3ce] BitFlags v0.1.9
  [62783981] BitTwiddlingConvenienceFunctions v0.1.6
  [8e7c35d0] BlockArrays v1.4.0
  [764a87c0] BoundaryValueDiffEq v5.15.0
  [7227322d] BoundaryValueDiffEqAscher v1.4.0
  [56b672f2] BoundaryValueDiffEqCore v1.7.0
  [85d9eb09] BoundaryValueDiffEqFIRK v1.5.0
  [1a22d4ce] BoundaryValueDiffEqMIRK v1.5.0
  [9255f1d6] BoundaryValueDiffEqMIRKN v1.4.0
  [ed55bfe0] BoundaryValueDiffEqShooting v1.5.0
  [70df07ce] BracketingNonlinearSolve v1.1.0
...
  [8e850b90] libblastrampoline_jll v5.11.0+0
  [8e850ede] nghttp2_jll v1.59.0+0
  [3f19e933] p7zip_jll v17.4.0+2
Info Packages marked with ⌃ and ⌅ have new versions available. Those with ⌃ may be upgradable, but those with ⌅ are restricted by compatibility constraints from upgrading. To see why use `status --outdated -m`

• Output of versioninfo()

Julia Version 1.11.3
Commit d63adeda50 (2025-01-21 19:42 UTC)
Build Info:
  Official https://julialang.org/ release
Platform Info:
  OS: Windows (x86_64-w64-mingw32)
  CPU: 12 × AMD Ryzen 5 4600H with Radeon Graphics         
  WORD_SIZE: 64
  LLVM: libLLVM-16.0.6 (ORCJIT, znver2)
Threads: 1 default, 0 interactive, 1 GC (on 12 virtual cores)
Environment:
  JULIA_NUM_THREADS = 

[levyung]

The same bug occurs for other discrete callbacks like SavingCallback and IntegratingCallback from https://docs.sciml.ai/DiffEqCallbacks/dev/integrating/
All of them skip the the time steps produced by continuous callbacks.

[levyung]

Closing it because I moved it to DifferentialEquations.jl