From 21b1e786e6dee042a8902c8c7f5cf0a4cc13e2dc Mon Sep 17 00:00:00 2001
From: mtfishman <mfishman@flatironinstitute.org>
Date: Sun, 19 Jan 2025 22:07:21 -0500
Subject: [PATCH] Reorganize and simplify code

---
 src/ITensorQuantumOperatorDefinitions.jl |  15 +-
 src/has_fermion_string.jl                |   1 +
 src/itensor/has_fermion_string.jl        |  25 +
 src/itensor/op.jl                        | 138 ++++
 src/itensor/siteinds.jl                  |  78 +++
 src/itensor/state.jl                     |  24 +
 src/itensor/val.jl                       |  20 +
 src/op.jl                                |  24 +
 src/sitetype.jl                          | 830 +----------------------
 src/sitetypes/boson.jl                   |  19 +-
 src/sitetypes/fermion.jl                 |  57 +-
 src/sitetypes/generic_sites.jl           | 102 +--
 src/sitetypes/qudit.jl                   |  26 +-
 src/sitetypes/spinone.jl                 |   6 +-
 src/sitetypes/tj.jl                      | 295 ++++----
 src/space.jl                             |   7 +
 src/state.jl                             |  13 +
 src/val.jl                               |  14 +
 18 files changed, 613 insertions(+), 1081 deletions(-)
 create mode 100644 src/has_fermion_string.jl
 create mode 100644 src/itensor/has_fermion_string.jl
 create mode 100644 src/itensor/op.jl
 create mode 100644 src/itensor/siteinds.jl
 create mode 100644 src/itensor/state.jl
 create mode 100644 src/itensor/val.jl
 create mode 100644 src/op.jl
 create mode 100644 src/space.jl
 create mode 100644 src/state.jl
 create mode 100644 src/val.jl

diff --git a/src/ITensorQuantumOperatorDefinitions.jl b/src/ITensorQuantumOperatorDefinitions.jl
index daa4441..11c65a8 100644
--- a/src/ITensorQuantumOperatorDefinitions.jl
+++ b/src/ITensorQuantumOperatorDefinitions.jl
@@ -1,7 +1,12 @@
 module ITensorQuantumOperatorDefinitions
 
 include("sitetype.jl")
-include("ITensorQuantumOperatorDefinitionsChainRulesCoreExt.jl")
+include("space.jl")
+include("val.jl")
+include("state.jl")
+include("op.jl")
+include("has_fermion_string.jl")
+
 include("sitetypes/aliases.jl")
 include("sitetypes/generic_sites.jl")
 include("sitetypes/qubit.jl")
@@ -13,4 +18,12 @@ include("sitetypes/tj.jl")
 include("sitetypes/qudit.jl")
 include("sitetypes/boson.jl")
 
+include("ITensorQuantumOperatorDefinitionsChainRulesCoreExt.jl")
+
+include("itensor/siteinds.jl")
+include("itensor/val.jl")
+include("itensor/state.jl")
+include("itensor/op.jl")
+include("itensor/has_fermion_string.jl")
+
 end
