WIP: Function interface specs

[mauro3]

Concering comment SciML/Roadmap#5 (comment):

I only saw the jac_exists function after I put this together. However, neither has_jac nor jac_exits passes the test in above comment.

However, when using a trait-function it would be fully inferred. Note however, that the traits would struggle with adding, say a jacobian, to the function. It may would still dispatch to the wrong one.

[ChrisRackauckas]

Note however, that the traits would struggle with adding, say a jacobian, to the function. It may would still dispatch to the wrong one.

I think this is a fundamental issue, not an implementation issue. If having a trait can be inferred at compile time, then it must be constant (or at least, the has_jac(f) on a function has to be constant and essentially declared when f is declared). If has_jac(f) can notice that you later added a Jacobian function, then it has to be finding that out at runtime and thus not have the information at compile time. From this understanding (if it's correct), I don't think there's an easy way around that.

What we could do is have a different has_jac(f) dispatch for ParameterizedFunction which uses a trait and has_jac(f) for Function which uses the methods table. Then the easy way of defining a Jacobian will work, even if it's not necessarily good, and the more advanced way of using call-overloading with subtyping could then be documented after that (similar to how out-of-place functions are allowed, though in-place functions are preferred). With this, solvers can be written without care as to whether users took the easy or the hard way?

[mauro3]

Actually, the way around it, is to call @traitimpl HasJac{F} after a jacobian is added. Then (barring issues with Julia issue 265) it should work.

(Note the constant function used for trait-dispatch is trait here. Thus the use of a generated function, otherwise it would lead to dynamic dispatch in trait-functions.)

[ChrisRackauckas]

Actually, the way around it, is to call @traitimpl HasJac{F} after a jacobian is added. Then (barring issues with Julia issue 265) it should work.

That seems fine to me. It would just go in the same place of the documentation as defining the Jacobian (I want to make a full example for that first, and then below have the full interface. That way most users will likely define just the Jacobian which is fine: the others probably can't be done by hand).

The second probably needs a better title: https://juliadiffeq.github.io/DiffEqDocs.jl/latest/features/performance_overloads.html#Declaring-Explicit-Jacobians-1

[ChrisRackauckas]

At least as stated before, the interface has a capital here:

https://juliadiffeq.github.io/DiffEqDocs.jl/latest/features/performance_overloads.html

f(Val{:Jac},t,u,J) # Call the explicit Jacobian function

[mauro3]

Yes, sorry, I'm not up to date with the exact specific. Thanks for the link.

[mauro3]

Actually, would you mind updating the readme with links to all relevant places in the docs? (I have no overview)

[ChrisRackauckas]

Sure thing.

[ChrisRackauckas]

Updated.

[ChrisRackauckas]

That is weird but not a show-stopper

[ChrisRackauckas]

Just need to change to :Jac

[ChrisRackauckas]

The tests should move to /test and I think the following test should be used to make sure the traits option is a compile-time check:

function testing(f)
  if has_jac(f)
    a = 2
  else
    a = 3.0
  end
  a
end
@inferred testing(f)

I think we should also just go forward with the SimpleTraits.jl dependency: it has a compelling enough use case. After this goes through I'll set ParameterizedFunctions.jl up with traits.

Could you also expand the interface to have the other checks as well?

[ChrisRackauckas]

The other way is that ParameterizedFunctions and the docs can switch to :jac. Is there a guide on how value types should be named? I know functions should be lower-case according to the guide, and types should be upper case, but I don't think they thought about a value-type for declaring a new function. It poses quite an interesting dilemma.

[mauro3]

I think it's fine either way as long as we're consistent. Actually, to be consistent with variable use, I'd think :jac is better (or :J but that is too short).

[ChrisRackauckas]

I took it to JuliaPraxis, and after some discussion agree it should change to :jac.

https://gitter.im/JuliaPraxis/Naming/archives/2016/11/16

[ChrisRackauckas]

ParameterizedFunctions.jl was changed to have everything de-captialized (:jac). The docs were changed to match. I'll implement the traits on all of the ParameterizedFunctions in a bit.

[ChrisRackauckas]

I looked into this. It seems that as I was saying before, any method which uses the methods list / table will not be able to compute this at compile time. We could go through the trouble of creating a workaround that is at runtime and etc. etc., but really, the user experience isn't bad at all:

f123(x,y) = 1
f123(::Val{:jac}, x) = 2
@traitimpl HasJac{typeof(f123)}

If we just have it documented to put that one line after the Jacobian declaration, then the three lines:

@traitdef HasJac{F}
@traitfn has_jac(x::::HasJac) = true
@traitfn has_jac(x::::(!HasJac)) = false

makes everything work at compile time. If someone is looking through the manual for how to declare Jacobians by hand (something which is already done automatically for many systems of ODEs), I think they can spare that one line of code. If we just accept that, then everything works perfectly.

I am willing to take that option, and if you want to work on a macro to make it easier:

@add_jac f(t,u,J) = ... # Makes f(Val{:jac},t,u,J) and does @traitimpl

then I think that's a nice and easy solution for both developers and users.

[ChrisRackauckas]

@traitfn testing(f::::HasJac)   =2
@traitfn testing(f::::!(HasJac))=3.0

@code_llvm testing(f123)
@code_llvm testing(g123)
@code_llvm testing(h123)

I am really dumb. After playing with it a bit, I see if you use the @traitfn correctly, everything is inferred. So this all works out.

Since this all works, I'll pull this in and add the other traits. Looks like there's a merge conflict to handle though.

[ChrisRackauckas]

I got it!

@traitdef HasJac{F}
__has_jac(f) = check_first_arg(f, Val{:jac})
@generated SimpleTraits.trait{F}(::Type{HasJac{F}}) = __has_jac(F) ? :(HasJac{F}) : :(Not{HasJac{F}})
has_jac{T}(f::T)   = istrait(HasJac{T})

Only check the methods table to generate the trait and you're good. Now it auto-finds the trait and has_jac is the right value. Thanks for the help!