Balanced gram_eigh_full convention, add one and operator-construction primitives#177
Merged
Conversation
Switch gram_eigh_full and gram_eigh_full_with_pinv from A ≈ X' * X to A ≈ X * X', i.e. X = V * sqrt(D) instead of sqrt(D) * V'. This matches the convention used by eigh and SVD (factor on the left, with codomain axes leading and a trailing rank axis), so the output of gram_eigh_full can be right-multiplied into a tensor in the operator's input space without an axis-order swap. The companion left inverse Y from gram_eigh_full_with_pinv now satisfies Y * X ≈ I on the rank subspace (it was a right inverse under the old convention).
Mirrors the existing factorization API (qr/svd/eigen/gram_eigh_full):
a layered dispatch chain accepting labels, perm-tuples, biperms, or
Val{ndims_codomain}, all funneling into a canonical worker that
matricizes, calls MatrixAlgebraKit.one!, and unmatricizes. The
result sits in canonical (codomain, domain) order.
Not exported, since exporting would clash with the implicit Base.one.
Qualify as `TensorAlgebra.one(A, ...)`.
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #177 +/- ##
==========================================
+ Coverage 80.61% 80.97% +0.36%
==========================================
Files 24 26 +2
Lines 1011 1041 +30
==========================================
+ Hits 815 843 +28
- Misses 196 198 +2
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
mtfishman
changed the title
Allows callers to supply codomain and domain axes in the same direction. Storage axes become `(codomain..., conj.(domain)...)`, baking the bra/ket flip into the primitive instead of pushing it onto every caller.
mtfishman
added a commit
to ITensor/NamedDimsArrays.jl
that referenced
this pull request
Jun 2, 2026
`TensorAlgebra.similar_map` now flips the domain direction itself when laying out storage (ITensor/TensorAlgebra.jl#177), so the wrapper passes both sides in the same direction and lets the primitive handle the bra/ket flip.
`projectto!(dest, src) -> dest` orthogonally projects `src` onto the representable subspace of `dest`, in place. Magnitude-blind and tolerance-free: the non-representable component is dropped unconditionally. `checked_projectto!(dest, src; atol, rtol) -> dest` runs the projection and then verifies via `isapprox(src, dest)`, throwing on tolerance violations. The default `projectto!` falls through to `copyto!`, which is correct for dense backends where the representable subspace is everything. Block-structured backends (`AbelianGradedArray`, `FusionTensor`) overload it to copy only the symmetry-allowed entries. This pairs with `similar_map` to give a uniform "allocate then project" idiom for operator constructors across backends, mirroring the standard `copyto!(similar(...), src)` pattern but with backend-aware projection in place of faithful copy. Co-authored-by: Claude <noreply@anthropic.com>
…ppers The data-bearing sibling of `similar_map` / `zeros_map`: allocate via `similar_map` and project `raw` in via `projectto!`. `checked_project_map` is the same on top of `checked_projectto!`, throwing if the discarded weight exceeds the tolerance. Co-authored-by: Claude <noreply@anthropic.com>
mtfishman
changed the title
Unmatricize `S` so all three SVD outputs share one backend. Previously `S` came back as the matricized form from `MatrixAlgebra.svd!!`, while `U` and `Vᴴ` went through `unmatricize`, leaving downstream consumers with mixed backends (`FusedGradedMatrix` for `S`, the unmatricized backend for the others). `S` now uses `tuplemortar(((axes(U, 2),), (axes(Vᴴ, 1),)))` as its 1+1 axis pair.
…gs to isapprox Replace the internal-design-doc framing in the four `projectto.jl` docstrings with high-level one-liners that describe the user-facing contract. `checked_projectto!` now takes `kwargs...` and forwards them to `isapprox`, instead of a fixed `atol` / `rtol` signature, leaving the tolerance policy to the caller (and to `isapprox`'s own defaults). `checked_project_map` forwards through unchanged. Co-authored-by: Claude <noreply@anthropic.com>
…hange The function forwards kwargs to `isapprox`, so when the caller doesn't pass any the behavior tracks `isapprox`'s own defaults. Flag that this is not a stability contract so we can revisit later (a tighter built-in default, or a per-backend hook). Co-authored-by: Claude <noreply@anthropic.com>
Replace the conj / graded-axes design exposition with a one-line default-implementation reference. Rename the return value from `O` (which is the convention for operators / square maps) to `M` for a general linear map. Compact the doctest example. Co-authored-by: Claude <noreply@anthropic.com>
Match the project convention of `using Module: Module [as Alias]` over `import Module [as Alias]`. No behavior change. Co-authored-by: Claude <noreply@anthropic.com>
mtfishman
added a commit
to ITensor/GradedArrays.jl
that referenced
this pull request
Jun 3, 2026
GradedArrays now references `TensorAlgebra.projectto!`, added in ITensor/TensorAlgebra.jl#177. Without the pin, CI resolves against the registered TensorAlgebra 0.9.4 which does not have `projectto!`, and precompilation fails. Drop the pin once that branch registers. Co-authored-by: Claude <noreply@anthropic.com>
Both helpers were threading `eltype(raw)` through to `similar_map`, which is exactly the default in the 3-arg `similar_map` method. Call the 3-arg form instead. Co-authored-by: Claude <noreply@anthropic.com>
Adds `TensorAlgebra.trivialrange(::Type{<:AbstractUnitRange})`, the identity range for `tensor_product_axis`. A one-dimensional range that leaves any other range of the same family unchanged when fused. Defaults to `Base.OneTo(1)`. Downstream packages overload the type-level method to return their own identity (for example a charge-0 one-dimensional sector for graded ranges).
Co-authored-by: Claude <noreply@anthropic.com>
Use `axes(S, 1)` and `axes(S, 2)` for the unmatricized S instead of `axes(U, 2)` and `axes(Vᴴ, 1)`. For non-graded backends the two are equal, but for graded operators U's rank axis and S's first axis are conj-related rather than identical, so pulling S's axes from U and Vᴴ gave it the conj of its own axes and `U * S * V` failed to contract.
Y is the pseudoinverse of X, so Y's domain contracts with X's codomain (`axes_codomain`) and must carry the conjugate axes, not `axes_codomain` itself. For non-graded backends `conj` is a no-op and the original form worked, but for graded operators it flips the isdual, and the previous version wrote forbidden-charge blocks into Y's unmatricize buffer once the operator had non-trivial sector structure (manifesting as `Block (i, j) is not stored` deep in `gram_eigh_full_with_pinv` on a 2nd-sweep BP gauge step in the U(1) AKLT validation testbed).
`trivialrange` and the `projectto!` / `checked_projectto!` / `project_map` / `checked_project_map` family already shipped without test coverage. The new `test_projectto.jl` wraps a `Rounded` mock array whose `projectto!` actually projects, so the tolerance check has something to reject. The default `projectto!` is `copyto!`, which can never fail that check on its own. Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>
Inverts the underlying flat permutation and rewraps with the same block structure. Needed by `LittleSet(invperm(s.values))` in NamedDimsArrays, which the broadcast-alignment path hits when a named array's `dimnames` come from a `BlockedTuple`. Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>
…trings Stops exporting the `!!` (may-destroy-input) variants of the MatrixAlgebra factorizations and drops them from the `gram_eigh_full` docstrings, treating them as internal. Nothing outside TensorAlgebra uses them, and they stay callable as qualified `MatrixAlgebra.svd!!`-style functions. `test_exports.jl` is updated to match.
mtfishman
added a commit
to ITensor/NamedDimsArrays.jl
that referenced
this pull request
Jun 15, 2026
…_full (#229) ## Summary - Named-operator constructors on `AbstractNamedDimsArray` / `AbstractNamedDimsOperator`: `Base.one(operator)` and `Base.one(array, codomain, domain)`, `similar_operator(prototype, [T,] axes, [codomain_names,] domain_names)` (routed through `TensorAlgebra.similar_map`), generic `Random.randn!` and `Random.rand!`. - `Base.conj` on `AbstractNamedUnitRange` and `AbstractNamedDimsArray`, forwarded through the underlying so that graded axes flip their sector arrows. Without the array-level method, `conj` falls through to an element-wise broadcast that leaves the named axes' duality flags untouched, silently producing a `(-1)` per dual bond in graded bra-ket contractions. - `domain` on the operator `Bijection` now returns the `keys` of the domain-keyed dict (a `Base.KeySet`) instead of the `values` of the codomain-keyed dict, so the result compares correctly with `==`. Both dicts are built from the same pairs, so the positional lock-step between `codomain` and `domain` is unchanged. - `LinearAlgebra.normalize`, `Base.fill!`, and `LinearAlgebra.promote_leaf_eltypes` on `AbstractNamedDimsArray` routed through `denamed`, so they hit the concrete storage instead of iterating or scalar-indexing block-structured backends. - `denamed(a::AbstractNamedDimsArray, inds)` now compares dim-name order via `Tuple` instead of `LittleSet`. The old check used `LittleSet ==`, which is set-equality, so a permutation looked identical to the identity order and the lazy-permute step got skipped, silently breaking permuted broadcast. - `gram_eigh_full` aligned with the balanced convention from ITensor/TensorAlgebra.jl#177. --------- Co-authored-by: Claude <noreply@anthropic.com>
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
gram_eigh_fullandgram_eigh_full_with_pinvflipped to the balancedA ≈ X * X'convention (was right-GramA ≈ X' * X). Breaking, but lands as a patch bump sincegram_eigh_fullis new in v0.9.4 with no known downstream users.TensorAlgebra.one(A, codomain, domain)for identity operator tensors. Not exported (clashes withBase.one).similar_map(prototype, [T,] codomain_axes, domain_axes)for allocating linear-map-shaped arrays.NamedDimsArrays.similar_operatorroutes through this (Named-operator constructors, projection verbs, and balanced gram_eigh_full NamedDimsArrays.jl#229).projectto!andchecked_projectto!for projecting into a restricted subspace.project_mapandchecked_project_mapassimilar_map+projectto!allocators.trivialrange(::Type{<:AbstractUnitRange})for the identity range undertensor_product_axis. Defaults toBase.OneTo(1), overloaded by downstream packages (for example a charge-0 graded range).!!(may-destroy-input) factorization variants inMatrixAlgebraare no longer exported and no longer appear in thegram_eigh_fulldocstrings, treating them as internal. They remain callable as qualifiedMatrixAlgebra.svd!!-style functions. Nothing outside TensorAlgebra used them.