diff --git a/src/has_fermion_string.jl b/src/has_fermion_string.jl
new file mode 100644
index 0000000..9620299
--- /dev/null
+++ b/src/has_fermion_string.jl
@@ -0,0 +1 @@
+has_fermion_string(::OpName, ::SiteType) = nothing
diff --git a/src/itensor/has_fermion_string.jl b/src/itensor/has_fermion_string.jl
new file mode 100644
index 0000000..22330b2
--- /dev/null
+++ b/src/itensor/has_fermion_string.jl
@@ -0,0 +1,25 @@
+using ITensorBase: Index
+
+has_fermion_string(operator::AbstractArray{<:Number}, s::Index; kwargs...)::Bool = false
+
+function has_fermion_string(opname::AbstractString, s::Index; kwargs...)::Bool
+  opname = strip(opname)
+
+  # Interpret operator names joined by *
+  # as acting sequentially on the same site
+  starpos = findfirst(isequal('*'), opname)
+  if !isnothing(starpos)
+    op1 = opname[1:prevind(opname, starpos)]
+    op2 = opname[nextind(opname, starpos):end]
+    return xor(has_fermion_string(op1, s; kwargs...), has_fermion_string(op2, s; kwargs...))
+  end
+
+  Ntags = length(tags(s))
+  stypes = _sitetypes(s)
+  opn = OpName(opname)
+  for st in stypes
+    res = has_fermion_string(opn, st)
+    !isnothing(res) && return res
+  end
+  return false
+end
diff --git a/src/itensor/op.jl b/src/itensor/op.jl
new file mode 100644
index 0000000..40a8074
--- /dev/null
+++ b/src/itensor/op.jl
@@ -0,0 +1,138 @@
+using ITensorBase: Index, ITensor, dag, prime, tags
+
+op(::OpName, ::SiteType, ::Index...; kwargs...) = nothing
+function op(
+  ::OpName, ::SiteType, ::SiteType, sitetypes_inds::Union{SiteType,Index}...; kwargs...
+)
+  return nothing
+end
+
+_sitetypes(i::Index) = _sitetypes(tags(i))
+
+function commontags(i::Index...)
+  return union(tags.(i)...)
+end
+
+op(X::AbstractArray, s::Index...) = ITensor(X, (prime.(s)..., dag.(s)...))
+
+# TODO: Delete these.
+op(opname, s::Vector{<:Index}; kwargs...) = op(opname, s...; kwargs...)
+op(s::Vector{<:Index}, opname; kwargs...) = op(opname, s...; kwargs...)
+
+function op(name::AbstractString, s::Index...; adjoint::Bool=false, kwargs...)
+  name = strip(name)
+  # TODO: filter out only commons tags
+  # if there are multiple indices
+  commontags_s = commontags(s...)
+  # first we handle the + and - algebra, which requires a space between ops to avoid clashing
+  name_split = nothing
+  @ignore_derivatives name_split = String.(split(name, " "))
+  oplocs = findall(x -> x ∈ ("+", "-"), name_split)
+  if !isempty(oplocs)
+    @ignore_derivatives !isempty(kwargs) &&
+      error("Lazy algebra on parametric gates not allowed")
+    # the string representation of algebra ops: ex ["+", "-", "+"]
+    labels = name_split[oplocs]
+    # assign coefficients to each term: ex [+1, -1, +1]
+    coeffs = [1, [(-1)^Int(label == "-") for label in labels]...]
+    # grad the name of each operator block separated by an algebra op, and do so by
+    # making sure blank spaces between opnames are kept when building the new block.
+    start, opnames = 0, String[]
+    for oploc in oplocs
+      finish = oploc
+      opnames = vcat(
+        opnames, [prod([name_split[k] * " " for k in (start + 1):(finish - 1)])]
+      )
+      start = oploc
+    end
+    opnames = vcat(
+      opnames, [prod([name_split[k] * " " for k in (start + 1):length(name_split)])]
+    )
+    # build the vector of blocks and sum
+    op_list = [
+      coeff * (op(opname, s...; kwargs...)) for (coeff, opname) in zip(coeffs, opnames)
+    ]
+    return sum(op_list)
+  end
+  # the the multiplication come after
+  oploc = findfirst("*", name)
+  if !isnothing(oploc)
+    op1, op2 = nothing, nothing
+    @ignore_derivatives begin
+      op1 = name[1:prevind(name, oploc.start)]
+      op2 = name[nextind(name, oploc.start):end]
+      if !(op1[end] == ' ' && op2[1] == ' ')
+        @warn "($op1*$op2) composite op definition `A*B` deprecated: please use `A * B` instead (with spaces)"
+      end
+    end
+    return product(op(op1, s...; kwargs...), op(op2, s...; kwargs...))
+  end
+  common_stypes = _sitetypes(commontags_s)
+  @ignore_derivatives push!(common_stypes, SiteType("Generic"))
+  opn = OpName(name)
+  for st in common_stypes
+    op_mat = op(opn, st; kwargs...)
+    if !isnothing(op_mat)
+      rs = reverse(s)
+      res = ITensor(op_mat, (prime.(rs)..., dag.(rs)...))
+      adjoint && return swapprime(dag(res), 0 => 1)
+      return res
+    end
+  end
+  return throw(
+    ArgumentError(
+      "Overload of \"op\" or \"op!\" functions not found for operator name \"$name\" and Index tags: $(tags.(s)).",
+    ),
+  )
+end
+
+function op(opname, s::Vector{<:Index}, ns::NTuple{N,Integer}; kwargs...) where {N}
+  return op(opname, ntuple(n -> s[ns[n]], Val(N))...; kwargs...)
+end
+
+function op(opname, s::Vector{<:Index}, ns::Vararg{Integer}; kwargs...)
+  return op(opname, s, ns; kwargs...)
+end
+
+function op(s::Vector{<:Index}, opname, ns::Tuple{Vararg{Integer}}; kwargs...)
+  return op(opname, s, ns...; kwargs...)
+end
+
+function op(s::Vector{<:Index}, opname, ns::Integer...; kwargs...)
+  return op(opname, s, ns; kwargs...)
+end
+
+function op(s::Vector{<:Index}, opname, ns::Tuple{Vararg{Integer}}, kwargs::NamedTuple)
+  return op(opname, s, ns; kwargs...)
+end
+
+function op(s::Vector{<:Index}, opname, ns::Integer, kwargs::NamedTuple)
+  return op(opname, s, (ns,); kwargs...)
+end
+
+op(s::Vector{<:Index}, o::Tuple) = op(s, o...)
+
+op(o::Tuple, s::Vector{<:Index}) = op(s, o...)
+
+function op(
+  s::Vector{<:Index},
+  f::Function,
+  opname::AbstractString,
+  ns::Tuple{Vararg{Integer}};
+  kwargs...,
+)
+  return f(op(opname, s, ns...; kwargs...))
+end
+
+function op(
+  s::Vector{<:Index}, f::Function, opname::AbstractString, ns::Integer...; kwargs...
+)
+  return f(op(opname, s, ns; kwargs...))
+end
+
+# Here, Ref is used to not broadcast over the vector of indices
+# TODO: consider overloading broadcast for `op` with the example
+# here: https://discourse.julialang.org/t/how-to-broadcast-over-only-certain-function-arguments/19274/5
+# so that `Ref` isn't needed.
+ops(s::Vector{<:Index}, os::AbstractArray) = [op(oₙ, s) for oₙ in os]
+ops(os::AbstractVector, s::Vector{<:Index}) = [op(oₙ, s) for oₙ in os]
diff --git a/src/itensor/siteinds.jl b/src/itensor/siteinds.jl
new file mode 100644
index 0000000..be81f46
--- /dev/null
+++ b/src/itensor/siteinds.jl
@@ -0,0 +1,78 @@
+using ITensorBase: Index, addtags
+
+function siteind(st::SiteType; addtags="", kwargs...)
+  sp = space(st; kwargs...)
+  isnothing(sp) && return nothing
+  return Index(sp, "Site, $(value(st)), $addtags")
+end
+
+function siteind(st::SiteType, n; kwargs...)
+  s = siteind(st; kwargs...)
+  !isnothing(s) && return addtags(s, "n=$n")
+  sp = space(st, n; kwargs...)
+  isnothing(sp) && error(space_error_message(st))
+  return Index(sp, "Site, $(value(st)), n=$n")
+end
+
+siteind(tag::String; kwargs...) = siteind(SiteType(tag); kwargs...)
+
+siteind(tag::String, n; kwargs...) = siteind(SiteType(tag), n; kwargs...)
+
+# Special case of `siteind` where integer (dim) provided
+# instead of a tag string
+#siteind(d::Integer, n::Integer; kwargs...) = Index(d, "Site,n=$n")
+function siteind(d::Integer, n::Integer; addtags="", kwargs...)
+  return Index(d, "Site,n=$n, $addtags")
+end
+
+siteinds(::SiteType, N; kwargs...) = nothing
+
+"""
+    siteinds(tag::String, N::Integer; kwargs...)
+
+Create an array of `N` physical site indices of type `tag`.
+Keyword arguments can be used to specify quantum number conservation,
+see the `space` function corresponding to the site type `tag` for
+supported keyword arguments.
+
+# Example
+
+```julia
+N = 10
+s = siteinds("S=1/2", N; conserve_qns=true)
+```
+"""
+function siteinds(tag::String, N::Integer; kwargs...)
+  st = SiteType(tag)
+
+  si = siteinds(st, N; kwargs...)
+  if !isnothing(si)
+    return si
+  end
+
+  return [siteind(st, j; kwargs...) for j in 1:N]
+end
+
+"""
+    siteinds(f::Function, N::Integer; kwargs...)
+
+Create an array of `N` physical site indices where the site type at site `n` is given
+by `f(n)` (`f` should return a string).
+"""
+function siteinds(f::Function, N::Integer; kwargs...)
+  return [siteind(f(n), n; kwargs...) for n in 1:N]
+end
+
+# Special case of `siteinds` where integer (dim)
+# provided instead of a tag string
+"""
+    siteinds(d::Integer, N::Integer; kwargs...)
+
+Create an array of `N` site indices, each of dimension `d`.
+
+# Keywords
+- `addtags::String`: additional tags to be added to all indices
+"""
+function siteinds(d::Integer, N::Integer; kwargs...)
+  return [siteind(d, n; kwargs...) for n in 1:N]
+end
diff --git a/src/itensor/state.jl b/src/itensor/state.jl
new file mode 100644
index 0000000..52b355f
--- /dev/null
+++ b/src/itensor/state.jl
@@ -0,0 +1,24 @@
+using ITensorBase: ITensor, Index, onehot
+
+state(::StateName, ::SiteType, ::Index; kwargs...) = nothing
+
+# Syntax `state("Up", Index(2, "S=1/2"))`
+state(sn::String, i::Index; kwargs...) = state(i, sn; kwargs...)
+
+function state(s::Index, name::AbstractString; kwargs...)::ITensor
+  stypes = _sitetypes(s)
+  sname = StateName(name)
+  for st in stypes
+    v = state(sname, st; kwargs...)
+    !isnothing(v) && return ITensor(v, (s,))
+  end
+  return throw(
+    ArgumentError(
+      "Overload of \"state\" functions not found for state name \"$name\" and Index tags $(tags(s))",
+    ),
+  )
+end
+
+state(s::Index, n::Integer) = onehot(s => n)
+
+state(sset::Vector{<:Index}, j::Integer, st; kwargs...) = state(sset[j], st; kwargs...)
diff --git a/src/itensor/val.jl b/src/itensor/val.jl
new file mode 100644
index 0000000..2bc938c
--- /dev/null
+++ b/src/itensor/val.jl
@@ -0,0 +1,20 @@
+using ITensorBase: Index
+
+function val(s::Index, name::AbstractString)::Int
+  stypes = _sitetypes(s)
+  sname = ValName(name)
+
+  # Try calling val(::StateName"Name",::SiteType"Tag",)
+  for st in stypes
+    res = val(sname, st)
+    !isnothing(res) && return res
+  end
+
+  return throw(
+    ArgumentError("Overload of \"val\" function not found for Index tags $(tags(s))")
+  )
+end
+
+val(s::Index, n::Integer) = n
+
+val(sset::Vector{<:Index}, j::Integer, st) = val(sset[j], st)
diff --git a/src/op.jl b/src/op.jl
new file mode 100644
index 0000000..529be1b
--- /dev/null
+++ b/src/op.jl
@@ -0,0 +1,24 @@
+# TODO: Add `params<:NamedTuple` field.
+struct OpName{Name} end
+OpName(s::AbstractString) = OpName{Symbol(s)}()
+OpName(s::Symbol) = OpName{s}()
+name(::OpName{N}) where {N} = N
+macro OpName_str(s)
+  return OpName{Symbol(s)}
+end
+
+# Default implementations of op
+op(::OpName; kwargs...) = nothing
+op(::OpName, ::SiteType; kwargs...) = nothing
+
+function _sitetypes(ts::Set)
+  return collect(SiteType, SiteType.(ts))
+end
+
+op(name::AbstractString; kwargs...) = error("Must input indices when creating an `op`.")
+
+# To ease calling of other op overloads,
+# allow passing a string as the op name
+op(opname::AbstractString, t::SiteType; kwargs...) = op(OpName(opname), t; kwargs...)
+
+op(f::Function, args...; kwargs...) = f(op(args...; kwargs...))
diff --git a/src/sitetype.jl b/src/sitetype.jl
index 50aecfd..d6459cd 100644
--- a/src/sitetype.jl
+++ b/src/sitetype.jl
@@ -1,835 +1,13 @@
 using ChainRulesCore: @ignore_derivatives
-using ITensorBase: Index, ITensor, itensor, addtags, dag, onehot, prime, tags
 
 # TODO: Need to define or replace.
 # using ITensorBase: product, swapprime
 
