Skip to content
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

re-export VectorInterface? #88

Closed
rveltz opened this issue May 31, 2024 · 10 comments
Closed

re-export VectorInterface? #88

rveltz opened this issue May 31, 2024 · 10 comments

Comments

@rveltz
Copy link
Contributor

rveltz commented May 31, 2024

Hi,

I rely on KK for BifurcationKit.jl. In particular, I implemented a type BorderedArray

It seems the only method I need to add is add. First, I thought having axpby! defined would not require to implement add since I have similar and copy. Can you confirm this ?

In the case I really need to add the method add, I need to import VectorInferface: add which requires to add VectorInterface as a dependency to BifurcationKit. Is this really necessary?

@lkdvos
Copy link
Collaborator

lkdvos commented May 31, 2024

In principle, the required interface to be fully compatible with KrylovKit is the full interface specified by VectorInterface. I think for now though, there are still fallback methods available that will automatically dispatch to LinearAlgebra's axpy! and axpby!, as well as rmul!, whenever these are defined.

I think the easiest is to indeed add VectorInterface as a dependency, and just implement these methods for your custom type. Note that VectorInterface is quite lightweight, specifically for this reason, and as KrylovKit already depends on this it would already be in the implicit dependencies anyways.

An alternative solution is to make use of the fact that KrylovKit already depends on VectorInterface, and do something like:

import KrylovKit: VectorInterface
function VectorInterface.add(x::MyType, y::MyType, a::Number, b::Number)
    # my implementation
end

At the very least in the REPL this seems to work, and does not require an explicit dependency.

Note that I do not know which of these is considered best practice, but I hope to get people's opinion on this matter soon, see this discourse post

@Jutho
Copy link
Owner

Jutho commented May 31, 2024

Up till version 0.6, we required axpby! and a few other methods, as you can find in the old docs:
https://jutho.github.io/KrylovKit.jl/v0.6/

Since v0.7, we switched to depending on VectorInterface, by which we can support a much larger range of types (such as nested arrays, tuples, …) out of the box. In particular, from VectorInterface we use both the non-mutating functions (scale, add, ...) as well as the possibly mutating methods (scale!!, add!!).

As @lkdvos mentioned, VectorInterface is a dependency of KrylovKit. We did not reexport it in order not to pollute the namespace, but you can indeed import VectorInterface from KrylovKit as mentioned there.

As for which methods to extend, it also depends a bit. For AbstractArray, there are default implementations that use broadcasting, or use BLAS methods (such as axpby!) for suitable types. However, for specific types, more optimised implementations might be preferred. For generic types, there are fallback that use + and *, or test (dynamically!) wheter axpby! and rmul! and friends can be applied. These fallbacks probably have some overhead, and I am not sure if we want to keep those in the future.

We might think about some way to automate generating these method definitions (since there are quite a few) for new types, e.g. using a macro where you just specify your type, whether the data can be mutated, and whether the implementation should be lowering to e.g. LinearAlgebra.axpby! and friends, or to broadcasting, ...). Any input from actual use cases like yours would be very welcome for that.

@rveltz
Copy link
Contributor Author

rveltz commented May 31, 2024

Thank you for your helpful comments. In the end I implemented the interface using import KrylovKit: VectorInterface.

@rveltz rveltz closed this as completed May 31, 2024
@rveltz rveltz reopened this Jun 1, 2024
@rveltz
Copy link
Contributor Author

rveltz commented Jun 1, 2024

Can you help me a bit with it please? It gives me NaN that were not there before and this block tagging a new version of BifurcationKit

I am doing this

VectorInterface.add!(y::BorderedArray, x::BorderedArray, α::Number = 1, β::Number = 1) = axpby!(β, x, α, y)
VectorInterface.add!!(y::BorderedArray, x::BorderedArray, α::Number = 1, β::Number = 1) = VectorInterface.add!(y, x, α, β)
@inline VectorInterface.scalartype(W::BorderedArray{vectype, Tv1}) where {vectype, Tv1} = eltype(BorderedArray{vectype, Tv1})
VectorInterface.zerovector(b::BorderedArray, S::Type{<:Number} = VectorInterface.scalartype(b)) = 0 * similar(b, S)
VectorInterface.scale(x::BorderedArray, α::Number) = mul!(similar(x), x, α)# α * x
VectorInterface.scale!(y::BorderedArray, x::BorderedArray, α::Number) = throw("")#mul!(y, x, α)
VectorInterface.scale!!(x::BorderedArray, α::Number) = α * x
VectorInterface.scale!!(y::BorderedArray, x::BorderedArray, α::Number) = mul!(y, x, α)
VectorInterface.inner(x::BorderedArray, y::BorderedArray) = dot(x, y)

@Jutho
Copy link
Owner

Jutho commented Jun 1, 2024

Just guessing, but if similar has nans in its unitialized memory, then 0 * … will preserve those. That’s a difference between Julia semantics and BLAS. Maybe use some version of fill!

@rveltz
Copy link
Contributor Author

rveltz commented Jun 1, 2024

OK will try, thank you

@lkdvos
Copy link
Collaborator

lkdvos commented Jun 1, 2024

I think the cleanest way is something like:
VectorInterface.zerovector(x::BorderedArray, ::Type{S}) where {S} = BorderedArray(zerovector(x.u, S), zerovector(x.p, S))

In principle you can use a similar approach for all other functions, by defining them recursively, as vectorinterface should work on most AbstractArray, Tuple, Number, ... out of the box

@rveltz
Copy link
Contributor Author

rveltz commented Jun 1, 2024

Ah interesting

@lkdvos
Copy link
Collaborator

lkdvos commented Jun 12, 2024

Did everything work out for you?
I think we will not be re-exporting VectorInterface, as it really is meant to be a lightweight dependency and there should be no cost associated to depending on it. According to the discussion on discourse, best practice would be to explicitly add VectorInterface to your dependencies, because now you risk that a breaking vectorinterface update would not be caught by your compat settings. I do admit that this is rather unlikely, and in the end it's a bit a matter of preference too.
If you are happy with this, can I close this issue?

@rveltz
Copy link
Contributor Author

rveltz commented Jun 12, 2024

Yeah, I agree, I should add VI to my project.toml which is not what I do at the moment

@rveltz rveltz closed this as completed Jun 12, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants