created coherentstate! for inplace operations

src/QuantumOptics.jl

@@ -16,7 +16,7 @@ export bases, Basis, GenericBasis, CompositeBasis, basis,
         superoperators, SuperOperator, DenseSuperOperator, SparseSuperOperator,
                 spre, spost, liouvillian,
         fock, FockBasis, number, destroy, create,
-                fockstate, coherentstate, displace,
+                fockstate, coherentstate, coherentstate!, displace,
         randstate, randoperator, thermalstate, coherentthermalstate, phase_average, passive_state,
         spin, SpinBasis, sigmax, sigmay, sigmaz, sigmap, sigmam, spinup, spindown,
         subspace, SubspaceBasis, projector,

src/fock.jl

@@ -1,6 +1,6 @@
 module fock
 
-export FockBasis, number, destroy, create, displace, fockstate, coherentstate
+export FockBasis, number, destroy, create, displace, fockstate, coherentstate, coherentstate!
 
 import Base: ==
 
@@ -83,14 +83,24 @@ end
 
 Coherent state ``|α⟩`` for the specified Fock space.
 """
-function coherentstate(b::FockBasis, alpha::Number, result=Ket(b, Vector{ComplexF64}(undef, length(b))))
+function coherentstate(b::FockBasis, alpha::Number)
+    result = Ket(b, Vector{ComplexF64}(undef, length(b)))
+    coherentstate!(result, b, alpha)
+    return result
+end
+
+"""
+    coherentstate!(ket::Ket, b::FockBasis, alpha)
+
+Inplace creation of coherent state ``|α⟩`` for the specified Fock space.
+"""
+function coherentstate!(ket::Ket, b::FockBasis, alpha::Number)
     alpha = complex(alpha)
-    data = result.data
+    data = ket.data
     data[1] = exp(-abs2(alpha)/2)
     @inbounds for n=1:b.N
         data[n+1] = data[n]*alpha/sqrt(n)
     end
-    return result
 end
 
 end # module

src/phasespace.jl

@@ -76,7 +76,7 @@ function qfunc(state::Union{Ket{B}, AbstractOperator{B,B}}, x::Number, y::Number
 end
 
 function _qfunc_operator(rho::AbstractOperator{B,B}, alpha::ComplexF64, tmp1::Ket, tmp2::Ket) where B<:FockBasis
-    coherentstate(basis(rho), alpha, tmp1)
+    coherentstate!(tmp1, basis(rho), alpha)
     operators.gemv!(complex(1.), rho, tmp1, complex(0.), tmp2)
     a = dot(tmp1.data, tmp2.data)
     return a/pi