Skip to content

Commit

Permalink
Apply runic in src/ and test/
Browse files Browse the repository at this point in the history
  • Loading branch information
sebastiendesignolle committed Nov 25, 2024
1 parent 0794dc8 commit 851abe0
Show file tree
Hide file tree
Showing 20 changed files with 475 additions and 476 deletions.
102 changes: 51 additions & 51 deletions src/basic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@
Produces a ket of dimension `d` with nonzero element `i`.
"""
function ket(::Type{T}, i::Integer, d::Integer = 2) where {T<:Number}
function ket(::Type{T}, i::Integer, d::Integer = 2) where {T <: Number}
psi = zeros(T, d)
psi[i] = 1
return psi
Expand All @@ -26,15 +26,15 @@ export ketbra
Produces a projector onto the basis state `i` in dimension `d`.
"""
function proj(::Type{T}, i::Integer, d::Integer = 2) where {T<:Number}
function proj(::Type{T}, i::Integer, d::Integer = 2) where {T <: Number}
p = Hermitian(zeros(T, d, d))
p[i, i] = 1
return p
end
proj(i::Integer, d::Integer = 2) = proj(Bool, i, d)
export proj

const Measurement{T} = Vector{Hermitian{T,Matrix{T}}}
const Measurement{T} = Vector{Hermitian{T, Matrix{T}}}
export Measurement

