[WIP] Name-based index-set algebra, whole-tensor index ops, and the id/uuid rename

Project.toml

@@ -1,6 +1,6 @@
 name = "ITensorBase"
 uuid = "4795dd04-0d67-49bb-8f44-b89c448a1dc7"
-version = "0.10.9"
+version = "0.11.0"
 authors = ["ITensor developers <support@itensor.org> and contributors"]
 
 [workspace]
@@ -32,6 +32,10 @@ OMEinsumContractionOrders = "6f22d1fd-8eed-4bb7-9776-e7d684900715"
 TensorKit = "07d1fe3e-3e46-537d-9eac-e9e13d0d4cec"
 TensorOperations = "6aa20fa7-93e2-5fca-9bc0-fbd0db3c71a2"
 
+[sources.TensorAlgebra]
+rev = "mf/nextgen-followups-round1"
+url = "https://github.com/ITensor/TensorAlgebra.jl"
+
 [extensions]
 ITensorBaseAdaptExt = "Adapt"
 ITensorBaseMooncakeExt = "Mooncake"

docs/Project.toml

@@ -12,7 +12,7 @@ path = ".."
 
 [compat]
 Documenter = "1"
-ITensorBase = "0.10"
+ITensorBase = "0.11"
 ITensorFormatter = "0.2.27"
 Literate = "2"
 MatrixAlgebraKit = "0.2, 0.3, 0.4, 0.5, 0.6"

examples/Project.toml

@@ -6,5 +6,5 @@ MatrixAlgebraKit = "6c742aac-3347-4629-af66-fc926824e5e4"
 path = ".."
 
 [compat]
-ITensorBase = "0.10"
+ITensorBase = "0.11"
 MatrixAlgebraKit = "0.2, 0.3, 0.4, 0.5, 0.6"

src/ITensorBase.jl

@@ -1,9 +1,10 @@
 module ITensorBase
 
 export AbstractNamedTensor, NamedTensor, AbstractITensor, ITensor, Index,
-    NamedUnitRange, aligndims, aligneddims, apply, codomainnames, dimnames,
-    dimnametype, domainnames, inds, named, nameddims, noprime, operator, prime,
-    similar_operator, state, uniquename
+    NamedUnitRange, aligndims, aligneddims, apply, codomainnames, commonind, commoninds,
+    dimnames, dimnametype, domainnames, hascommoninds, id, inds, mapinds, named, nameddims,
+    noncommoninds, noprime, operator, prime, replaceinds, sim, similar_operator, state,
+    trycommonind, trynoncommonind, uniqueind, uniqueinds, unioninds, uniquename
 using Compat: @compat
 @compat public @names
 @compat public IndexName, name, nametype, replacedimnames, setname, unnamed, unnamedtype

src/abstractnamedtensor.jl

@@ -454,6 +454,24 @@ function Base.similar(
     )
     return similar_nameddims(a, elt, inds)
 end
+
+# Rank-0 (empty named axes): a scalar tensor on `a`'s backend, e.g. a backend-matched unit
+# for a product accumulator. Spelled out separately because the tuple forms above require at
+# least one `NamedUnitRange`.
+Base.similar(a::AbstractNamedTensor, inds::Tuple{}) = similar(a, eltype(a), inds)
+function Base.similar(a::AbstractNamedTensor, elt::Type, inds::Tuple{})
+    return similar_nameddims(a, elt, inds)
+end
+
+"""
+    one(a::AbstractNamedTensor) -> AbstractNamedTensor
+
+Return a rank-0 (scalar) tensor holding `one(scalartype(a))` on `a`'s backend. This is the
+multiplicative unit matching `a`'s element type and backend (dense, graded, `TensorMap`, …),
+useful as the seed of a product accumulator.
+"""
+Base.one(a::AbstractNamedTensor) = fill!(similar(a, ()), one(scalartype(a)))
+
 function setdimnames(a::AbstractNamedTensor, dimnames)
     return nameddims(unnamed(a), dimnames)
 end
