Conversation
…mber of cotangent columns
| """ | ||
| antihermitian_columns!(X, K) | ||
|
|
||
| Given the columns `K` of a square matrix that is nonzero only in its columns `K`, overwrite `X` | ||
| with the same columns of the antihermitian part of that matrix. | ||
| """ | ||
| function antihermitian_columns!(X, K) | ||
| # NOTE: all columns in order (e.g. `ind = Colon()`): the original in-place projection | ||
| is_leading_index(K, size(X, 1)) && return project_antihermitian!(X) | ||
| XKK = project_antihermitian!(X[K, :]) | ||
| X ./= 2 | ||
| X[K, :] .= XKK | ||
| return X | ||
| end |
There was a problem hiding this comment.
Any reason to use K instead of ind as variable name? I find this quite confusing as I often denote certain matrices coming up in the pullbacks as K (i.e. the antihermitian matrix associated with the infinitesimal in-space rotation of an isometry).
There was a problem hiding this comment.
Also, I don't understand what is happening here. Is X square, or is X already the K columns of a larger square matrix. And then X[K, :] is square, corresponding to the diagonal block of that matrix? I think it is the latter, but the doc string could be a bit clearer.
There was a problem hiding this comment.
No good reason at all, switched it back to use ind. I also updated the docstring and variable names to make it more clear what's actually happening. X holds the columns ind of a larger square matrix M that is nonzero only in those columns, and X[ind, :] = M[ind, ind] is its diagonal block, the only part of X to which M' contributes.
Codecov Report✅ All modified and coverable lines are covered by tests.
🚀 New features to boost your workflow:
|
| indD = select_indices(axes(D, 1), ind) | ||
| indV = select_indices(axes(V, 2), ind) | ||
| K = select_indices(axes(D, 1), ind) | ||
| k = length(K) |
There was a problem hiding this comment.
Or ind′ in case just want a single variable, and since ind is already in use.
There was a problem hiding this comment.
I renamed to ind′ for consistency.
This branch has not been deployed
2 participants
When the cotangents of
svd_pullback!are given on only k of the r singular vectors (throughind), the pullback pads them with zeros to all r columns and continues with r × r matrices:check_and_prepare_svd_cotangentsformsU₁' * ΔU₁andΔU₁ - U₁ * (U₁' * ΔU₁), the same forV, and the result is applied asU₁ * (UᴴΔAV * V₁ᴴ). The cost is therefore O(m n r) whatever k is.eigh_pullback!does the same withV' * ΔV₁andV * VᴴΔAV * V', at cost O(n³) for an n × n matrix. This is the common case for a full pullback of a truncated decomposition, whereindholds the kept indices.With cotangents on the columns
Konly,U₁ᴴΔU₁andV₁ᴴΔV₁are nonzero only in the columnsK. ThenUᴴΔAVis nonzero only in the rows and columnsK, and its rows follow from its columns by antihermiticity. In this PR,check_and_prepare_svd_cotangentstherefore no longer pads the cotangents, but computes only the r × k block of columnsKofUᴴΔAVand the corresponding block of its rows. For 2k ≤ r,svd_pullback!applies these two blocks directly as rank-k updates, at cost O(m n k); otherwise it assemblesUᴴΔAVfrom them and applies it as before.check_and_prepare_eigh_cotangentsandeigh_pullback!are changed in the same way, so that for 2k ≤ n the cost ofeigh_pullback!is O(n² k) instead of O(n³). The gauge check covers the same entries as before, since the columnsKcontain every nonzero entry up to conjugation. Cotangents on columns beyond a full rank have components along all ofU₁orV₁ᴴ, so in that case the block still spans all r columns.svd_trunc_pullback!andeigh_trunc_pullback!call the same functions with all columns and receive the same full matrices as before.Bug fix. This PR also fixes
svd_pullback!for a nonzeroΔSwith anindother than1:k. Since #232,check_and_prepare_svd_cotangentsonmainindexesΔSby column number instead of by position withinind, so that for exampleind = [3, 1, 7, 2]throws aBoundsError. The new code adds each entry ofΔSto the diagonal entry of its own column; the comparison below includes this case.On random real and complex, square and rectangular, full-rank and rank-deficient matrices, with
ind = 1:5,[3, 1, 7, 2]and1:6, and with zeroΔU,ΔVᴴorΔS, the result agrees withmain(called with the same cotangents zero-padded to all columns) to 7.9e-16 relative. On square matrices with exponentially decaying singular values or eigenvalues and k between n / 36 and n / 9 (the spectra and sizes of a CTMRG step in PEPSKit.jl, with n = χD² and k = χ, which motivated this investigation), it is 4–22× faster:main(s)Minimum times on a laptop with 4 BLAS threads.
main'ssrc/pullbacks/svd.jlandsrc/pullbacks/eigh.jlwere loaded into the same process, and the two versions were timed alternately with the same arguments.On random n × n matrices (n = 400 and 1000, real and complex; same timing method), the speedup over
mainfor cotangents on the first k columns (for eigh, the k eigenvalues of largest magnitude) is:ind = Colon())With 2k ≤ r the new path is faster in every case. Above that, the r × r matrix is formed from the k computed columns and applied as before, which is slower than
mainin five of the 24 cases with 2k > r: by up to 1.3× in two SVD cases at n = 1000, and by at most 6% in the other three. Withind = Colon()the result is identical to that ofmain.Benchmark
On
main:With this PR: