vprusso/toqito

Add unambiguous discrimination strategy to state distinguishability and state exclusion

vprusso opened this issue ยท 7 comments

Enable the ability for the user to unambiguously exclude or discriminate amongst of quantum states in the state_distinguishability.py and state_exclusion.py files, respectively.

Refer to Section 2.4 in arXiv:1707.02571 for further information on the constraints and definition of the unambiguous case for the SDP and refer to ppt_distinguishability.py for an example of how the unambiguous strategy is supported within the SDP and function arguments.

@vprusso I have a code that works with solver="mosek", but for a reason I can't figure out, not with solver="cvxopt". Of course, cvxopt being the default one, it's not quite convenient.

Thus, I see two options here:

  • Try another formulation (that is, not the immediate one), for which cvxopt may succeed. The problem with this approach is that the implementation would be quite different, that would be a whole new function altogether, rater than just adding constraints on the existing code.
  • Implement it as is, and open an issue on picos to report the issue.

What is your preferred solution?

Of course, it might be that my implementation is wrong altogether! But the fact that I obtain the correct results with Mosek seems weird to me. If you want, I don't mind sharing my code here so that we can discuss about it ๐Ÿ™‚

Implement it as is, and open an issue on picos to report the issue.

Not Vincent but I think this is the way to go if it is indeed a bug in picos. mosek giving you the expected results makes me think something in picos might not be working correctly.

I don't mind sharing my code here so that we can discuss about it

Sharing your code for further discussion would be a better option as well!

@purva-thakre Sure, here it is:

def state_distinguishability(
    vectors: list[np.ndarray],
    probs: list[float] = None,
    strategy: str = "min_error",
    solver: str = "cvxopt",
    primal_dual: str = "dual",
) -> float:
    r"""[Docstring]
    """
    if not has_same_dimension(vectors):
        raise ValueError("Vectors for state distinguishability must all have the same dimension.")

    # Assumes a uniform probabilities distribution among the states if one is not explicitly provided.
    n = len(vectors)
    probs = [1 / n] * n if probs is None else probs
    dim = calculate_vector_matrix_dimension(vectors[0])

    if primal_dual == "primal":
        return _min_error_primal(vectors=vectors, dim=dim, probs=probs, solver=solver, strategy=strategy)
    return _min_error_dual(vectors=vectors, dim=dim, probs=probs, solver=solver, strategy=strategy)


def _min_error_primal(
    vectors: list[np.ndarray],
    dim: int,
    probs: list[float] = None,
    solver: str = "cvxopt",
    strategy: str = "min_error",
) -> tuple[float, list[picos.HermitianVariable]]:
    """Find the primal problem for minimum-error quantum state distinguishability SDP."""
    n = len(vectors)

    problem = picos.Problem()
    num_measurements = len(vectors)
    measurements = [picos.HermitianVariable(f"M[{i}]", (dim, dim)) for i in range(num_measurements)]

    if strategy == "unambiguous":
        measurements.append(picos.HermitianVariable(f"M[{num_measurements}]", (dim, dim)))

    problem.add_list_of_constraints([meas >> 0 for meas in measurements])
    problem.add_constraint(picos.sum(measurements) == picos.I(dim))

    dms = [vector_to_density_matrix(vector) for vector in vectors]

    if strategy == "unambiguous":
        problem.add_list_of_constraints(measurements[i] | dms[j] == 0 for i in range(n) for j in range(n) if i != j)

    problem.set_objective("max", np.real(picos.sum([(probs[i] * dms[i] | measurements[i]) for i in range(n)])))
    solution = problem.solve(solver=solver)
    return solution.value, measurements


def _min_error_dual(
    vectors: list[np.ndarray],
    dim: int,
    probs: list[float] = None,
    solver: str = "cvxopt",
    strategy: str = "min_error",
) -> tuple[float, list[picos.HermitianVariable]]:
    """Find the dual problem for minimum-error quantum state distinguishability SDP."""
    if strategy != "min_error":
        raise ValueError("The dual method only supports minimal-error distinguishability at this time.")

    n = len(vectors)
    problem = picos.Problem()

    # Set up variables and constraints for SDP:
    y_var = picos.HermitianVariable("Y", (dim, dim))
    problem.add_list_of_constraints(
        [y_var >> probs[i] * vector_to_density_matrix(vector) for i, vector in enumerate(vectors)]
    )

    # Objective function:
    problem.set_objective("min", picos.trace(y_var))
    solution = problem.solve(solver=solver, primals=None)

    measurements = [problem.get_constraint(k).dual for k in range(n)]

    return solution.value, measurements

