Example: Standard POVM Tomography on 1-Qubit with Cvxpy
[1]:
import numpy as np
from quara.objects.qoperation_typical import (
generate_qoperation,
generate_qoperation_object,
)
from quara.objects.composite_system_typical import generate_composite_system
from quara.protocol.qtomography.standard.standard_povmt import (
StandardPovmt,
)
[2]:
# System
num_qubits = 1
c_sys = generate_composite_system(mode="qubit", num=num_qubits)
[3]:
# Tester Objects
state_names = ["x0", "y0", "z0", "z1"]
tester_states = [
generate_qoperation_object(
mode="state", object_name="state", name=name, c_sys=c_sys
)
for name in state_names
]
[4]:
# True Object
true_object_name = "z"
true_object = generate_qoperation(
mode="povm", name=true_object_name, c_sys=c_sys
)
[5]:
# Povmt
on_para_eq_constraint = True
povmt = StandardPovmt(
states=tester_states,
num_outcomes=true_object.num_outcomes,
on_para_eq_constraint=on_para_eq_constraint,
schedules="all",
)
[6]:
# empi_dists
prob_dists = povmt.calc_prob_dists(true_object)
empi_dists = [(10, prob_dist) for prob_dist in prob_dists]
Povmt with Cvxpy
[7]:
from quara.interface.cvxpy.qtomography.standard.loss_function import (
CvxpyLossFunctionOption,
CvxpyRelativeEntropy,
CvxpyUniformSquaredError,
CvxpyApproximateRelativeEntropyWithZeroProbabilityTerm,
)
from quara.interface.cvxpy.qtomography.standard.minimization_algorithm import (
CvxpyMinimizationAlgorithm,
CvxpyMinimizationAlgorithmOption,
)
from quara.interface.cvxpy.qtomography.standard.estimator import (
CvxpyLossMinimizationEstimator,
)
[8]:
# Loss and Algorithm
mode_constraint = "physical"
#mode_constraint = "unconstraint"
name_solver = "mosek"
#name_solver = "scs"
#name_solver = "cvxopt"
loss = CvxpyRelativeEntropy()
#loss = CvxpyUniformSquaredError()
#loss = CvxpyApproximateRelativeEntropyWithZeroProbabilityTerm()
loss_option = CvxpyLossFunctionOption()
algo = CvxpyMinimizationAlgorithm()
algo_option = CvxpyMinimizationAlgorithmOption(
name_solver=name_solver, mode_constraint=mode_constraint
)
estimator = CvxpyLossMinimizationEstimator()
[9]:
import numpy.testing as npt
[10]:
# Estimation
sqt = povmt
result = estimator.calc_estimate(
qtomography=sqt,
empi_dists=empi_dists,
loss=loss,
loss_option=loss_option,
algo=algo,
algo_option=algo_option,
is_computation_time_required=True,
)
var_estimate = result.estimated_var
# Test
actual = var_estimate
expected = true_object.to_var()
decimal = 1e-8
npt.assert_almost_equal(actual, expected, decimal=decimal)
print("\n")
print("mode_constraint =", mode_constraint)
print("actual =", actual)
print("expected =", expected)
print("squared_error =", np.inner(actual-expected, actual-expected))
mode_constraint = physical
actual = [7.07106781e-01 1.77604451e-06 9.16673177e-08 7.07106774e-01]
expected = [0.70710678 0. 0. 0.70710678]
squared_error = 3.1627875890750366e-12