-@eval struct SiteType{T}
-  (f::Type{<:SiteType})() = $(Expr(:new, :f))
-end
-
-# Note that the complicated definition of
-# SiteType above is a workaround for performance
-# issues when creating parameterized types
-# in Julia 1.4 and 1.5-beta. Ideally we
-# can just use the following in the future:
-# struct SiteType{T}
-# end
-
-"""
-SiteType is a parameterized type which allows
-making Index tags into Julia types. Use cases
-include overloading functions such as `op`,
-`siteinds`, and `state` which generate custom
-operators, Index arrays, and IndexVals associated
-with Index objects having a certain tag.
-
-To make a SiteType type, you can use the string
-macro notation: `SiteType"MyTag"`
-
-To make a SiteType value or object, you can use
-the notation: `SiteType("MyTag")`
-
-There are currently a few built-in site types
-recognized by `jl`. The system is easily extensible
-by users. To add new operators to an existing site type,
-or to create new site types, you can follow the instructions
-[here](https://itensor.github.io/jl/stable/examples/Physics.html).
-
-The current built-in site types are:
-
-- `SiteType"S=1/2"` (or `SiteType"S=½"`)
-- `SiteType"S=1"`
-- `SiteType"Qubit"`
-- `SiteType"Qudit"`
-- `SiteType"Boson"`
-- `SiteType"Fermion"`
-- `SiteType"tJ"`
-- `SiteType"Electron"`
-
-# Examples
-
-Tags on indices get turned into SiteTypes internally, and then
-we search for overloads of functions like `op` and `siteind`.
-For example:
-
-```julia
-julia> s = siteind("S=1/2")
-(dim=2|id=862|"S=1/2,Site")
-
-julia> @show op("Sz", s);
-op(s, "Sz") = ITensor ord=2
-Dim 1: (dim=2|id=862|"S=1/2,Site")'
-Dim 2: (dim=2|id=862|"S=1/2,Site")
-NDTensors.Dense{Float64,Array{Float64,1}}
- 2×2
- 0.5   0.0
- 0.0  -0.5
-
-julia> @show op("Sx", s);
-op(s, "Sx") = ITensor ord=2
-Dim 1: (dim=2|id=862|"S=1/2,Site")'
-Dim 2: (dim=2|id=862|"S=1/2,Site")
-NDTensors.Dense{Float64,Array{Float64,1}}
- 2×2
- 0.0  0.5
- 0.5  0.0
-
-julia> @show op("Sy", s);
-op(s, "Sy") = ITensor ord=2
-Dim 1: (dim=2|id=862|"S=1/2,Site")'
-Dim 2: (dim=2|id=862|"S=1/2,Site")
-NDTensors.Dense{Complex{Float64},Array{Complex{Float64},1}}
- 2×2
- 0.0 + 0.0im  -0.0 - 0.5im
- 0.0 + 0.5im   0.0 + 0.0im
-
-julia> s = siteind("Electron")
-(dim=4|id=734|"Electron,Site")
-
-julia> @show op("Nup", s);
-op(s, "Nup") = ITensor ord=2
-Dim 1: (dim=4|id=734|"Electron,Site")'
-Dim 2: (dim=4|id=734|"Electron,Site")
-NDTensors.Dense{Float64,Array{Float64,1}}
- 4×4
- 0.0  0.0  0.0  0.0
- 0.0  1.0  0.0  0.0
- 0.0  0.0  0.0  0.0
- 0.0  0.0  0.0  1.0
-```
-
-Many operators are available, for example:
-
-- `SiteType"S=1/2"`: `"Sz"`, `"Sx"`, `"Sy"`, `"S+"`, `"S-"`, ...
-- `SiteType"Electron"`: `"Nup"`, `"Ndn"`, `"Nupdn"`, `"Ntot"`, `"Cup"`,
-   `"Cdagup"`, `"Cdn"`, `"Cdagdn"`, `"Sz"`, `"Sx"`, `"Sy"`, `"S+"`, `"S-"`, ...
-- ...
-
-You can view the source code for the internal SiteType definitions
-and operators that are defined [here](https://github.com/ITensor/jl/tree/main/src/physics/site_types).
-"""
+# TODO: Add `params<:NamedTuple` field.
+struct SiteType{T} end
 SiteType(s::AbstractString) = SiteType{Symbol(s)}()
-
-SiteType(t::Integer) = SiteType{Symbol(t)}()
-
-tag(::SiteType{T}) where {T} = T
-
+SiteType(i::Integer) = SiteType{Symbol(i)}()
+value(::SiteType{T}) where {T} = T
 macro SiteType_str(s)
   return SiteType{Symbol(s)}
 end
