-
-
Notifications
You must be signed in to change notification settings - Fork 65
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Update implicit methods to use lazy W #96
Conversation
The version for constant caches isn't used at the moment, but is needed for exponential solvers in the future.
StochasticDiffEq doesn't seem to use transformed W operator at the moment.
We should try to be consistent where possible, so if it's not to difficult we should change this.
That's fine. |
http://docs.juliadiffeq.org/latest/solvers/sde_solve.html#Stiff-Methods-1 Set |
A lot of those don't need the same solution, for example just checking the W operator. But when you do need to check it, you can use a noise process generated by a pre-determined noise grid. For an example see this: https://github.com/JuliaDiffEq/DiffEqDevTools.jl/blob/master/src/benchmark.jl#L236-L241 |
Codecov Report
@@ Coverage Diff @@
## master #96 +/- ##
=========================================
Coverage ? 88.73%
=========================================
Files ? 33
Lines ? 1447
Branches ? 0
=========================================
Hits ? 1284
Misses ? 163
Partials ? 0
Continue to review full report at Codecov.
|
For the integration, I figured the simplest solution is to reset the random seed every time Currently all tests except for split_tests.jl pass. The split tests error on master as well so I don't think it has anything to do with the current PR. Stack trace:
|
Also about using StochasticDiffEq
println("Scalar g")
A = [-1.0 0.0; 0.0 -0.5]
u0 = [1.0, 1.0]; tspan = (0.0,1.0)
_f = (u,p,t) -> t*(A*u)
_g = (u,p,t) -> 1.0
prob = SDEProblem(SDEFunction(_f, _g), _g, u0, tspan)
integrator = init(prob, SKenCarp(); adaptive=false, dt=0.01)
step!(integrator) # fine
println("Vector g")
_g = (u,p,t) -> [1.0, 1.0]
prob = SDEProblem(SDEFunction(_f, _g), _g, u0, tspan)
println("Implicit EM")
integrator = init(prob, ImplicitEM(); adaptive=false, dt=0.01)
step!(integrator) # fine
println("SKenCarp")
integrator = init(prob, SKenCarp(); adaptive=false, dt=0.01)
step!(integrator) # error Gives the following output:
@ChrisRackauckas Can you confirm these errors indeed exist and are not a result of my environment/configuration? |
That only holds for some regression testing due to the fact that adaptivity can make the same seed on different problems generate a different time-sequence of random numbers |
Those two issues were issues on master and the fixes were merged into the PR |
Mimics what is done to OrdinaryDiffEq SciML/OrdinaryDiffEq.jl#460 in StochasticDiffEq.
Notes:
I'm not sure if this is a mistake or deliberate choice:
https://github.com/JuliaDiffEq/StochasticDiffEq.jl/compare/master...MSeeker1340:implicit-update?expand=1#diff-35ab888388acf2e7239e17f33e5eefe7R20
In OrdinaryDiffEq
integrator.eigen_est
is updated regardless of which branchcalc_J!
takes, but this is not the case for StochasticDiffEq.In StochasticDiffEq the return values of
calc_W!
for constant caches is(J, W)
, whereas in OrdinaryDiffEq onlyW
is returned. This seems a bit inconsistent (but on the other hand they are separate packages, so maybe it is ok?).alg_cache
in StochasticDiffEq does not have thedt
parameter (again an inconsistency), so I've chosen to usezero(t)
as the initial value ofdtgamma
in the construction ofWOperator
. The initial value should not matter though as it is always reset incalc_W!
.When I try to use custom mass matrices I got the error message "Algorithm must be set as symplectic or theta=1 for mass matrices". @ChrisRackauckas can you explain what this means? (Currently all mass matrices in the test scripts are
I
).It occurred to me that the integration tests for OrdinaryDiffEq (in utility_tests.jl) do not work for StochasticDiffEq because the solution is indeterministic, thus cannot be compared. I'm not sure how to modify it though. Maybe do a Monte Carlo simulation?