Merge pull request #33 from JuliaQuantum/general-ct

Generalized CTranspose + conversion methods and promotion rules. Closes #31

src/arrays/arraymath.jl

@@ -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
 
@@ -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")
@@ -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")
@@ -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))
@@ -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))
@@ -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)
 
 ###############
 # Commutators #

src/arrays/constructors.jl

@@ -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))

src/arrays/quarray.jl

@@ -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))
@@ -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 #
@@ -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)'
@@ -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 #
 ####################

src/bases/finitebasis.jl

@@ -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) 
 
     ##########################