M.E. Irizarry-Gelpí

Physics impostor. Mathematics interloper. Husband. Father.

CHSH Game 3


Here is the Python code that I used to show that there is a strategy with a 2-qubit state share by Alice and Bob to win the CHSH game with a probability that is larger than if they had classical bits instead.

I used SymPy to do the symbolic calculations. In the future, I would love to do a numerical simulation of the quantum game. First you import SymPy and enable pretty printing:

import sympy as sym
sym.init_printing()

You create some symbols for using in expressions.

alpha, beta = sym.symbols('alpha beta')
alpha0, alpha1 = sym.symbols('alpha0 alpha1')
beta0, beta1 = sym.symbols('beta0 beta1')

The 1-qubit state vectors are just 2-dimensional vectors. Here are some functions that return the two ortho-normal basis states:

def bra0(theta):
    return sym.Matrix([[sym.cos(theta), sym.sin(theta)]])

def bra1(theta):
    return sym.Matrix([[sym.sin(theta), -sym.cos(theta)]])

def ket0(theta):
    return bra0(theta).T

def ket1(theta):
    return bra1(theta).T

You can use the sym.tensor.tensorproduct function to compute the tensor product and get the 2-qubit basis states:

def ketA0B0(a, b):
    return sym.Matrix(sym.tensor.tensorproduct(bra0(a), bra0(b)))

def ketA0B1(a, b):
    return sym.Matrix(sym.tensor.tensorproduct(bra0(a), bra1(b)))

def ketA1B0(a, b):
    return sym.Matrix(sym.tensor.tensorproduct(bra1(a), bra0(b)))

def ketA1B1(a, b):
    return sym.Matrix(sym.tensor.tensorproduct(bra1(a), bra1(b)))

def braA0B0(a, b):
    return ketA0B0(a,b).T

def braA0B1(a, b):
    return ketA0B1(a,b).T

def braA1B0(a, b):
    return ketA1B0(a,b).T

def braA1B1(a, b):
    return ketA1B1(a,b).T

It is also convenient to introduce the four Bell states:

braPhiP = sym.Matrix([[1,0,0,1]]) / sym.sqrt(2)
braPhiN = sym.Matrix([[1,0,0,-1]]) / sym.sqrt(2)
braPsiP = sym.Matrix([[0,1,1,0]]) / sym.sqrt(2)
braPsiN = sym.Matrix([[0,1,-1,0]]) / sym.sqrt(2)

ketPhiP = braPhiP.T
ketPhiN = braPhiN.T
ketPsiP = braPsiP.T
ketPsiN = braPsiN.T

Here ketPhiP stands for \(| \Phi^{+} \rangle\) and ketPhiN stands for \(| \Phi^{-} \rangle\). There are to distinct contributions to probabilities given by the following functions:

def p1(a, b, ket):
    return ((braA0B0(a, b) @ ket)[0,0])**2 + ((braA1B1(a, b) @ ket)[0,0])**2

def p2(a, b, ket):
    return ((braA0B1(a, b) @ ket)[0,0])**2 + ((braA0B1(a, b) @ ket)[0,0])**2

Here @ is the matrix multiplication operator. The probability of winning the game as a function of the two sets of angles is calculated by the following function:

def p_win(a0, a1, b0, b1, ket):
    return (p1(a0, b0, ket) + p1(a0, b1, ket) + p1(a1, b0, ket) + p2(a1, b1, ket)) / 4

You can check that certain values for the set of angles give a maximal probability. You can get a numerical result:

sym.simplify(p_win(alpha0, alpha1, beta0, beta1, ketPhiP)).subs([
    (alpha0, 0),
    (alpha1, sym.pi/4),
    (beta0, sym.pi/8),
    (beta1, -sym.pi/8),
]).evalf()

Or you can get a symbolic result:

sym.simplify(p_win(alpha0, alpha1, beta0, beta1, ketPhiN)).subs([
    (alpha0, 0),
    (alpha1, sym.pi/4),
    (beta0, -sym.pi/8),
    (beta1, sym.pi/8),
])

Both of these give the maximal amount of probability:

$$ p_{\text{win}} = \frac{1}{2} + \frac{1}{\sqrt{8}} \approx 0.853553390593274 $$