seeded random problems and solutions

SciML · Jun 26, 2017 · d2f221b · d2f221b
1 parent 2ec6fa1
commit d2f221b
Show file tree

Hide file tree

Showing 4 changed files with 28 additions and 14 deletions.
diff --git a/src/problems/noise_problems.jl b/src/problems/noise_problems.jl
@@ -1,5 +1,8 @@
 type NoiseProblem{N<:AbstractNoiseProcess,T} <: AbstractNoiseProblem
   noise::N
   tspan::Tuple{T,T}
+  seed::UInt64
 end
-NoiseProblem(noise,tspan) = NoiseProblem{typeof(noise),promote_type(map(typeof,tspan)...)}(noise,tspan)
+NoiseProblem(noise,tspan;seed=UInt64(0)) = NoiseProblem{typeof(noise),
+                                                      promote_type(map(typeof,tspan)...)
+                                                      }(noise,tspan,seed)
diff --git a/src/problems/rode_problems.jl b/src/problems/rode_problems.jl
@@ -6,11 +6,15 @@ type RODEProblem{uType,tType,isinplace,NP,F,C,MM,ND} <: AbstractRODEProblem{uTyp
   callback::C
   mass_matrix::MM
   rand_prototype::ND
+  seed::UInt64
 end
 
 function RODEProblem(f,u0,tspan; iip = isinplace(f,4),
                       rand_prototype = nothing,
-                      noise= nothing,
+                      noise= nothing, seed = UInt64(0),
                      callback=nothing,mass_matrix=I)
-  RODEProblem{typeof(u0),promote_type(map(typeof,tspan)...),iip,typeof(noise),typeof(f),typeof(callback),typeof(mass_matrix),typeof(rand_prototype)}(f,u0,tspan,noise,callback,mass_matrix,rand_prototype)
+  RODEProblem{typeof(u0),promote_type(map(typeof,tspan)...),
+              iip,typeof(noise),typeof(f),typeof(callback),
+              typeof(mass_matrix),typeof(rand_prototype)}(
+              f,u0,tspan,noise,callback,mass_matrix,rand_prototype,seed)
 end
diff --git a/src/problems/sde_problems.jl b/src/problems/sde_problems.jl
@@ -7,15 +7,18 @@ type SDEProblem{uType,tType,isinplace,NP,F,G,C,MM,ND} <: AbstractSDEProblem{uTyp
   callback::C
   mass_matrix::MM
   noise_rate_prototype::ND
+  seed::UInt64
 end
 
 function SDEProblem(f,g,u0,tspan;iip = isinplace(f,3),
         noise_rate_prototype = nothing,
-        noise= nothing,
+        noise= nothing, seed = UInt64(0),
         callback = nothing,mass_matrix=I)
 
-  SDEProblem{typeof(u0),promote_type(map(typeof,tspan)...),iip,typeof(noise),typeof(f),typeof(g),typeof(callback),typeof(mass_matrix),typeof(noise_rate_prototype)}(
-  f,g,u0,tspan,noise,callback,mass_matrix,noise_rate_prototype)
+  SDEProblem{typeof(u0),promote_type(map(typeof,tspan)...),iip,typeof(noise),
+             typeof(f),typeof(g),typeof(callback),typeof(mass_matrix),
+             typeof(noise_rate_prototype)}(
+  f,g,u0,tspan,noise,callback,mass_matrix,noise_rate_prototype,seed)
 end
 
 # Mu' = f[1] + f[2] + ...
@@ -28,16 +31,17 @@ type SplitSDEProblem{uType,tType,isinplace,NP,F,G,C,MM,ND} <: AbstractSDEProblem
   callback::C
   mass_matrix::MM
   noise_rate_prototype::ND
+  seed::UInt64
 end
 
 function SplitSDEProblem(f,g,u0,tspan; iip = isinplace(f[2],3),
-                         noise = nothing,
+                         noise = nothing, seed = UInt64(0),
                          noise_rate_prototype = nothing,
                          callback=nothing,mass_matrix=I)
   @assert typeof(f) <: Tuple
   SplitSDEProblem{typeof(u0),promote_type(map(typeof,tspan)...),iip,typeof(noise),typeof(f),
                   typeof(g),typeof(callback),typeof(mass_matrix),
-                  typeof(noise_rate_prototype)}(f,g,u0,tspan,noise,callback,mass_matrix,noise_rate_prototype)
+                  typeof(noise_rate_prototype)}(f,g,u0,tspan,noise,callback,mass_matrix,noise_rate_prototype,seed)
 end
 
 # M[i]*u[i]' = f[i]
@@ -50,10 +54,11 @@ type PartitionedSDEProblem{uType,tType,isinplace,NP,F,G,C,MM,ND} <: AbstractSDEP
   callback::C
   mass_matrix::MM
   noise_rate_prototype::ND
+  seed::UInt64
 end
 
 function PartitionedSDEProblem(f,g,u0,tspan; iip = isinplace(f[1],3),
-                               noise = nothing,
+                               noise = nothing, seed = UInt64(0),
                                noise_rate_prototype = nothing,
                                callback=nothing,mass_matrix=nothing)
   @assert typeof(f) <: Tuple
@@ -64,5 +69,5 @@ function PartitionedSDEProblem(f,g,u0,tspan; iip = isinplace(f[1],3),
   end
   PartitionedSDEProblem{typeof(u0),promote_type(map(typeof,tspan)...),iip,typeof(noise),typeof(f),typeof(g),
                         typeof(callback),typeof(_mm),typeof(noise_rate_prototype)}(
-                        f,g,u0,tspan,noise,callback,_mm,noise_rate_prototype)
+                        f,g,u0,tspan,noise,callback,_mm,noise_rate_prototype,seed)
 end
diff --git a/src/solutions/rode_solutions.jl b/src/solutions/rode_solutions.jl
@@ -12,6 +12,7 @@ type RODESolution{T,N,uType,uType2,DType,tType,randType,P,A,IType} <: AbstractRO
   dense::Bool
   tslocation::Int
   retcode::Symbol
+  seed::UInt64
 end
 (sol::RODESolution)(t,deriv::Type=Val{0};idxs=nothing) = sol.interp(t,idxs,deriv)
 (sol::RODESolution)(v,t,deriv::Type=Val{0};idxs=nothing) = sol.interp(v,t,idxs,deriv)
@@ -21,7 +22,8 @@ function build_solution(
         alg,t,u;W=[],timeseries_errors = length(u) > 2,
         dense = false,dense_errors=dense,calculate_error=true,
         interp = LinearInterpolation(t,u),
-        retcode = :Default, kwargs...)
+        retcode = :Default,
+        seed = UInt64(0), kwargs...)
 
   T = eltype(eltype(u))
   if typeof(prob.u0) <: Tuple
