Merge pull request #25 from amitjamadagni/ops

Add operators: squeezingop, displaceop + tests. Add function coherentstatevec + test. Add functions for QuArray: expm, full. Fix issue with multiplication of two DualMatrices.

src/arrays/arraymath.jl

@@ -59,7 +59,7 @@ function *{B<:OrthonormalBasis}(qm1::QuMatrix{B}, qm2::QuMatrix{B})
                    bases(qm2,2))
 end
 
-*(dm1::DualMatrix, dm2::DualMatrix) = (dm1.qarr*dm2.qarr)'
+*(dm1::DualMatrix, dm2::DualMatrix) = (dm2.qarr*dm1.qarr)'
 
 function +{B<:OrthonormalBasis,N}(qarr1::AbstractQuArray{B,N}, qarr2::AbstractQuArray{B,N})
     if bases(qarr1) == bases(qarr2)
@@ -90,6 +90,12 @@ Base.scale(qarr::Union(QuArray,CTranspose), num::Number) = scale(num, qarr)
 *(qarr::Union(QuArray,CTranspose), num::Number) = scale(qarr, num)
 /(qarr::Union(QuArray,CTranspose), num::Number) = scale(1/num, qarr)
 
+# sparse to dense
+Base.full(qarr::AbstractQuMatrix) = QuArray(full(coeffs(qarr)),bases(qarr))
+
+# exponential of dense matrix
+Base.expm(qarr::AbstractQuMatrix) = QuArray(expm(full(coeffs(qarr))),bases(qarr))
+
 # normalization
 Base.norm(qarr::Union(QuArray,CTranspose)) = vecnorm(rawcoeffs(qarr))
 

src/arrays/constructors.jl

@@ -16,4 +16,21 @@
     one_at_ind!(arr, i) = setindex!(arr, one(eltype(arr)), i)
     single_coeff(i, lens...) = one_at_ind!(zeros(lens), i)
 
-export statevec
+# Reference :
+# Section - 1.1, http://cds.cern.ch/record/331607/files/9708012.pdf
+function coherentstatevec_inf(n::Int,alpha::Number)
+    s = zeros(typeof(float(alpha)),n)
+    s[1] = one(alpha)
+    for i in [2:n]
+        s[i] = alpha/sqrt(i-1)*s[i-1]
+    end
+    z = QuArray(s)
+    return scale!(exp(-abs2(alpha)/2),z)
+end
+
+# Reference :
+# 2nd Defintion : http://en.wikipedia.org/wiki/Coherent_states#Mathematical_features_of_the_canonical_coherent_states
+coherentstatevec(n::Int, alpha::Number) = displaceop(n,alpha)*statevec(1,FiniteBasis(n))
+
+export statevec,
+    coherentstatevec

src/arrays/ladderops.jl

@@ -37,7 +37,15 @@ function momentumop(n::Int)
     return scale(im/sqrt(2.), cop-cop')
 end
 
+squeeze_construct(a::AbstractQuMatrix,b::AbstractQuMatrix, z::Number) = scale!(0.5,(scale!(z',a*b)-scale!(z,a'*b')))
+squeezingop(a::AbstractQuMatrix, b::AbstractQuMatrix, z::Number) = expm(squeeze_construct(a,b,z))
+squeezingop(n::Int,z::Number) = squeezingop(lowerop(n),lowerop(n),z)
+
+displaceop(n::Int,alpha::Number) = expm(scale(alpha,lowerop(n)')-scale(alpha',lowerop(n)))
+
 export raiseop,
     lowerop,
     positionop,
-    momentumop
+    momentumop,
+    squeezingop,
+    displaceop

test/multest.jl

@@ -30,4 +30,4 @@ qv = QuArray(v)
 @test_approx_eq rawcoeffs(normalize(qv)) rawcoeffs(qv)/norm(v)
 
 # a simple test of the `==` operator
-@assert qm*3im == scale!(3im, copy(qm))
+@assert qm*3im == scale!(3im, copy(qm))

test/operatortest.jl

@@ -10,11 +10,22 @@
 
 @assert coeffs(commutator(sigma_x, sigma_x)) == spzeros(2,2)
 
-######################################
-# Position & Momentum Operators Test #
-######################################
+####################################################
+# Position, Displacement & Momentum Operators Test #
+####################################################
 
 p = positionop(2)
 m = momentumop(2)
 @assert coeffs(commutator(sigma_x, p)) == spzeros(2,2)
 @assert coeffs(commutator(sigma_y, m)) == spzeros(2,2)
+
+coherentstate_inf = QuBase.coherentstatevec_inf(20,1)
+coherentstate = displaceop(20,1)*statevec(1,FiniteBasis(20))
+@test_approx_eq_eps coeffs(coherentstate_inf)[1] coeffs(coherentstate)[1] 1e-8
+
+############################
+# Squeezing Operators Test #
+############################
+
+@assert squeezingop(2,1.0)== QuArray(eye(2))
+@test_approx_eq coeffs(squeezingop(lowerop(2), QuArray(eye(2)), 2.0)') coeffs(displaceop(2,1.0))