From 3f71173138017b2c186ccf6333b361e1cabd6683 Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Sun, 1 Mar 2015 15:50:16 -0500 Subject: [PATCH 1/8] rm'd Conjugate/Transpose types, moved remaining CTranspose type into quarray.jl --- src/arrays/dual.jl | 88 ------------------------------------------- src/arrays/quarray.jl | 62 +++++++++++++++++++++++++----- 2 files changed, 53 insertions(+), 97 deletions(-) delete mode 100644 src/arrays/dual.jl diff --git a/src/arrays/dual.jl b/src/arrays/dual.jl deleted file mode 100644 index 4feb269..0000000 --- a/src/arrays/dual.jl +++ /dev/null @@ -1,88 +0,0 @@ -##################################### -# Conjugate/Transpose Wrapper types # -##################################### - abstract DualWrapper{B,T,N,A} <: AbstractQuArray{B,T,N} - - immutable Transpose{B,T,N,A} <: DualWrapper{B,T,N,A} - data::QuArray{B,T,N,A} - Transpose(data::QuVector{B,T,A}) = new(data) - Transpose(data::QuMatrix{B,T,A}) = new(data) - Transpose(data::QuArray) = error("Transposition is unsupported for QuArrays of dimension $N") - end - - immutable Conjugate{B,T,N,A} <: DualWrapper{B,T,N,A} - data::QuArray{B,T,N,A} - end - - immutable CTranspose{B,T,N,A} <: DualWrapper{B,T,N,A} - data::QuArray{B,T,N,A} - CTranspose(data::QuVector{B,T,A}) = new(data) - CTranspose(data::QuMatrix{B,T,A}) = new(data) - CTranspose(data::QuArray) = error("Conjugate transpose is unsupported for QuArrays of dimension $N") - end - - Transpose{B,T,N,A}(qa::QuArray{B,T,N,A}) = Transpose{B,T,N,A}(qa) - Conjugate{B,T,N,A}(qa::QuArray{B,T,N,A}) = Conjugate{B,T,N,A}(qa) - CTranspose{B,T,N,A}(qa::QuArray{B,T,N,A}) = CTranspose{B,T,N,A}(qa) - - typealias TranVector{B,T,A} Transpose{B,T,1,A} - typealias TranMatrix{B,T,A} Transpose{B,T,2,A} - - typealias ConjVector{B,T,A} Conjugate{B,T,1,A} - typealias ConjMatrix{B,T,A} Conjugate{B,T,2,A} - - typealias CTranVector{B,T,A} CTranspose{B,T,1,A} - typealias CTranMatrix{B,T,A} CTranspose{B,T,2,A} - - ###################### - # Accessor functions # - ###################### - data(tarr::Transpose) = tarr.data - data(carr::Conjugate) = carr.data - data(ctarr::CTranspose) = ctarr.data - - rawcoeffs(dw::DualWrapper) = rawcoeffs(data(dw)) - coefftype{B,T,N,A}(::DualWrapper{B,T,N,A}) = A - - revind(len, i) = len - (i-1) - bases(tmat::Union(CTranMatrix, TranMatrix), i) = bases(data(tmat), revind(ndims(tmat), i)) - - isconj(::Union(Conjugate, CTranspose)) = true - isconj(x) = false - - istran(::Union(Transpose, CTranspose)) = true - istran(x) = false - - ######################## - # Array-like functions # - ######################## - Base.size(dw::DualWrapper, i...) = size(data(dw), i...) - Base.size(tmat::Union(TranMatrix,CTranMatrix)) = reverse(size(data(tmat))) - Base.size(tmat::Union(TranMatrix,CTranMatrix), i) = size(tmat, revind(ndims(tmat), i)) - - Base.ndims(dw::DualWrapper) = ndims(data(dw)) - Base.length(dw::DualWrapper) = length(data(dw)) - - Base.getindex(tarr::Transpose, i, j) = getindex(data(tarr), j, i).' - Base.getindex(tvec::TranVector, i) = getindex(data(tvec), i).' - Base.getindex(ctarr::CTranspose, i, j) = getindex(data(ctarr), j, i)' - Base.getindex(ctvec::CTranVector, i) = getindex(data(ctvec), i)' - Base.getindex(carr::Conjugate, i...) = conj(getindex(data(carr), i...)) - - Base.setindex!(tarr::Transpose, x, i, j) = setindex!(data(tarr), x.', j, i) - Base.setindex!(tvec::TranVector, x, i) = setindex!(data(tvec), x.', i) - Base.setindex!(ctarr::CTranspose, x, i, j) = setindex!(data(ctarr), x', j, i) - Base.setindex!(ctvec::CTranVector, x, i) = setindex!(data(ctvec), x', i) - Base.setindex!(carr::Conjugate, x, i...) = setindex!(data(carr), conj(x), i...) - - Base.conj(qarr::QuArray) = Conjugate(qarr) - Base.conj(ctarr::CTranspose) = Transpose(data(ctarr)) - Base.conj(tarr::Transpose) = CTranspose(data(tarr)) - Base.conj(carr::Conjugate) = data(carr) - - Base.transpose(qarr::QuArray) = Transpose(qarr) - Base.transpose(ctarr::CTranspose) = Conjugate(data(ctarr)) - Base.transpose(tarr::Transpose) = data(tarr) - Base.transpose(carr::Conjugate) = CTranspose(data(carr)) - - Base.ctranspose(qarr::AbstractQuArray) = transpose(conj(qarr)) diff --git a/src/arrays/quarray.jl b/src/arrays/quarray.jl index 73405fd..9bf4478 100644 --- a/src/arrays/quarray.jl +++ b/src/arrays/quarray.jl @@ -46,24 +46,68 @@ Base.in(c, qarr::QuArray) = in(c, rawcoeffs(qarr)) +############## +# CTranspose # +############## + immutable CTranspose{B,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") + end + + CTranspose{B,T,N,A}(qa::QuArray{B,T,N,A}) = CTranspose{B,T,N,A}(qa) + + typealias DualVector{B,T,A} CTranspose{B,T,1,A} + typealias DualMatrix{B,T,A} CTranspose{B,T,2,A} + ###################### - # Printing Functions # + # Accessor functions # ###################### - Base.summary{B}(qarr::AbstractQuArray{B}) = "$(sizenotation(size(qarr))) $(typenotation(qarr)) in $B" + rawcoeffs(ct::CTranspose) = rawcoeffs(ct.qarr) + coeffs(ct::CTranspose) = rawcoeffs(ct)' + coefftype{B,T,N,A}(::CTranspose{B,T,N,A}) = A + + revind(len, i) = len - (i-1) + bases(ct::CTranspose, i) = bases(ct.qarr, revind(ndims(ct), i)) + + ######################## + # Array-like functions # + ######################## + 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.getindex(ct::CTranspose, i, j) = getindex(ct.qarr, j, i)' + Base.getindex(dv::DualVector, i) = getindex(dv.qarr, i)' + + 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 + +###################### +# Printing Functions # +###################### + Base.summary{B}(qarr::AbstractQuArray{B}) = "$(sizenotation(size(qarr))) $(typerepr(qarr)) in $B" function Base.show(io::IO, qarr::AbstractQuArray) println(io, summary(qarr)*":") - println(io, "...transposed/conjugated?: $(istran(qarr))/$(isconj(qarr))") println(io, "...original coefficients: $(coefftype(qarr))") print(io, repr(rawcoeffs(qarr))) end - #################### - # Helper Functions # - #################### - typenotation(::AbstractQuVector) = "QuVector" - typenotation(::AbstractQuMatrix) = "QuMatrix" - typenotation(::AbstractQuArray) = "QuArray" +#################### +# Helper Functions # +#################### + typerepr(::QuVector) = "QuVector" + typerepr(::QuMatrix) = "QuMatrix" + typerepr(::DualVector) = "DualVector" + typerepr(::DualMatrix) = "DualMatrix" + typerepr(::QuArray) = "QuArray" sizenotation(tup::(Int,)) = "$(first(tup))-element" sizenotation(tup::(Int...)) = reduce(*, map(s->"$(s)x", tup))[1:end-1] From 52b0fb8dcc5155c00164c45d030797a30816c82a Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Sun, 1 Mar 2015 16:21:56 -0500 Subject: [PATCH 2/8] rm'd overly complicated default basis generating function --- src/arrays/quarray.jl | 8 +------- 1 file changed, 1 insertion(+), 7 deletions(-) diff --git a/src/arrays/quarray.jl b/src/arrays/quarray.jl index 9bf4478..9edbd07 100644 --- a/src/arrays/quarray.jl +++ b/src/arrays/quarray.jl @@ -20,7 +20,7 @@ QuArray{B<:AbstractBasis,T,N}(coeffs::AbstractArray{T,N}, bases::NTuple{N,B}) = QuArray{B,T,N,typeof(coeffs)}(coeffs, bases) QuArray(coeffs, bases::AbstractBasis...) = QuArray(coeffs, bases) - QuArray(coeffs) = QuArray(coeffs, basesfordims(size(coeffs))) + QuArray(coeffs) = QuArray(coeffs, map(FiniteBasis, size(coeffs))) typealias QuVector{B<:AbstractBasis,T,A} QuArray{B,T,1,A} typealias QuMatrix{B<:AbstractBasis,T,A} QuArray{B,T,2,A} @@ -119,10 +119,4 @@ return false end - # Assumes that every basis type passed in - # has a constructor B(::eltype(lens)) - function basesfordims(lens::Tuple, B=ntuple(length(lens), x->FiniteBasis)) - return ntuple(length(lens), n->B[n](lens[n])) - end - export QuArray From bd64c3e204a99628b73a257079f54d186f7534cd Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Sun, 1 Mar 2015 16:25:50 -0500 Subject: [PATCH 3/8] defined quarray multiplication methods --- src/arrays/arraymath.jl | 55 +++++++++++++++++++++++++++++++++++++++++ src/arrays/arrays.jl | 2 +- 2 files changed, 56 insertions(+), 1 deletion(-) create mode 100644 src/arrays/arraymath.jl diff --git a/src/arrays/arraymath.jl b/src/arrays/arraymath.jl new file mode 100644 index 0000000..1458532 --- /dev/null +++ b/src/arrays/arraymath.jl @@ -0,0 +1,55 @@ +import Base: * + +################## +# Multiplication # +################## +# bra * ket -> scalar +function *{B}(bra::DualVector{B}, ket::QuVector{B}) + return first(Ac_mul_B(rawcoeffs(bra), rawcoeffs(ket))) +end + +# bra * operator -> bra +function *{B}(bra::DualVector{B}, op::QuMatrix{B}) + return QuArray(Ac_mul_B(rawcoeffs(op), rawcoeffs(bra)), bases(op,2))' +end + +function *{B}(bra::DualVector{B}, op::DualMatrix{B}) + return QuArray(rawcoeffs(op)*rawcoeffs(bra), bases(op,2))' +end + +# operator * ket -> ket +function *{B}(op::QuMatrix{B}, ket::QuVector{B}) + return QuArray(rawcoeffs(op)*rawcoeffs(ket), bases(op,1)) +end + +function *{B}(op::DualMatrix{B}, ket::QuVector{B}) + return QuArray(Ac_mul_B(rawcoeffs(op), rawcoeffs(ket)), bases(op,1)) +end + +# ket * bra -> operator +function *{B}(ket::QuVector{B}, bra::DualVector{B}) + return QuArray(A_mul_Bc(rawcoeffs(ket), rawcoeffs(bra)), + bases(ket,1), + bases(bra,1)) +end + +# operator * operator -> operator +function *{B}(qm::QuMatrix{B}, dm::DualMatrix{B}) + return QuArray(A_mul_Bc(rawcoeffs(qm), rawcoeffs(dm)), + bases(qm,1), + bases(dm,2)) +end + +function *{B}(dm::DualMatrix{B}, qm::QuMatrix{B}) + return QuArray(Ac_mul_B(rawcoeffs(dm), rawcoeffs(qm)), + bases(dm,1), + bases(qm,2)) +end + +function *{B}(qm1::QuMatrix{B}, qm2::QuMatrix{B}) + return QuArray(rawcoeffs(qm1)*rawcoeffs(qm2), + bases(qm1,1), + bases(qm2,2)) +end + +*(dm1::DualMatrix, dm2::DualMatrix) = (dm1.qarr*dm2.qarr)' diff --git a/src/arrays/arrays.jl b/src/arrays/arrays.jl index a96f40d..d568143 100644 --- a/src/arrays/arrays.jl +++ b/src/arrays/arrays.jl @@ -2,7 +2,7 @@ # Include Statements # ###################### include("quarray.jl") - include("dual.jl") + include("arraymath.jl") include("constructors.jl") include("ladderops.jl") \ No newline at end of file From 1afc417c8c8f8ccd09f2e109b3721badd4ea7255 Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Sun, 1 Mar 2015 16:26:43 -0500 Subject: [PATCH 4/8] added quarray multiplication tests --- test/multest.jl | 23 +++++++++++++++++++++++ test/runtests.jl | 2 ++ 2 files changed, 25 insertions(+) create mode 100644 test/multest.jl diff --git a/test/multest.jl b/test/multest.jl new file mode 100644 index 0000000..558cf90 --- /dev/null +++ b/test/multest.jl @@ -0,0 +1,23 @@ +m = [1+1im 2+2im 3+3im; + 4+4im 5+5im 6+6im; + 7+7im 8+8im 9+9im] + +v = [1.0,2.0,3.0] + +qm = QuArray(m) +qv = QuArray(v) + +@assert rawcoeffs(qm*qm) == m*m +@assert rawcoeffs(qm'*qm') == m*m +@assert rawcoeffs(qm'*qm) == m'*m +@assert rawcoeffs(qm*qm') == m*m' + +@assert rawcoeffs(qm*qv) == m*v +@assert rawcoeffs(qv'*qm') == m*v +@assert rawcoeffs(qm'*qv) == m'*v +@assert rawcoeffs(qv'*qm) == m'*v + +@assert qv'*qv == first(v'*v) +@assert rawcoeffs(qv*qv') == v*v' + +@assert qv'*qm*qv == first(v'*m*v) \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index 4479748..7c227f6 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,3 +1,5 @@ using QuBase +using QuBase.rawcoeffs using Base.Test +include("multest.jl") From f52b14bd93a557f749da6a77a4a8c5a611e3906c Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Sun, 1 Mar 2015 16:39:02 -0500 Subject: [PATCH 5/8] removed lingering coeffs() definition --- src/arrays/quarray.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/src/arrays/quarray.jl b/src/arrays/quarray.jl index 9edbd07..99c41d9 100644 --- a/src/arrays/quarray.jl +++ b/src/arrays/quarray.jl @@ -65,7 +65,6 @@ # Accessor functions # ###################### rawcoeffs(ct::CTranspose) = rawcoeffs(ct.qarr) - coeffs(ct::CTranspose) = rawcoeffs(ct)' coefftype{B,T,N,A}(::CTranspose{B,T,N,A}) = A revind(len, i) = len - (i-1) From 268e448257b5dd7ef7a1dfd933d7ebec73a3bac0 Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Mon, 2 Mar 2015 11:31:02 -0500 Subject: [PATCH 6/8] use dot for inner product of vectors --- src/arrays/arraymath.jl | 2 +- test/multest.jl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/arrays/arraymath.jl b/src/arrays/arraymath.jl index 1458532..117e664 100644 --- a/src/arrays/arraymath.jl +++ b/src/arrays/arraymath.jl @@ -5,7 +5,7 @@ import Base: * ################## # bra * ket -> scalar function *{B}(bra::DualVector{B}, ket::QuVector{B}) - return first(Ac_mul_B(rawcoeffs(bra), rawcoeffs(ket))) + return dot(rawcoeffs(bra), rawcoeffs(ket)) end # bra * operator -> bra diff --git a/test/multest.jl b/test/multest.jl index 558cf90..61735c1 100644 --- a/test/multest.jl +++ b/test/multest.jl @@ -17,7 +17,7 @@ qv = QuArray(v) @assert rawcoeffs(qm'*qv) == m'*v @assert rawcoeffs(qv'*qm) == m'*v -@assert qv'*qv == first(v'*v) +@assert qv'*qv == dot(v,v) @assert rawcoeffs(qv*qv') == v*v' @assert qv'*qm*qv == first(v'*m*v) \ No newline at end of file From c1224d36c0235881d4832b2f68fd9c1503c3c9a2 Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Mon, 2 Mar 2015 11:54:03 -0500 Subject: [PATCH 7/8] restricted basis type for specific functions --- src/arrays/arraymath.jl | 25 ++++++++++++++++--------- 1 file changed, 16 insertions(+), 9 deletions(-) diff --git a/src/arrays/arraymath.jl b/src/arrays/arraymath.jl index 117e664..bff1106 100644 --- a/src/arrays/arraymath.jl +++ b/src/arrays/arraymath.jl @@ -3,50 +3,57 @@ import Base: * ################## # Multiplication # ################## + +# The below needs to be defined more generically; +# the below only defines multiplication on orthonormal, +# N<3 objects + +typealias OrthonormalBasis{S<:OrthonormalBasis} AbstractBasis{S} + # bra * ket -> scalar -function *{B}(bra::DualVector{B}, ket::QuVector{B}) +function *{B<:OrthonormalBasis}(bra::DualVector{B}, ket::QuVector{B}) return dot(rawcoeffs(bra), rawcoeffs(ket)) end # bra * operator -> bra -function *{B}(bra::DualVector{B}, op::QuMatrix{B}) +function *{B<:OrthonormalBasis}(bra::DualVector{B}, op::QuMatrix{B}) return QuArray(Ac_mul_B(rawcoeffs(op), rawcoeffs(bra)), bases(op,2))' end -function *{B}(bra::DualVector{B}, op::DualMatrix{B}) +function *{B<:OrthonormalBasis}(bra::DualVector{B}, op::DualMatrix{B}) return QuArray(rawcoeffs(op)*rawcoeffs(bra), bases(op,2))' end # operator * ket -> ket -function *{B}(op::QuMatrix{B}, ket::QuVector{B}) +function *{B<:OrthonormalBasis}(op::QuMatrix{B}, ket::QuVector{B}) return QuArray(rawcoeffs(op)*rawcoeffs(ket), bases(op,1)) end -function *{B}(op::DualMatrix{B}, ket::QuVector{B}) +function *{B<:OrthonormalBasis}(op::DualMatrix{B}, ket::QuVector{B}) return QuArray(Ac_mul_B(rawcoeffs(op), rawcoeffs(ket)), bases(op,1)) end # ket * bra -> operator -function *{B}(ket::QuVector{B}, bra::DualVector{B}) +function *{B<:OrthonormalBasis}(ket::QuVector{B}, bra::DualVector{B}) return QuArray(A_mul_Bc(rawcoeffs(ket), rawcoeffs(bra)), bases(ket,1), bases(bra,1)) end # operator * operator -> operator -function *{B}(qm::QuMatrix{B}, dm::DualMatrix{B}) +function *{B<:OrthonormalBasis}(qm::QuMatrix{B}, dm::DualMatrix{B}) return QuArray(A_mul_Bc(rawcoeffs(qm), rawcoeffs(dm)), bases(qm,1), bases(dm,2)) end -function *{B}(dm::DualMatrix{B}, qm::QuMatrix{B}) +function *{B<:OrthonormalBasis}(dm::DualMatrix{B}, qm::QuMatrix{B}) return QuArray(Ac_mul_B(rawcoeffs(dm), rawcoeffs(qm)), bases(dm,1), bases(qm,2)) end -function *{B}(qm1::QuMatrix{B}, qm2::QuMatrix{B}) +function *{B<:OrthonormalBasis}(qm1::QuMatrix{B}, qm2::QuMatrix{B}) return QuArray(rawcoeffs(qm1)*rawcoeffs(qm2), bases(qm1,1), bases(qm2,2)) From aacf45509a37f458aba7adaade36c4d88674129f Mon Sep 17 00:00:00 2001 From: Jarrett Revels Date: Mon, 2 Mar 2015 11:56:16 -0500 Subject: [PATCH 8/8] fixed typo/set default structure to orthonormal --- src/arrays/arraymath.jl | 2 +- src/bases/finitebasis.jl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/arrays/arraymath.jl b/src/arrays/arraymath.jl index bff1106..0bd5413 100644 --- a/src/arrays/arraymath.jl +++ b/src/arrays/arraymath.jl @@ -8,7 +8,7 @@ import Base: * # the below only defines multiplication on orthonormal, # N<3 objects -typealias OrthonormalBasis{S<:OrthonormalBasis} AbstractBasis{S} +typealias OrthonormalBasis{S<:Orthonormal} AbstractBasis{S} # bra * ket -> scalar function *{B<:OrthonormalBasis}(bra::DualVector{B}, ket::QuVector{B}) diff --git a/src/bases/finitebasis.jl b/src/bases/finitebasis.jl index 344f059..1196929 100644 --- a/src/bases/finitebasis.jl +++ b/src/bases/finitebasis.jl @@ -24,7 +24,7 @@ FiniteBasis(lens::Int...) = new(lens) end - FiniteBasis{N, S<:AbstractStructure}(lens::NTuple{N,Int}, ::Type{S}=AbstractStructure) = FiniteBasis{S,N}(lens) + FiniteBasis{N, S<:AbstractStructure}(lens::NTuple{N,Int}, ::Type{S}=Orthonormal) = FiniteBasis{S,N}(lens) FiniteBasis(lens::Int...) = FiniteBasis(lens) Base.convert{S,N}(::Type{FiniteBasis{S,N}}, f::FiniteBasis) = FiniteBasis(f.lens, S)