Essentially, the problem of unambiguous is the same as the minimal error one, with the additional constraint that $\mathrm{Tr}\left[M_i\rho_j\right]=0$ for $i\neq j$, so the implementation is essentially just: add a measurement operator, and add the constraints.

If you want to test by yourself, this part seems to be the problematic one:

if strategy == "unambiguous":
        measurements.append(picos.HermitianVariable(f"M[{num_measurements}]", (dim, dim)))

Indeed, even if you comment the other if case, cvxopt fails with an ArithmeticError.

But maybe there's something trivial I missed, that's definitely a possibility ๐Ÿ˜„

Hey @tnemoz,

Thank you for your interest and contribution!

Hmm, so I was thinking about this a bit yesterday. Here is another approach to unambiguous distinguishability that makes use of the associated Gram matrix of the set of states. The nice thing about this is that the problem scales very well (it's a function of the number of states as opposed to the dimension of states):

Primal:

def unambig_primal(gram_matrix: np.ndarray, probs: list[float] | None = None, solver="mosek") -> float:
    """SDP from Equation (5) of arXiv:2402.06365 for unambiguous state discrimination."""
    problem = picos.Problem()
    n = gram_matrix.shape[0]
    probs = [1 / n] * n if probs is None else probs    

    x_var = picos.RealVariable("x_var", n, lower=0)
    problem.set_objective("max", probs * x_var)

    x_mat = picos.diag(x_var)
    problem.add_constraint(gram_matrix - x_mat >> 0)

    problem.solve(solver=solver)
    print(np.around(x_mat, decimals=5))

    return x_var.value[0]

Dual:

def unambig_dual(gram_matrix: np.ndarray, probs: list[float] | None = None, solver="mosek") -> float:
    """SDP for the dual problem for unambiguous state discrimination."""
    problem = picos.Problem()
    n = gram_matrix.shape[0]
    
    probs = [1 / n] * n if probs is None else probs
    m = len(probs)

    # Variables
    Z = picos.SymmetricVariable("Z", n)
    z = picos.RealVariable("z", m, lower=0)
    
    # Objective
    problem.set_objective("min", picos.trace(gram_matrix * Z))
    
    # Constraints
    for i in range(m):
        # Define F_i with -1 at position (i, i)
        Fi = np.zeros((n, n))
        Fi[i, i] = -1
        problem.add_constraint(z[i] + probs[i] + picos.trace(Fi * Z) == 0)
    
    problem.add_constraint(Z >> 0)
    problem.add_constraint(z >= 0)
    
    problem.solve(solver=solver)

    print(z.value)
    print(np.around(Z.value, decimals=5))
    
    return problem.value

Not sure if providing this is helpful, but these might be a bit more cooperative with cvxopt. I prefer mosek as well, but sadly, it's a paid option. I definitely want to support that and encourage users to make use of it if available, but failing that, I certainly do want to ensure that support for the opensource option is accessible.

Let me know if that helps, and thanks again!

@vprusso Thanks! I've managed to implement the unambiguous state discrimination and some tests with it.
However, for now I'm stuck on the quantum exclusion part. I thought it would be as easy as turning the "min" into "max" and vice-versa, but not only does that not work, I also don't understand why the SDP provided by the paper you provided me with is correct.

I mean, it certainly is (it holds numerically at least :D), but I don't understand why. I'll keep working on it and will keep you posted if I find my error!

I thought it would be as easy as turning the "min" into "max" and vice-versa

The way I understand it, I also think this should work for the state exclusion case.

I also don't understand why the SDP provided by the paper you provided me with is correct

I found a paper that explicitly defines the unambiguous SDPs for state exclusion. See Appendix E. https://arxiv.org/abs/1306.4683
Maybe this helps in finding your error?

@vprusso Thanks! I've managed to implement the unambiguous state discrimination and some tests with it. However, for now I'm stuck on the quantum exclusion part. I thought it would be as easy as turning the "min" into "max" and vice-versa, but not only does that not work, I also don't understand why the SDP provided by the paper you provided me with is correct.

I mean, it certainly is (it holds numerically at least :D), but I don't understand why. I'll keep working on it and will keep you posted if I find my error!

One thing to note here is that the max and min are indeed flipped for state exclusion (as you mentioned). Depending on where the issue is that you're encountering, another thing to be aware of is that the PSD condition in the dual also has a flipped inequality. That is, the condition should be Y <= ... in the dual for state exclusion, where it is Y >= in the dual for state distinguishability.

Not sure if that helps out or not. Also, please do feel free to ping here in the chat thread with any additional details if that's helpful. As @purva-thakre mentioned, the paper she linked does indeed contain the SDP for the unambiguous state exclusion problem which may prove useful.