@@ -41,7 +43,7 @@ function build_solution(
     errors = Dict{Symbol,eltype(prob.u0)}()
     sol = RODESolution{T,N,typeof(u),typeof(u_analytic),typeof(errors),typeof(t),typeof(W),
                        typeof(prob),typeof(alg),typeof(interp)}(
-                       u,u_analytic,errors,t,W,prob,alg,interp,dense,0,retcode)
+                       u,u_analytic,errors,t,W,prob,alg,interp,dense,0,retcode,seed)
 
     if calculate_error
       calculate_solution_errors!(sol;timeseries_errors=timeseries_errors,dense_errors=dense_errors)
@@ -51,7 +53,7 @@ function build_solution(
   else
     return RODESolution{T,N,typeof(u),Void,Void,typeof(t),
                         typeof(W),typeof(prob),typeof(alg),typeof(interp)}(
-                        u,nothing,nothing,t,W,prob,alg,interp,dense,0,retcode)
+                        u,nothing,nothing,t,W,prob,alg,interp,dense,0,retcode,seed)
   end
 end
 
@@ -98,5 +100,5 @@ function build_solution(sol::AbstractRODESolution,u_analytic,errors)
   RODESolution{T,N,typeof(sol.u),typeof(u_analytic),typeof(errors),typeof(sol.t),
                typeof(sol.W),typeof(sol.prob),typeof(sol.alg),typeof(sol.interp)}(
                sol.u,u_analytic,errors,sol.t,sol.W,sol.prob,sol.alg,sol.interp,
-               sol.dense,sol.tslocation,sol.retcode)
+               sol.dense,sol.tslocation,sol.retcode,sol.seed)
 end