optimize_measurment()

docs/src/BellScenario/quantum_opt.md

@@ -0,0 +1,8 @@
+```@meta
+CurrentModule = BellScenario
+```
+# Quantum Optimization
+
+```@docs
+optimize_measurement
+```

src/BellScenario.jl

@@ -1,6 +1,6 @@
 module BellScenario
 
-using QBase
+using QBase, LinearAlgebra
 
 import Base: *, convert
 
@@ -15,6 +15,7 @@ include("./games.jl")
 
 # quantum scenarios
 include("./quantum_strategies.jl")
+include("./quantum_opt.jl")
 
 # read/write and printing
 include("./file_io.jl")

src/quantum_opt.jl

@@ -0,0 +1,60 @@
+using Convex, SCS
+
+export optimize_measurement
+
+"""
+    optimize_measurement(
+        game :: BellGame,
+        ρ_states :: Vector{<:States.AbstractDensityMatrix},
+        PM :: PrepareAndMeasure
+    )
+
+Perform semi-definite programming to find the optimal quantum measurement which
+maximizes the score for the specified set of quantum  states `ρ_states` and `BellGame`.
+"""
+function optimize_measurement(
+    game::BellGame,
+    ρ_states::Vector{<:States.AbstractDensityMatrix},
+    PM::PrepareAndMeasure,
+) :: Dict
+    if PM.X != length(ρ_states)
+        throw(DomainError(PM, "expexted length of ρ_states is $(PM.d), but got $(length(ρ_states)) instead"))
+    end
+
+    if size(ρ_states[1]) != (PM.d,PM.d)
+        throw(DomainError(ρ_states, "dimension of ρ_states is not $(PM.d)"))
+    end
+
+    norm_game_vector = convert(Vector{Int64}, game)
+    norm_bound = norm_game_vector[end]
+    norm_game = reshape(norm_game_vector[1:(end-1)], (PM.B-1, PM.X))
+
+
+
+    # add povm variables and constraints
+    Π_vars = map(i -> HermitianSemidefinite(PM.d), 1:PM.B)
+    constraints = (sum(map(Π_b -> real(Π_b), Π_vars)) == Matrix{Float64}(I, PM.d, PM.d))
+    constraints += (sum(map(Π_b -> imag(Π_b), Π_vars)) == zeros(Float64, PM.d, PM.d))
+
+    # sum up the state contibutions for each row
+    H_b = map(row_id -> sum(norm_game[row_id,:] .* ρ_states), 1:PM.B-1)
+
+    # add the objective
+    objective = maximize(real(tr(sum(Π_vars[1:end-1] .* H_b))), constraints)
+
+    # optimize model
+    solve!(objective, SCS.Optimizer(verbose=0))
+
+    # parse/return results
+    score = objective.optval
+    violation = score - norm_bound
+    Π_opt = map(Π_b -> Π_b.value, Π_vars)
+
+    Dict(
+        "violation" => violation,
+        "povm" => Π_opt,
+        "game" => game,
+        "scenario" => PM,
+        "states" => ρ_states
+    )
+end

src/quantum_strategies.jl

@@ -36,9 +36,8 @@ function quantum_strategy(
     ρ_states :: Vector{<:States.AbstractDensityMatrix},
     PM :: PrepareAndMeasure
 ) :: Strategy
-    dits = PM.d
-    if (size(Π[1]) != (dits, dits)) || (size(ρ_states[1]) != (dits, dits))
-        throw(DomainError((Π, ρ_states), "POVM or States are not dimension `d=$dits`."))
+    if (size(Π[1]) != (PM.d, PM.d)) || (size(ρ_states[1]) != (PM.d, PM.d))
+        throw(DomainError((Π, ρ_states), "POVM or States are not dimension `d=$(PM.d)`."))
     end
     conditionals = measurement_probs(Π, ρ_states)
     Strategy(conditionals, PM)

test/unit/quantum_opt.jl

@@ -0,0 +1,65 @@
+using Test, QBase
+
+@testset "./src/quantum_opt.jl" begin
+
+using BellScenario
+
+@testset "optimize_measurement(PrepareAndMeasure)" begin
+    @testset "trine states" begin
+        PM = PrepareAndMeasure(3,3,2)
+        game = BellGame([1 0 0;0 1 0;0 0 1], 2)
+        ρ_states = States.trine_qubits
+
+        dict = optimize_measurement(game, ρ_states, PM)
+
+        @test isapprox(dict["violation"], 0.0, atol=1e-6)
+        @test all(isapprox.(dict["povm"][1], 2/3*ρ_states[1], atol=1e-6))
+        @test all(isapprox.(dict["povm"][2], 2/3*ρ_states[2], atol=1e-6))
+        @test all(isapprox.(dict["povm"][3], 2/3*ρ_states[3], atol=1e-6))
+
+        @test dict["states"] == ρ_states
+        @test dict["game"] == game
+        @test dict["scenario"] == PM
+    end
+
+    @testset "bb84 states" begin
+        PM = PrepareAndMeasure(4,4,2)
+        game = BellGame([1 0 0 0;0 1 0 0;0 0 1 0;0 0 0 1], 2)
+        ρ_states = States.bb84_qubits
+
+        dict = optimize_measurement(game, ρ_states, PM)
+
+        @test isapprox(dict["violation"], 0.0, atol=1e-6)
+        @test all(isapprox.(dict["povm"][1], 1/2*ρ_states[1], atol=1e-3))
+        @test all(isapprox.(dict["povm"][2], 1/2*ρ_states[2], atol=1e-3))
+        @test all(isapprox.(dict["povm"][3], 1/2*ρ_states[3], atol=1e-3))
+        @test all(isapprox.(dict["povm"][4], 1/2*ρ_states[4], atol=1e-3))
+    end
+
+    @testset "sic qubit states" begin
+        PM = PrepareAndMeasure(4,4,2)
+        game = BellGame([1 0 0 0;0 1 0 0;0 0 1 0;0 0 0 1], 2)
+        ρ_states = States.sic_qubits
+
+        dict = optimize_measurement(game, ρ_states, PM)
+
+        @test isapprox(dict["violation"], 0.0, atol=1e-5)
+        @test all(isapprox.(dict["povm"][1], 1/2*ρ_states[1], atol=1e-5))
+        @test all(isapprox.(dict["povm"][2], 1/2*ρ_states[2], atol=1e-5))
+        @test all(isapprox.(dict["povm"][3], 1/2*ρ_states[3], atol=1e-5))
+        @test all(isapprox.(dict["povm"][4], 1/2*ρ_states[4], atol=1e-5))
+    end
+
+    @testset "Errors" begin
+        PM = PrepareAndMeasure(3,3,2)
+        states = States.bb84_qubits
+        game = BellGame([1 0 0 0;0 1 0 0;0 0 1 0;0 0 0 1], 2)
+
+        @test_throws DomainError optimize_measurement(game, states, PM)
+
+        PM = PrepareAndMeasure(4,4,3)
+        @test_throws DomainError  optimize_measurement(game, states, PM)
+    end
+end
+
+end