Skip to content

Commit 4a11ce8

Browse files
mtfishmanmtfishman
committed
Make more generic operators
1 parent af105ed commit 4a11ce8

File tree

13 files changed

+500
-539
lines changed

13 files changed

+500
-539
lines changed

‎README.md‎

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -73,8 +73,9 @@ TODO: Delete.
7373
@test Vector{Float32}(StateName("0"), SiteType("Qubit")) == Float32[1, 0]
7474
@test Vector{Float32}(StateName("1"), SiteType("Qubit")) == Float32[0, 1]
7575

76-
@test AbstractArray(StateName("Up"), SiteType("Qubit")) == [1, 0]
77-
@test AbstractArray(StateName("Dn"), SiteType("Qubit")) == [0, 1]
76+
# TODO: Add this back.
77+
# @test AbstractArray(StateName("Up"), SiteType("Qubit")) == [1, 0]
78+
# @test AbstractArray(StateName("Dn"), SiteType("Qubit")) == [0, 1]
7879

7980
@test Matrix(OpName("X"), SiteType("Qubit")) == [0 1; 1 0]
8081
@test Matrix(OpName("σx"), SiteType("Qubit")) == [0 1; 1 0]

‎examples/README.jl‎

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -75,8 +75,9 @@ using Test: @test
7575
@test Vector{Float32}(StateName("0"), SiteType("Qubit")) == Float32[1, 0]
7676
@test Vector{Float32}(StateName("1"), SiteType("Qubit")) == Float32[0, 1]
7777

78-
@test AbstractArray(StateName("Up"), SiteType("Qubit")) == [1, 0]
79-
@test AbstractArray(StateName("Dn"), SiteType("Qubit")) == [0, 1]
78+
## TODO: Add this back.
79+
## @test AbstractArray(StateName("Up"), SiteType("Qubit")) == [1, 0]
80+
## @test AbstractArray(StateName("Dn"), SiteType("Qubit")) == [0, 1]
8081

8182
@test Matrix(OpName("X"), SiteType("Qubit")) == [0 1; 1 0]
8283
@test Matrix(OpName("σx"), SiteType("Qubit")) == [0 1; 1 0]

‎src/ITensorQuantumOperatorDefinitions.jl‎

Lines changed: 0 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -7,9 +7,6 @@ include("state.jl")
77
include("op.jl")
88
include("has_fermion_string.jl")
99

10-
include("sitetypes/generic.jl")
11-
include("sitetypes/aliases.jl")
12-
include("sitetypes/generic_sites.jl")
1310
include("sitetypes/qubit.jl")
1411
include("sitetypes/spinhalf.jl")
1512
include("sitetypes/spinone.jl")

‎src/op.jl‎

Lines changed: 327 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -15,6 +15,84 @@ macro OpName_str(s)
1515
return OpName{Symbol(s)}
1616
end
1717

