Draft: ImplicitDiscreteSolve
#534
Conversation
|
@ChrisRackauckas This is barebones but working. Any important features I'm missing here? Also should this actually go in OrdinaryDiffEq? |
|
Can this to use SimpleNonlinearSolve.jl's Newton and make this be SimpleImplicitDiscreteSolve.jl? I think the complete version will want to use the Jacobian reuse aspect and such, so for now having a static version would be a good way to handle this. |
|
Did that refactor, just needs more tests now I think |
lib/SimpleImplicitDiscreteSolve/src/SimpleImplicitDiscreteSolve.jl
Outdated
Show resolved
Hide resolved
|
@ChrisRackauckas I think this is basically ready, initialization for problems without initializationproblems needs work but I think it's fine for the MTK purposes. How do I register it? |
|
Fixed these |
|
Okay so this implementation finds itself in a somewhat odd spot. It's not a "Simple" alg because it's not statically compile-able, which is a property that all others have. It would need to avoid the OrdinaryDiffEqCore part. But it's not a fully complete version either? But it's a rough start towards it. I guess seeing it at this stage, how about moving this to OrdinaryDiffEq and making this be ImplicitDiscreteSolve.jl, and we can use this as a base to do the required optimizations to make the full version. We can do the SimpleImplicitDiscreteSolve.jl as a somewhat separate version for static support. |
|
I see, didn't realize I couldn't use OrdinaryDiffEqCore, will move. I'll keep this open too to keep iterating it toward a SimpleImplicitDiscreteSolve. |
|
See for example https://github.com/SciML/SimpleDiffEq.jl/blob/master/src/rk4/looprk4.jl which is the most basic RK4 implementation. For IDS, it would just be a loop over nonlinear solve calls. With SimpleNonlinearSolve being similar, the whole thing would statically compile since no mutable structs would be in the middle. |
|
@ChrisRackauckas tried to refactor this to a simple solve |
| for i in 2:length(ts) | ||
| uprev = u | ||
| t = ts[i] | ||
| nlf = isinplace(f) ? (out, u, p) -> f(out, u, uprev, p, t) : (u, p) -> f(u, uprev, p, t) |
There was a problem hiding this comment.
is the compiler always smart enough to optimize this closure? It probably is these days. I would've just made a callable type, but I think this got better. I guess we can just patch this if we find a case where it does allocate.
|
This folder will need a license file in order to register, but it should be good to go now. |
|
Alright I think this is now done, fingers crossed. But, if ImplicitDiscreteSolve is in OrdinaryDiffEq.jl, it may make sense to move this there before registering. But let's just make sure tests pass and moving code around is then simple. |
| if save_everystep && save_start | ||
| us = Vector{typeof(u0)}(undef, length(ts)) | ||
| us[1] = u0 | ||
| elseif save_everystep | ||
| us = Vector{typeof(u0)}(undef, length(ts) - 1) | ||
| elseif save_start | ||
| us = Vector{typeof(u0)}(undef, 2) | ||
| us[1] = u0 | ||
| else | ||
| us = Vector{typeof(u0)}(undef, 1) # for interface compatibility |
There was a problem hiding this comment.
It would be good to make this possibly be an SArray to skip the allocation when possible. This is done for example in https://github.com/SciML/SimpleDiffEq.jl/blob/master/src/rk4/gpurk4.jl#L40
|
Alright it looks like it's passing. Let's move this version over and register. |
|
Oops I just added this staticarrays thing |
| l = save_everystep ? length(ts) - 1 : 1 | ||
| save_start && (l = l + 1) | ||
| u0type = typeof(u0) | ||
| us = u0type <: StaticArray ? MVector{l, u0type}(undef) : Vector{u0type}(undef, l) |
There was a problem hiding this comment.
This won't infer though
There was a problem hiding this comment.
Hm how would I do this then?
There was a problem hiding this comment.
Maybe a switch in the alg that is type based? We can follow up with that though. Let's get a first version and then get that non-allocating.
There was a problem hiding this comment.
OK I'll open it in OrdinaryDiffEq
Requires:
AbstractDiscreteProbleminstead ofDiscreteProbleminsolveOrdinaryDiffEq.jl#2607