Closed Gertian closed 1 month ago
These test failures seem unrelated, we should probably check if the master even runs...
As a side-comment, this is not really a bug. The SDD algorithm should be a lot faster for large matrices, and this change really does have some performance implications. Nevertheless, I agree that in our case, as the stability of SDD turns out to be insufficient for the cases at hand, we should probably default to the (slower) more stable SVD algorithm
We should probably adopt a SafeSVD
version that does
try
# use SDD
catch
# fall back to SVD
end
A singular value decomposition is probably sufficiently costly in most cases that the try ... catch
overhead is negligible.
@lkdvos I'm quite certain that the failing tests are indeed unrelated since I've ran many hours of Gradient Descent with this change without any issues.
@Jutho This was my original patch, the problem is that SDD already changes the input tensor so we should really do :
Make hardcopy(tensor) try catch
but then there is a substantial overhead which is why I finally just accepted to SDD.
In this case without the exclamation mark, tsvd
should not change the input parameter. I agree that it is more tricky with the mutating method. So maybe SafeSVD
should only be valid for tsvd
but not for tsvd!
.
Hi all,
I updated the PR.
From the TensorKitManifolds perspective the only changes compared to master are now that all methods (I think) now take the alg keyword which can be SVD or SDD.
I simultaneously made a PR in MPSKit which ads an extra argument svd_alg to the GradientGrassman struct. This can be set to SVD or SDD and is passed to the underlying TensorKitManifolds algorithms.
I added some changes to make the tests run again:
getproperty
is treated specially in the Julia compiler, in the sense that it has a default aggressive constprop. This ensures that it is type stable with the different symbols. Adding a keyword argument breaks this, so I had to revert that change. In principle, we could consider manually adding a Base.@constprop :aggressive
, but somehow I don't think getproperty
with a keyword is very Julian.This is probably a ridiculous edge case but I'm wondering what might happen if a user sets the TensorKitManifolds default svd algorithm to SDD and then uses the dispatching in MPSKit to use SVD everywhere.
In this case the SDD call inside getproperty might still lead to a crash ?
Or do you think users should simply not use SDD anymore ?
The current approach does not actually have an option to set the default algorithm, it is a compile time constant.
Note that your approach also did not really expose the alg keyword in the getproperty
function, this is the function that gets called when you access fields of a tangent: y.U
is completely equivalent to getproperty(y, :U)
. (This is also why this function is treated differently).
Thus, the current implementation has the stable one as default, and if you would want to change that, it would require quite a bit of rewrite. I am more or less okay with anything, so I would wait for @Jutho 's opinion
Ok, that makes sense thanks !
My apologies for not getting to this sooner, which will probably result in double work. I have the following requests / suggestions:
default_svd_alg(A::TensorMap) = SVD()
in the main module file ("src/TensorKitManifolds.jl")alg
keyword of projectisometric!
(also in the TensorKitManifolds.jl file) to defualt_svd_alg(A)
. Also use this value in any other tsvd
and tsvd!
call that you can find, I think they are only in "src/grassman.jl".alg
value (reset it to nothing
) in Grassmann.retract
and friends (invretract
, relativegauge
). This alg
argument is actually aimed at different retraction schemes (which exist for example for the Stiefel case). This makes all svd calls in TensorKitManifolds.jl completely hard coded, but you can still change it by redefining TensorKitManifolds.default_svd_alg(::TensorMap) = SDD()
.
Ultimately, we should have algorithms structures at the level of retract
and transport
and friends, where then the algorithm to be used for the svd
would be a field in that structure. But that should be part of a larger overhaul, possibly where we consider to buy in to the Manifolds.jl ecosystem.
Thanks @Jutho for the suggestions and @lkdvos for implementing them faster then I could even read them :sweat:
Does this mean that
using TensorKitManifolds
TensorKitManifolds.default_svd_alg(::TensorMap) = SDD()
using MPSKit
rest of my code
will result in using the stable version of svd inside all MPSKit
algorithms ?
No, SDD()
is the unstable/fast algorithm, so by default the stable version will already be selected, you won't have to do anything.
Oh yes, indeed. Thanks for pointing that out !
Great. I will merge when lights turn green. There seems one svd!
call hidden in "src/auxiliary.jl" in _stiefellog
, but that is unused I assume so let's keep that for another time.
Also, we should probably switch to the LinearAlgebra
name of these algorithms in TensorKit.jl at some point. default_svd_alg
is also defined in LinearAlgebra
, but it might actually be useful to have independent default_svd_alg
routines in the different packages for more granular control, until we have a batter strategy implemented.
I saw the talk about ScopedValues on JuliaCon last week, it seems like these kinds of defaults and configurations are precisely what is being targeted with that, and from what I understand, this would solve our issues. Just something to keep in the back of our heads when we revisit this.
fixed tiny bug