Not correctly matching the SciML interface

[ChrisRackauckas]

The value of variable 2 at time 0.0005

sol(2,0.0005)

That is incorrect. It should be matching the interface, i.e. sol(0.0005, idxs = 2)

[ChrisRackauckas]

# The total number of steps to end the simulation
sol.totalSteps
498
# The number of simultaneous steps during the simulation
sol.simulStepCount
132
# The total number of events during the simulation
sol.evCount
53
# The actual data is stored in two vectors:
sol.savedTimes
sol.savedVars

Those would be designed to be in the sol.stats

[ChrisRackauckas]

Is there a reason why this isn't able to take a standard ODEProblem? Could the codegen to QSS be handled automatically in ModelingToolkit?

[mongibellili]

It should be matching the interface, i.e. sol(0.0005, idxs = 2) Those would be designed to be in the sol.stats

Noted. I will add these.

Is there a reason why this isn't able to take a standard ODEProblem?

At a time step, the QSS method does not update all variables (it needs to decide which variables to update). So, it uses some dependency matrices (Vector{Vector{Int}} in code) and a personalized function that holds the differential equations and events in a special format**. These datastructures are stored in a NLODEProblem{PRTYPE,T,Z,D,CS} wich is different from a standard ODEProblem. However, we can have it as a subtype then make the solve function dispatch on the QSS algorithm and NLODEProblem.

** The QSS solver uses variables of type "Taylor" from the TaylorSeries.jl for approximation to make calculation easier and cheaper. But all operations from the TaylorSeries.jl allocate, so I redid internally this package with added cache to every operation. To exploit this, I had to change the AST of the differential equations upfront to obtain the personalized function. Example:

du[1] = t
for k in 2:5 
        du[k]=discrete[1]*(u[k]-u[k-1]) ;
end
if t-5.0>0.0
        discrete[1]=0.0
end

Becomes

RuntimeGeneratedFunction{(:i, :zc, :ev, :q, :d, :t, :cache)
    if i == 1
        createT(t, cache[1])  
   elseif 2 <= i <= 5
        mulT(d[1], subT(q[i], q[i - 1], cache[2]), cache[1])
   end
   if zc == 1
        subT(t, 5.0, cache[1])
   end
   if ev == 1
        d[1] = 0.0
   end
  return nothing

Could the codegen to QSS be handled automatically in ModelingToolkit?

As mentioned above, QSS needs the datastructures stored in the NLODEProblem. So if ModelingToolkit can transform a set of differential equations with events (if-statements) to a NLODEProblem{PRTYPE,T,Z,D,CS} then YES.