Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Generalized CTranspose + conversion methods and promotion rules #33

Merged
merged 9 commits into from
May 20, 2015
49 changes: 21 additions & 28 deletions src/arrays/arraymath.jl
Original file line number Diff line number Diff line change
Expand Up @@ -27,13 +27,13 @@ end
# operator * ket -> ket
function *{B<:OrthonormalBasis}(op::AbstractQuMatrix{B}, ket::AbstractQuVector{B})
vc = rawcoeffs(op)*rawcoeffs(ket)
QAT = similar_type(typeof(ket))
QAT = similar_type(ket)
return QAT(vc, bases(op,1))
end

function *{B<:OrthonormalBasis}(op::DualMatrix{B}, ket::AbstractQuVector{B})
vc = Ac_mul_B(rawcoeffs(op), rawcoeffs(ket))
QAT = similar_type(typeof(ket))
QAT = similar_type(ket)
return QAT(vc, bases(op,1))
end

Expand All @@ -50,29 +50,27 @@ end

function *{B<:OrthonormalBasis}(dm::DualMatrix{B}, qm::AbstractQuMatrix{B})
mc = Ac_mul_B(rawcoeffs(dm), rawcoeffs(qm))
QAT = similar_type(typeof(qm))
QAT = similar_type(qm)
return QAT(mc, bases(dm,1), bases(qm,2))
end

function *{B<:OrthonormalBasis}(qm::AbstractQuMatrix{B}, dm::DualMatrix{B})
mc = A_mul_Bc(rawcoeffs(dm), rawcoeffs(qm))
QAT = similar_type(typeof(qm))
QAT = similar_type(qm)
return QAT(mc, bases(qm,1), bases(dm,2))
end

function *{B<:OrthonormalBasis}(qm1::AbstractQuMatrix{B}, qm2::AbstractQuMatrix{B})
mc = rawcoeffs(qm1)*rawcoeffs(qm2)
QAT = similar_type(promote_type(typeof(qm1), typeof(qm2)))
QAT = similar_type(qm1, qm2)
return QAT(mc, bases(qm1,1), bases(qm2,2))
end



# addition and subtraction
function +{B<:AbstractBasis}(qarr1::AbstractQuArray{B}, qarr2::AbstractQuArray{B})
if bases(qarr1) == bases(qarr2)
sc = coeffs(qarr1) + coeffs(qarr2)
QAT = similar_type(promote_type(typeof(qarr1), typeof(qarr2)))
QAT = similar_type(qarr1, qarr2)
return QAT(sc, bases(qarr1))
else
error("Bases not compatible")
Expand All @@ -82,7 +80,7 @@ end
function -{B<:AbstractBasis}(qarr1::AbstractQuArray{B}, qarr2::AbstractQuArray{B})
if bases(qarr1) == bases(qarr2)
sc = coeffs(qarr1) - coeffs(qarr2)
QAT = similar_type(promote_type(typeof(qarr1), typeof(qarr2)))
QAT = similar_type(qarr1, qarr2)
return QAT(sc, bases(qarr1))
else
error("Bases not compatible")
Expand All @@ -92,15 +90,14 @@ end
+(ct1::CTranspose, ct2::CTranspose) = (ct1.qarr+ct2.qarr)'
-(ct1::CTranspose, ct2::CTranspose) = (ct1.qarr-ct2.qarr)'


