From fb60ee54efdf21332ba244f9fcca44499dfbf76a Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Designolle?=
 <sebastien.designolle@gmail.com>
Date: Thu, 16 Jan 2025 10:35:42 +0100
Subject: [PATCH] Add Dicke states

---
 TODO            |  1 -
 docs/src/api.md |  2 ++
 src/states.jl   | 38 +++++++++++++++++++++++++++++++++++++-
 test/states.jl  |  3 +++
 4 files changed, 42 insertions(+), 2 deletions(-)

diff --git a/TODO b/TODO
index 274f4f5..8051d91 100644
--- a/TODO
+++ b/TODO
@@ -8,7 +8,6 @@ SD - Properly add measurements/states from https://gitlab.com/tqo/sqma/ (add gen
 CG - In _dps_constraints! bypass new variables when n = 1, as then only PPT is needed (minor technical issue)
 CG - GME detection with DPS
 CG - Basic state discrimination stuff
-CG - state_dicke
 LP - seesaw with fixed state
 MA - signalling and non-signalling bounds for Bell inequalities
 MA - canonical representation of Bell inequalities
diff --git a/docs/src/api.md b/docs/src/api.md
index 7b2f3ca..d1ca4d4 100644
--- a/docs/src/api.md
+++ b/docs/src/api.md
@@ -117,6 +117,8 @@ state_ghz_ket
 state_ghz
 state_w_ket
 state_w
+state_dicke_ket
+state_dicke
 white_noise
 white_noise!
 ```
diff --git a/src/states.jl b/src/states.jl
index 3e1046b..cbef786 100644
--- a/src/states.jl
+++ b/src/states.jl
@@ -204,7 +204,7 @@ state_supersinglet_ket(N::Integer = 3; kwargs...) = state_supersinglet_ket(Compl
 export state_supersinglet_ket
 
 """
-    state_supersinglet([T=ComplexF64,] N::Integer = 3; v::Real = 1, coeff = 1/√d)
+    state_supersinglet([T=ComplexF64,] N::Integer = 3; v::Real = 1)
 
 Produces the `N`-partite `N`-level singlet state with visibility `v`.
 This state is invariant under simultaneous rotations on all parties: `(U⊗…⊗U)ρ(U⊗…⊗U)'=ρ`.
@@ -218,3 +218,39 @@ function state_supersinglet(::Type{T}, N::Integer = 3; v::Real = 1) where {T<:Nu
 end
 state_supersinglet(N::Integer = 3; kwargs...) = state_supersinglet(ComplexF64, N; kwargs...)
 export state_supersinglet
+
+"""
+    state_dicke_ket([T=ComplexF64,] k::Integer, N::Integer; coeff = 1/√Cᴺₖ)
+
+Produces the ket of the `N`-partite Dicke state with `k` excitations.
+
+Reference: Robert H. Dicke [doi:10.1103/PhysRev.93.99](https://doi.org/10.1103/PhysRev.93.99)
+"""
+function state_dicke_ket(::Type{T}, k::Integer, N::Integer; coeff = inv(_sqrt(T, binomial(N, k)))) where {T<:Number}
+    psi = zeros(T, 2^N)
+    ind = zeros(Int8, N)
+    @inbounds for i in eachindex(psi)
+        if sum(ind) == k
+            psi[i] = coeff
+        end
+        _update_odometer!(ind, 2)
+    end
+    return psi
+end
+state_dicke_ket(k::Integer, N::Integer; kwargs...) = state_dicke_ket(ComplexF64, k, N; kwargs...)
+export state_dicke_ket
+
+"""
+    state_dicke([T=ComplexF64,] k::Integer, N::Integer; v::Real = 1)
+
+Produces the `N`-partite Dicke state with `k` excitations.
+
+Reference: Robert H. Dicke [doi:10.1103/PhysRev.93.99](https://doi.org/10.1103/PhysRev.93.99)
+"""
+function state_dicke(::Type{T}, k::Integer, N::Integer; v::Real = 1) where {T<:Number}
+    rho = ketbra(state_dicke_ket(T, k, N; coeff = one(T)))
+    parent(rho) ./= binomial(N, k)
+    return white_noise!(rho, v)
+end
+state_dicke(k::Integer, N::Integer; kwargs...) = state_dicke(ComplexF64, k, N; kwargs...)
+export state_dicke
diff --git a/test/states.jl b/test/states.jl
index fc4892d..eddf2f9 100644
--- a/test/states.jl
+++ b/test/states.jl
@@ -12,6 +12,9 @@
         @test ketbra(ψ) ≈ state_w(T)
         coeff = T.([1, 2, 2]) / 9
         @test state_w_ket(T; coeff) == T.(ket(2, 8) + 2 * ket(3, 8) + 2 * ket(5, 8)) / 9
+        ψ = state_dicke_ket(T, 2, 3)
+        @test ψ == inv(sqrt(R(3))) * (ket(4, 8) + ket(6, 8) + ket(7, 8))
+        @test ketbra(ψ) ≈ state_dicke(T, 2, 3)
         ψ = state_ghz_ket(T)
         @test ψ == inv(sqrt(R(2))) * (ket(1, 8) + ket(8, 8))
         @test ketbra(ψ) ≈ state_ghz(T)