Merge pull request #17 from JuliaDiffEq/post_interp

Post-Solution interpolation
SciML · May 12, 2017 · 396480a · 396480a
2 parents b5a70bd + ee8f429
commit 396480a
Show file tree

Hide file tree

Showing 2 changed files with 20 additions and 9 deletions.
diff --git a/src/common.jl b/src/common.jl
@@ -17,14 +17,13 @@ export daskr
 function solve{uType,duType,tType,isinplace,LinearSolver}(
     prob::AbstractDAEProblem{uType,duType,tType,isinplace},
     alg::DASKRDAEAlgorithm{LinearSolver},
-    timeseries = [], ts = [], ks = [];
-
+    timeseries = [], ts = [], ks = []; dense = true,
     verbose=true,
     callback = nothing, abstol = 1/10^6, reltol = 1/10^3,
     saveat = Float64[], adaptive = true, maxiter = Int(1e5),
     timeseries_errors = true, save_everystep = isempty(saveat),
     save_start = true, save_timeseries = nothing,
-    userdata = nothing, isdiff = fill(true, length(prob.u0)), kwargs...)
+    userdata = nothing, kwargs...)
 
     verbose && !isempty(kwargs) && check_keywords(alg, kwargs, warnlist)
 
@@ -92,7 +91,12 @@ function solve{uType,duType,tType,isinplace,LinearSolver}(
                           u = vec(u); du=vec(du); 0)
     end
 
-    id = Int32[x ? 1 : -1 for x in isdiff]
+    if prob.differential_vars == nothing
+      id = ones(Int32,length(u0))
+    else
+      id = Int32[x ? 1 : -1 for x in prob.differential_vars]
+    end
+
     tstart = 0.0
     tstop = 50.0
     Nsteps = 500
@@ -126,9 +130,11 @@ function solve{uType,duType,tType,isinplace,LinearSolver}(
     psol = Int32[0]
 
     ures = Vector{Vector{Float64}}()
+    dures = Vector{Vector{Float64}}()
     save_start ? ts = [t0] : ts = Float64[]
     save_start ? start_idx = 1 : start_idx = 2
     save_start && push!(ures, copy(u0))
+    save_start && dense && push!(dures, copy(du0))
 
     u = copy(u0)
     du = copy(du0)
@@ -139,6 +145,9 @@ function solve{uType,duType,tType,isinplace,LinearSolver}(
             DASKR.unsafe_solve(res, N, t, u, du, tout, info, rtol, atol, idid, rwork, lrw, iwork, liw, rpar, ipar, jac, psol, rt, nrt, jroot)
             push!(ures,copy(u))
             push!(ts, t[1])
+            if dense
+              push!(dures,copy(du))
+            end
         end
     end
     ### Finishing Routine
@@ -157,6 +166,8 @@ function solve{uType,duType,tType,isinplace,LinearSolver}(
     end
 
     build_solution(prob,alg,ts,timeseries,
+                      du = dures,
+                      dense = dense,
                       timeseries_errors = timeseries_errors,
                       retcode = :Success)
 end
diff --git a/test/runtests.jl b/test/runtests.jl
@@ -70,13 +70,13 @@ let
     @test length(sol.t) > 2
     sol = solve(prob, daskr(),save_everystep=false)
     @test length(sol.u) == length(sol.t) == 2
-    sol = solve(prob, daskr(), saveat = saveat, isdiff = [true, true, false])
+    prob2 = DAEProblem(resrob,u0,du0,(0.0,100000.0),differential_vars = [true, true, false])
+    sol = solve(prob2, daskr(), saveat = saveat)
     @test sol.t == saveat
-    sol = solve(prob, daskr(), saveat = dt, isdiff = [true, true, false])
+    sol = solve(prob2, daskr(), saveat = dt)
     @test sol.t == saveat
-    sol = solve(prob, daskr(), saveat = saveat,
-                      save_everystep = true,
-                      isdiff = [true, true, false])
+    sol = solve(prob2, daskr(), saveat = saveat,
+                      save_everystep = true)
     @test minimum([t ∈ sol.t for t in saveat])
     sol = solve(prob, daskr(), saveat = saveat, save_everystep = true)
     @test intersect(sol.t, saveat) == saveat