18+
function op_alias_expr(name1, name2, pars...)
19+
return :(function alias(n::OpName{Symbol($name1)})
20+
return OpName{Symbol($name2)}(; params(n)..., $(esc.(pars)...))
21+
end)
22+
end
23+
macro op_alias(name1, name2, pars...)
24+
return op_alias_expr(name1, name2, pars...)
25+
end
26+
27+
alias(n::OpName) = n
28+
function op_convert(
29+
arrtype::Type{<:AbstractArray{<:Any,N}},
30+
domain_size::Tuple{Vararg{Integer}},
31+
a::AbstractArray{<:Any,N},
32+
) where {N}
33+
# TODO: Check the dimensions.
34+
return convert(arrtype, a)
35+
end
36+
function op_convert(
37+
arrtype::Type{<:AbstractArray}, domain_size::Tuple{Vararg{Integer}}, a::AbstractArray
38+
)
39+
# TODO: Check the dimensions.
40+
return convert(arrtype, a)
41+
end
42+
function op_convert(
43+
arrtype::Type{<:AbstractArray{<:Any,N}},
44+
domain_size::Tuple{Vararg{Integer}},
45+
a::AbstractArray,
46+
) where {N}
47+
size = (domain_size..., domain_size...)
48+
@assert length(size) == N
49+
return convert(arrtype, reshape(a, size))
50+
end
51+
function (arrtype::Type{<:AbstractArray})(n::OpName, ts::Tuple{Vararg{SiteType}})
52+
return op_convert(arrtype, length.(ts), AbstractArray(n, ts))
53+
end
54+
function (arrtype::Type{<:AbstractArray})(n::OpName, domain_size::Tuple{Vararg{Integer}})
55+
return op_convert(arrtype, domain_size, AbstractArray(n, Int.(domain_size)))
56+
end
57+
function (arrtype::Type{<:AbstractArray})(n::OpName, domain_size::Integer...)
58+
return arrtype(n, domain_size)
59+
end
60+
(arrtype::Type{<:AbstractArray})(n::OpName, ts::SiteType...) = arrtype(n, ts)
61+
Base.AbstractArray(n::OpName, ts::SiteType...) = AbstractArray(n, ts)
62+
function Base.AbstractArray(n::OpName, ts::Tuple{Vararg{SiteType}})
63+
n′ = alias(n)
64+
ts′ = alias.(ts)
65+
if n′ == n && ts′ == ts
66+
return AbstractArray(n′, length.(ts′))
67+
end
68+
return AbstractArray(n′, ts′)
69+
end
70+
function Base.AbstractArray(n::OpName, domain_size::Tuple{Vararg{Int}})
71+
n′ = alias(n)
72+
if n′ == n
73+
error("Not implemented.")
74+
end
75+
return AbstractArray(n′, domain_size)
76+
end
77+
78+
# TODO: Decide on this.
79+
function Base.AbstractArray(n::OpName)
80+
return AbstractArray(n, ntuple(Returns(default_sitetype()), nsites(n)))
81+
end
82+
function (arrtype::Type{<:AbstractArray})(n::OpName)
83+
return arrtype(n, ntuple(Returns(default_sitetype()), nsites(n)))
84+
end
85+
86+
# `Int(ndims // 2)`, i.e. `ndims_domain`/`ndims_codomain`.
87+
function nsites(n::Union{StateName,OpName})
88+
n′ = alias(n)
89+
if n′ == n
90+
# Default value, assume 1-site operator/state.
91+
return 1
92+
end
93+
return nsites(n′)
94+
end
95+
1896
## TODO: Delete.
1997
## # Default implementations of op
2098
## op(::OpName; kwargs...) = nothing
@@ -34,3 +112,252 @@ end
34112