@@ -482,7 +500,12 @@ julia> dimnames(replacedimnames(a, :i => :k))
 See also [`dimnames`](@ref).
 """
 function replacedimnames end
+# `name` strips an `Index`/`NamedUnitRange` to its dimension name and passes a bare name
+# through unchanged, so an index-keyed pair (`i => j`) relabels like the name-keyed pair
+# (`name(i) => name(j)`). `dimnames(a)` holds names, so a raw-index key would never match and
+# silently no-op.
 function replacedimnames(a::AbstractNamedTensor, replacements::Pair...)
+    replacements = map(p -> name(first(p)) => name(last(p)), replacements)
     new_dimnames = replace(dimnames(a), replacements...)
     return nameddims(unnamed(a), new_dimnames)
 end
@@ -492,20 +515,161 @@ function replacedimnames(f, a::AbstractNamedTensor)
 end
 mapdimnames(f, a::AbstractNamedTensor) = replacedimnames(f, a)
 
-# Replace over `axes` (a `Tuple`) rather than `inds` (a `Vector`): `replace` on a `Vector`
-# is homogeneous and would fail to convert a replacement index backed by a different range
-# type (e.g. `UnitRange` into a `OneTo`-backed vector), whereas a `Tuple` admits the mixed
-# element types. The result is splatted into `getindex`, so only the order matters.
+"""
+    replaceinds(a::AbstractNamedTensor, replacements::Pair...)
+    replaceinds(f, a::AbstractNamedTensor)
+
+Return a tensor with the same data as `a` but with its indices relabeled. Unlike
+[`replacedimnames`](@ref), whose function `f` receives a dimension name, `replaceinds` works
+at the index level: the pair form takes `old => new` index pairs, and the function form
+relabels each index `i` using `f(i)`. Either way this is a name-only relabel, taking just the
+name of the replacement and leaving the underlying space untouched (replacing the space
+instead would scalar-index a graded axis). Pairs whose index is absent are ignored.
+
+See also [`mapinds`](@ref), [`replacedimnames`](@ref).
+"""
 function replaceinds(a::AbstractNamedTensor, replacements::Pair...)
-    new_inds = replace(axes(a), replacements...)
-    return unnamed(a)[new_inds...]
-end
-function replaceinds(f, a::AbstractNamedTensor)
-    new_inds = replace(f, axes(a))
-    return unnamed(a)[new_inds...]
+    return replacedimnames(a, replacements...)
 end
+replaceinds(f, a::AbstractNamedTensor) = replaceinds(a, map(i -> i => f(i), inds(a))...)
+
+"""
+    mapinds(f, a::AbstractNamedTensor)
+
+Return a tensor with the same data as `a` but with each index `i` relabeled using `f(i)`, a
+name-only relabel that leaves the underlying space untouched. This is the whole-tensor
+index-map primitive behind [`prime`](@ref), [`noprime`](@ref), and [`sim`](@ref).
+
+See also [`replaceinds`](@ref).
+"""
 mapinds(f, a::AbstractNamedTensor) = replaceinds(f, a)
 
+# Name-based index-set algebra (the `commoninds`/etc. surface). Layered in three:
+# small order-preserving set ops on arbitrary collections, the same ops keyed by index
+# name, and the tensor-family functions built on top.
+
+# Small-collection set operations keyed by a transform `by` (elements compare equal when
+# `by(x) == by(y)`). These scan linearly rather than building `Set`s: Base's `Set`-based ops
+# hash, and hashing a whole `Index` can be expensive or fall back to iterating a graded axis.
+# The intersect/setdiff/union/symdiff forms return elements of the first argument as `Vector`s.
+function smallintersect(a, b; by = identity)
+    return (kb = Iterators.map(by, b); [x for x in a if by(x) ∈ kb])
+end
+function smallsetdiff(a, b; by = identity)
+    return (kb = Iterators.map(by, b); [x for x in a if by(x) ∉ kb])
+end
+smallunion(a, b; by = identity) = vcat(collect(a), smallsetdiff(b, a; by))
+smallsymdiff(a, b; by = identity) = vcat(smallsetdiff(a, b; by), smallsetdiff(b, a; by))
+smallisdisjoint(a, b; by = identity) = (kb = Iterators.map(by, b); !any(x -> by(x) ∈ kb, a))
+smallissubset(a, b; by = identity) = (kb = Iterators.map(by, b); all(x -> by(x) ∈ kb, a))
+smallissetequal(a, b; by = identity) = smallissubset(a, b; by) && smallissubset(b, a; by)
+
+# The small ops keyed by index name (`by = name`) rather than full `Index` equality. On a
+# graded axis a shared bond appears as an index on one tensor and its dual (`conj`) on the
+# other (same name, opposite arrow), so full-`Index` `==` misses it while the names match.
+# On the dense backend the two coincide.
+nameintersect(a, b) = smallintersect(a, b; by = name)
+namesetdiff(a, b) = smallsetdiff(a, b; by = name)
+nameunion(a, b) = smallunion(a, b; by = name)
+namesymdiff(a, b) = smallsymdiff(a, b; by = name)
+nameisdisjoint(a, b) = smallisdisjoint(a, b; by = name)
+nameissubset(a, b) = smallissubset(a, b; by = name)
+nameissetequal(a, b) = smallissetequal(a, b; by = name)
+
+"""
+    commoninds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The indices shared by name between `a` and `b`, as a `Vector` in the order they appear in `a`.
+
+See also [`commonind`](@ref), [`uniqueinds`](@ref), [`hascommoninds`](@ref).
+"""
+commoninds(a::AbstractNamedTensor, b::AbstractNamedTensor) = nameintersect(inds(a), inds(b))
+
+"""
+    uniqueinds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The indices of `a` that do not appear by name in `b`, as a `Vector` in the order they appear
+in `a`.
+
+See also [`uniqueind`](@ref), [`commoninds`](@ref), [`noncommoninds`](@ref).
+"""
+uniqueinds(a::AbstractNamedTensor, b::AbstractNamedTensor) = namesetdiff(inds(a), inds(b))
+
+"""
+    unioninds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The union by name of the indices of `a` and `b`, as a `Vector`: the indices of `a` followed
+by the indices of `b` not already present in `a`.
+
+See also [`commoninds`](@ref), [`noncommoninds`](@ref).
+"""
+unioninds(a::AbstractNamedTensor, b::AbstractNamedTensor) = nameunion(inds(a), inds(b))
+
+"""
+    noncommoninds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The indices not shared by name between `a` and `b` (the symmetric difference), as a `Vector`:
+the indices unique to `a` followed by those unique to `b`.
+
+See also [`uniqueinds`](@ref), [`commoninds`](@ref).
+"""
+function noncommoninds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+    return namesymdiff(inds(a), inds(b))
+end
+
+"""
+    hascommoninds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+Whether `a` and `b` share any index by name.
+
+See also [`commoninds`](@ref).
+"""
+function hascommoninds(a::AbstractNamedTensor, b::AbstractNamedTensor)
+    return !nameisdisjoint(inds(a), inds(b))
+end
+
+"""
+    commonind(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The single index shared by name between `a` and `b`. Errors unless there is exactly one
+shared index. Use [`trycommonind`](@ref) to get `nothing` instead of an error.
+
+See also [`commoninds`](@ref), [`uniqueind`](@ref).
+"""
+commonind(a::AbstractNamedTensor, b::AbstractNamedTensor) = only(commoninds(a, b))
+
+"""
+    uniqueind(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The single index of `a` that does not appear by name in `b`. Errors unless there is exactly
+one such index. Use [`trynoncommonind`](@ref) to get `nothing` instead of an error.
+
+See also [`uniqueinds`](@ref), [`commonind`](@ref).
+"""
+uniqueind(a::AbstractNamedTensor, b::AbstractNamedTensor) = only(uniqueinds(a, b))
+
+"""
+    trycommonind(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The single index shared by name between `a` and `b`, or `nothing` if they share no index or
+more than one. The non-erroring counterpart of [`commonind`](@ref).
+"""
+function trycommonind(a::AbstractNamedTensor, b::AbstractNamedTensor)
+    cs = commoninds(a, b)
+    return length(cs) == 1 ? only(cs) : nothing
+end
+
+"""
+    trynoncommonind(a::AbstractNamedTensor, b::AbstractNamedTensor)
+
+The single index of `a` that does not appear by name in `b`, or `nothing` if there is no such
+index or more than one. The non-erroring counterpart of [`uniqueind`](@ref).
+"""
+function trynoncommonind(a::AbstractNamedTensor, b::AbstractNamedTensor)
+    us = uniqueinds(a, b)
+    return length(us) == 1 ? only(us) : nothing
+end
+
 # `Base.isempty(a::AbstractArray)` is defined as `length(a) == 0`,
 # which involves comparing a named integer to an unnamed integer
 # which isn't well defined.