-
-# Keep TagType defined for backwards
-# compatibility; will be deprecated later
-const TagType = SiteType
-macro TagType_str(s)
-  return TagType{Symbol(s)}
-end
-
-#---------------------------------------
-#
-# op system
-#
-#---------------------------------------
-
-@eval struct OpName{Name}
-  (f::Type{<:OpName})() = $(Expr(:new, :f))
-end
-
-# Note that the complicated definition of
-# OpName above is a workaround for performance
-# issues when creating parameterized types
-# in Julia 1.4 and 1.5-beta. Ideally we
-# can just use the following in the future:
-# struct OpName{Name}
-# end
-
-"""
-OpName is a parameterized type which allows
-making strings into Julia types for the purpose
-of representing operator names.
-The main use of OpName is overloading the
-`op!` method which generates operators
-for indices with certain tags such as "S=1/2".
-
-To make a OpName type, you can use the string
-macro notation: `OpName"MyTag"`.
-
-To make an OpName value or object, you can use
-the notation: `OpName("myop")`
-"""
-OpName(s::AbstractString) = OpName{Symbol(s)}()
-OpName(s::Symbol) = OpName{s}()
-name(::OpName{N}) where {N} = N
-
-macro OpName_str(s)
-  return OpName{Symbol(s)}
-end
-
-# Default implementations of op and op!
-op(::OpName; kwargs...) = nothing
-op(::OpName, ::SiteType; kwargs...) = nothing
-op(::OpName, ::SiteType, ::Index...; kwargs...) = nothing
-function op(
-  ::OpName, ::SiteType, ::SiteType, sitetypes_inds::Union{SiteType,Index}...; kwargs...
-)
-  return nothing
-end
-op!(::ITensor, ::OpName, ::SiteType, ::Index...; kwargs...) = nothing
-function op!(
-  ::ITensor,
-  ::OpName,
-  ::SiteType,
-  ::SiteType,
-  sitetypes_inds::Union{SiteType,Index}...;
-  kwargs...,
-)
-  return nothing
-end
-
-# Deprecated version, for backwards compatibility
-op(::SiteType, ::Index, ::AbstractString; kwargs...) = nothing
-
-function _sitetypes(ts::Set)
-  return collect(SiteType, SiteType.(ts))
-end
-
-_sitetypes(i::Index) = _sitetypes(tags(i))
-
-function commontags(i::Index...)
-  return union(tags.(i)...)
-end
-
-"""
-    op(opname::String, s::Index; kwargs...)
-
-Return an ITensor corresponding to the operator
-named `opname` for the Index `s`. The operator
-is constructed by calling an overload of either
-the `op` or `op!` methods which take a `SiteType`
-argument that corresponds to one of the tags of
-the Index `s` and an `OpName"opname"` argument
-that corresponds to the input operator name.
-
-Operator names can be combined using the `"*"`
-symbol, for example `"S+*S-"` or `"Sz*Sz*Sz"`.
-The result is an ITensor made by forming each operator
-then contracting them together in a way corresponding
-to the usual operator product or matrix multiplication.
-
-The `op` system is used by the OpSum
-system to convert operator names into ITensors,
-and can be used directly such as for applying
-operators to MPS.
-
-# Example
-
-```julia
-s = Index(2, "Site,S=1/2")
-Sz = op("Sz", s)
-```
-
-To see all of the operator names defined for the site types included with
-ITensor, please view the [source code](https://github.com/ITensor/jl/tree/main/src/physics/site_types)
-for each site type. Note that some site types such as "S=1/2" and "Qubit"
-are aliases for each other and share operator definitions.
-"""
-function op(name::AbstractString, s::Index...; adjoint::Bool=false, kwargs...)
-  name = strip(name)
-  # TODO: filter out only commons tags
-  # if there are multiple indices
-  commontags_s = commontags(s...)
-
-  # first we handle the + and - algebra, which requires a space between ops to avoid clashing
-  name_split = nothing
-  @ignore_derivatives name_split = String.(split(name, " "))
-  oplocs = findall(x -> x ∈ ("+", "-"), name_split)
-
-  if !isempty(oplocs)
-    @ignore_derivatives !isempty(kwargs) &&
-      error("Lazy algebra on parametric gates not allowed")
-
-    # the string representation of algebra ops: ex ["+", "-", "+"]
-    labels = name_split[oplocs]
-    # assign coefficients to each term: ex [+1, -1, +1]
-    coeffs = [1, [(-1)^Int(label == "-") for label in labels]...]
-
-    # grad the name of each operator block separated by an algebra op, and do so by
-    # making sure blank spaces between opnames are kept when building the new block.
-    start, opnames = 0, String[]
-    for oploc in oplocs
-      finish = oploc
-      opnames = vcat(
-        opnames, [prod([name_split[k] * " " for k in (start + 1):(finish - 1)])]
-      )
-      start = oploc
-    end
-    opnames = vcat(
-      opnames, [prod([name_split[k] * " " for k in (start + 1):length(name_split)])]
-    )
-
-    # build the vector of blocks and sum
-    op_list = [
-      coeff * (op(opname, s...; kwargs...)) for (coeff, opname) in zip(coeffs, opnames)
-    ]
-    return sum(op_list)
-  end
-
-  # the the multiplication come after
-  oploc = findfirst("*", name)
-  if !isnothing(oploc)
-    op1, op2 = nothing, nothing
-    @ignore_derivatives begin
-      op1 = name[1:prevind(name, oploc.start)]
-      op2 = name[nextind(name, oploc.start):end]
-      if !(op1[end] == ' ' && op2[1] == ' ')
-        @warn "($op1*$op2) composite op definition `A*B` deprecated: please use `A * B` instead (with spaces)"
-      end
-    end
-    return product(op(op1, s...; kwargs...), op(op2, s...; kwargs...))
-  end
-
-  common_stypes = _sitetypes(commontags_s)
-  @ignore_derivatives push!(common_stypes, SiteType("Generic"))
-  opn = OpName(name)
-
-  #
-  # Try calling a function of the form:
-  #    op(::OpName, ::SiteType, ::Index...; kwargs...)
-  #
-  for st in common_stypes
-    res = op(opn, st, s...; kwargs...)
-    if !isnothing(res)
-      adjoint && return swapprime(dag(res), 0 => 1)
-      return res
-    end
-  end
-
-  #
-  # Try calling a function of the form:
-  #    op(::OpName; kwargs...)
-  #  for backward compatibility with previous
-  #  gate system in PastaQ.jl
-  #
-  op_mat = op(opn; kwargs...)
-  if !isnothing(op_mat)
-    rs = reverse(s)
-    res = itensor(op_mat, prime.(rs)..., dag.(rs)...)
-    adjoint && return swapprime(dag(res), 0 => 1)
-    return res
-  end
-  #
-  # otherwise try calling a function of the form:
-  #    op(::OpName, ::SiteType; kwargs...)
-  # which returns a Julia matrix
-  #
-  for st in common_stypes
-    op_mat = op(opn, st; kwargs...)
-    if !isnothing(op_mat)
-      rs = reverse(s)
-      #return itensor(op_mat, prime.(rs)..., dag.(rs)...)
-      res = itensor(op_mat, prime.(rs)..., dag.(rs)...)
-      adjoint && return swapprime(dag(res), 0 => 1)
-      return res
-    end
-  end
-
-  # otherwise try calling a function of the form:
-  #    op!(::ITensor, ::OpName, ::SiteType, ::Index...; kwargs...)
-  #
-  ## Op = ITensor(prime.(s)..., dag.(s)...)
-  ## for st in common_stypes
-  ##   op!(Op, opn, st, s...; kwargs...)
-  ##   if !isempty(Op)
-  ##     adjoint && return swapprime(dag(Op), 0 => 1)
-  ##     return Op
-  ##   end
-  ## end
-
-  if length(s) > 1
-    # No overloads for common tags found. It might be a
-    # case of making an operator with mixed site types,
-    # searching for overloads like:
-    #   op(::OpName,
-    #      ::SiteType...,
-    #      ::Index...;
-    #      kwargs...)
-    #   op!(::ITensor, ::OpName,
-    #       ::SiteType...,
-    #       ::Index...;
-    #       kwargs...)
-    stypes = _sitetypes.(s)
-
-    for st in Iterators.product(stypes...)
-      res = op(opn, st..., s...; kwargs...)
-      if !isnothing(res)
-        adjoint && return swapprime(dag(res), 0 => 1)
-        return res
-      end
-    end
-
-    Op = ITensor(prime.(s)..., dag.(s)...)
-    for st in Iterators.product(stypes...)
-      op!(Op, opn, st..., s...; kwargs...)
-      if !isempty(Op)
-        adjoint && return swapprime(dag(Op), 0 => 1)
-        return Op
-      end
-    end
-
-    throw(
-      ArgumentError(
-        "Overload of \"op\" or \"op!\" functions not found for operator name \"$name\" and Index tags: $(tags.(s)).",
-      ),
-    )
-  end
-
-  #
-  # otherwise try calling a function of the form:
-  #   op(::SiteType, ::Index, ::AbstractString)
-  #
-  # (Note: this version is for backwards compatibility
-  #  after version 0.1.10, and may be eventually
-  #  deprecated)
-  #
-  for st in common_stypes
-    res = op(st, s[1], name; kwargs...)
-    if !isnothing(res)
-      adjoint && return dag(res)
-      return res
-    end
-  end
-
-  return throw(
-    ArgumentError(
-      "Overload of \"op\" or \"op!\" functions not found for operator name \"$name\" and Index tags: $(tags.(s)).",
-    ),
-  )
-end
-
-op(name::AbstractString; kwargs...) = error("Must input indices when creating an `op`.")
-
-"""
-    op(X::AbstractArray, s::Index...)
-    op(M::Matrix, s::Index...)
-
-Given a matrix M and a set of indices s,t,...
-return an operator ITensor with matrix elements given
-by M and indices s, s', t, t'
-
-# Example
-
-```julia
-julia> s = siteind("S=1/2")
-(dim=2|id=575|"S=1/2,Site")
-
-julia> Sz = op([1/2 0; 0 -1/2],s)
-ITensor ord=2 (dim=2|id=575|"S=1/2,Site")' (dim=2|id=575|"S=1/2,Site")
-NDTensors.Dense{Float64, Vector{Float64}}
-
-julia> @show Sz
-Sz = ITensor ord=2
-Dim 1: (dim=2|id=575|"S=1/2,Site")'
-Dim 2: (dim=2|id=575|"S=1/2,Site")
-NDTensors.Dense{Float64, Vector{Float64}}
- 2×2
- 0.5   0.0
- 0.0  -0.5
-ITensor ord=2 (dim=2|id=575|"S=1/2,Site")' (dim=2|id=575|"S=1/2,Site")
-NDTensors.Dense{Float64, Vector{Float64}}
-```
-"""
-op(X::AbstractArray, s::Index...) = itensor(X, prime.([s...]), dag.([s...]))
-
-op(opname, s::Vector{<:Index}; kwargs...) = op(opname, s...; kwargs...)
-
-op(s::Vector{<:Index}, opname; kwargs...) = op(opname, s...; kwargs...)
-
-# For backwards compatibility, version of `op`
-# taking the arguments in the other order:
-op(s::Index, opname; kwargs...) = op(opname, s; kwargs...)
-
-# To ease calling of other op overloads,
-# allow passing a string as the op name
-op(opname::AbstractString, t::SiteType; kwargs...) = op(OpName(opname), t; kwargs...)
-
-"""
-    op(opname::String,sites::Vector{<:Index},n::Int; kwargs...)
-
-Return an ITensor corresponding to the operator
-named `opname` for the n'th Index in the array
-`sites`.
-
-# Example
-
-```julia
-s = siteinds("S=1/2", 4)
-Sz2 = op("Sz", s, 2)
-```
-"""
-function op(opname, s::Vector{<:Index}, ns::NTuple{N,Integer}; kwargs...) where {N}
-  return op(opname, ntuple(n -> s[ns[n]], Val(N))...; kwargs...)
-end
-
-function op(opname, s::Vector{<:Index}, ns::Vararg{Integer}; kwargs...)
-  return op(opname, s, ns; kwargs...)
-end
-
-function op(s::Vector{<:Index}, opname, ns::Tuple{Vararg{Integer}}; kwargs...)
-  return op(opname, s, ns...; kwargs...)
-end
-
-function op(s::Vector{<:Index}, opname, ns::Integer...; kwargs...)
-  return op(opname, s, ns; kwargs...)
-end
-
-function op(s::Vector{<:Index}, opname, ns::Tuple{Vararg{Integer}}, kwargs::NamedTuple)
-  return op(opname, s, ns; kwargs...)
-end
-
-function op(s::Vector{<:Index}, opname, ns::Integer, kwargs::NamedTuple)
-  return op(opname, s, (ns,); kwargs...)
-end
-
-op(s::Vector{<:Index}, o::Tuple) = op(s, o...)
-
-op(o::Tuple, s::Vector{<:Index}) = op(s, o...)
-
-op(f::Function, args...; kwargs...) = f(op(args...; kwargs...))
-
-function op(
-  s::Vector{<:Index},
-  f::Function,
-  opname::AbstractString,
-  ns::Tuple{Vararg{Integer}};
-  kwargs...,
-)
-  return f(op(opname, s, ns...; kwargs...))
-end
-
-function op(
-  s::Vector{<:Index}, f::Function, opname::AbstractString, ns::Integer...; kwargs...
-)
-  return f(op(opname, s, ns; kwargs...))
-end
-
-# Here, Ref is used to not broadcast over the vector of indices
-# TODO: consider overloading broadcast for `op` with the example
-# here: https://discourse.julialang.org/t/how-to-broadcast-over-only-certain-function-arguments/19274/5
-# so that `Ref` isn't needed.
-ops(s::Vector{<:Index}, os::AbstractArray) = [op(oₙ, s) for oₙ in os]
-ops(os::AbstractVector, s::Vector{<:Index}) = [op(oₙ, s) for oₙ in os]
-
-@doc """
-    ops(s::Vector{<:Index}, os::Vector)
-    ops(os::Vector, s::Vector{<:Index})
-
-Given a list of operators, create ITensors using the collection
-of indices.
-
-# Examples
-
-```julia
-s = siteinds("Qubit", 4)
-os = [("H", 1), ("X", 2), ("CX", 2, 4)]
-# gates = ops(s, os)
-gates = ops(os, s)
-```
-""" ops(::Vector{<:Index}, ::AbstractArray)
-
-#---------------------------------------
-#
-# state system
-#
-#---------------------------------------
-
-@eval struct StateName{Name}
-  (f::Type{<:StateName})() = $(Expr(:new, :f))
-end
-
-StateName(s::AbstractString) = StateName{Symbol(s)}()
-StateName(s::Symbol) = StateName{s}()
-name(::StateName{N}) where {N} = N
-
-macro StateName_str(s)
-  return StateName{Symbol(s)}
-end
-
-state(::StateName, ::SiteType; kwargs...) = nothing
-state(::StateName, ::SiteType, ::Index; kwargs...) = nothing
-state!(::ITensor, ::StateName, ::SiteType, ::Index; kwargs...) = nothing
-
-# Syntax `state("Up", Index(2, "S=1/2"))`
-state(sn::String, i::Index; kwargs...) = state(i, sn; kwargs...)
-
-"""
-    state(s::Index, name::String; kwargs...)
-
-Return an ITensor corresponding to the state
-named `name` for the Index `s`. The returned
-ITensor will have `s` as its only index.
-
-The terminology here is based on the idea of a
-single-site state or wavefunction in physics.
-
-The `state` function is implemented for various
-Index tags by overloading either the
-`state` or `state!` methods which take a `SiteType`
-argument corresponding to one of the tags of
-the Index `s` and an `StateName"name"` argument
-that corresponds to the input state name.
-
-The `state` system is used by the MPS type
-to construct product-state MPS and for other purposes.
-
-# Example
-
-```julia
-s = Index(2, "Site,S=1/2")
-sup = state(s,"Up")
-sdn = state(s,"Dn")
-sxp = state(s,"X+")
-sxm = state(s,"X-")
-```
-"""
-function state(s::Index, name::AbstractString; kwargs...)::ITensor
-  stypes = _sitetypes(s)
-  sname = StateName(name)
-
-  # Try calling state(::StateName"Name",::SiteType"Tag",s::Index; kwargs...)
-  for st in stypes
-    v = state(sname, st, s; kwargs...)
-    if !isnothing(v)
-      if v isa ITensor
-        return v
-      else
-        # TODO: deprecate, only for backwards compatibility.
-        return itensor(v, s)
-      end
-    end
-  end
-
-  ## TODO: Bring this back, it is broken right now.
-  ## `isempty` isn't the correct function, it should use `isallocated`.
-  ## # Try calling state!(::ITensor,::StateName"Name",::SiteType"Tag",s::Index;kwargs...)
-  ## T = ITensor(s)
-  ## for st in stypes
-  ##   state!(T, sname, st, s; kwargs...)
-  ##   !isempty(T) && return T
-  ## end
-
-  #
-  # otherwise try calling a function of the form:
-  #    state(::StateName"Name", ::SiteType"Tag"; kwargs...)
-  # which returns a Julia vector
-  #
-  for st in stypes
-    v = state(sname, st; kwargs...)
-    !isnothing(v) && return itensor(v, s)
-  end
-
-  return throw(
-    ArgumentError(
-      "Overload of \"state\" or \"state!\" functions not found for state name \"$name\" and Index tags $(tags(s))",
-    ),
-  )
-end
-
-state(s::Index, n::Integer) = onehot(s => n)
-
-state(sset::Vector{<:Index}, j::Integer, st; kwargs...) = state(sset[j], st; kwargs...)
-
-#---------------------------------------
-#
-# val system
-#
-#---------------------------------------
-
-@eval struct ValName{Name}
-  (f::Type{<:ValName})() = $(Expr(:new, :f))
-end
-
-ValName(s::AbstractString) = ValName{Symbol(s)}()
-ValName(s::Symbol) = ValName{s}()
-name(::ValName{N}) where {N} = N
-
-macro ValName_str(s)
-  return ValName{Symbol(s)}
-end
-
-val(::ValName, ::SiteType) = nothing
-val(::AbstractString, ::SiteType) = nothing
-
-"""
-    val(s::Index, name::String)
-
-Return an integer corresponding to the `name`
-of a certain value the Index `s` can take.
-In other words, the `val` function maps strings
-to specific integer values within the range `1:dim(s)`.
-
-The `val` function is implemented for various
-Index tags by overloading methods named `val`
-which take a `SiteType` argument corresponding to
-one of the tags of the Index `s` and an `ValName"name"`
-argument that corresponds to the input name.
-
-# Example
-
-```julia
-s = Index(2, "Site,S=1/2")
-val(s,"Up") == 1
-val(s,"Dn") == 2
-
-s = Index(2, "Site,Fermion")
-val(s,"Emp") == 1
-val(s,"Occ") == 2
-```
-"""
-function val(s::Index, name::AbstractString)::Int
-  stypes = _sitetypes(s)
-  sname = ValName(name)
-
-  # Try calling val(::StateName"Name",::SiteType"Tag",)
-  for st in stypes
-    res = val(sname, st)
-    !isnothing(res) && return res
-  end
-
-  return throw(
-    ArgumentError("Overload of \"val\" function not found for Index tags $(tags(s))")
-  )
-end
-
-val(s::Index, n::Integer) = n
-
-val(sset::Vector{<:Index}, j::Integer, st) = val(sset[j], st)
-
-#---------------------------------------
-#
-# siteind system
-#
-#---------------------------------------
-
-space(st::SiteType; kwargs...) = nothing
-
-space(st::SiteType, n::Int; kwargs...) = space(st; kwargs...)
-
-function space_error_message(st::SiteType)
-  return "Overload of \"space\",\"siteind\", or \"siteinds\" functions not found for Index tag: $(tag(st))"
-end
-
-function siteind(st::SiteType; addtags="", kwargs...)
-  sp = space(st; kwargs...)
-  isnothing(sp) && return nothing
-  return Index(sp, "Site, $(tag(st)), $addtags")
-end
-
-function siteind(st::SiteType, n; kwargs...)
-  s = siteind(st; kwargs...)
-  !isnothing(s) && return addtags(s, "n=$n")
-  sp = space(st, n; kwargs...)
-  isnothing(sp) && error(space_error_message(st))
-  return Index(sp, "Site, $(tag(st)), n=$n")
-end
-
-siteind(tag::String; kwargs...) = siteind(SiteType(tag); kwargs...)
-
-siteind(tag::String, n; kwargs...) = siteind(SiteType(tag), n; kwargs...)
-
-# Special case of `siteind` where integer (dim) provided
-# instead of a tag string
-#siteind(d::Integer, n::Integer; kwargs...) = Index(d, "Site,n=$n")
-function siteind(d::Integer, n::Integer; addtags="", kwargs...)
-  return Index(d, "Site,n=$n, $addtags")
-end
-
-#---------------------------------------
-#
-# siteinds system
-#
-#---------------------------------------
-
-siteinds(::SiteType, N; kwargs...) = nothing
-
-"""
-    siteinds(tag::String, N::Integer; kwargs...)
-
-Create an array of `N` physical site indices of type `tag`.
-Keyword arguments can be used to specify quantum number conservation,
-see the `space` function corresponding to the site type `tag` for
-supported keyword arguments.
-
-# Example
-
-```julia
-N = 10
-s = siteinds("S=1/2", N; conserve_qns=true)
-```
-"""
-function siteinds(tag::String, N::Integer; kwargs...)
-  st = SiteType(tag)
-
-  si = siteinds(st, N; kwargs...)
-  if !isnothing(si)
-    return si
-  end
-
-  return [siteind(st, j; kwargs...) for j in 1:N]
-end
-
-"""
-    siteinds(f::Function, N::Integer; kwargs...)
-
-Create an array of `N` physical site indices where the site type at site `n` is given
-by `f(n)` (`f` should return a string).
-"""
-function siteinds(f::Function, N::Integer; kwargs...)
-  return [siteind(f(n), n; kwargs...) for n in 1:N]
-end
-
-# Special case of `siteinds` where integer (dim)
-# provided instead of a tag string
-"""
-    siteinds(d::Integer, N::Integer; kwargs...)
-
-Create an array of `N` site indices, each of dimension `d`.
-
-# Keywords
-- `addtags::String`: additional tags to be added to all indices
-"""
-function siteinds(d::Integer, N::Integer; kwargs...)
-  return [siteind(d, n; kwargs...) for n in 1:N]
-end
-
-#---------------------------------------
-#
-# has_fermion_string system
-#
-#---------------------------------------
-
-has_fermion_string(operator::AbstractArray{<:Number}, s::Index; kwargs...)::Bool = false
-
-has_fermion_string(::OpName, ::SiteType) = nothing
-
-function has_fermion_string(opname::AbstractString, s::Index; kwargs...)::Bool
-  opname = strip(opname)
-
-  # Interpret operator names joined by *
-  # as acting sequentially on the same site
-  starpos = findfirst(isequal('*'), opname)
-  if !isnothing(starpos)
-    op1 = opname[1:prevind(opname, starpos)]
-    op2 = opname[nextind(opname, starpos):end]
-    return xor(has_fermion_string(op1, s; kwargs...), has_fermion_string(op2, s; kwargs...))
-  end
-
-  Ntags = length(tags(s))
-  stypes = _sitetypes(s)
-  opn = OpName(opname)
-  for st in stypes
-    res = has_fermion_string(opn, st)
-    !isnothing(res) && return res
-  end
-  return false
-end
diff --git a/src/sitetypes/boson.jl b/src/sitetypes/boson.jl
index 387edad..6a011a4 100644
--- a/src/sitetypes/boson.jl
+++ b/src/sitetypes/boson.jl
@@ -1,4 +1,3 @@
-
 alias(::SiteType"Boson") = SiteType"Qudit"()
 
 """
@@ -16,17 +15,17 @@ space(st::SiteType"Boson"; kwargs...) = space(alias(st); kwargs...)
 
 val(vn::ValName, st::SiteType"Boson") = val(vn, alias(st))
 
-function state(sn::StateName, st::SiteType"Boson", s::Index; kwargs...)
-  return state(sn, alias(st), s; kwargs...)
-end
+## function state(sn::StateName, st::SiteType"Boson", s::Index; kwargs...)
+##   return state(sn, alias(st), s; kwargs...)
+## end
 
 function op(on::OpName, st::SiteType"Boson", ds::Int...; kwargs...)
   return op(on, alias(st), ds...; kwargs...)
 end
 
-function op(on::OpName, st::SiteType"Boson", s1::Index, s_tail::Index...; kwargs...)
-  rs = reverse((s1, s_tail...))
-  ds = dim.(rs)
-  opmat = op(on, st, ds...; kwargs...)
-  return itensor(opmat, prime.(rs)..., dag.(rs)...)
-end
+## function op(on::OpName, st::SiteType"Boson", s1::Index, s_tail::Index...; kwargs...)
+##   rs = reverse((s1, s_tail...))
+##   ds = dim.(rs)
+##   opmat = op(on, st, ds...; kwargs...)
+##   return itensor(opmat, prime.(rs)..., dag.(rs)...)
+## end
diff --git a/src/sitetypes/fermion.jl b/src/sitetypes/fermion.jl
index 008bdf0..f9735c2 100644
--- a/src/sitetypes/fermion.jl
+++ b/src/sitetypes/fermion.jl
@@ -1,4 +1,3 @@
-
 """
     space(::SiteType"Fermion";
           conserve_qns=false,
@@ -69,34 +68,34 @@ state(::StateName"Occ", ::SiteType"Fermion") = [0.0 1.0]
 state(::StateName"0", st::SiteType"Fermion") = state(StateName("Emp"), st)
 state(::StateName"1", st::SiteType"Fermion") = state(StateName("Occ"), st)
 
-function op!(Op::ITensor, ::OpName"N", ::SiteType"Fermion", s::Index)
-  return Op[s' => 2, s => 2] = 1.0
-end
-function op!(Op::ITensor, on::OpName"n", st::SiteType"Fermion", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"C", ::SiteType"Fermion", s::Index)
-  return Op[s' => 1, s => 2] = 1.0
-end
-function op!(Op::ITensor, on::OpName"c", st::SiteType"Fermion", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Cdag", ::SiteType"Fermion", s::Index)
-  return Op[s' => 2, s => 1] = 1.0
-end
-function op!(Op::ITensor, on::OpName"c†", st::SiteType"Fermion", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-function op!(Op::ITensor, on::OpName"cdag", st::SiteType"Fermion", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"F", ::SiteType"Fermion", s::Index)
-  Op[s' => 1, s => 1] = +1.0
-  return Op[s' => 2, s => 2] = -1.0
-end
+## function op!(Op::ITensor, ::OpName"N", ::SiteType"Fermion", s::Index)
+##   return Op[s' => 2, s => 2] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"n", st::SiteType"Fermion", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"C", ::SiteType"Fermion", s::Index)
+##   return Op[s' => 1, s => 2] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"c", st::SiteType"Fermion", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Cdag", ::SiteType"Fermion", s::Index)
+##   return Op[s' => 2, s => 1] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"c†", st::SiteType"Fermion", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## function op!(Op::ITensor, on::OpName"cdag", st::SiteType"Fermion", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"F", ::SiteType"Fermion", s::Index)
+##   Op[s' => 1, s => 1] = +1.0
+##   return Op[s' => 2, s => 2] = -1.0
+## end
 
 has_fermion_string(::OpName"C", ::SiteType"Fermion") = true
 function has_fermion_string(on::OpName"c", st::SiteType"Fermion")
diff --git a/src/sitetypes/generic_sites.jl b/src/sitetypes/generic_sites.jl
index b91bba9..9f06047 100644
--- a/src/sitetypes/generic_sites.jl
+++ b/src/sitetypes/generic_sites.jl
@@ -1,51 +1,51 @@
-using ITensorBase: ITensor, itensor
-using LinearAlgebra: I
-
-# TODO: Need to define or replace.
-# using ITensorBase: settensor!
-
-function op!(
-  o::ITensor, ::OpName"Id", ::SiteType"Generic", s1::Index, sn::Index...; eltype=Float64
-)
-  s = (s1, sn...)
-  n = prod(dim.(s))
-  t = itensor(Matrix(one(eltype) * I, n, n), prime.(s)..., dag.(s)...)
-  return settensor!(o, tensor(t))
-end
-
-function op!(o::ITensor, on::OpName"I", st::SiteType"Generic", s::Index...; kwargs...)
-  return op!(o, alias(on), st, s...; kwargs...)
-end
-
-function op!(o::ITensor, ::OpName"F", st::SiteType"Generic", s::Index; kwargs...)
-  return op!(o, OpName("Id"), st, s; kwargs...)
-end
-
-function default_random_matrix(eltype::Type, s::Index...)
-  n = prod(dim.(s))
-  return randn(eltype, n, n)
-end
-
-# Haar-random unitary
-#
-# Reference:
-# Section 4.6
-# http://math.mit.edu/~edelman/publications/random_matrix_theory.pdf
-function op!(
-  o::ITensor,
-  ::OpName"RandomUnitary",
-  ::SiteType"Generic",
-  s1::Index,
-  sn::Index...;
-  eltype=ComplexF64,
-  random_matrix=default_random_matrix(eltype, s1, sn...),
-)
-  s = (s1, sn...)
-  Q, _ = NDTensors.qr_positive(random_matrix)
-  t = itensor(Q, prime.(s)..., dag.(s)...)
-  return settensor!(o, tensor(t))
-end
-
-function op!(o::ITensor, ::OpName"randU", st::SiteType"Generic", s::Index...; kwargs...)
-  return op!(o, OpName("RandomUnitary"), st, s...; kwargs...)
-end
+## using ITensorBase: ITensor, itensor
+## using LinearAlgebra: I
+## 
+## # TODO: Need to define or replace.
+## # using ITensorBase: settensor!
+## 
+## function op!(
+##   o::ITensor, ::OpName"Id", ::SiteType"Generic", s1::Index, sn::Index...; eltype=Float64
+## )
+##   s = (s1, sn...)
+##   n = prod(dim.(s))
+##   t = itensor(Matrix(one(eltype) * I, n, n), prime.(s)..., dag.(s)...)
+##   return settensor!(o, tensor(t))
+## end
+## 
+## function op!(o::ITensor, on::OpName"I", st::SiteType"Generic", s::Index...; kwargs...)
+##   return op!(o, alias(on), st, s...; kwargs...)
+## end
+## 
+## function op!(o::ITensor, ::OpName"F", st::SiteType"Generic", s::Index; kwargs...)
+##   return op!(o, OpName("Id"), st, s; kwargs...)
+## end
+## 
+## function default_random_matrix(eltype::Type, s::Index...)
+##   n = prod(dim.(s))
+##   return randn(eltype, n, n)
+## end
+## 
+## # Haar-random unitary
+## #
+## # Reference:
+## # Section 4.6
+## # http://math.mit.edu/~edelman/publications/random_matrix_theory.pdf
+## function op!(
+##   o::ITensor,
+##   ::OpName"RandomUnitary",
+##   ::SiteType"Generic",
+##   s1::Index,
+##   sn::Index...;
+##   eltype=ComplexF64,
+##   random_matrix=default_random_matrix(eltype, s1, sn...),
+## )
+##   s = (s1, sn...)
+##   Q, _ = NDTensors.qr_positive(random_matrix)
+##   t = itensor(Q, prime.(s)..., dag.(s)...)
+##   return settensor!(o, tensor(t))
+## end
+## 
+## function op!(o::ITensor, ::OpName"randU", st::SiteType"Generic", s::Index...; kwargs...)
+##   return op!(o, OpName("RandomUnitary"), st, s...; kwargs...)
+## end
diff --git a/src/sitetypes/qudit.jl b/src/sitetypes/qudit.jl
index 6644cf0..93b3330 100644
--- a/src/sitetypes/qudit.jl
+++ b/src/sitetypes/qudit.jl
@@ -28,12 +28,12 @@ function val(::ValName{N}, ::SiteType"Qudit") where {N}
   return parse(Int, String(N)) + 1
 end
 
-function state(::StateName{N}, ::SiteType"Qudit", s::Index) where {N}
-  n = parse(Int, String(N))
-  st = zeros(dim(s))
-  st[n + 1] = 1.0
-  return itensor(st, s)
-end
+## function state(::StateName{N}, ::SiteType"Qudit", s::Index) where {N}
+##   n = parse(Int, String(N))
+##   st = zeros(dim(s))
+##   st[n + 1] = 1.0
+##   return itensor(st, s)
+## end
 
 # one-body operators
 function op(::OpName"Id", ::SiteType"Qudit", ds::Int...)
@@ -88,13 +88,13 @@ function op(::OpName"a†b†", st::SiteType"Qudit", d1::Int, d2::Int)
   return kron(op(OpName("a†"), st, d1), op(OpName("a†"), st, d2))
 end
 
-# interface
-function op(on::OpName, st::SiteType"Qudit", s1::Index, s_tail::Index...; kwargs...)
-  rs = reverse((s1, s_tail...))
-  ds = dim.(rs)
-  opmat = op(on, st, ds...; kwargs...)
-  return itensor(opmat, prime.(rs)..., dag.(rs)...)
-end
+## # interface
+## function op(on::OpName, st::SiteType"Qudit", s1::Index, s_tail::Index...; kwargs...)
+##   rs = reverse((s1, s_tail...))
+##   ds = dim.(rs)
+##   opmat = op(on, st, ds...; kwargs...)
+##   return itensor(opmat, prime.(rs)..., dag.(rs)...)
+## end
 
 function op(on::OpName, st::SiteType"Qudit"; kwargs...)
   return error("`op` can't be called without indices or dimensions.")
diff --git a/src/sitetypes/spinone.jl b/src/sitetypes/spinone.jl
index 269fda7..91766f1 100644
--- a/src/sitetypes/spinone.jl
+++ b/src/sitetypes/spinone.jl
@@ -143,8 +143,8 @@ space(st::SiteType"SpinOne"; kwargs...) = space(alias(st); kwargs...)
 state(name::StateName, st::SiteType"SpinOne") = state(name, alias(st))
 val(name::ValName, st::SiteType"SpinOne") = val(name, alias(st))
 
-function op!(Op::ITensor, o::OpName, st::SiteType"SpinOne", s::Index)
-  return op!(Op, o, alias(st), s)
-end
+## function op!(Op::ITensor, o::OpName, st::SiteType"SpinOne", s::Index)
+##   return op!(Op, o, alias(st), s)
+## end
 
 op(o::OpName, st::SiteType"SpinOne") = op(o, alias(st))
diff --git a/src/sitetypes/tj.jl b/src/sitetypes/tj.jl
index 063f46d..e603f19 100644
--- a/src/sitetypes/tj.jl
+++ b/src/sitetypes/tj.jl
@@ -1,4 +1,3 @@
-
 """
     space(::SiteType"tJ";
           conserve_qns = false,
@@ -68,153 +67,153 @@ state(::StateName"0", st::SiteType"tJ") = state(StateName("Emp"), st)
 state(::StateName"↑", st::SiteType"tJ") = state(StateName("Up"), st)
 state(::StateName"↓", st::SiteType"tJ") = state(StateName("Dn"), st)
 
-function op!(Op::ITensor, ::OpName"Nup", ::SiteType"tJ", s::Index)
-  return Op[s' => 2, s => 2] = 1.0
-end
-function op!(Op::ITensor, on::OpName"n↑", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Ndn", ::SiteType"tJ", s::Index)
-  return Op[s' => 3, s => 3] = 1.0
-end
-function op!(Op::ITensor, on::OpName"n↓", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Ntot", ::SiteType"tJ", s::Index)
-  Op[s' => 2, s => 2] = 1.0
-  return Op[s' => 3, s => 3] = 1.0
-end
-function op!(Op::ITensor, on::OpName"ntot", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Cup", ::SiteType"tJ", s::Index)
-  return Op[s' => 1, s => 2] = 1.0
-end
-function op!(Op::ITensor, on::OpName"c↑", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Cdagup", ::SiteType"tJ", s::Index)
-  return Op[s' => 2, s => 1] = 1.0
-end
-function op!(Op::ITensor, on::OpName"c†↑", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Cdn", ::SiteType"tJ", s::Index)
-  return Op[s' => 1, s => 3] = 1.0
-end
-function op!(Op::ITensor, on::OpName"c↓", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Cdagdn", ::SiteType"tJ", s::Index)
-  return Op[s' => 3, s => 1] = 1.0
-end
-function op!(Op::ITensor, on::OpName"c†↓", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Aup", ::SiteType"tJ", s::Index)
-  return Op[s' => 1, s => 2] = 1.0
-end
-function op!(Op::ITensor, on::OpName"a↑", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Adagup", ::SiteType"tJ", s::Index)
-  return Op[s' => 2, s => 1] = 1.0
-end
-function op!(Op::ITensor, on::OpName"a†↑", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Adn", ::SiteType"tJ", s::Index)
-  return Op[s' => 1, s => 3] = 1.0
-end
-function op!(Op::ITensor, on::OpName"a↓", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Adagdn", ::SiteType"tJ", s::Index)
-  return Op[s' => 3, s => 1] = 1.0
-end
-function op!(Op::ITensor, on::OpName"a†↓", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"F", ::SiteType"tJ", s::Index)
-  Op[s' => 1, s => 1] = +1.0
-  Op[s' => 2, s => 2] = -1.0
-  return Op[s' => 3, s => 3] = -1.0
-end
-
-function op!(Op::ITensor, ::OpName"Fup", ::SiteType"tJ", s::Index)
-  Op[s' => 1, s => 1] = +1.0
-  Op[s' => 2, s => 2] = -1.0
-  return Op[s' => 3, s => 3] = +1.0
-end
-function op!(Op::ITensor, on::OpName"F↑", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Fdn", ::SiteType"tJ", s::Index)
-  Op[s' => 1, s => 1] = +1.0
-  Op[s' => 2, s => 2] = +1.0
-  return Op[s' => 3, s => 3] = -1.0
-end
-function op!(Op::ITensor, on::OpName"F↓", st::SiteType"tJ", s::Index)
-  return op!(Op, alias(on), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Sz", ::SiteType"tJ", s::Index)
-  Op[s' => 2, s => 2] = +0.5
-  return Op[s' => 3, s => 3] = -0.5
-end
-
-function op!(Op::ITensor, ::OpName"Sᶻ", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("Sz"), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"Sx", ::SiteType"tJ", s::Index)
-  Op[s' => 2, s => 3] = 0.5
-  return Op[s' => 3, s => 2] = 0.5
-end
-
-function op!(Op::ITensor, ::OpName"Sˣ", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("Sx"), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"S+", ::SiteType"tJ", s::Index)
-  return Op[s' => 2, s => 3] = 1.0
-end
-
-function op!(Op::ITensor, ::OpName"S⁺", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("S+"), st, s)
-end
-function op!(Op::ITensor, ::OpName"Sp", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("S+"), st, s)
-end
-function op!(Op::ITensor, ::OpName"Splus", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("S+"), st, s)
-end
-
-function op!(Op::ITensor, ::OpName"S-", ::SiteType"tJ", s::Index)
-  return Op[s' => 3, s => 2] = 1.0
-end
-
-function op!(Op::ITensor, ::OpName"S⁻", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("S-"), st, s)
-end
-function op!(Op::ITensor, ::OpName"Sm", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("S-"), st, s)
-end
-function op!(Op::ITensor, ::OpName"Sminus", st::SiteType"tJ", s::Index)
-  return op!(Op, OpName("S-"), st, s)
-end
+## function op!(Op::ITensor, ::OpName"Nup", ::SiteType"tJ", s::Index)
+##   return Op[s' => 2, s => 2] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"n↑", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Ndn", ::SiteType"tJ", s::Index)
+##   return Op[s' => 3, s => 3] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"n↓", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Ntot", ::SiteType"tJ", s::Index)
+##   Op[s' => 2, s => 2] = 1.0
+##   return Op[s' => 3, s => 3] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"ntot", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Cup", ::SiteType"tJ", s::Index)
+##   return Op[s' => 1, s => 2] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"c↑", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Cdagup", ::SiteType"tJ", s::Index)
+##   return Op[s' => 2, s => 1] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"c†↑", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Cdn", ::SiteType"tJ", s::Index)
+##   return Op[s' => 1, s => 3] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"c↓", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Cdagdn", ::SiteType"tJ", s::Index)
+##   return Op[s' => 3, s => 1] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"c†↓", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Aup", ::SiteType"tJ", s::Index)
+##   return Op[s' => 1, s => 2] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"a↑", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Adagup", ::SiteType"tJ", s::Index)
+##   return Op[s' => 2, s => 1] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"a†↑", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Adn", ::SiteType"tJ", s::Index)
+##   return Op[s' => 1, s => 3] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"a↓", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Adagdn", ::SiteType"tJ", s::Index)
+##   return Op[s' => 3, s => 1] = 1.0
+## end
+## function op!(Op::ITensor, on::OpName"a†↓", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"F", ::SiteType"tJ", s::Index)
+##   Op[s' => 1, s => 1] = +1.0
+##   Op[s' => 2, s => 2] = -1.0
+##   return Op[s' => 3, s => 3] = -1.0
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Fup", ::SiteType"tJ", s::Index)
+##   Op[s' => 1, s => 1] = +1.0
+##   Op[s' => 2, s => 2] = -1.0
+##   return Op[s' => 3, s => 3] = +1.0
+## end
+## function op!(Op::ITensor, on::OpName"F↑", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Fdn", ::SiteType"tJ", s::Index)
+##   Op[s' => 1, s => 1] = +1.0
+##   Op[s' => 2, s => 2] = +1.0
+##   return Op[s' => 3, s => 3] = -1.0
+## end
+## function op!(Op::ITensor, on::OpName"F↓", st::SiteType"tJ", s::Index)
+##   return op!(Op, alias(on), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Sz", ::SiteType"tJ", s::Index)
+##   Op[s' => 2, s => 2] = +0.5
+##   return Op[s' => 3, s => 3] = -0.5
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Sᶻ", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("Sz"), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Sx", ::SiteType"tJ", s::Index)
+##   Op[s' => 2, s => 3] = 0.5
+##   return Op[s' => 3, s => 2] = 0.5
+## end
+## 
+## function op!(Op::ITensor, ::OpName"Sˣ", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("Sx"), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"S+", ::SiteType"tJ", s::Index)
+##   return Op[s' => 2, s => 3] = 1.0
+## end
+## 
+## function op!(Op::ITensor, ::OpName"S⁺", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("S+"), st, s)
+## end
+## function op!(Op::ITensor, ::OpName"Sp", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("S+"), st, s)
+## end
+## function op!(Op::ITensor, ::OpName"Splus", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("S+"), st, s)
+## end
+## 
+## function op!(Op::ITensor, ::OpName"S-", ::SiteType"tJ", s::Index)
+##   return Op[s' => 3, s => 2] = 1.0
+## end
+## 
+## function op!(Op::ITensor, ::OpName"S⁻", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("S-"), st, s)
+## end
+## function op!(Op::ITensor, ::OpName"Sm", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("S-"), st, s)
+## end
+## function op!(Op::ITensor, ::OpName"Sminus", st::SiteType"tJ", s::Index)
+##   return op!(Op, OpName("S-"), st, s)
+## end
 
 has_fermion_string(::OpName"Cup", ::SiteType"tJ") = true
 function has_fermion_string(on::OpName"c↑", st::SiteType"tJ")
diff --git a/src/space.jl b/src/space.jl
new file mode 100644
index 0000000..78f8fad
--- /dev/null
+++ b/src/space.jl
@@ -0,0 +1,7 @@
+space(st::SiteType; kwargs...) = nothing
+
+space(st::SiteType, n::Int; kwargs...) = space(st; kwargs...)
+
+function space_error_message(st::SiteType)
+  return "Overload of \"space\",\"siteind\", or \"siteinds\" functions not found for Index tag: $(value(st))"
+end
diff --git a/src/state.jl b/src/state.jl
new file mode 100644
index 0000000..1e2ccb8
--- /dev/null
+++ b/src/state.jl
@@ -0,0 +1,13 @@
+@eval struct StateName{Name}
+  (f::Type{<:StateName})() = $(Expr(:new, :f))
+end
+
+StateName(s::AbstractString) = StateName{Symbol(s)}()
+StateName(s::Symbol) = StateName{s}()
+name(::StateName{N}) where {N} = N
+
+macro StateName_str(s)
+  return StateName{Symbol(s)}
+end
+
+state(::StateName, ::SiteType; kwargs...) = nothing
diff --git a/src/val.jl b/src/val.jl
new file mode 100644
index 0000000..f45634e
--- /dev/null
+++ b/src/val.jl
@@ -0,0 +1,14 @@
+@eval struct ValName{Name}
+  (f::Type{<:ValName})() = $(Expr(:new, :f))
+end
+
+ValName(s::AbstractString) = ValName{Symbol(s)}()
+ValName(s::Symbol) = ValName{s}()
+name(::ValName{N}) where {N} = N
+
+macro ValName_str(s)
+  return ValName{Symbol(s)}
+end
+
+val(::ValName, ::SiteType) = nothing
+val(::AbstractString, ::SiteType) = nothing