diff --git a/Project.toml b/Project.toml index 8db9c7a..6499f65 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ITensorBase" uuid = "4795dd04-0d67-49bb-8f44-b89c448a1dc7" -version = "0.12.4" +version = "0.13.0" authors = ["ITensor developers and contributors"] [workspace] @@ -14,7 +14,6 @@ Combinatorics = "861a8166-3701-5b0c-9a16-15d98fcdc6aa" ConstructionBase = "187b0558-2788-49d3-abe0-74a17ed4e7c9" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MatrixAlgebraKit = "6c742aac-3347-4629-af66-fc926824e5e4" -OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" SimpleTraits = "699a6c99-e7fa-54fc-8d76-47d257e15c1d" TensorAlgebra = "68bd88dc-f39d-4e12-b2ca-f046b68fcc6a" @@ -47,7 +46,6 @@ LinearAlgebra = "1.10" MatrixAlgebraKit = "0.2, 0.3, 0.4, 0.5, 0.6" Mooncake = "0.4.202, 0.5" OMEinsumContractionOrders = "1.3" -OrderedCollections = "1.6" Random = "1.10" SimpleTraits = "0.9.4" TensorAlgebra = "0.17.5" diff --git a/docs/Project.toml b/docs/Project.toml index b5ba24f..7273ead 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -12,7 +12,7 @@ path = ".." [compat] Documenter = "1" -ITensorBase = "0.12" +ITensorBase = "0.13" ITensorFormatter = "0.2.27" Literate = "2" MatrixAlgebraKit = "0.2, 0.3, 0.4, 0.5, 0.6" diff --git a/docs/src/dev_interface.md b/docs/src/dev_interface.md index 7346fb7..56a1135 100644 --- a/docs/src/dev_interface.md +++ b/docs/src/dev_interface.md @@ -50,17 +50,17 @@ aligneddims ## Experimental These features support building and applying operators, where an operator is a tensor whose -dimension names are split into a codomain (output) set and a domain (input) set. The API is +dimension names are split into an output set and an input set. The API is still being refined and is subject to change. Build an operator with [`operator`](@ref) or allocate one with [`similar_operator`](@ref), apply it to a tensor with [`apply`](@ref), and -recover its underlying tensor and name sets with [`state`](@ref), [`codomainnames`](@ref), -and [`domainnames`](@ref). +recover its underlying tensor and name sets with [`state`](@ref), [`outputnames`](@ref), +and [`inputnames`](@ref). ```@docs; canonical=false operator similar_operator apply state -codomainnames -domainnames +outputnames +inputnames ``` diff --git a/examples/Project.toml b/examples/Project.toml index 1b8dfcb..998e0a1 100644 --- a/examples/Project.toml +++ b/examples/Project.toml @@ -6,5 +6,5 @@ MatrixAlgebraKit = "6c742aac-3347-4629-af66-fc926824e5e4" path = ".." [compat] -ITensorBase = "0.12" +ITensorBase = "0.13" MatrixAlgebraKit = "0.2, 0.3, 0.4, 0.5, 0.6" diff --git a/src/ITensorBase.jl b/src/ITensorBase.jl index 0e06acf..95c71fd 100644 --- a/src/ITensorBase.jl +++ b/src/ITensorBase.jl @@ -1,9 +1,9 @@ module ITensorBase export AbstractNamedTensor, NamedTensor, AbstractITensor, ITensor, Index, - NamedUnitRange, aligndims, aligneddims, apply, codomainnames, commonind, commoninds, - dimnames, dimnametype, domainnames, hascommoninds, id, inds, mapinds, named, nameddims, - noncommonind, noncommoninds, noprime, operator, prime, replaceinds, sim, + NamedUnitRange, aligndims, aligneddims, apply, commonind, commoninds, + dimnames, dimnametype, hascommoninds, id, inds, inputnames, mapinds, named, nameddims, + noncommonind, noncommoninds, noprime, operator, outputnames, prime, replaceinds, sim, similar_operator, state, trycommonind, trynoncommonind, tryuniqueind, uniqueind, uniqueinds, unioninds, uniquename diff --git a/src/namedtensoroperator.jl b/src/namedtensoroperator.jl index 77cc09b..6759bae 100644 --- a/src/namedtensoroperator.jl +++ b/src/namedtensoroperator.jl @@ -1,6 +1,5 @@ using Base.Broadcast: Broadcast as BC, Broadcasted using LinearAlgebra: LinearAlgebra as LA -using OrderedCollections: OrderedDict using Random: Random # Named dimension operator minimal interface. @@ -10,9 +9,9 @@ using Random: Random """ state(a) -The underlying tensor of a named operator, with its codomain/domain structure -forgotten. An operator carries a tensor together with a pairing of its codomain and -domain dimension names (its Choi, or state, representation). `state` returns that tensor +The underlying tensor of a named operator, with its output/input structure +forgotten. An operator carries a tensor together with a pairing of its output and +input dimension names (its Choi, or state, representation). `state` returns that tensor on its own. For a plain tensor that is not an operator, `state` returns it unchanged. # Examples @@ -24,69 +23,104 @@ julia> state(a) == a true ``` -See also [`operator`](@ref), [`codomainnames`](@ref), [`domainnames`](@ref). +See also [`operator`](@ref), [`outputnames`](@ref), [`inputnames`](@ref). """ state(a) = throw(MethodError(state, (a,))) -# Operator representation of the named state given pairs of named codomain and domain indices. -operator(a, codomain, domain) = throw(MethodError(operator, (a, codomain, domain))) +# Operator representation of the named state given pairs of named output and input indices. +operator(a, output, input) = throw(MethodError(operator, (a, output, input))) -# Get the codomain dimension names of the operator. +# Get the output dimension names of the operator. """ - codomainnames(a) + outputnames(a) -The codomain (output) dimension names of an operator `a`. An operator pairs each of its -codomain names with a domain name. Applying the operator contracts over the domain and -leaves the codomain. +The output dimension names of an operator `a`. An operator pairs each of its +output names with an input name. Applying the operator contracts over the input and +leaves the output. A plain tensor is a trivial operator with no pairing, so its output +names are empty. # Examples ```jldoctest julia> op = operator(zeros(2, 2), ("i",), ("j",)); -julia> collect(codomainnames(op)) +julia> outputnames(op) 1-element Vector{String}: "i" ``` -See also [`domainnames`](@ref), [`operator`](@ref), [`apply`](@ref). +See also [`inputnames`](@ref), [`operator`](@ref), [`apply`](@ref). """ -codomainnames(a) = throw(MethodError(codomainnames, (a,))) +outputnames(a::AbstractNamedTensor{DimName}) where {DimName} = DimName[] -# Get the domain dimension names of the operator. +# Get the input dimension names of the operator. """ - domainnames(a) + inputnames(a) -The domain (input) dimension names of an operator `a`. These are the names contracted -over when the operator is applied to a tensor. +The input dimension names of an operator `a`. These are the names contracted +over when the operator is applied to a tensor. A plain tensor is a trivial operator with +no pairing, so its input names are empty. # Examples ```jldoctest julia> op = operator(zeros(2, 2), ("i",), ("j",)); -julia> collect(domainnames(op)) +julia> inputnames(op) 1-element Vector{String}: "j" ``` -See also [`codomainnames`](@ref), [`operator`](@ref), [`apply`](@ref). +See also [`outputnames`](@ref), [`operator`](@ref), [`apply`](@ref). """ -domainnames(a) = throw(MethodError(domainnames, (a,))) +inputnames(a::AbstractNamedTensor{DimName}) where {DimName} = DimName[] + +# `outputname(a, i, default)` returns the output name paired with input name `i`, and +# `inputname(a, i, default)` returns the input name paired with output name `i`, each +# returning `default` when `i` is not one of the operator's paired names (mirroring +# `Base.get`). A tensor with no operator pairing returns `default` for every name. +outputname(a::AbstractNamedTensor, i, default) = default +inputname(a::AbstractNamedTensor, i, default) = default +# The two-argument forms return the paired name and throw when `i` is unpaired (a name on +# no wire, or any name of a plain tensor), like indexing versus `get`. `nothing` is never a +# valid name, so it is the "unpaired" sentinel; the `@something` fallback only builds the +# error on a miss. +function outputname(a::AbstractNamedTensor, i) + return @something outputname(a, i, nothing) throw( + ArgumentError("`$i` is not a paired input name") + ) +end +function inputname(a::AbstractNamedTensor, i) + return @something inputname(a, i, nothing) throw( + ArgumentError("`$i` is not a paired output name") + ) +end -# Given a domain dimension name, return the corresponding codomain dimension name. -# If it doesn't exist, return the index itself. -get_codomain_name(a, i) = throw(MethodError(get_codomain_name, (a, i))) -# Given a codomain dimension name, return the corresponding domain dimension name. -# If it doesn't exist, return the index itself. -get_domain_name(a, i) = throw(MethodError(get_domain_name, (a, i))) +# `apply(x, y)` lands each input of `x` on a matching output (or dangling leg) of `y`, so +# every name shared between `x` and `y` must be an input of `x` that is not also an input of +# `y`. Landing an input of `x` on an input of `y`, or sharing any other name, is rejected. +function check_apply(x::AbstractNamedTensor, y::AbstractNamedTensor) + for s in intersect(dimnames(x), dimnames(y)) + if !(s in inputnames(x)) || (s in inputnames(y)) + throw( + ArgumentError( + "`apply` contracts an input of the first operator with a matching " * + "output of the second; the shared name `$s` is not such a pair." + ) + ) + end + end + return nothing +end """ apply(x::AbstractNamedTensor, y::AbstractNamedTensor) -Apply the operator `x` to `y`. This contracts the state tensors of `x` and `y` over -their shared names, then renames each surviving codomain name of `x` back to its paired -domain name, so the result carries the same names `y` would map to. Applying the -identity operator leaves `y` unchanged. +Apply the operator `x` to `y`, contracting each input of `x` with the matching output (or +dangling leg) of `y` and renaming each consumed output of `x` back to its paired input, so +the result sits on `x`'s input space. Uncontracted structure passes through: `y`'s remaining +input wires stay wires, and a part of `x` disjoint from `y` is tensored in. Applying an +operator to a bare state gives a bare state, so applying the identity operator leaves `y` +unchanged; applying it to another operator gives an operator. # Examples @@ -99,89 +133,100 @@ julia> apply(op, v) == v true ``` -See also [`operator`](@ref), [`state`](@ref), [`codomainnames`](@ref), -[`domainnames`](@ref). +See also [`operator`](@ref), [`state`](@ref), [`outputnames`](@ref), +[`inputnames`](@ref). """ function apply(x::AbstractNamedTensor, y::AbstractNamedTensor) - xy = state(x) * state(y) - return mapdimnames(xy) do i - return get_domain_name(x, i) - end -end - -function apply_dag(x::AbstractNamedTensor, y::AbstractNamedTensor) - xy = state(x) * state(y) - return mapdimnames(xy) do i - return get_codomain_name(y, i) - end -end - -# TODO: Define versions that accept codomain and domain names, -# i.e. `transpose(a, codomain, domain)` and `adjoint(a, codomain, domain)` (?). + check_apply(x, y) + xy = x * y + relabels = [ + ox => ix for (ox, ix) in zip(outputnames(x), inputnames(x)) if ix in dimnames(y) + ] + result = replacedimnames(xy, relabels...) + # A result with no surviving pairing is a plain state (e.g. an operator applied to a + # bare state), so return it unwrapped. + return isempty(outputnames(result)) ? state(result) : result +end + +# TODO: Define versions that accept output and input names, +# i.e. `transpose(a, output, input)` and `adjoint(a, output, input)` (?). function Base.transpose(a::AbstractNamedTensor) - c = codomainnames(a) - d = domainnames(a) - a_map = merge(Dict(c .=> d), Dict(d .=> c)) + out = outputnames(a) + inp = inputnames(a) + a_map = merge(Dict(out .=> inp), Dict(inp .=> out)) a′ = mapdimnames(state(a)) do i return get(a_map, i, i) end - return operator(a′, c, d) + return operator(a′, out, inp) end function Base.adjoint(a::AbstractNamedTensor) return transpose(conj(a)) end -function product(x::AbstractNamedTensor, y::AbstractNamedTensor) - c = codomainnames(x) - d = domainnames(x) - c′ = uniquename.(c) - x′_map = merge(Dict(c .=> c′), Dict(d .=> c)) - x′ = mapdimnames(parent(x)) do i - return get(x′_map, i, i) +# A wire of `a` and a wire of `b` may either be the same wire (composed by `product`) or +# share no name (independent). Any other overlap is rejected: two different wires ending on +# the same input site or leaving the same output site is ambiguous, and connecting two +# operators end-to-end through a shared bond is a job for `*`, not `product`. +function check_product(a::AbstractNamedTensor, b::AbstractNamedTensor) + for (out_a, in_a) in zip(outputnames(a), inputnames(a)), + (out_b, in_b) in zip(outputnames(b), inputnames(b)) + + overlap = !isdisjoint((out_a, in_a), (out_b, in_b)) + match = (out_a == out_b) && (in_a == in_b) + if overlap && !match + throw( + ArgumentError( + "`product` operands have overlapping but mismatched wires " * + "($out_a←$in_a vs $out_b←$in_b); wires must match to compose or share " * + "no name. Use `*` to connect operators end-to-end." + ) + ) + end end - x′y = x′ * parent(y) - x′y_map = Dict(c′ .=> c) - xy = mapdimnames(x′y) do i - return get(x′y_map, i, i) + return nothing +end + +# `product(a, b)` composes operators on matching sites. A name that is an input of both +# operands marks a shared site: its (matching) wire is welded — `a`'s input and `b`'s +# output are renamed to a fresh bond — so the site composes instead of contracting +# input↔input. Every other shared name — dangling Kraus/batch dimensions, a plain tensor's +# legs — contracts like `*`. See `check_product` for the wires that are rejected. +function product(a::AbstractNamedTensor, b::AbstractNamedTensor) + check_product(a, b) + a′, b′ = a, b + for s in intersect(inputnames(a), inputnames(b)) + bond = uniquename(s) + a′ = replacedimnames(a′, s => bond) # a's input site → bond + b′ = replacedimnames(b′, outputname(b, s, s) => bond) # b's matching output → bond end - return operator(xy, c, d) + return operator_product(a′, b′) end -struct Bijection{Codomain, Domain} <: AbstractDict{Domain, Codomain} - domain_to_codomain::OrderedDict{Domain, Codomain} - codomain_to_domain::OrderedDict{Codomain, Domain} -end -function Bijection(domain, codomain) - pairs = domain .=> codomain - domain_to_codomain = OrderedDict(pairs) - codomain_to_domain = OrderedDict(reverse(kv) for kv in pairs) - return Bijection(domain_to_codomain, codomain_to_domain) -end -function Base.get(b::Bijection, k, default) - return get(b.domain_to_codomain, k, default) -end -function inverse(b::Bijection) - return Bijection(b.codomain_to_domain, b.domain_to_codomain) -end -# Both accessors return the `keys(::OrderedDict)` of the dict that has the -# requested side as its key type, so the result is a `Base.KeySet` that -# compares correctly with `==`. The two dicts are constructed from the same -# pairs in the constructor, so `codomain(b)[i]` and `domain(b)[i]` remain in -# lock-step positional order. -function codomain(b::Bijection) - return keys(b.codomain_to_domain) -end -function domain(b::Bijection) - return keys(b.domain_to_codomain) -end -Base.iterate(b::Bijection) = iterate(b.domain_to_codomain) -Base.iterate(b::Bijection, state) = iterate(b.domain_to_codomain, state) -Base.length(b::Bijection) = length(b.domain_to_codomain) - -struct NamedTensorOperator{DimName, P <: AbstractNamedTensor{DimName}, D, C} <: +struct NamedTensorOperator{DimName, P <: AbstractNamedTensor{DimName}} <: AbstractNamedTensor{DimName} parent::P - dimnames_bijection::Bijection{D, C} + # `outputnames[i]` is paired with `inputnames[i]` (positional pairing), so the two + # vectors always have equal length. Names in the parent that appear in neither vector + # are dangling (not part of the operator pairing). + outputnames::Vector{DimName} + inputnames::Vector{DimName} +end + +function NamedTensorOperator( + parent::AbstractNamedTensor{DimName}, outputnames, inputnames + ) where {DimName} + if length(outputnames) != length(inputnames) + throw( + ArgumentError( + "Operator `outputnames` and `inputnames` must have equal length " * + "(positional pairing), got $(length(outputnames)) and " * + "$(length(inputnames))." + ) + ) + end + return NamedTensorOperator{DimName, typeof(parent)}( + parent, collect(DimName, outputnames), collect(DimName, inputnames) + ) end state(a::AbstractNamedTensor) = a @@ -190,106 +235,125 @@ Base.parent(a::NamedTensorOperator) = state(a) unnamed(a::NamedTensorOperator) = unnamed(state(a)) dimnames(a::NamedTensorOperator) = dimnames(state(a)) -function NamedTensorOperator(a::AbstractNamedTensor, codomainnames, domainnames) - return NamedTensorOperator(a, Bijection(domainnames, codomainnames)) -end - parenttype(type::Type{<:NamedTensorOperator}) = fieldtype(type, :parent) statetype(type::Type{<:NamedTensorOperator}) = parenttype(type) function nameddimsof(a::NamedTensorOperator, b::AbstractArray) - return NamedTensorOperator(nameddimsof(state(a), b), a.dimnames_bijection) + return NamedTensorOperator(nameddimsof(state(a), b), a.outputnames, a.inputnames) +end + +outputnames(a::NamedTensorOperator) = a.outputnames +inputnames(a::NamedTensorOperator) = a.inputnames + +# Relabeling an operator's dimension names updates both its state and its pairing (the +# generic `AbstractNamedTensor` methods reconstruct via `nameddims` and would drop the +# pairing). `mapdimnames(f, op)` routes through the function form via the generic +# `mapdimnames(f, ::AbstractNamedTensor) = replacedimnames(f, ...)`. +function replacedimnames(op::NamedTensorOperator, replacements::Pair...) + isempty(replacements) && return op + ps = map(p -> name(first(p)) => name(last(p)), replacements) + return operator( + replacedimnames(state(op), ps...), + replace(outputnames(op), ps...), + replace(inputnames(op), ps...) + ) +end +function replacedimnames(f, op::NamedTensorOperator) + return operator( + replacedimnames(f, state(op)), map(f, outputnames(op)), map(f, inputnames(op)) + ) end -codomainnames(a::NamedTensorOperator) = codomain(a.dimnames_bijection) -domainnames(a::NamedTensorOperator) = domain(a.dimnames_bijection) - -function get_codomain_name(a::NamedTensorOperator, i) - return get(a.dimnames_bijection, i, i) +# Linear scan over the paired name vectors: for the small number of paired legs an +# operator carries, a `findfirst` is faster than maintaining a hashed lookup, and it +# keeps the operator representation to just the two vectors. +function outputname(a::NamedTensorOperator, i, default) + k = findfirst(==(i), a.inputnames) + return isnothing(k) ? default : a.outputnames[k] end -function get_domain_name(a::NamedTensorOperator, i) - return get(inverse(a.dimnames_bijection), i, i) +function inputname(a::NamedTensorOperator, i, default) + k = findfirst(==(i), a.outputnames) + return isnothing(k) ? default : a.inputnames[k] end """ - operator(a, codomain, domain) + operator(a, output, input) Build a named operator from a tensor (or plain array) `a` by partitioning its dimension -names into a `codomain` (output) set and a `domain` (input) set. The operator pairs each -codomain name with a domain name, so it can be applied to a tensor with -[`apply`](@ref), contracting over the domain. `codomain` and `domain` may be given as -dimension names or as named ranges such as `Index`es. Recover the underlying tensor -with [`state`](@ref) and the name sets with [`codomainnames`](@ref) and -[`domainnames`](@ref). +names into an `output` set and an `input` set. The operator pairs each output name with +an input name, so it can be applied to a tensor with [`apply`](@ref), contracting over +the input. `output` and `input` may be given as dimension names or as named ranges such +as `Index`es. Recover the underlying tensor with [`state`](@ref) and the name sets with +[`outputnames`](@ref) and [`inputnames`](@ref). # Examples ```jldoctest julia> op = operator(zeros(2, 2), ("i",), ("j",)); -julia> collect(codomainnames(op)) +julia> outputnames(op) 1-element Vector{String}: "i" -julia> collect(domainnames(op)) +julia> inputnames(op) 1-element Vector{String}: "j" ``` -See also [`state`](@ref), [`codomainnames`](@ref), [`domainnames`](@ref), +See also [`state`](@ref), [`outputnames`](@ref), [`inputnames`](@ref), [`apply`](@ref), [`similar_operator`](@ref). """ function operator end -# `codomain` and `domain` may be given as dimension names or as named ranges +# `output` and `input` may be given as dimension names or as named ranges # (such as `Index`es); `name` maps the latter to their names and leaves names as-is. # TODO: Unify these two functions. -function operator(a::AbstractArray, codomain, domain) - codomain, domain = name.(codomain), name.(domain) - na = nameddims(a, (codomain..., domain...)) - return operator(na, codomain, domain) +function operator(a::AbstractArray, output, input) + output, input = name.(output), name.(input) + na = nameddims(a, (output..., input...)) + return operator(na, output, input) end -function operator(a::AbstractNamedTensor, codomain, domain) - return NamedTensorOperator(a, name.(codomain), name.(domain)) +function operator(a::AbstractNamedTensor, output, input) + return NamedTensorOperator(a, name.(output), name.(input)) end # Operator-preserving contraction. Contracting two named arrays sums over their -# shared names, so the result keeps each operand's surviving codomain/domain +# shared names, so the result keeps each operand's surviving output/input # structure. A non-operator tensor contributes no pairs (all its names are # dangling from the operator point of view). The result is always an -# `NamedTensorOperator`, even when its codomain and domain both come out empty, so +# `NamedTensorOperator`, even when its output and input both come out empty, so # the product type does not depend on the runtime names being contracted. -operator_pairs(a::NamedTensorOperator) = a.dimnames_bijection.domain_to_codomain +operator_pairs(a::NamedTensorOperator) = a.inputnames .=> a.outputnames operator_pairs(a::AbstractNamedTensor) = () -# Compose the codomain/domain of `a * b`. The `domain => codomain` pairs of both +# Compose the output/input of `a * b`. The `input => output` pairs of both # operands form a graph of maximum degree two (each name is paired at most once -# per operand), so it is a disjoint union of simple paths. Each surviving domain -# name reaches its surviving codomain partner by following the pairing through +# per operand), so it is a disjoint union of simple paths. Each surviving input +# name reaches its surviving output partner by following the pairing through # any contracted (shared) names. A name whose chain dead-ends on a contracted # index is left dangling, so the result is well defined for any contraction. -function product_codomain_domain(a::AbstractNamedTensor, b::AbstractNamedTensor) +function product_output_input(a::AbstractNamedTensor, b::AbstractNamedTensor) shared = intersect(dimnames(a), dimnames(b)) pairs = collect(Iterators.flatten((operator_pairs(a), operator_pairs(b)))) forward = Dict(pairs) - domain = eltype(keys(forward))[] - codomain = eltype(values(forward))[] - for (d, c) in pairs - d in shared && continue - while c in shared && haskey(forward, c) - c = forward[c] + input = eltype(keys(forward))[] + output = eltype(values(forward))[] + for (inp, out) in pairs + inp in shared && continue + while out in shared && haskey(forward, out) + out = forward[out] end - c in shared && continue - push!(domain, d) - push!(codomain, c) + out in shared && continue + push!(input, inp) + push!(output, out) end - return codomain, domain + return output, input end function operator_product(a::AbstractNamedTensor, b::AbstractNamedTensor) ab = state(a) * state(b) - codomain, domain = product_codomain_domain(a, b) - return operator(ab, codomain, domain) + output, input = product_output_input(a, b) + return operator(ab, output, input) end Base.:*(a::NamedTensorOperator, b::NamedTensorOperator) = operator_product(a, b) @@ -302,12 +366,12 @@ Base.:*(a::AbstractNamedTensor, b::NamedTensorOperator) = operator_product(a, b) # preservation since they lower to broadcasting. `copy` / `similar` unwrap each # operator operand to its `state` (the shared `NamedTensorStyle` machinery does this via # `unnamed`), build the `NamedTensor` result, then rewrap as an operator using the -# codomain/domain split recovered from the operands. +# output/input split recovered from the operands. function BC.BroadcastStyle(arraytype::Type{<:NamedTensorOperator}) return NamedTensorOperatorStyle{ndims(arraytype)}() end -# Recover the codomain/domain split shared by all operator operands of `bc`, +# Recover the output/input split shared by all operator operands of `bc`, # erroring if any two operators disagree. operator_operands(bc::Broadcasted) = operator_operands(bc.args...) function operator_operands(arg::NamedTensorOperator, args...) @@ -319,36 +383,36 @@ end operator_operands(arg, args...) = operator_operands(args...) operator_operands() = () -function broadcast_operator_codomain_domain(bc::Broadcasted) +function broadcast_operator_output_input(bc::Broadcasted) ops = operator_operands(bc) op1 = first(ops) - cod1 = codomainnames(op1) - dom1 = domainnames(op1) + out1 = outputnames(op1) + inp1 = inputnames(op1) for op in Base.tail(ops) - (issetequal(codomainnames(op), cod1) && issetequal(domainnames(op), dom1)) || + (issetequal(outputnames(op), out1) && issetequal(inputnames(op), inp1)) || throw( ArgumentError( - "Operator operands disagree on their codomain/domain split: " * - "$((cod1, dom1)) vs $((codomainnames(op), domainnames(op))). " * + "Operator operands disagree on their output/input split: " * + "$((out1, inp1)) vs $((outputnames(op), inputnames(op))). " * "Broadcasting operators requires a matching split." ) ) end - return cod1, dom1 + return out1, inp1 end function Base.copy(bc::Broadcasted{<:NamedTensorOperatorStyle}) - cod, dom = broadcast_operator_codomain_domain(bc) + out, inp = broadcast_operator_output_input(bc) result = copy(statebroadcasted(bc)) - return operator(result, cod, dom) + return operator(result, out, inp) end for f in MATRIX_FUNCTIONS @eval begin function Base.$f(a::NamedTensorOperator) - c = codomainnames(a) - d = domainnames(a) - return operator($f(state(a), c, d), c, d) + out = outputnames(a) + inp = inputnames(a) + return operator($f(state(a), out, inp), out, inp) end end end @@ -367,10 +431,10 @@ const _gram_eigh_full_operator_docstring = """ Gram factorization of a Hermitian positive semi-definite named operator `a`, returning `x` such that `x * x_cod ≈ state(a)`, where `x_cod` is -`conj(x)` with its domain dimension names replaced by the corresponding -codomain names of `a`. `x` carries `a`'s domain dimension names and a -fresh trailing rank name. The codomain and domain partition is taken from -`codomainnames(a)` and `domainnames(a)`. +`conj(x)` with its input dimension names replaced by the corresponding +output names of `a`. `x` carries `a`'s input dimension names and a +fresh trailing rank name. The output and input partition is taken from +`outputnames(a)` and `inputnames(a)`. `kwargs` are forwarded to `TensorAlgebra.MatrixAlgebra.gram_eigh_full` on the underlying named array (e.g. `atol`, `rtol`). @@ -401,8 +465,8 @@ const _gram_eigh_full_with_pinv_operator_docstring = """ Like `TensorAlgebra.MatrixAlgebra.gram_eigh_full`, but additionally returns a named array `y` that is a left inverse of `x`: `y * x ≈ I` on the rank subspace (equal to the identity when `a` is full rank). The -codomain and domain partition is taken from `codomainnames(a)` and -`domainnames(a)`. +output and input partition is taken from `outputnames(a)` and +`inputnames(a)`. # Examples @@ -433,36 +497,36 @@ for f in (:gram_eigh_full, :gram_eigh_full_with_pinv) doc_sym = Symbol("_", f, "_operator_docstring") @eval begin @doc $doc_sym function MA.$f(a::NamedTensorOperator; kwargs...) - return MA.$f(state(a), codomainnames(a), domainnames(a); kwargs...) + return MA.$f(state(a), outputnames(a), inputnames(a); kwargs...) end end end function MAK.project_hermitian(a::NamedTensorOperator; kwargs...) - h = MAK.project_hermitian(state(a), codomainnames(a), domainnames(a); kwargs...) - return operator(h, codomainnames(a), domainnames(a)) + h = MAK.project_hermitian(state(a), outputnames(a), inputnames(a); kwargs...) + return operator(h, outputnames(a), inputnames(a)) end for f in (:sqrth_safe, :invsqrth_safe) @eval function MA.$f(a::NamedTensorOperator; kwargs...) - x = MA.$f(state(a), codomainnames(a), domainnames(a); kwargs...) - return operator(x, codomainnames(a), domainnames(a)) + x = MA.$f(state(a), outputnames(a), inputnames(a); kwargs...) + return operator(x, outputnames(a), inputnames(a)) end end function MA.sqrth_invsqrth_safe(a::NamedTensorOperator; kwargs...) - x, y = MA.sqrth_invsqrth_safe(state(a), codomainnames(a), domainnames(a); kwargs...) - return operator(x, codomainnames(a), domainnames(a)), - operator(y, codomainnames(a), domainnames(a)) + x, y = MA.sqrth_invsqrth_safe(state(a), outputnames(a), inputnames(a); kwargs...) + return operator(x, outputnames(a), inputnames(a)), + operator(y, outputnames(a), inputnames(a)) end """ Base.one(op::NamedTensorOperator) -> Id -Return the identity operator with the same codomain/domain names and shape as +Return the identity operator with the same output/input names and shape as `op`. `op` is treated as a shape prototype and is not mutated. The identity acts as the multiplicative identity for `ITensorBase.apply`: it -contracts on the domain names and renames the resulting codomain names back to -the domain names, leaving the input unchanged. +contracts on the input names and renames the resulting output names back to +the input names, leaving the input unchanged. # Examples @@ -482,42 +546,41 @@ true ``` """ function Base.one(op::NamedTensorOperator) - co, dom = codomainnames(op), domainnames(op) - return operator(one(state(op), co, dom), co, dom) + out, inp = outputnames(op), inputnames(op) + return operator(one(state(op), out, inp), out, inp) end """ LinearAlgebra.tr(op::NamedTensorOperator) -> scalar -Trace of a named operator: contracts each codomain index with its paired domain index and -sums the diagonal, over the operator's intrinsic codomain/domain split. +Trace of a named operator: contracts each output index with its paired input index and +sums the diagonal, over the operator's intrinsic output/input split. """ function LA.tr(op::NamedTensorOperator) a = state(op) byname = Dict(name(i) => i for i in inds(a)) - codomain = map(n -> byname[n], Tuple(codomainnames(op))) - domain = map(n -> byname[n], Tuple(domainnames(op))) - return LA.tr(a, codomain, domain) + output = map(n -> byname[n], Tuple(outputnames(op))) + input = map(n -> byname[n], Tuple(inputnames(op))) + return LA.tr(a, output, input) end # === similar_operator === # -# Allocate an operator with the user-supplied axes as the domain (input). The -# codomain (output) shares the domain direction and either takes -# explicitly-supplied names or fresh `uniquename` outputs. The 5-arg form is -# canonical, the others fill in defaults. The bra/ket flip on the storage side -# is handled inside `TA.similar_map`. +# Allocate an operator with the user-supplied axes as the input. The output +# shares the input direction and either takes explicitly-supplied names or fresh +# `uniquename` outputs. The 5-arg form is canonical, the others fill in defaults. +# The bra/ket flip on the storage side is handled inside `TA.similar_map`. """ - similar_operator(prototype, [T,] unnamed_domain_axes, [codomain_names,] domain_names) -> op - similar_operator(prototype, [T,] named_domain_axes) -> op + similar_operator(prototype, [T,] unnamed_input_axes, [outputnames,] inputnames) -> op + similar_operator(prototype, [T,] named_input_axes) -> op Allocate an operator-shaped named array with undefined data, with the -user-supplied side as the domain (input) and a matching codomain (output). -Element type defaults to `eltype(prototype)`. Codomain names default to fresh +user-supplied side as the input and a matching output. +Element type defaults to `eltype(prototype)`. Output names default to fresh `uniquename`-generated names. The first form takes unnamed (raw) axes and explicit names, the second takes already-named axes and reuses their names as -the domain. Storage layout (including the bra/ket flip on the domain side for +the input. Storage layout (including the bra/ket flip on the input side for graded axes) is delegated to `TensorAlgebra.similar_map`. # Examples @@ -525,7 +588,7 @@ graded axes) is delegated to `TensorAlgebra.similar_map`. ```jldoctest julia> op = similar_operator(zeros(2, 2), (Base.OneTo(2),), (:i,), (:j,)); -julia> collect(domainnames(op)) +julia> inputnames(op) 1-element Vector{Symbol}: :j ``` @@ -533,36 +596,36 @@ julia> collect(domainnames(op)) See also [`operator`](@ref), [`uniquename`](@ref). """ function similar_operator( - prototype, ::Type{T}, unnamed_domain_axes, codomain_names, domain_names + prototype, ::Type{T}, unnamed_input_axes, outputnames, inputnames ) where {T} - codomain_axes = named.(unnamed_domain_axes, codomain_names) - domain_axes = named.(unnamed_domain_axes, domain_names) - raw = TA.similar_map(prototype, T, codomain_axes, domain_axes) - return operator(raw, codomain_names, domain_names) + output_axes = named.(unnamed_input_axes, outputnames) + input_axes = named.(unnamed_input_axes, inputnames) + raw = TA.similar_map(prototype, T, output_axes, input_axes) + return operator(raw, outputnames, inputnames) end function similar_operator( - prototype, ::Type{T}, unnamed_domain_axes, domain_names + prototype, ::Type{T}, unnamed_input_axes, inputnames ) where {T} - codomain_names = uniquename.(domain_names) + outputnames = uniquename.(inputnames) return similar_operator( - prototype, T, unnamed_domain_axes, codomain_names, domain_names + prototype, T, unnamed_input_axes, outputnames, inputnames ) end -function similar_operator(prototype, ::Type{T}, named_domain_axes) where {T} +function similar_operator(prototype, ::Type{T}, named_input_axes) where {T} return similar_operator( - prototype, T, unnamed.(named_domain_axes), name.(named_domain_axes) + prototype, T, unnamed.(named_input_axes), name.(named_input_axes) ) end -function similar_operator(prototype, unnamed_domain_axes, codomain_names, domain_names) +function similar_operator(prototype, unnamed_input_axes, outputnames, inputnames) return similar_operator( - prototype, eltype(prototype), unnamed_domain_axes, codomain_names, domain_names + prototype, eltype(prototype), unnamed_input_axes, outputnames, inputnames ) end -function similar_operator(prototype, unnamed_domain_axes, domain_names) - return similar_operator(prototype, eltype(prototype), unnamed_domain_axes, domain_names) +function similar_operator(prototype, unnamed_input_axes, inputnames) + return similar_operator(prototype, eltype(prototype), unnamed_input_axes, inputnames) end -function similar_operator(prototype, named_domain_axes) - return similar_operator(prototype, eltype(prototype), named_domain_axes) +function similar_operator(prototype, named_input_axes) + return similar_operator(prototype, eltype(prototype), named_input_axes) end # Forward `Random.randn!` / `Random.rand!` to the operator's state, which diff --git a/test/Project.toml b/test/Project.toml index 2093879..af44b3a 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -30,7 +30,7 @@ AbstractTrees = "0.4.5" Adapt = "4" Aqua = "0.8.9" Combinatorics = "1" -ITensorBase = "0.12" +ITensorBase = "0.13" ITensorPkgSkeleton = "0.3.42" JLArrays = "0.2, 0.3" LinearAlgebra = "1.10" diff --git a/test/test_exports.jl b/test/test_exports.jl index 5ccd5a4..e575ae4 100644 --- a/test/test_exports.jl +++ b/test/test_exports.jl @@ -4,10 +4,10 @@ using Test: @test, @testset exports = [ :ITensorBase, :AbstractNamedTensor, :NamedTensor, :AbstractITensor, :ITensor, :Index, :NamedUnitRange, - :aligndims, :aligneddims, :apply, :codomainnames, :commonind, :commoninds, - :dimnames, :dimnametype, :domainnames, :hascommoninds, :id, - :inds, :mapinds, :named, :nameddims, :noncommonind, :noncommoninds, :noprime, - :operator, + :aligndims, :aligneddims, :apply, :commonind, :commoninds, + :dimnames, :dimnametype, :hascommoninds, :id, + :inds, :inputnames, :mapinds, :named, :nameddims, :noncommonind, :noncommoninds, + :noprime, :operator, :outputnames, :prime, :replaceinds, :sim, :similar_operator, :state, :trycommonind, :trynoncommonind, :tryuniqueind, :uniqueind, :uniqueinds, :unioninds, :uniquename, diff --git a/test/test_operator.jl b/test/test_operator.jl index b3d3bba..8a98b2a 100644 --- a/test/test_operator.jl +++ b/test/test_operator.jl @@ -1,6 +1,6 @@ -using ITensorBase: ITensorBase as NDA, NamedTensor, NamedTensorOperator, apply, - codomainnames, dimnames, domainnames, id, nameddims, namedoneto, operator, product, - replacedimnames, similar_operator, state, unname, unnamed +using ITensorBase: ITensorBase as NDA, NamedTensor, NamedTensorOperator, apply, dimnames, + id, inputname, inputnames, nameddims, namedoneto, operator, outputname, outputnames, + product, replacedimnames, similar_operator, state, unname, unnamed using LinearAlgebra: I, norm using MatrixAlgebraKit: project_hermitian using Random: Random @@ -14,22 +14,22 @@ using Test: @test, @test_throws, @testset o = operator(randn(2, 2, 2, 2), ("i'", "j'"), ("i", "j")) @test o isa NamedTensorOperator{String} @test eltype(o) ≡ Float64 - @test issetequal(NDA.codomainnames(o), ("i'", "j'")) - @test issetequal(NDA.domainnames(o), ("i", "j")) + @test issetequal(NDA.outputnames(o), ("i'", "j'")) + @test issetequal(NDA.inputnames(o), ("i", "j")) o = operator(randn(2, 2, 2, 2), ("i'", "j'"), ("i", "j")) õ = similar(o) @test õ isa NamedTensorOperator{String} @test eltype(õ) ≡ Float64 - @test issetequal(NDA.codomainnames(õ), ("i'", "j'")) - @test issetequal(NDA.domainnames(õ), ("i", "j")) + @test issetequal(NDA.outputnames(õ), ("i'", "j'")) + @test issetequal(NDA.inputnames(õ), ("i", "j")) o = operator(randn(2, 2, 2, 2), ("i'", "j'"), ("i", "j")) õ = similar(o, Float32) @test õ isa NamedTensorOperator{String} @test eltype(õ) ≡ Float32 - @test issetequal(NDA.codomainnames(õ), ("i'", "j'")) - @test issetequal(NDA.domainnames(õ), ("i", "j")) + @test issetequal(NDA.outputnames(õ), ("i'", "j'")) + @test issetequal(NDA.inputnames(õ), ("i", "j")) o = operator(randn(2, 2, 2, 2), ("i'", "j'"), ("i", "j")) @test o isa NamedTensorOperator @@ -48,13 +48,148 @@ using Test: @test, @test_throws, @testset @test ov ≈ replacedimnames(o * v, "i'" => "i", "j'" => "j") end +@testset "wire lookups" begin + o = operator(randn(2, 2), ("i'",), ("i",)) + # Three-argument form returns the default for an unpaired name; two-argument form throws. + @test outputname(o, "i", "x") == "i'" + @test inputname(o, "i'", "x") == "i" + @test outputname(o, "z", "x") == "x" + @test outputname(o, "i") == "i'" + @test inputname(o, "i'") == "i" + @test_throws ArgumentError outputname(o, "z") + # A plain tensor has no pairing, so every two-argument lookup throws. + @test_throws ArgumentError inputname(NamedTensor(randn(2), ("i",)), "i") +end + +@testset "apply composition" begin + # Operator applied to another operator: `x` lands on `y`'s output, and the result is an + # operator on `x`'s input space keeping `y`'s input. + A = operator(randn(2, 2), ("k",), ("j",)) + B = operator(randn(2, 2), ("j",), ("i",)) + AB = apply(A, B) + @test AB isa NamedTensorOperator + @test issetequal(outputnames(AB), ("j",)) + @test issetequal(inputnames(AB), ("i",)) + Am = unname(state(A), ("k", "j")) + Bm = unname(state(B), ("j", "i")) + @test unname(state(AB), ("j", "i")) ≈ Am * Bm + + # Applying to a bare state gives a bare state. + v = NamedTensor(randn(2, 2), ("j", "i")) + Av = apply(A, v) + @test Av isa NamedTensor + @test issetequal(dimnames(Av), ("j", "i")) + @test unname(Av, ("j", "i")) ≈ Am * unname(v, ("j", "i")) + + # Disjoint apply tensors, like `*`/product. + C = operator(randn(2, 2), ("m",), ("n",)) + D = operator(randn(2, 2), ("j",), ("i",)) + CD = apply(C, D) + @test CD isa NamedTensorOperator + @test issetequal(outputnames(CD), ("m", "j")) + @test issetequal(inputnames(CD), ("n", "i")) + + # Applying into `y`'s input side (not its codomain) errors. + @test_throws ArgumentError apply( + operator(randn(2, 2), ("k",), ("j",)), operator(randn(2, 2), ("i",), ("j",)) + ) +end + +@testset "product composition" begin + # Same-space single site: product is matrix multiplication, and it is + # bidirectional (product(B, A) is the other matrix order). + A = operator(randn(2, 2), ("a'",), ("a",)) + B = operator(randn(2, 2), ("a'",), ("a",)) + Am = unname(state(A), ("a'", "a")) + Bm = unname(state(B), ("a'", "a")) + + AB = product(A, B) + @test AB isa NamedTensorOperator + @test issetequal(outputnames(AB), ("a'",)) + @test issetequal(inputnames(AB), ("a",)) + @test unname(state(AB), ("a'", "a")) ≈ Am * Bm + + BA = product(B, A) + @test unname(state(BA), ("a'", "a")) ≈ Bm * Am + + # Disjoint operators tensor into a two-site operator. + C = operator(randn(2, 2), ("i'",), ("i",)) + D = operator(randn(3, 3), ("j'",), ("j",)) + CD = product(C, D) + @test issetequal(outputnames(CD), ("i'", "j'")) + @test issetequal(inputnames(CD), ("i", "j")) + + # Partial overlap: A on sites (1, 2), B on sites (2, 3), sharing site 2 → + # a three-site operator composed on 2 and tensored on 1 and 3. + A3 = operator( + nameddims(randn(2, 2, 2, 2), ("1'", "2'", "1", "2")), + ("1'", "2'"), + ("1", "2") + ) + B3 = operator( + nameddims(randn(2, 2, 2, 2), ("2'", "3'", "2", "3")), + ("2'", "3'"), + ("2", "3") + ) + AB3 = product(A3, B3) + @test issetequal(outputnames(AB3), ("1'", "2'", "3'")) + @test issetequal(inputnames(AB3), ("1", "2", "3")) + # Value: weld A's input 2 to B's output 2' and contract. + manual = + replacedimnames(state(A3), "2" => "bond") * + replacedimnames(state(B3), "2'" => "bond") + order = ("1'", "2'", "3'", "1", "2", "3") + @test unname(state(AB3), order) ≈ unname(manual, order) + + # Batched: an operator applied to a state with a spectator index. The shared `i` + # contracts (a's input meets b's leg), `j` rides along, and the output stays on `i'`. + Aop = operator(randn(2, 2), ("i'",), ("i",)) + v = NamedTensor(randn(2, 3), ("i", "j")) + Av = product(Aop, v) + @test issetequal(dimnames(Av), ("i'", "j")) + @test isempty(outputnames(Av)) + @test isempty(inputnames(Av)) + @test unname(state(Av), ("i'", "j")) ≈ + unname(state(Aop), ("i'", "i")) * unname(v, ("i", "j")) + + # Kraus: composing two operators that share both the i-wire and a dangling Kraus + # index `j` welds the wire and sums over `j`. + K = operator(nameddims(randn(2, 2, 3), ("i'", "i", "j")), ("i'",), ("i",)) + KK = product(K, K) + @test issetequal(outputnames(KK), ("i'",)) + @test issetequal(inputnames(KK), ("i",)) + @test issetequal(dimnames(KK), ("i'", "i")) # j contracted away + Km = unname(state(K), ("i'", "i", "j")) + manual_kraus = sum(Km[:, :, jj] * Km[:, :, jj] for jj in axes(Km, 3)) + @test unname(state(KK), ("i'", "i")) ≈ manual_kraus + + # Overlapping but mismatched wires are rejected. Shared input, differing output: + @test_throws ArgumentError product( + operator(randn(2, 2), ("j",), ("i",)), operator(randn(2, 2), ("k",), ("i",)) + ) + # ...and, symmetrically, shared output, differing input: + @test_throws ArgumentError product( + operator(randn(2, 2), ("j",), ("i",)), operator(randn(2, 2), ("j",), ("k",)) + ) + + # Chaining (`a`'s output is `b`'s input) is not composition on a shared site; it is + # end-to-end contraction, so `product` rejects it and `*` handles it, returning the + # composed operator `a'←b`. + chain_a = operator(randn(2, 2), ("a'",), ("m",)) + chain_b = operator(randn(2, 2), ("m",), ("b",)) + @test_throws ArgumentError product(chain_a, chain_b) + chain = chain_a * chain_b + @test issetequal(outputnames(chain), ("a'",)) + @test issetequal(inputnames(chain), ("b",)) +end + @testset "operator from named ranges" begin - # Codomain/domain may be given as named ranges, not just names. + # Output/input may be given as named ranges, not just names. i, ip = namedoneto(2, "i"), namedoneto(2, "i'") o = operator(randn(2, 2), [ip], [i]) @test o isa NamedTensorOperator{String} - @test issetequal(codomainnames(o), ("i'",)) - @test issetequal(domainnames(o), ("i",)) + @test issetequal(outputnames(o), ("i'",)) + @test issetequal(inputnames(o), ("i",)) end @testset "one(::NamedTensorOperator)" begin @@ -63,8 +198,8 @@ end op = operator(randn(i, j, k, l), ("i", "j"), ("k", "l")) Id = one(op) @test Id isa NamedTensorOperator{String} - @test codomainnames(Id) == codomainnames(op) - @test domainnames(Id) == domainnames(op) + @test outputnames(Id) == outputnames(op) + @test inputnames(Id) == inputnames(op) Id_mat = matricize(state(Id), (i, j) => "row", (k, l) => "col") @test unname(Id_mat, ("row", "col")) ≈ I(6) end @@ -100,23 +235,23 @@ end end @testset "similar_operator" begin - # Five-arg canonical: explicit element type, axes, codomain, domain names. + # Five-arg canonical: explicit element type, axes, output, input names. op = similar_operator(randn(3, 3), Float32, (Base.OneTo(3),), ("i'",), ("i",)) @test op isa NamedTensorOperator{String} - @test issetequal(codomainnames(op), ("i'",)) - @test issetequal(domainnames(op), ("i",)) + @test issetequal(outputnames(op), ("i'",)) + @test issetequal(inputnames(op), ("i",)) - # Codomain names default to fresh `uniquename` outputs. + # Output names default to fresh `uniquename` outputs. op = similar_operator(randn(3, 3), Float64, (Base.OneTo(3),), ("i",)) @test op isa NamedTensorOperator{String} - @test issetequal(domainnames(op), ("i",)) - @test only(codomainnames(op)) != "i" + @test issetequal(inputnames(op), ("i",)) + @test only(outputnames(op)) != "i" - # Named-axes form reuses each axis's name as the domain. + # Named-axes form reuses each axis's name as the input. i = namedoneto(3, "i") op = similar_operator(randn(3, 3), Float64, (i,)) - @test issetequal(domainnames(op), ("i",)) - @test only(codomainnames(op)) != "i" + @test issetequal(inputnames(op), ("i",)) + @test only(outputnames(op)) != "i" # Element type defaults to `eltype(prototype)`. op = similar_operator(randn(ComplexF32, 3, 3), (Base.OneTo(3),), ("i'",), ("i",)) @@ -137,7 +272,7 @@ end @testset "operator-preserving broadcasting" begin # `+`, `-`, and scalar multiplication lower to broadcasting. An operator # broadcasts as itself (it is not peeled to its `state`), so these operations - # preserve the `NamedTensorOperator` wrapper and its codomain/domain bijection. + # preserve the `NamedTensorOperator` wrapper and its output/input pairing. # (Contraction `*` is operator-preserving too, in its own testset below.) o = operator(randn(2, 2), ("i'",), ("i",)) s = state(o) @@ -145,8 +280,8 @@ end for r in (o + o, o - o, -o, 2 * o, o * 2, 2 .* o, o .* 2, o ./ 2) @test r isa NamedTensorOperator - @test issetequal(codomainnames(r), ("i'",)) - @test issetequal(domainnames(r), ("i",)) + @test issetequal(outputnames(r), ("i'",)) + @test issetequal(inputnames(r), ("i",)) end @test unname(state(o + o), nms) ≈ 2 .* unname(s, nms) @@ -159,54 +294,54 @@ end @test unname(state(o ./ 2), nms) ≈ unname(s, nms) ./ 2 # `o` shares both its names with itself, so `o * o` fully contracts to a - # scalar with no surviving codomain/domain. It is still an `NamedTensorOperator` - # (with empty codomain/domain), so the product type does not depend on which + # scalar with no surviving output/input. It is still an `NamedTensorOperator` + # (with empty output/input), so the product type does not depend on which # names happen to contract. oo = o * o @test oo isa NamedTensorOperator - @test isempty(codomainnames(oo)) - @test isempty(domainnames(oo)) + @test isempty(outputnames(oo)) + @test isempty(inputnames(oo)) # Operator combined with a non-operator tensor is rejected. plain = NamedTensor(randn(2, 2), ("i'", "i")) @test_throws ArgumentError o .+ plain - # Two operators whose name sets match but whose codomain/domain split differs + # Two operators whose name sets match but whose output/input split differs # are rejected (the split would otherwise be ambiguous). o_swapped = operator(randn(2, 2), ("i",), ("i'",)) @test_throws ArgumentError o .+ o_swapped end @testset "operator-preserving contraction" begin - # A shared *dangling* leg (in neither bijection) is summed away, and the - # surviving codomain/domain of each operand combine. This is the `c† * c` + # A shared *dangling* leg (in neither pairing) is summed away, and the + # surviving output/input of each operand combine. This is the `c† * c` # hopping pattern: two operators paired over an auxiliary link. a = operator(nameddims(randn(2, 2, 3), ("i'", "i", "aux")), ["i'"], ["i"]) b = operator(nameddims(randn(2, 2, 3), ("j'", "j", "aux")), ["j'"], ["j"]) ab = a * b @test ab isa NamedTensorOperator - @test issetequal(codomainnames(ab), ("i'", "j'")) - @test issetequal(domainnames(ab), ("i", "j")) + @test issetequal(outputnames(ab), ("i'", "j'")) + @test issetequal(inputnames(ab), ("i", "j")) @test !("aux" in dimnames(ab)) @test state(ab) ≈ state(a) * state(b) - # A shared *paired* index (a's domain equals b's codomain) chains through the - # contraction: a's codomain partner pairs with b's domain partner. + # A shared *paired* index (a's input equals b's output) chains through the + # contraction: a's output partner pairs with b's input partner. A = operator(randn(2, 2), ("a'",), ("m",)) B = operator(randn(2, 2), ("m",), ("b",)) AB = A * B @test AB isa NamedTensorOperator - @test issetequal(codomainnames(AB), ("a'",)) - @test issetequal(domainnames(AB), ("b",)) + @test issetequal(outputnames(AB), ("a'",)) + @test issetequal(inputnames(AB), ("b",)) - # Applying an operator to a plain state contracts the operator's domain and + # Applying an operator to a plain state contracts the operator's input and # leaves its output dangling. The result stays an `NamedTensorOperator` with empty - # codomain/domain (the surviving `a'` leg is dangling, in neither). + # output/input (the surviving `a'` leg is dangling, in neither). v = nameddims(randn(2), ("m",)) Av = operator(randn(2, 2), ("a'",), ("m",)) * v @test Av isa NamedTensorOperator - @test isempty(codomainnames(Av)) - @test isempty(domainnames(Av)) + @test isempty(outputnames(Av)) + @test isempty(inputnames(Av)) @test issetequal(dimnames(Av), ("a'",)) end @@ -241,8 +376,8 @@ end # Hermitian part. H_op = project_hermitian(M_op) @test H_op isa NamedTensorOperator - @test codomainnames(H_op) == codomainnames(M_op) - @test domainnames(H_op) == domainnames(M_op) + @test outputnames(H_op) == outputnames(M_op) + @test inputnames(H_op) == inputnames(M_op) @test H_op ≈ M_op @test project_hermitian(operator(B, ["ket"], ["bra"])) ≈ operator((B + B') / 2, ["ket"], ["bra"]) @@ -250,8 +385,8 @@ end # The roots are again bond operators, with the same codomain/domain as the input. for X in (sqrth_safe(M_op), invsqrth_safe(M_op), sqrth_invsqrth_safe(M_op)...) @test X isa NamedTensorOperator - @test codomainnames(X) == codomainnames(M_op) - @test domainnames(X) == domainnames(M_op) + @test outputnames(X) == outputnames(M_op) + @test inputnames(X) == inputnames(M_op) end P = unnamed(state(sqrth_safe(M_op)))