# scaling
Base.scale!(num::Number, qarr::AbstractQuArray) = (scale!(num, rawcoeffs(qarr)); return qarr)
Base.scale!(num::Number, ct::CTranspose) = CTranspose(scale!(num', ct.qarr))
Base.scale!(qarr::AbstractQuArray, num::Number) = scale!(num, qarr)

function Base.scale(num::Number, qarr::AbstractQuArray)
fc = scale(num, rawcoeffs(qarr))
QAT = QuBase.similar_type(typeof(qarr))
QAT = similar_type(qarr)
return QAT(fc, bases(qarr))
end
Base.scale(num::Number, ct::CTranspose) = CTranspose(scale(num', ct.qarr))
Expand All @@ -110,7 +107,6 @@ Base.scale(qarr::AbstractQuArray, num::Number) = scale(num, qarr)
*(qarr::AbstractQuArray, num::Number) = scale(qarr, num)
/(qarr::AbstractQuArray, num::Number) = scale(1/num, qarr)


# matrix operations returning a scalar
# normalization
Base.norm(qarr::AbstractQuArray) = vecnorm(rawcoeffs(qarr))
Expand All @@ -122,48 +118,45 @@ end

normalize(qarr::AbstractQuArray) = normalize!(copy(qarr))


# matrix operations returning an array
# sparse to dense
function Base.full(qarr::AbstractQuMatrix)
fc = full(rawcoeffs(qarr))
QAT = similar_type(typeof(qarr))
QAT = similar_type(qarr)
return QAT(fc, bases(qarr))
end
#Base.full(ct::CTranspose) = full(ct.qarr)'

# exponential of dense matrix
function Base.expm(qarr::AbstractQuMatrix)
fc = expm(full(rawcoeffs(qarr)))
QAT = similar_type(typeof(qarr))
QAT = similar_type(qarr)
return QAT(fc, bases(qarr))
end
#Base.expm(ct::CTranspose) = expm(ct.qarr)'


##################
# Tensor Product #
##################

# General tensor product definitions for orthonormal bases
function tensor{B<:OrthonormalBasis,T1,T2,N}(qarr1::AbstractQuArray{B,T1,N}, qarr2::AbstractQuArray{B,T2,N})
# General tensor product definitions
function tensor{B1,B2,T1,T2,N}(qarr1::AbstractQuArray{B1,T1,N}, qarr2::AbstractQuArray{B2,T2,N})
tc = kron(coeffs(qarr1), coeffs(qarr2))
QAT = similar_type(promote_type(typeof(qarr1), typeof(qarr2)))
QAT = similar_type(qarr1, qarr2)
return QAT(tc, map(tensor, bases(qarr1), bases(qarr2)))
end

function tensor{B<:OrthonormalBasis,T1,T2,N}(qarr1::CTranspose{B,T1,N}, qarr2::CTranspose{B,T2,N})
return QuArray(kron(rawcoeffs(qarr1), rawcoeffs(qarr2)),
map(tensor, rawbases(qarr1), rawbases(qarr2)))'
end
# defined to resolve ambiguity warnings
tensor(ct1::CTranspose, ct2::CTranspose) = tensor(ct1.qarr, ct2.qarr)'
tensor(ct1::DualVector, ct2::DualVector) = tensor(ct1.qarr, ct2.qarr)'

function tensor{B<:OrthonormalBasis}(ket::QuVector{B}, bra::DualVector{B})
return QuArray(kron(coeffs(ket), coeffs(bra)),
bases(ket,1),
bases(bra,1))
function tensor(ket::AbstractQuVector, bra::DualVector)
tc = kron(coeffs(ket), coeffs(bra))
QAT = similar_type(ket)
return QAT(tc, bases(ket,1), bases(bra,1))
end

tensor{B<:OrthonormalBasis}(bra::DualVector{B}, ket::QuVector{B}) = tensor(ket, bra)
tensor(bra::DualVector, ket::AbstractQuVector) = tensor(ket, bra)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am probably just missing it, but does that make sense?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess it's a matter of intent...I mainly defined that thinking that if one takes the tensor product of Bra and a Ket, they would expect the outer product of the states? For example, isn't the following is an accepted use of the ⊗ notation:

⟨ ψ | ⊗ | ϕ ⟩ = | ϕ ⟩ ⊗ ⟨ ψ | = | ϕ ⟩⟨ ψ |

Not a peer-reviewed publication, obviously, but: http://physics.stackexchange.com/questions/134575/tensor-product-of-a-bra-and-a-ket

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see. Never thought about it that way, but it makes perfect sense of course!


###############
# Commutators #
Expand Down
1 change: 1 addition & 0 deletions src/arrays/constructors.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@
QAT = similar_type(typeof(qa))
return QAT(fc, bases(qa))
end

function Base.eye(qa::AbstractQuArray)
fc = eye(coeffs(qa))
QAT = similar_type(typeof(qa))
Expand Down
77 changes: 53 additions & 24 deletions src/arrays/quarray.jl
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,9 @@
########################
# Array-like functions #
########################
Base.eltype{B,T,N}(::AbstractQuArray{B,T,N}) = T
Base.eltype{B,T,N}(::Type{AbstractQuArray{B,T,N}}) = T

Base.size(qarr::AbstractQuArray, i...) = size(rawcoeffs(qarr), i...)
Base.ndims(qarr::AbstractQuArray) = ndims(rawcoeffs(qarr))
Base.length(qarr::AbstractQuArray) = length(rawcoeffs(qarr))
Expand All @@ -33,6 +36,8 @@
typealias AbstractQuVector{B<:AbstractBasis,T} AbstractQuArray{B,T,1}
typealias AbstractQuMatrix{B<:AbstractBasis,T} AbstractQuArray{B,T,2}

similar_type{Q<:AbstractQuArray}(::Q) = similar_type(Q)
similar_type{A<:AbstractQuArray,B<:AbstractQuArray}(::A,::B) = similar_type(promote_type(A,B))

###########
# QuArray #
Expand Down Expand Up @@ -61,43 +66,44 @@
# and a tuple of bases (see QuArray(coeffs, bases) above)
similar_type{Q<:QuArray}(::Type{Q}) = QuArray


######################
# Accessor functions #
######################
coefftype{B,T,N,A}(::QuArray{B,T,N,A}) = A
coefftype{B,T,N,A}(::Type{QuArray{B,T,N,A}}) = A

rawcoeffs(qarr::QuArray) = qarr.coeffs

rawbases(qarr::QuArray) = qarr.bases
rawbases(qarr::QuArray, i) = qarr.bases[i]

########################
# Array-like functions #
########################
Base.copy(qa::QuArray) = QuArray(copy(qa.coeffs), copy(qa.bases))

Base.eltype{B,T,N,A}(::Type{QuArray{B,T,N,A}}) = T

##############
# CTranspose #
##############
immutable CTranspose{B<:AbstractBasis,T,N,A} <: AbstractQuArray{B,T,N}
qarr::QuArray{B,T,N,A}
CTranspose(qarr::QuVector{B,T,A}) = new(qarr)
CTranspose(qarr::QuMatrix{B,T,A}) = new(qarr)
CTranspose(qarr::QuArray) = error("Conjugate transpose is unsupported for QuArrays of dimension $N")
immutable CTranspose{B<:AbstractBasis,T,N,Q} <: AbstractQuArray{B,T,N}
qarr::Q
CTranspose(qarr::AbstractQuVector{B,T}) = new(qarr)
CTranspose(qarr::AbstractQuMatrix{B,T}) = new(qarr)
CTranspose(qarr::AbstractQuArray) = error("Conjugate transpose is unsupported for AbstractQuArrays of dimension $N")
end

CTranspose{B<:AbstractBasis,T,N,A}(qa::QuArray{B,T,N,A}) = CTranspose{B,T,N,A}(qa)
CTranspose{B<:AbstractBasis,T,N}(qa::AbstractQuArray{B,T,N}) = CTranspose{B,T,N,typeof(qa)}(qa)
CTranspose{B<:AbstractBasis,T,N}(coeffs::AbstractArray{T,N}, bases::NTuple{N,B}) = CTranspose(QuArray(coeffs, bases))

typealias DualVector{B,T,A} CTranspose{B,T,1,A}
typealias DualMatrix{B,T,A} CTranspose{B,T,2,A}
typealias DualVector{B<:AbstractBasis,T,Q} CTranspose{B,T,1,Q}
typealias DualMatrix{B<:AbstractBasis,T,Q} CTranspose{B,T,2,Q}

similar_type{QC<:CTranspose}(::Type{QC}) = CTranspose


######################
# Accessor functions #
######################
coefftype{B,T,N,A}(::CTranspose{B,T,N,A}) = A
coefftype(ct::CTranspose) = coefftype(ct.qarr)
coefftype{B,T,N,Q}(::Type{CTranspose{B,T,N,Q}}) = coefftype(Q)

rawcoeffs(ct::CTranspose) = rawcoeffs(ct.qarr)
coeffs(ct::CTranspose) = rawcoeffs(ct)'
Expand All @@ -109,57 +115,80 @@

bases(ct::CTranspose, i) = rawbases(ct, revind(ndims(ct), i))

qarr_type{B,T,N,Q}(::Type{CTranspose{B,T,N,Q}}) = Q
qarr_type(ct::CTranspose) = qarr_type(typeof(ct))

similar_type{QC<:CTranspose}(::Type{QC}) = CTranspose

########################
# Array-like functions #
########################
Base.eltype{B,T,N,Q}(::Type{CTranspose{B,T,N,Q}}) = T

Base.copy(ct::CTranspose) = CTranspose(copy(ct.qarr))

Base.ndims(ct::CTranspose) = ndims(ct.qarr)
Base.length(ct::CTranspose) = length(ct.qarr)

Base.size(ct::CTranspose) = reverse(size(ct.qarr))
Base.size(ct::CTranspose, i) = size(ct, revind(ndims(ct), i))
Base.size(ct::CTranspose, i) = size(ct.qarr, revind(ndims(ct), i))

Base.getindex(ct::CTranspose, i, j) = getindex(ct.qarr, j, i)'
Base.getindex(dv::DualVector, i) = getindex(dv.qarr, i)'

Base.in(c, ct::CTranspose) = in(c', ct.qarr)

Base.setindex!(ct::CTranspose, x, i, j) = setindex!(ct.qarr, x', j, i)
Base.setindex!(dv::DualVector, x, i) = setindex!(dv.qarr, x', i)

Base.ctranspose(qarr::QuArray) = CTranspose(qarr)
Base.ctranspose(ct::CTranspose) = ct.qarr

eager_ctranspose(qarr::AbstractQuArray) = similar_type(qarr)(coeffs(qarr)', reverse(bases(qarr)))

################
# LabelQuArray #
################
typealias LabelQuArray{B<:LabelBasis,T,N,A} QuArray{B,T,N,A}
typealias LabeledQuArray{B<:LabelBasis,T,N} AbstractQuArray{B,T,N}
typealias TupleArray{T<:Tuple,N} Array{T,N}

Base.getindex(larr::LabelQuArray, tups::Union(Tuple,TupleArray)...) = getindex(larr, map(getindex, bases(larr), tups)...)
Base.setindex!(larr::LabelQuArray, x, tups::Union(Tuple,TupleArray)...) = setindex!(larr, x, map(getindex, bases(larr), tups)...)
# redudant definitions added to resolve ambiguity warnings
Base.getindex{B<:LabelBasis}(larr::DualVector{B}, tups::Union(Tuple,TupleArray)) = getindex(larr, bases(larr,1)[tups])
Base.getindex{B<:LabelBasis}(larr::CTranspose{B}, tups::Union(Tuple,TupleArray)...) = getindex(larr, map(getindex, bases(larr), tups)...)
Base.getindex(larr::LabeledQuArray, tups::Union(Tuple,TupleArray)...) = getindex(larr, map(getindex, bases(larr), tups)...)

# redudant definitions added to resolve ambiguity warnings
Base.setindex!{B<:LabelBasis}(larr::DualVector{B}, x, tups::Union(Tuple,TupleArray)) = setindex(larr, x, bases(larr,1)[tups])
Base.setindex!{B<:LabelBasis}(larr::CTranspose{B}, x, tups::Union(Tuple,TupleArray)...) = setindex!(larr, x, map(getindex, bases(larr), tups)...)
Base.setindex!(larr::LabeledQuArray, x, tups::Union(Tuple,TupleArray)...) = setindex!(larr, x, map(getindex, bases(larr), tups)...)

###################
# Promotion rules #
###################
########################
# Conversion/Promotion #
########################
Base.promote_rule{B1,B2,T1,T2,N,A1,A2}(::Type{QuArray{B1,T1,N,A1}}, ::Type{QuArray{B2,T2,N,A2}}) = QuArray{promote_type(B1,B2),promote_type(T1,T2),N,promote_type(A1,A2)}
Base.promote_rule{C<:CTranspose,B,T,N,A}(::Type{C}, ::Type{QuArray{B,T,N,A}}) = promote_type(qarr_type(C), QuArray{B,T,N,A})
Base.promote_rule{DV<:DualVector, QV<:QuVector}(::Type{DV}, ::Type{QV}) = typejoin(DV,QV)

Base.promote_rule{B1,B2,T1,T2,N,A1,A2}(::Type{CTranspose{B1,T1,N,A1}}, ::Type{QuArray{B2,T2,N,A2}}) = QuArray{promote_type(B1,B2),promote_type(T1,T2),N,promote_type(A1,A2)}
Base.convert{B,T,N,A}(::Type{QuArray{B,T,N,A}}, qa::QuArray{B}) = QuArray(convert(A, coeffs(qa)), bases(qa))
Base.convert{B,T,N,A}(::Type{QuArray{B,T,N,A}}, qa::QuArray) = QuArray(convert(A, coeffs(qa)), map(i -> convert(B,i), bases(qa)))

Base.convert{C<:CTranspose,B,T,N,Q<:AbstractQuArray}(::Type{C}, ct::CTranspose{B,T,N,Q}) = CTranspose(convert(qarr_type(C),ct.qarr))
Base.convert{B,T,N,Q<:AbstractQuArray}(::Type{Q}, ct::CTranspose{B,T,N,Q}) = eager_ctranspose(ct.qarr)
Base.convert{QV<:AbstractQuVector, DV<:DualVector}(::Type{QV}, ::DV) = error("Cannot convert Bra to Ket. Use ctranspose if you would like to take the conjugate transpose.")
Base.convert{DV<:DualVector,B,T,QV<:AbstractQuVector}(::Type{DV}, ct::DualVector{B,T,QV}) = CTranspose(convert(qarr_type(DV),ct.qarr))

######################
# Printing Functions #
######################
Base.summary{B}(qarr::AbstractQuArray{B}) = "$(sizenotation(size(qarr))) $(typerepr(qarr)) in $B"
Base.summary(larr::LabelQuArray) = "$(sizenotation(size(larr))) $(typerepr(larr)) in labeled bases $(bases(larr))"
Base.summary(larr::LabeledQuArray) = "$(sizenotation(size(larr))) $(typerepr(larr)) in labeled bases $(bases(larr))"

function Base.show(io::IO, qarr::AbstractQuArray)
println(io, summary(qarr)*":")
println(io, "...coefficients: $(coefftype(qarr))")
print(io, repr(coeffs(qarr)))
end


####################
# Helper Functions #
####################
Expand Down
13 changes: 6 additions & 7 deletions src/bases/finitebasis.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,28 +18,27 @@
# This tensor product basis, F_4 x F_3 x F_2, could be represented
# by calling FiniteBasis(4, 3, 2).

immutable FiniteBasis{S,N} <: AbstractFiniteBasis{S}
lens::NTuple{N,Int}
FiniteBasis(lens::NTuple{N,Int}) = new(lens)
immutable FiniteBasis{S} <: AbstractFiniteBasis{S}
lens::@compat(Tuple{Vararg{Int}})
FiniteBasis(lens::@compat(Tuple{Vararg{Int}})) = new(lens)
FiniteBasis(lens::Int...) = new(lens)
end

FiniteBasis{N, S<:AbstractStructure}(lens::NTuple{N,Int}, ::Type{S}=Orthonormal) = FiniteBasis{S,N}(lens)
FiniteBasis{S<:AbstractStructure}(lens::@compat(Tuple{Vararg{Int}}), ::Type{S}=Orthonormal) = FiniteBasis{S}(lens)
FiniteBasis(lens::Int...) = FiniteBasis(lens)

Base.convert{A,B,N}(::Type{FiniteBasis{A,N}}, f::FiniteBasis{B,N}) = FiniteBasis(f.lens, A)
Base.convert{A,B}(::Type{FiniteBasis{A}}, f::FiniteBasis{B}) = FiniteBasis(f.lens, A)

######################
# Property Functions #
######################
structure{S}(::Type{FiniteBasis{S}}) = S
structure{S,N}(::Type{FiniteBasis{S,N}}) = S

Base.length(basis::FiniteBasis) = prod(basis.lens)
Base.size(basis::FiniteBasis) = basis.lens
Base.size(basis::FiniteBasis, i) = basis.lens[i]

nfactors{S,N}(::FiniteBasis{S,N}) = N
nfactors(basis::FiniteBasis) = length(basis.lens)
checkcoeffs(coeffs, dim::Int, basis::FiniteBasis) = size(coeffs, dim) == length(basis)

##########################
Expand Down
Loading