From 9ebfc08478395d788f167a992b2e66790f65d437 Mon Sep 17 00:00:00 2001 From: ScottPJones Date: Fri, 12 May 2017 11:35:58 -0400 Subject: [PATCH] Add checks for passing explicit t-gradient or Jacobian to stiff solver --- src/common.jl | 96 ++++++++++++++++++++++++++++++--------------------- 1 file changed, 56 insertions(+), 40 deletions(-) diff --git a/src/common.jl b/src/common.jl index 6120c29..9153d42 100644 --- a/src/common.jl +++ b/src/common.jl @@ -8,14 +8,29 @@ function solve{uType,tType,isinplace}( verbose=true, abstol=1/10^6,reltol=1/10^3, tstops=Float64[], - saveat=Float64[],maxiter=Int(1e5), + saveat=Float64[], maxiter=Int(1e5), save_start=true, - callback = nothing, - timeseries_errors=true,save_everystep= isempty(saveat), - save_timeseries = nothing, - userdata=nothing,kwargs...) - - verbose && !isempty(kwargs) && check_keywords(alg, kwargs, warnlist) + callback=nothing, + timeseries_errors=true, + save_everystep=isempty(saveat), + save_timeseries=nothing, + userdata=nothing, + kwargs...) + + if verbose + warned = false + isempty(kwargs) || check_keywords(alg, kwargs, warnlist) + if has_tgrad(prob.f) + warn("Explicit t-gradient given to this stiff solver is ignored.") + warned = true + end + if has_jac(prob.f) + warn("Explicit Jacobian given to this stiff solver is ignored.") + warned = true + end + warned && + warn("See http://docs.juliadiffeq.org/latest/basics/compatibility_chart.html") + end if save_timeseries != nothing warn("save_timeseries is deprecated. Use save_everystep instead") @@ -131,44 +146,45 @@ function solve{uType,tType,isinplace}( lsoda_prepare(ctx,opt) for k in 2:length(save_ts) - ttmp[1] = save_ts[k] - if t[1]ttmp[1] # overstepd, interpolate back + ttmp[1] = save_ts[k] + if t[1] < ttmp[1] + while t[1] < ttmp[1] + lsoda(ctx, utmp, t, ttmp[1]) + if t[1] > ttmp[1] # overstepd, interpolate back + t2[1] = t[1] # save step values + copy!(utmp2,utmp) # save step values + opt.itask = 1 # change to interpolating + lsoda(ctx, utmp, t, ttmp[1]) + opt.itask = itask_tmp + push!(ures, copy(utmp)) + push!(ts, t[1]) + # don't overstep the last timestep + if k != length(save_ts) && save_ts[k+1] > t2[1] + push!(ures, copy(utmp2)) + push!(ts, t2[1]) + end + copy!(utmp, utmp2) + t[1] = t2[1] + else + push!(ures, copy(utmp)) + push!(ts,t[1]) + end + end + else t2[1] = t[1] # save step values - copy!(utmp2,utmp) # save step values + copy!(utmp2, utmp) # save step values opt.itask = 1 # change to interpolating - lsoda(ctx,utmp,t,ttmp[1]) + lsoda(ctx, utmp, t, ttmp[1]) opt.itask = itask_tmp - push!(ures,copy(utmp)) - push!(ts,t[1]) + push!(ures, copy(utmp)) + push!(ts, t[1]) if k != length(save_ts) && save_ts[k+1] > t2[1] # don't overstep the last timestep - push!(ures,copy(utmp2)) - push!(ts,t2[1]) + push!(ures,copy(utmp2)) + push!(ts,t2[1]) end copy!(utmp,utmp2) t[1] = t2[1] - else - push!(ures,copy(utmp)) - push!(ts,t[1]) - end - end - else - t2[1] = t[1] # save step values - copy!(utmp2,utmp) # save step values - opt.itask = 1 # change to interpolating - lsoda(ctx,utmp,t,ttmp[1]) - opt.itask = itask_tmp - push!(ures,copy(utmp)) - push!(ts,t[1]) - if k != length(save_ts) && save_ts[k+1] > t2[1] # don't overstep the last timestep - push!(ures,copy(utmp2)) - push!(ts,t2[1]) end - copy!(utmp,utmp2) - t[1] = t2[1] - end end ### Finishing Routine @@ -185,7 +201,7 @@ function solve{uType,tType,isinplace}( end end - build_solution(prob,alg,ts,timeseries, - timeseries_errors = timeseries_errors, - retcode = :Success) + build_solution(prob, alg, ts, timeseries, + timeseries_errors = timeseries_errors, + retcode = :Success) end