35113
# TODO: Bring this back?
36114
# op(f::Function, args...; kwargs...) = f(op(args...; kwargs...))
115+
116+
using LinearAlgebra: Diagonal
117+
function Base.AbstractArray(::OpName"Id", domain_size::Tuple{Int})
118+
return Diagonal(trues(only(domain_size)))
119+
end
120+
function Base.AbstractArray(n::OpName"Id", domain_size::Tuple{Int,Vararg{Int}})
121+
return Base.AbstractArray(kron(ntuple(Returns(n), length(domain_size))...), domain_size)
122+
end
123+
@op_alias "I" "Id"
124+
@op_alias "σ0" "Id"
125+
@op_alias "σ⁰" "Id"
126+
@op_alias "σ₀" "Id"
127+
# TODO: Is this a good definition?
128+
@op_alias "F" "Id"
129+
130+
# TODO: Define as `::Tuple{Int}`.
131+
function Base.AbstractArray(n::OpName"a†", domain_size::Tuple{Int})
132+
d = only(domain_size)
133+
a = zeros(d, d)
134+
for k in 1:(d - 1)
135+
a[k + 1, k] = √k
136+
end
137+
return a
138+
end
139+
@op_alias "Adag" "a†"
140+
@op_alias "adag" "a†"
141+
alias(::OpName"a") = OpName"a†"()'
142+
@op_alias "A" "a"
143+
alias(::OpName"n") = OpName"a†"() * OpName"a"()
144+
@op_alias "N" "n"
145+
alias(::OpName"aa") = OpName("a") ⊗ OpName("a")
146+
alias(::OpName"a†a") = OpName("a†") ⊗ OpName("a")
147+
alias(::OpName"aa†") = OpName("a") ⊗ OpName("a†")
148+
alias(::OpName"a†a†") = OpName("a†") ⊗ OpName("a†")
149+
150+
δ(x, y) = (x == y) ? 1 : 0
151+
152+
function Base.AbstractArray(n::OpName"σ⁺", domain_size::Tuple{Int})
153+
d = only(domain_size)
154+
s = (d - 1) / 2
155+
return [2 * δ(i + 1, j) * √((s + 1) * (i + j - 1) - i * j) for i in 1:d, j in 1:d]
156+
end
157+
alias(::OpName"S⁺") = OpName("σ⁺") / 2
158+
@op_alias "S+" "S⁺"
159+
@op_alias "Splus" "S+"
160+
@op_alias "Sp" "S+"
161+
162+
alias(::OpName"σ⁻") = OpName"σ⁺"()'
163+
alias(::OpName"S⁻") = OpName("σ⁻") / 2
164+
@op_alias "S-" "S⁻"
165+
@op_alias "Sminus" "S-"
166+
@op_alias "Sm" "S-"
167+
168+
alias(::OpName"X") = (OpName"σ⁺"() + OpName"σ⁻"()) / 2
169+
@op_alias "σx" "X"
170+
@op_alias "σˣ" "X"
171+
@op_alias "σₓ" "X"
172+
@op_alias "σ1" "X"
173+
@op_alias "σ¹" "X"
174+
@op_alias "σ₁" "X"
175+
@op_alias "σy" "Y"
176+
@op_alias "iX" "im" op = OpName"X"()
177+
@op_alias "√X" "√" op = OpName"X"()
178+
@op_alias "√NOT" "√" op = OpName"X"()
179+
alias(n::OpName"Sx") = OpName("X") / 2
180+
@op_alias "Sˣ" "Sx"
181+
@op_alias "Sₓ" "Sx"
182+
183+
alias(::OpName"Y") = -im * (OpName"σ⁺"() - OpName"σ⁻"()) / 2
184+
# TODO: No subsript `\_y` available
185+
# in unicode.
186+
@op_alias "σʸ" "Y"
187+
@op_alias "σ2" "Y"
188+
@op_alias "σ²" "Y"
189+
@op_alias "σ₂" "Y"
190+
function alias(::OpName"iY")
191+
return real(OpName"Y"()im)
192+
end
193+
@op_alias "iσy" "iY"
194+
# TODO: No subsript `\_y` available
195+
# in unicode.
196+
@op_alias "iσʸ" "iY"
197+
@op_alias "iσ2" "iY"
198+
@op_alias "iσ²" "iY"
199+
@op_alias "iσ₂" "iY"
200+
alias(n::OpName"Sy") = OpName("Y") / 2
201+
@op_alias "Sʸ" "Sy"
202+
alias(n::OpName"iSy") = OpName("iY") / 2
203+
@op_alias "iSʸ" "iSy"
204+
205+
function Base.AbstractArray(n::OpName"σᶻ", domain_size::Tuple{Int})
206+
d = only(domain_size)
207+
s = (d - 1) / 2
208+
return Diagonal([2 * (s + 1 - i) for i in 1:d])
209+
end
210+
# TODO: No subsript `\_z` available
211+
# in unicode.
212+
@op_alias "Z" "σᶻ"
213+
@op_alias "σ3" "Z"
214+
@op_alias "σ³" "Z"
215+
@op_alias "σ₃" "Z"
216+
@op_alias "σz" "Z"
217+
@op_alias "iZ" "im" op = OpName"Z"()
218+
alias(n::OpName"Sz") = OpName("Z") / 2
219+
@op_alias "Sᶻ" "Sz"
220+
# TODO: Make sure it ends up real, using `S⁺` and `S⁻`,
221+
# define `Sx²`, `Sy²`, and `Sz²`, etc.
222+
# Should be equal to:
223+
# ```julia
224+
# α * I = s * (s + 1) * I
225+
# = (d - 1) * (d + 1) / 4 * I
226+
# ```
227+
# where `s = (d - 1) / 2`. See https://en.wikipedia.org/wiki/Spin_(physics)#Higher_spins.
228+
alias(n::OpName"S2") = (OpName("Sˣ")^2 + OpName("Sʸ")^2 + OpName("Sᶻ")^2)
229+
@op_alias "S²" "S2"
230+
231+
using LinearAlgebra: eigen
232+
function Base.AbstractArray(n::OpName"H", domain_size::Tuple{Int})
233+
Λ, H = eigen(AbstractArray(OpName("X"), domain_size))
234+
p = sortperm(Λ; rev=true)
235+
return H[:, p]
236+
end
237+
238+
## TODO: Bring back these definitions.
239+
## function default_random_matrix(eltype::Type, s::Index...)
240+
## n = prod(dim.(s))
241+
## return randn(eltype, n, n)
242+
## end
243+
##
244+
## # Haar-random unitary
245+
## #
246+
## # Reference:
247+
## # Section 4.6
248+
## # http://math.mit.edu/~edelman/publications/random_matrix_theory.pdf
249+
## function op!(
250+
## o::ITensor,
251+
## ::OpName"RandomUnitary",
252+
## ::SiteType"Generic",
253+
## s1::Index,
254+
## sn::Index...;
255+
## eltype=ComplexF64,
256+
## random_matrix=default_random_matrix(eltype, s1, sn...),
257+
## )
258+
## s = (s1, sn...)
259+
## Q, _ = NDTensors.qr_positive(random_matrix)
260+
## t = itensor(Q, prime.(s)..., dag.(s)...)
261+
## return settensor!(o, tensor(t))
262+
## end
263+
##
264+
## function op!(o::ITensor, ::OpName"randU", st::SiteType"Generic", s::Index...; kwargs...)
265+
## return op!(o, OpName("RandomUnitary"), st, s...; kwargs...)
266+
## end
267+
268+
# Unary operations
269+
nsites(n::OpName"f") = nsites(n.op)
270+
function Base.AbstractArray(n::OpName"f", domain_size::Tuple{Vararg{Int}})
271+
return n.f(AbstractArray(n.op, domain_size))
272+
end
273+
274+
for f in (
275+
:(Base.sqrt),
276+
:(Base.real),
277+
:(Base.imag),
278+
:(Base.complex),
279+
:(Base.exp),
280+
:(Base.cos),
281+
:(Base.sin),
282+
:(Base.adjoint),
283+
:(Base.:+),
284+
:(Base.:-),
285+
)
286+
@eval begin
287+
$f(n::OpName) = OpName"f"(; f=$f, op=n)
288+
end
289+
end
290+
291+
nsites(n::OpName"^") = nsites(n.op)
292+
function Base.AbstractArray(n::OpName"^", domain_size::Tuple{Vararg{Int}})
293+
return AbstractArray(n.op, domain_size)^n.exponent
294+
end
295+
Base.:^(n::OpName, exponent) = OpName"^"(; op=n, exponent)
296+
297+
nsites(n::OpName"kron") = nsites(n.op1) + nsites(n.op2)
298+
function Base.AbstractArray(n::OpName"kron", domain_size::Tuple{Vararg{Int}})
299+
domain_size1 = domain_size[1:nsites(n.op1)]
300+
domain_size2 = domain_size[(nsites(n.op1) + 1):end]
301+
@assert length(domain_size2) == nsites(n.op2)
302+
return kron(AbstractArray(n.op1, domain_size1), AbstractArray(n.op2, domain_size2))
303+
end
304+
Base.kron(n1::OpName, n2::OpName) = OpName"kron"(; op1=n1, op2=n2)
305+
⊗(n1::OpName, n2::OpName) = kron(n1, n2)
306+
307+
function nsites(n::OpName"+")
308+
@assert nsites(n.op1) == nsites(n.op2)
309+
return nsites(n.op1)
310+
end
311+
function Base.AbstractArray(n::OpName"+", domain_size::Tuple{Vararg{Int}})
312+
return AbstractArray(n.op1, domain_size) + AbstractArray(n.op2, domain_size)
313+
end
314+
Base.:+(n1::OpName, n2::OpName) = OpName"+"(; op1=n1, op2=n2)
315+
Base.:-(n1::OpName, n2::OpName) = n1 + (-n2)
316+
317+
function nsites(n::OpName"*")
318+
@assert nsites(n.op1) == nsites(n.op2)
319+
return nsites(n.op1)
320+
end
321+
function Base.AbstractArray(n::OpName"*", domain_size::Tuple{Vararg{Int}})
322+
return AbstractArray(n.op1, domain_size) * AbstractArray(n.op2, domain_size)
323+
end
324+
Base.:*(n1::OpName, n2::OpName) = OpName"*"(; op1=n1, op2=n2)
325+
326+
nsites(n::OpName"scaled") = nsites(n.op)
327+
function Base.AbstractArray(n::OpName"scaled", domain_size::Tuple{Vararg{Int}})
328+
return AbstractArray(n.op, domain_size) * n.c
329+
end
330+
function Base.:*(c::Number, n::OpName)
331+
return OpName"scaled"(; op=n, c)
332+
end
333+
function Base.:*(n::OpName, c::Number)
334+
return OpName"scaled"(; op=n, c)
335+
end
336+
function Base.:/(n::OpName, c::Number)
337+
return OpName"scaled"(; op=n, c=inv(c))
338+
end
339+
340+
function Base.:*(c::Number, n::OpName"scaled")
341+
return OpName"scaled"(; op=n.op, c=(c * n.c))
342+
end
343+
function Base.:*(n::OpName"scaled", c::Number)
344+
return OpName"scaled"(; op=n.op, c=(n.c * c))
345+
end
346+
function Base.:/(n::OpName"scaled", c::Number)
347+
return OpName"scaled"(; op=n.op, c=(n.c / c))
348+
end
349+
350+
alias(n::OpName"im") = OpName"scaled"(; op=n.op, c=im)
351+
352+
controlled(n::OpName; ncontrol=1) = OpName"Control"(; ncontrol, op=n)
353+
354+
# Expand the operator in a new basis.
355+
using LinearAlgebra: ⋅
356+
function expand(v::AbstractArray, basis)
357+
gramian = [basisᵢ ⋅ basisⱼ for basisᵢ in basis, basisⱼ in basis]
358+
vbasis = [basisᵢ ⋅ v for basisᵢ in basis]
359+
return gramian \ vbasis
360+
end
361+
function expand(n::OpName, basis, ts...)
362+
return expand(AbstractArray(n, ts...), map(basisᵢ -> AbstractArray(basisᵢ, ts...), basis))
363+
end

‎src/sitetype.jl‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,3 +14,10 @@ value(::SiteType{T}) where {T} = T
1414
macro SiteType_str(s)
1515
return SiteType{Symbol(s)}
1616
end
17+
18+
alias(t::SiteType) = t
19+
Base.AbstractUnitRange(t::SiteType) = Base.OneTo(length(t))
20+
Base.size(t::SiteType) = (length(t),)
21+
Base.size(t::SiteType, dim::Integer) = size(t)[dim]
22+
Base.axes(t::SiteType) = (AbstractUnitRange(t),)
23+
Base.axes(t::SiteType, dim::Integer) = axes(t)[dim]

0 commit comments

Comments
 (0)