Skip to content

Commit 138f5f2

Browse files
committed
Fix minor bug if DAE integrator is used
1 parent c060428 commit 138f5f2

File tree

1 file changed

+25
-24
lines changed

1 file changed

+25
-24
lines changed

src/SimulateAndPlot.jl

Lines changed: 25 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -372,32 +372,33 @@ end
372372

373373
function simulateSegment!(m::SimulationModel{FloatType,TimeType}, algorithm=missing; kwargs...) where {FloatType,TimeType}
374374
solution = nothing
375+
options = m.options
375376

376377
sizesOfLinearEquationSystems = Int[length(leq.b) for leq in m.linearEquations]
377378

378379
# Define problem and callbacks based on algorithm and model type
379-
interval = m.options.interval
380-
if abs(m.options.stopTime - m.options.startTime) <= 0
380+
interval = options.interval
381+
if abs(options.stopTime - options.startTime) <= 0
381382
interval = 1.0
382-
tspan2 = [m.options.startTime]
383-
elseif abs(m.options.interval) < abs(m.options.stopTime-m.options.startTime)
383+
tspan2 = [options.startTime]
384+
elseif abs(options.interval) < abs(options.stopTime-options.startTime)
384385
if m.nsegments == 1
385-
tspan2 = m.options.startTime:interval:m.options.stopTime
386+
tspan2 = options.startTime:interval:options.stopTime
386387
else
387-
i = ceil( (m.options.startTime - m.options.startTimeFirstSegment)/interval )
388-
tnext = m.options.startTimeFirstSegment + i*interval
389-
tspan2 = tnext:interval:m.options.stopTime
390-
if tspan2[1] > m.options.startTime
391-
tspan2 = [m.options.startTime, tspan2...]
388+
i = ceil( (options.startTime - options.startTimeFirstSegment)/interval )
389+
tnext = options.startTimeFirstSegment + i*interval
390+
tspan2 = tnext:interval:options.stopTime
391+
if tspan2[1] > options.startTime
392+
tspan2 = [options.startTime, tspan2...]
392393
end
393394
end
394-
if tspan2[end] < m.options.stopTime
395-
tspan2 = [tspan2..., m.options.stopTime]
395+
if tspan2[end] < options.stopTime
396+
tspan2 = [tspan2..., options.stopTime]
396397
end
397398
else
398-
tspan2 = [m.options.startTime, m.options.stopTime]
399+
tspan2 = [options.startTime, options.stopTime]
399400
end
400-
tspan = (m.options.startTime, m.options.stopTime)
401+
tspan = (options.startTime, options.stopTime)
401402

402403
eh = m.eventHandler
403404
m.odeMode = true
@@ -450,9 +451,9 @@ function simulateSegment!(m::SimulationModel{FloatType,TimeType}, algorithm=miss
450451
# FunctionalCallingCallback(outputs!, ...) is not correctly called when zero crossings are present.
451452
# The fix is to call outputs!(..) from the previous to the current event, when an event occurs.
452453
# (alternativey: callback4 = DifferentialEquations.PresetTimeCallback(tspan2, affect_outputs!) )
453-
callback1 = DifferentialEquations.FunctionCallingCallback(outputs!, funcat=[m.options.startTime]) # call outputs!(..) at startTime
454+
callback1 = DifferentialEquations.FunctionCallingCallback(outputs!, funcat=[options.startTime]) # call outputs!(..) at startTime
454455
callback3 = DifferentialEquations.VectorContinuousCallback(zeroCrossings!,
455-
affectStateEvent!, eh.nz, interp_points=m.options.interp_points, rootfind=DifferentialEquations.SciMLBase.RightRootFind)
456+
affectStateEvent!, eh.nz, interp_points=options.interp_points, rootfind=DifferentialEquations.SciMLBase.RightRootFind)
456457
#callback4 = DifferentialEquations.PresetTimeCallback(tspan2, affect_outputs!)
457458
callbacks = DifferentialEquations.CallbackSet(callback1, callback2, callback3) #, callback4)
458459
else
@@ -463,24 +464,24 @@ function simulateSegment!(m::SimulationModel{FloatType,TimeType}, algorithm=miss
463464

464465
# Initial step size (the default of DifferentialEquations integrators is too large) + step-size of fixed-step algorithm
465466
if !m.sundials
466-
dt = m.options.adaptive ? m.options.interval/10 : m.options.interval # initial step-size
467+
dt = options.adaptive ? options.interval/10 : options.interval # initial step-size
467468
end
468469

469470
# Compute solution
470-
abstol = 0.1*m.options.tolerance
471+
abstol = 0.1*options.tolerance
471472
tstops = (m.eventHandler.nextEventTime,)
472473
maxiters = Int(typemax(Int32)) # switch off maximum number of iterations (typemax(Int) gives an inexact error for Sundials)
473474
if ismissing(algorithm)
474-
TimerOutputs.@timeit m.timer "DifferentialEquations.solve" solution = DifferentialEquations.solve(problem, reltol=m.options.tolerance, abstol=abstol, save_everystep=false,
475-
callback=callbacks, adaptive=m.options.adaptive, saveat=tspan2, dt=dt, dtmax=m.options.dtmax, maxiters=maxiters, tstops = tstops,
475+
TimerOutputs.@timeit m.timer "DifferentialEquations.solve" solution = DifferentialEquations.solve(problem, reltol=options.tolerance, abstol=abstol, save_everystep=false,
476+
callback=callbacks, adaptive=options.adaptive, saveat=tspan2, dt=dt, dtmax=options.dtmax, maxiters=maxiters, tstops = tstops,
476477
initializealg = DifferentialEquations.NoInit())
477478
elseif m.sundials
478-
TimerOutputs.@timeit m.timer "DifferentialEquations.solve" solution = DifferentialEquations.solve(problem, algorithm, reltol=m.options.tolerance, abstol=abstol, save_everystep=false,
479-
callback=callbacks, adaptive=m.options.adaptive, saveat=tspan2, dtmax=m.options.dtmax, maxiters=maxiters, tstops = tstops,
479+
TimerOutputs.@timeit m.timer "DifferentialEquations.solve" solution = DifferentialEquations.solve(problem, algorithm, reltol=options.tolerance, abstol=abstol, save_everystep=false,
480+
callback=callbacks, adaptive=options.adaptive, saveat=tspan2, dtmax=options.dtmax, maxiters=maxiters, tstops = tstops,
480481
initializealg = DifferentialEquations.NoInit())
481482
else
482-
TimerOutputs.@timeit m.timer "DifferentialEquations.solve" solution = DifferentialEquations.solve(problem, algorithm, reltol=m.options.tolerance, abstol=abstol, save_everystep=false,
483-
callback=callbacks, adaptive=m.options.adaptive, saveat=tspan2, dt=dt, dtmax=m.options.dtmax, maxiters=maxiters, tstops = tstops,
483+
TimerOutputs.@timeit m.timer "DifferentialEquations.solve" solution = DifferentialEquations.solve(problem, algorithm, reltol=options.tolerance, abstol=abstol, save_everystep=false,
484+
callback=callbacks, adaptive=options.adaptive, saveat=tspan2, dt=dt, dtmax=options.dtmax, maxiters=maxiters, tstops = tstops,
484485
initializealg = DifferentialEquations.NoInit())
485486
end
486487

0 commit comments

Comments
 (0)