fix up silly change and add more tests

SciML · Apr 8, 2018 · a306b52 · a306b52
1 parent f5c1cd2
commit a306b52
Show file tree

Hide file tree

Showing 9 changed files with 56 additions and 14 deletions.
diff --git a/src/caches/basic_method_caches.jl b/src/caches/basic_method_caches.jl
@@ -14,7 +14,7 @@ du_cache(c::EMCache) = (c.rtmp1,c.rtmp2)
 alg_cache(alg::EM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{false}}) = EMConstantCache()
 
 function alg_cache(alg::EM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
-  tmp = zeros(ΔW); rtmp1 = zeros(rate_prototype);
+  tmp = zeros(u); rtmp1 = zeros(rate_prototype);
   rtmp2 = zeros(noise_rate_prototype)
   if is_diagonal_noise(prob)
     rtmp3 = rtmp2
@@ -40,7 +40,7 @@ du_cache(c::SplitEMCache) = (c.rtmp1,c.rtmp2)
 alg_cache(alg::SplitEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{false}}) = SplitEMConstantCache()
 
 function alg_cache(alg::SplitEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
-  tmp = zeros(ΔW); rtmp1 = zeros(rate_prototype);
+  tmp = zeros(u); rtmp1 = zeros(rate_prototype);
   rtmp2 = zeros(noise_rate_prototype)
   if is_diagonal_noise(prob)
     rtmp3 = rtmp2
@@ -68,7 +68,7 @@ du_cache(c::EulerHeunCache) = (c.rtmp1,c.rtmp2,c.rtmp3,c.rtmp4)
 alg_cache(alg::EulerHeun,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{false}}) = EulerHeunConstantCache()
 
 function alg_cache(alg::EulerHeun,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
-  tmp = zeros(ΔW); ftmp1 = zeros(rate_prototype); ftmp2 = zeros(rate_prototype)
+  tmp = zeros(u); ftmp1 = zeros(rate_prototype); ftmp2 = zeros(rate_prototype)
   nrtmp = zeros(rate_prototype)
   gtmp1 = zeros(noise_rate_prototype); gtmp2 = zeros(noise_rate_prototype)
   EulerHeunCache(u,uprev,tmp,ftmp1,ftmp2,nrtmp,gtmp1,gtmp2)
@@ -88,7 +88,7 @@ du_cache(c::RandomEMCache) = (c.rtmp1,c.rtmp2)
 alg_cache(alg::RandomEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{false}}) = RandomEMConstantCache()
 
 function alg_cache(alg::RandomEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
-  tmp = zeros(ΔW); rtmp = zeros(rate_prototype)
+  tmp = zeros(u); rtmp = zeros(rate_prototype)
   RandomEMCache(u,uprev,tmp,rtmp)
 end
 
@@ -110,7 +110,7 @@ alg_cache(alg::RKMil,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltyp
 
 function alg_cache(alg::RKMil,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
   du1 = zeros(rate_prototype); du2 = zeros(rate_prototype)
-  K = zeros(rate_prototype); tmp = zeros(ΔW); L = zeros(rate_prototype)
+  K = zeros(rate_prototype); tmp = zeros(u); L = zeros(rate_prototype)
   RKMilCache(u,uprev,du1,du2,K,tmp,L)
 end
 

diff --git a/src/caches/implicit_split_step_caches.jl b/src/caches/implicit_split_step_caches.jl
@@ -32,7 +32,7 @@ function alg_cache(alg::ISSEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototy
   J = zeros(uEltypeNoUnits,length(u),length(u)) # uEltype?
   W = zeros(J)
   z = zeros(u)
-  dz = zeros(u); tmp = zeros(ΔW); gtmp = zeros(noise_rate_prototype)
+  dz = zeros(u); tmp = zeros(u); gtmp = zeros(noise_rate_prototype)
   fsalfirst = zeros(rate_prototype)
   k = zeros(rate_prototype)
 
@@ -128,7 +128,7 @@ function alg_cache(alg::ISSEulerHeun,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_
   J = zeros(uEltypeNoUnits,length(u),length(u)) # uEltype?
   W = zeros(J)
   z = zeros(u)
-  dz = zeros(u); tmp = zeros(ΔW); gtmp = zeros(noise_rate_prototype)
+  dz = zeros(u); tmp = zeros(u); gtmp = zeros(noise_rate_prototype)
   fsalfirst = zeros(rate_prototype)
   k = zeros(rate_prototype)
 

diff --git a/src/caches/kencarp_caches.jl b/src/caches/kencarp_caches.jl
@@ -78,7 +78,7 @@ function alg_cache(alg::SKenCarp,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prot
   dz = similar(u,indices(u))
   fsalfirst = zeros(rate_prototype)
   k = zeros(rate_prototype)
-  tmp = zeros(ΔW); b = similar(u,indices(u));
+  tmp = zeros(u); b = similar(u,indices(u));
   atmp = zeros(u,uEltypeNoUnits,indices(u))
 
   if typeof(f) <: SplitFunction

diff --git a/src/caches/lamba_caches.jl b/src/caches/lamba_caches.jl
@@ -18,7 +18,7 @@ alg_cache(alg::LambaEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uElt
 
 function alg_cache(alg::LambaEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
   du1 = zeros(rate_prototype); du2 = zeros(rate_prototype)
-  K = zeros(rate_prototype); tmp = zeros(ΔW);
+  K = zeros(rate_prototype); tmp = zeros(u);
   L = zeros(noise_rate_prototype)
   gtmp = zeros(noise_rate_prototype)
   if is_diagonal_noise(prob)
@@ -49,7 +49,7 @@ alg_cache(alg::LambaEulerHeun,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototy
 
 function alg_cache(alg::LambaEulerHeun,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
   du1 = zeros(rate_prototype); du2 = zeros(rate_prototype)
-  K = zeros(rate_prototype); tmp = zeros(ΔW);
+  K = zeros(rate_prototype); tmp = zeros(u);
   L = zeros(noise_rate_prototype)
   gtmp = zeros(noise_rate_prototype)
   if is_diagonal_noise(prob)

diff --git a/src/caches/predcorr_caches.jl b/src/caches/predcorr_caches.jl
@@ -11,7 +11,7 @@ struct PCEulerCache{uType,rateType,rateNoiseType,rateNoiseCollectionType} <: Sto
 end
 
 function alg_cache(alg::PCEuler,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_prototype,uEltypeNoUnits,uBottomEltype,tTypeNoUnits,uprev,f,t,::Type{Val{true}})
-  utmp = zeros(ΔW); ftmp = zeros(rate_prototype);
+  utmp = zeros(u); ftmp = zeros(rate_prototype);
   gtmp = zeros(noise_rate_prototype)
   bbprimetmp = zeros(ftmp)
   if is_diagonal_noise(prob)

diff --git a/src/caches/sdirk_caches.jl b/src/caches/sdirk_caches.jl
@@ -30,7 +30,7 @@ function alg_cache(alg::ImplicitEM,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate_pr
   J = zeros(uEltypeNoUnits,length(u),length(u)) # uEltype?
   W = zeros(J)
   z = zeros(u)
-  dz = zeros(u); tmp = zeros(ΔW); gtmp = zeros(noise_rate_prototype)
+  dz = zeros(u); tmp = zeros(u); gtmp = zeros(noise_rate_prototype)
   fsalfirst = zeros(rate_prototype)
   k = zeros(rate_prototype)
 
@@ -124,7 +124,7 @@ function alg_cache(alg::ImplicitEulerHeun,prob,u,ΔW,ΔZ,p,rate_prototype,noise_
   J = zeros(uEltypeNoUnits,length(u),length(u)) # uEltype?
   W = zeros(J)
   z = zeros(u)
-  dz = zeros(u); tmp = zeros(ΔW); gtmp = zeros(noise_rate_prototype)
+  dz = zeros(u); tmp = zeros(u); gtmp = zeros(noise_rate_prototype)
   fsalfirst = zeros(rate_prototype)
   k = zeros(rate_prototype)
 
@@ -219,7 +219,7 @@ function alg_cache(alg::ImplicitRKMil,prob,u,ΔW,ΔZ,p,rate_prototype,noise_rate
   J = zeros(uEltypeNoUnits,length(u),length(u)) # uEltype?
   W = zeros(J)
   z = zeros(u)
-  dz = zeros(u); tmp = zeros(ΔW); gtmp = zeros(noise_rate_prototype)
+  dz = zeros(u); tmp = zeros(u); gtmp = zeros(noise_rate_prototype)
   fsalfirst = zeros(rate_prototype)
   k = zeros(rate_prototype)
 

diff --git a/test/nondiagonal_tests.jl b/test/nondiagonal_tests.jl
@@ -93,3 +93,24 @@ simPCEuler = test_convergence(dts,probiip,PCEuler(ggprime),numMonte=numMonte)
 # ylabel("l2 error")
 # xlabel("dt")
 # legend()
+
+f_morenoise = (du,u,p,t) -> du.=1.01u
+g_morenoise = function (du,u,p,t)
+  du[1,1] = 0.3u[1]
+  du[1,2] = 0.6u[1]
+  du[1,3] = 0.9u[1]
+  du[1,4] = 0.12u[2]
+  du[2,1] = 1.2u[1]
+  du[2,2] = 0.2u[2]
+  du[2,3] = 0.3u[2]
+  du[2,4] = 1.8u[2]
+end
+prob = SDEProblem(f_morenoise,g_morenoise,ones(2),(0.0,1.0),
+                  noise_rate_prototype=zeros(2,4))
+
+sol =solve(prob,dt=1/2^(3),EM())
+sol =solve(prob,dt=1/2^(3),ISSEM())
+sol =solve(prob,dt=1/2^(3),ImplicitEM())
+sol =solve(prob,dt=1/2^(3),EulerHeun())
+sol =solve(prob,dt=1/2^(3),ImplicitEulerHeun())
+sol =solve(prob,dt=1/2^(3),ISSEulerHeun())
diff --git a/test/runtests.jl b/test/runtests.jl
@@ -31,6 +31,7 @@ LONGER_TESTS && @time @testset "Weak Convergence Tests" begin include("weak_conv
 @time @testset "Events Tests" begin include("events_test.jl") end
 @time @testset "Cache Tests" begin include("cache_test.jl") end
 @time @testset "Element-wise Tolerances Tests" begin include("tolerances_tests.jl") end
+@time @testset "Zero'd Noise Tests" begin include("zerod_noise_test.jl") end
 #@time @testset "Scalar Tests" begin include("scalar_noise.jl") end # Fails because of bounds checks
 
 #Adaptive SDE

diff --git a/test/zerod_noise_test.jl b/test/zerod_noise_test.jl
@@ -0,0 +1,20 @@
+using StochasticDiffEq, Base.Test
+
+function f(du, u, p, t)
+    du[1] = u[1]
+end
+function g(du, u, p, t)
+    0.0
+end
+
+u0 = [1.0]
+prob = SDEProblem{true}(f, g, u0, (0.0, 0.1))
+
+sol_ito = solve(prob, RKMil{:Ito}())
+@test length(sol_ito) < 100
+
+sol_strato = solve(prob, RKMil{:Stratonovich}(); dt=1e-2)
+@test length(sol_strato) < 100
+
+sol_leh = solve(prob, LambaEulerHeun())
+@test length(sol_leh) < 100