"""
Expand All @@ -54,7 +54,7 @@ Converts a set of measurements in the common tensor format into a matrix of (her
By default, the second argument is fixed by the size of `A`.
It can also contain custom number of outcomes if there are measurements with less outcomes.
"""
function tensor_to_povm(Aax::Array{T,4}, o::Vector{Int64} = fill(size(Aax, 3), size(Aax, 4))) where {T}
function tensor_to_povm(Aax::Array{T, 4}, o::Vector{Int64} = fill(size(Aax, 3), size(Aax, 4))) where {T}
return [[Hermitian(Aax[:, :, a, x]) for a in 1:o[x]] for x in axes(Aax, 4)]
end
export tensor_to_povm
Expand All @@ -64,7 +64,7 @@ export tensor_to_povm
Converts a matrix of (hermitian) matrices into a set of measurements in the common tensor format.
"""
function povm_to_tensor(Axa::Vector{Measurement{T}}) where {T<:Number}
function povm_to_tensor(Axa::Vector{Measurement{T}}) where {T <: Number}
d, o, m = _measurements_parameters(Axa)
Aax = zeros(T, d, d, maximum(o), m)
for x in eachindex(Axa)
Expand All @@ -76,7 +76,7 @@ function povm_to_tensor(Axa::Vector{Measurement{T}}) where {T<:Number}
end
export povm_to_tensor

function _measurements_parameters(Axa::Vector{Measurement{T}}) where {T<:Number}
function _measurements_parameters(Axa::Vector{Measurement{T}}) where {T <: Number}
@assert !isempty(Axa)
# dimension on which the measurements act
d = size(Axa[1][1], 1)
Expand All @@ -93,11 +93,11 @@ _measurements_parameters(Aa::Measurement) = _measurements_parameters([Aa])
Checks if the measurement defined by A is valid (hermitian, semi-definite positive, and normalized).
"""
function test_povm(E::Vector{<:AbstractMatrix{T}}) where {T<:Number}
function test_povm(E::Vector{<:AbstractMatrix{T}}) where {T <: Number}
!all(ishermitian.(E)) && return false
d = size(E[1], 1)
!(sum(E) I(d)) && return false
for i = 1:length(E)
for i in 1:length(E)
minimum(eigvals(E[i])) < -_rtol(T) && return false
end
return true
Expand All @@ -111,10 +111,10 @@ Constructs the shift operator X of dimension `d` to the power `p`.
Reference: [Generalized Clifford algebra](https://en.wikipedia.org/wiki/Generalized_Clifford_algebra)
"""
function shift(::Type{T}, d::Integer, p::Integer = 1) where {T<:Number}
function shift(::Type{T}, d::Integer, p::Integer = 1) where {T <: Number}
X = zeros(T, d, d)
for i in 0:d-1
X[mod(i + p, d)+1, i+1] = 1
for i in 0:(d - 1)
X[mod(i + p, d) + 1, i + 1] = 1
end
return X
end
Expand All @@ -128,21 +128,21 @@ Constructs the clock operator Z of dimension `d` to the power `p`.
Reference: [Generalized Clifford algebra](https://en.wikipedia.org/wiki/Generalized_Clifford_algebra)
"""
function clock(::Type{T}, d::Integer, p::Integer = 1) where {T<:Number}
function clock(::Type{T}, d::Integer, p::Integer = 1) where {T <: Number}
z = zeros(T, d)
ω = _root_unity(T, d)
for i in 0:d-1
for i in 0:(d - 1)
exponent = mod(i * p, d)
if exponent == 0
z[i+1] = 1
z[i + 1] = 1
elseif 4exponent == d
z[i+1] = im
z[i + 1] = im
elseif 2exponent == d
z[i+1] = -1
z[i + 1] = -1
elseif 4exponent == 3d
z[i+1] = -im
z[i + 1] = -im
else
z[i+1] = ω^exponent
z[i + 1] = ω^exponent
end
end
return Diagonal(z)
Expand All @@ -160,17 +160,17 @@ Vectors of integers between 0 and 3 or strings of I, X, Y, Z
automatically generate Kronecker products of the corresponding
operators.
"""
function pauli(::Type{T}, i::Integer) where {T<:Number}
function pauli(::Type{T}, i::Integer) where {T <: Number}
return gell_mann(T, i ÷ 2 + 1, i % 2 + 1, 2)
end
function pauli(::Type{T}, ind::Vector{<:Integer}) where {T<:Number}
function pauli(::Type{T}, ind::Vector{<:Integer}) where {T <: Number}
if length(ind) == 1
return pauli(T, ind[1])
else
return kron([pauli(T, i) for i in ind]...)
end
end
function pauli(::Type{T}, str::String) where {T<:Number}
function pauli(::Type{T}, str::String) where {T <: Number}
ind = Int64[]
for c in str
if c in ['I', 'i', '1']
Expand All @@ -187,7 +187,7 @@ function pauli(::Type{T}, str::String) where {T<:Number}
end
return pauli(T, ind)
end
pauli(::Type{T}, c::Char) where {T<:Number} = pauli(T, string(c))
pauli(::Type{T}, c::Char) where {T <: Number} = pauli(T, string(c))
pauli(i::Integer) = pauli(ComplexF64, i)
pauli(ind::Vector{<:Integer}) = pauli(ComplexF64, ind)
pauli(str::String) = pauli(ComplexF64, str)
Expand All @@ -202,7 +202,7 @@ Constructs the set `G` of generalized Gell-Mann matrices in dimension `d` such t
Reference: [Generalizations of Pauli matrices](https://en.wikipedia.org/wiki/Generalizations_of_Pauli_matrices)
"""
function gell_mann(::Type{T}, d::Integer = 3) where {T<:Number}
function gell_mann(::Type{T}, d::Integer = 3) where {T <: Number}
return [gell_mann(T, i, j, d) for j in 1:d, i in 1:d][:]
# d=2, ρ = 1/2(σ0 + n*σ)
# d=3, ρ = 1/3(I + √3 n*λ)
Expand All @@ -221,7 +221,7 @@ Constructs the set `i`,`j`th Gell-Mann matrix of dimension `d`.
Reference: [Generalizations of Pauli matrices](https://en.wikipedia.org/wiki/Generalizations_of_Pauli_matrices)
"""
function gell_mann(::Type{T}, i::Integer, j::Integer, d::Integer = 3) where {T<:Number}
function gell_mann(::Type{T}, i::Integer, j::Integer, d::Integer = 3) where {T <: Number}
return gell_mann!(zeros(T, d, d), i, j, d)
end
gell_mann(i::Integer, j::Integer, d::Integer = 3) = gell_mann(ComplexF64, i, j, d)
Expand All @@ -231,7 +231,7 @@ gell_mann(i::Integer, j::Integer, d::Integer = 3) = gell_mann(ComplexF64, i, j,
In-place version of `gell_mann`.
"""
function gell_mann!(res::AbstractMatrix{T}, i::Integer, j::Integer, d::Integer = 3) where {T<:Number}
function gell_mann!(res::AbstractMatrix{T}, i::Integer, j::Integer, d::Integer = 3) where {T <: Number}
if i < j
res[i, j] = 1
res[j, i] = 1
Expand All @@ -244,35 +244,35 @@ function gell_mann!(res::AbstractMatrix{T}, i::Integer, j::Integer, d::Integer =
end
elseif i == d
tmp = _sqrt(T, 2) / _sqrt(T, d * (d - 1))
for k in 1:d-1
for k in 1:(d - 1)
res[k, k] = tmp
end
res[d, d] = -(d - 1) * tmp
else
gell_mann!(view(res, 1:d-1, 1:d-1), i, j, d - 1)
gell_mann!(view(res, 1:(d - 1), 1:(d - 1)), i, j, d - 1)
end
return res
end
export gell_mann!

_rtol(::Type{T}) where {T<:Number} = Base.rtoldefault(real(T))
_rtol(::Type{T}) where {T <: Number} = Base.rtoldefault(real(T))

_eps(::Type{T}) where {T<:Number} = _realeps(real(T))
_realeps(::Type{T}) where {T<:AbstractFloat} = eps(T)
_eps(::Type{T}) where {T <: Number} = _realeps(real(T))
_realeps(::Type{T}) where {T <: AbstractFloat} = eps(T)
_realeps(::Type{<:Real}) = 0

"""
cleanup!(M::AbstractArray{T}; tol = Base.rtoldefault(real(T)))
Zeroes out real or imaginary parts of `M` that are smaller than `tol`.
"""
function cleanup!(M::AbstractArray{T}; tol = _eps(T)) where {T<:Number}
function cleanup!(M::AbstractArray{T}; tol = _eps(T)) where {T <: Number}
wrapper = Base.typename(typeof(M)).wrapper
cleanup!(parent(M); tol)
return wrapper(M)
end

function cleanup!(M::Array{T}; tol = _eps(T)) where {T<:Number}
function cleanup!(M::Array{T}; tol = _eps(T)) where {T <: Number}
if isbitstype(T)
M2 = reinterpret(real(T), M) #this is a no-op when T<:Real
_cleanup!(M2; tol)
Expand All @@ -288,7 +288,7 @@ end
export cleanup!

function _cleanup!(M; tol)
return M[abs.(M).<tol] .= 0
return M[abs.(M) .< tol] .= 0
end

function applykraus(K, M)
Expand All @@ -297,20 +297,20 @@ end
export applykraus

function _orthonormal_range_svd!(
A::AbstractMatrix{T};
tol::Union{Real,Nothing} = nothing,
alg = LinearAlgebra.default_svd_alg(A)
) where {T<:Number}
A::AbstractMatrix{T};
tol::Union{Real, Nothing} = nothing,
alg = LinearAlgebra.default_svd_alg(A)
) where {T <: Number}
dec = svd!(A; alg = alg)
tol = isnothing(tol) ? maximum(dec.S) * _eps(T) * minimum(size(A)) : tol
rank = sum(dec.S .> tol)
dec.U[:, 1:rank]
end

_orthonormal_range_svd(A::AbstractMatrix; tol::Union{Real,Nothing} = nothing) =
_orthonormal_range_svd(A::AbstractMatrix; tol::Union{Real, Nothing} = nothing) =
_orthonormal_range_svd!(deepcopy(A); tol = tol)

function _orthonormal_range_qr(A::SA.AbstractSparseMatrix{T,M}; tol::Union{Real,Nothing} = nothing) where {T<:Number,M}
function _orthonormal_range_qr(A::SA.AbstractSparseMatrix{T, M}; tol::Union{Real, Nothing} = nothing) where {T <: Number, M}
dec = qr(A)
tol = isnothing(tol) ? maximum(abs.(dec.R)) * _eps(T) : tol
rank = sum(abs.(Diagonal(dec.R)) .> tol)
Expand All @@ -325,10 +325,10 @@ otherwise uses an SVD and returns a dense matrix (`mode = 1`). Input `A` will be
Tolerance `tol` is used to compute the rank and is automatically set if not provided.
"""
function orthonormal_range(
A::SA.AbstractMatrix{T};
mode::Integer = -1,
tol::Union{Real,Nothing} = nothing
) where {T<:Number}
A::SA.AbstractMatrix{T};
mode::Integer = -1,
tol::Union{Real, Nothing} = nothing
) where {T <: Number}
mode == 1 && SA.issparse(A) && throw(ArgumentError("SVD does not work with sparse matrices, use a dense matrix."))
mode == -1 && (mode = SA.issparse(A) ? 0 : 1)

Expand Down Expand Up @@ -383,12 +383,12 @@ This function returns a generator, which can then be used e.g. in for loops with
entire basis at once. If you need a vector, call `collect` on it.
"""
function n_body_basis(
n::Integer,
n_parties::Integer;
sb::AbstractVector{<:AbstractMatrix} = [pauli(1), pauli(2), pauli(3)],
sparse::Bool = true,
eye::AbstractMatrix = I(size(sb[1], 1))
)
n::Integer,
n_parties::Integer;
sb::AbstractVector{<:AbstractMatrix} = [pauli(1), pauli(2), pauli(3)],
sparse::Bool = true,
eye::AbstractMatrix = I(size(sb[1], 1))
)
(n >= 0 && n_parties >= 2) || throw(ArgumentError("Number of parties must be ≥ 2 and n ≥ 0."))
n <= n_parties || throw(ArgumentError("Number of parties cannot be larger than n."))

Expand All @@ -397,8 +397,8 @@ function n_body_basis(
idx_eye = length(sb)
(
Base.kron(sb[t]...)
for p in Combinatorics.with_replacement_combinations(1:nb, n)
for t in Combinatorics.multiset_permutations([p; repeat([idx_eye], n_parties - n)], n_parties)
for p in Combinatorics.with_replacement_combinations(1:nb, n)
for t in Combinatorics.multiset_permutations([p; repeat([idx_eye], n_parties - n)], n_parties)
)
end
export n_body_basis
Loading

0 comments on commit 851abe0

Please sign in to comment.