Conversation
add QR gauge projection
Codecov Report❌ Patch coverage is
🚀 New features to boost your workflow:
|
Co-authored-by: Jutho <Jutho@users.noreply.github.com>
|
I guess I really don't see the point of removing all the |
|
Once this merges, I'll update the other PR |
| Q2 = view(Q, (p + 1):size(Q, 1), :) | ||
| ΔQ2 = view(ΔQ, (p + 1):size(Q, 1), :) | ||
| ΔQ2Q1ᴴ = ΔQ2 * Q1' | ||
| check_lq_full_cotangents(Q1, ΔQ2, ΔQ2Q1ᴴ; gauge_atol) |
There was a problem hiding this comment.
For a rank-deficient matrix p < minmn, check_lq_full_cotangents can also be called by lq_compact (i.e. size(Q,1) == minmn for lq_compact and size(Q,1) == n for lq_full). Since this is not a change from the current PR, I will try to look into when this was introduced, and whether everything is still fully correct.
There was a problem hiding this comment.
Does that mean we can go ahead with this one? That seems somewhat unrelated to the testsuite refactor anyways, so might be clearer in a separate PR anyways
| remove_eig_gauge_dependence!(ΔV, D, V) | ||
|
|
||
| Remove the gauge-dependent part from the cotangent `ΔV` of the eigenvector matrix `V`. The | ||
| eigenvectors are only determined up to complex phase (and unitary mixing for degenerate |
There was a problem hiding this comment.
In the eig case there can also be non-unitary mixing. I think the only convention respected by LAPACK is that all the eigenvectors have norm 1, but that condition cannot easily be transformed into a restriction onto a restriction on the matrices. In fact, even for non-degenerate eigenvalues we actually require v' * Δv (the diagonal element of V' * ΔV) to be completely zero, not just the imaginary part of it (the imaginary part corresponds to phase rotations, the real part to norm changes).
| ΔU[:, (minmn + 1):end] .= 0 | ||
| ΔVᴴ[(minmn + 1):end, :] .= 0 |
There was a problem hiding this comment.
Do we also want to support rank deficient matrices, where we should use some p instead of minmn?
| Q₁ = @view Q[:, 1:r] | ||
| ΔQ₂ = @view ΔQ[:, (r + 1):end] |
There was a problem hiding this comment.
For consistency, since you do use regular view in the last line:
| Q₁ = @view Q[:, 1:r] | |
| ΔQ₂ = @view ΔQ[:, (r + 1):end] | |
| Q₁ = view(Q, :, 1:r) | |
| ΔQ₂ = view(ΔQ, :, (r + 1):end) |
| ambiguity. Additionally, rows of `ΔR` beyond the rank are zeroed out. | ||
| """ | ||
| function remove_qr_gauge_dependence!(ΔQ, ΔR, A, Q, R) | ||
| r = MatrixAlgebraKit.qr_rank(R) |
There was a problem hiding this comment.
Do we want to support a rank_atol keyword to be passed along from remove_qr_gauge_dependence! to qr_rank?
| Q₁ = @view Q[1:r, :] | ||
| ΔQ₂ = @view ΔQ[(r + 1):end, :] |
There was a problem hiding this comment.
| Q₁ = @view Q[1:r, :] | |
| ΔQ₂ = @view ΔQ[(r + 1):end, :] | |
| Q₁ = view(Q, 1:r, :) | |
| ΔQ₂ = view(ΔQ, (r + 1):end, :) |
| Additionally, columns of `ΔL` beyond the rank are zeroed out. | ||
| """ | ||
| function remove_lq_gauge_dependence!(ΔL, ΔQ, A, L, Q) | ||
| r = MatrixAlgebraKit.lq_rank(L) |
| Q, _ = qr_compact(A) | ||
| mul!(ΔN, Q, Q' * ΔN) | ||
| return ΔN |
There was a problem hiding this comment.
Any reason to not simply pass this on to remove_qr_null_gauge_dependence in order to avoid code duplication?
| null space basis is only determined up to a unitary rotation, so `ΔNᴴ` is projected onto the | ||
| row span of the compact LQ factor `Q₁` of `A`. | ||
| """ | ||
| function remove_right_null_gauge_dependence!(ΔNᴴ, A, Nᴴ) |
There was a problem hiding this comment.
Same question with remove_lq_null_gauge_dependence!?
| ΔD2 = Diagonal(randn!(similar(A, complex(T), m))) | ||
| return DV, (ΔD, ΔV), (ΔD2, ΔV) | ||
| ΔV = remove_eig_gauge_dependence!(ΔV, D, V) | ||
| ΔD = Diagonal(randn!(similar(A, complex(T), m))) |
There was a problem hiding this comment.
| ΔD = Diagonal(randn!(similar(A, complex(T), m))) | |
| ΔD = Diagonal(randn!(similar(diagview(D)))) |
| ΔV = randn!(similar(A.diag, T, m, m)) | ||
| ΔV = remove_eiggauge_dependence!(ΔV, D, V) | ||
| ΔV = remove_eig_gauge_dependence!(ΔV, D, V) | ||
| ΔD = Diagonal(randn!(similar(A.diag, T, m))) |
There was a problem hiding this comment.
| ΔD = Diagonal(randn!(similar(A.diag, T, m))) | |
| ΔD = Diagonal(randn!(similar(D.diag))) |
or also using diagview(D) if we don't want to access .diag.
| ΔD2 = Diagonal(randn!(similar(A, real(T), m))) | ||
| return DV, (ΔD, ΔV), (ΔD2, ΔV) | ||
| ΔV = remove_eigh_gauge_dependence!(ΔV, D, V) | ||
| ΔD = Diagonal(randn!(similar(A, real(T), m))) |
There was a problem hiding this comment.
| ΔD = Diagonal(randn!(similar(A, real(T), m))) | |
| ΔD = Diagonal(randn!(similar(diagview(D)))) |
| m, n = size(A) | ||
| T = complex(eltype(A)) | ||
| D = eig_vals(A) | ||
| ΔD = randn!(similar(A, complex(T), m)) |
There was a problem hiding this comment.
| ΔD = randn!(similar(D)) |
| m, n = size(A) | ||
| T = complex(eltype(A)) | ||
| D = eig_vals(A) | ||
| ΔD = randn!(similar(A.diag, T, m)) |
There was a problem hiding this comment.
| ΔD = randn!(similar(D)) |
| m, n = size(A) | ||
| T = eltype(A) | ||
| D = eigh_vals(A) | ||
| ΔD = randn!(similar(A, real(T), m)) |
There was a problem hiding this comment.
| ΔD = randn!(similar(D)) |
| function ad_svd_compact_setup(A) | ||
| m, n = size(A) | ||
| T = eltype(A) | ||
| minmn = min(m, n) | ||
| ΔU = randn!(similar(A, T, m, minmn)) | ||
| ΔS = randn!(similar(A, real(T), minmn, minmn)) | ||
| ΔS2 = Diagonal(randn!(similar(A, real(T), minmn))) | ||
| ΔS = Diagonal(randn!(similar(A, real(T), minmn))) | ||
| ΔVᴴ = randn!(similar(A, T, minmn, n)) | ||
| U, S, Vᴴ = svd_compact(A) | ||
| ΔU, ΔVᴴ = remove_svdgauge_dependence!(ΔU, ΔVᴴ, U, S, Vᴴ) | ||
| return (U, S, Vᴴ), (ΔU, ΔS, ΔVᴴ), (ΔU, ΔS2, ΔVᴴ) | ||
| ΔU, ΔVᴴ = remove_svd_gauge_dependence!(ΔU, ΔVᴴ, U, S, Vᴴ) | ||
| return (U, S, Vᴴ), (ΔU, ΔS, ΔVᴴ) | ||
| end |
There was a problem hiding this comment.
Not sure why I cannot make a suggestion here, but anyway
function ad_svd_compact_setup(A)
U, S, Vᴴ = svd_compact(A)
ΔU = randn!(similar(U))
ΔVᴴ = randn!(similar(Vᴴ))
ΔS = Diagonal(randn!(similar(diagview(S))))
ΔU, ΔVᴴ = remove_svd_gauge_dependence!(ΔU, ΔVᴴ, U, S, Vᴴ)
return (U, S, Vᴴ), (ΔU, ΔS, ΔVᴴ)
endThere was a problem hiding this comment.
That would probably even remove the need to special case A::Diagonal below.
| @@ -324,7 +440,7 @@ | |||
| U, S, Vᴴ = svd_full(A) | |||
There was a problem hiding this comment.
It seems a bit wasteful to do svd_compact in ad_svd_compact_setup, and the again svd_full here.
| diagview(ΔSfull)[1:minmn] .= diagview(ΔS2) | ||
| diagview(ΔSfull)[1:minmn] .= diagview(ΔS) | ||
| return (U, S, Vᴴ), (ΔUfull, ΔSfull, ΔVᴴfull) | ||
| end |
There was a problem hiding this comment.
function ad_svd_full_setup(A)
U, S, Vᴴ = svd_full(A)
ΔU = randn!(similar(U))
ΔVᴴ = randn!(similar(Vᴴ))
ΔS = zero(S)
rand!(diagview(ΔS)) # although I think nonzero random contributions in all of ΔS would also just work fine
ΔU, ΔVᴴ = remove_svd_gauge_dependence!(ΔU, ΔVᴴ, U, S, Vᴴ)
return (U, S, Vᴴ), (ΔU, ΔS, ΔVᴴ)
end
| m, n = size(A) | ||
| T = eltype(A) | ||
| WP = left_polar(A) | ||
| ΔWP = (randn!(similar(A, T, m, n)), randn!(similar(A, T, n, n))) |
There was a problem hiding this comment.
| ΔWP = randn!.(similar.(WP)) |
There was a problem hiding this comment.
Since this is not about this PR, I should probably make some changes to ad_utils.jl in a separate PR.
Jutho
left a comment
Jutho
left a comment
There was a problem hiding this comment.
I left some suggestions, but I already approve!
|
Yay congrats! |
|
I didn't realize automerge was in effect. |
|
I will try to make a PR with the |
3 participants
This is a somewhat large refactor of the Mooncake test suite that hinges on the idea to wrap the in-place functions in such a way that they become admissible to the internal Mooncake testing framework. (Based on ideas suggested by @Jutho, thanks!)
The basic idea is that for a function
f!(A, input, alg)we need to:Aspace to not count finite-differences resultsinputvariablesI also somewhat reorganized the tests to get a slightly better overview of the mooncake tests.
Finally, this also allowed me to find some small mistakes in the mutating rules, which I have here fixed.
This is also relevant for #174 and #173.