%autosave 0
%matplotlib inline
from __future__ import division
from sglib import *
from sympy import sqrt, E, Rational, latex, Symbol
L = latex # Shorthand
# Zero and One states (they're also + and - z)
s0 = col(1,0); s1 = col(0,1)
pz = s0; mz = s1
# Create a formatter for stuff like |00>
fmt = sg_format_state(['0','1'], '').format
zfmt = sg_format_state(['+z','-z'], '').format
xfmt = sg_format_state(['+x','-x'], '').format
# Create two-bit states and basis
s00=TP(s0,s0); s01=TP(s0,s1); s10=TP(s1,s0); s11=TP(s1,s1)
b2bit = [s00, s01, s10, s11]
# Some misc definitions
half = Rational(1,2) # Define 1/2
alpha = Symbol('alpha') # Define some symbols
beta = Symbol('beta')
# The Bell states
bell_1 = 1/sqrt(2)*s00 + 1/sqrt(2)*s11 # phi+
bell_2 = 1/sqrt(2)*s00 - 1/sqrt(2)*s11 # phi-
bell_3 = 1/sqrt(2)*s01 + 1/sqrt(2)*s10 # psi+ (Triplet)
bell_4 = 1/sqrt(2)*s01 - 1/sqrt(2)*s10 # psi- (Singlet)
bell = [bell_1, bell_2, bell_3, bell_4]
The density operator of a pure state $|\psi\rangle$ is:
$$ \rho = |\psi\rangle\langle\psi| $$
The off-diagonal terms are the "interference" terms, for the basis in which the expression is written.
A mixed state can be viewed as a state chosen from a classical ensemble of states, consisting of states $\{\psi_i\}$ for which you know the relative frequency of each state in the ensemble.
The density matrix for a mixed state is simply the weighted sum of the density matrices for the individual states:
$$ (2.20) \;\; \rho = \sum_i p_i|\psi_i\rangle\langle\psi_i| $$ where $p_i$ is the (classical) probability of drawing $\psi_i$ from the ensemble.
For example, say you have a 50% chance of getting $|0\rangle$ and a 50% chance of getting $|1\rangle$. Then we have: $$ \rho = \frac{1}{2}|0\rangle\langle 0| + \frac{1}{2}|1\rangle\langle 1| $$
The trace operation is:
$$ \mathrm{Tr}(\rho) = \sum_i\langle\phi_i|\rho|\phi_i\rangle $$
where $\{\phi_i\}$ are a complete set of orthonormal basis vectors. The expression $\langle\phi_i|A|\phi_i\rangle$ picks out the $i$th row and column of $A$. So we're actually just adding up the main diagonal of the matrix.
Notice that trace work in any basis. As long as you use a set of orthonormal basis vectors, you'll trace out the same number.
from sympy import var
a_11, a_12, a_21, a_22 = var('a_11, a_12, a_21, a_22')
A = mat(a_11, a_12, a_21, a_22)
MOP = matrix_as_outer_product(A) # Just to get access to the "demo" Tr()
from sglib import pZ, mZ, pX, mX, pY, mY
trace = MOP.Tr(A, [pZ, mZ], V=True)
trace = MOP.Tr(A, [pX, mX], V=True)
trace = MOP.Tr(A, [pY, mY], V=True)
Say you have a state $|\psi\rangle$, and it's density matrix $\rho$. You want to know the probability that a measurement on $|\psi\rangle$ will yield some arbitrary state $|\phi\rangle$. This can be calculated in the following way: $$ \mathrm{Prob}(\phi)=\langle \phi|\rho|\phi\rangle =\mathrm{Tr}(|\phi\rangle\langle\phi|\rho) $$
Purity is a measure of how mixed the density matrix is.
$\;\;\;\;\mathrm{purity}(\rho) = \mathrm{Tr}(\rho^2)$
The maximum purity is $1$ (which will be achieved for any pure state). The minimum purity, will be one divided by the dimension of the matrix. $1/2$ for a one-bit state, $1/4$ for a two-bit state, and so on.
The Entropy $S$ of of a system can be calculated by the following formula, where the $\lambda$ are the eigenvalues of the system's density (or reduced density) matrix:
$\;\;\;\;S = - \displaystyle\sum_{n} \lambda_n \mathrm{log}(\lambda_n) $
The entropy of any pure state will be 0. The maximum entropy (assuming you are taking the logs in base 2) will be the number of bits in the state.
Example: A one bit pure state
psi = 1/sqrt(2)*s0 + 1/sqrt(2)*s1
rho = OP(psi,psi)
examine_dm(rho)
Example: A one bit even mixture
half = Rational(1,2)
rho = half*OP(s0,s0) + half*OP(s1,s1)
examine_dm(rho)
Example: One bit, somewhat mixed
rho = .15*OP(s0,s0) + .85*OP(s1,s1)
examine_dm(rho)
Example: Two bit pure state (Bell 1)
rho = OP(bell_1, bell_1)
examine_dm(rho)
Example: Even mixture of all four Bell states.
rho = .25*OP(bell_1, bell_1) + .25*OP(bell_2, bell_2)\
+ .25*OP(bell_3, bell_3) + .25*OP(bell_4, bell_4)
examine_dm(rho)
Example: Partially mixed two bit state
rho = .25*OP(s00,s00)+.15*OP(s01,s01)+.05*OP(s10,s10)+.55*OP(s11,s11)
print 'trace is %s'%rho.trace() # Make sure it adds up!
examine_dm(rho)
Example: A two bit state with an even mixture of two states
has the same purity and entropy as an evenly mixed one bit state.
rho = .50*OP(s00,s00)+.50*OP(s10,s10)
print 'trace is %s'%rho.trace() # Make sure it adds up!
examine_dm(rho)
Example: Calculating the probability of an arbitrary state
# The state we are given
psi = 1/sqrt(2)*pz + 1/sqrt(2)*mz
rho = OP(psi,psi)
phi = col(1/sqrt(2), 1/sqrt(2)) # This is +x
phi = col(1,0) # This is +z
Print(r'Given: $%s=%s$' %( zfmt(psi), L(psi) ))
Print(r'Find the probability of measuring $%s=%s$' %(zfmt(phi),L(phi) ))
Print(r'$\langle\phi|\rho|\phi\rangle=%s%s%s=%s$'
%(L(phi.adjoint()),L(rho),L(phi), (phi.adjoint()*rho*phi)[0] ))
foo = phi*phi.adjoint()*rho
string = r'$\mathrm{Tr}(|\phi\rangle\langle\phi|\rho)'
string += r'=\mathrm{Tr}\left(%s%s%s\right)=\mathrm{Tr}\left(%s\right)=%s$'
Print(string
%(L(phi),L(phi.adjoint()),L(rho),L(foo),foo.trace()),font_size=3)
(1) $\rho_{AB} = \rho_A \otimes \rho_B$.
$\;\;\;\;$ In other words: $|a_1 b_1\rangle \langle a_2 b_2| = |a_1\rangle \langle a_2| \otimes |b_1\rangle \langle b_2|$
(2) The trace of an outer product equals the (reversed) inner product
$\;\;\;\;\;\mathrm{Tr}(|x\rangle\langle y|) = \langle y|x\rangle$
(3) The Partial Trace operation, where we "trace out" bit B is:
$\;\;\;\;\;\mathrm{Tr_B}(|a_1 b_1\rangle\langle a_2 b_2|) \\\;\;\;\;\;\;\; = \mathrm{Tr_B}(|a_1\rangle\langle a_2|\otimes|b_1\rangle\langle b_2|) \\\;\;\;\;\;\;\; = |a_1\rangle\langle a_2|\mathrm{Tr_B}(|b_1\rangle\langle b_2|) \\\;\;\;\;\;\;\; = |a_1\rangle\langle a_2|\langle b_2|b_1\rangle $
(4) The Reduced Density Matrix, $\rho_A$ is the matrix where we have "traced out" bit B:
$\;\;\;\;\;\rho_A = \mathrm{Tr_B}(\rho_{AB})$
one = 1/sqrt(2)*s0 + 1/sqrt(2)*s1
two = sqrt(3)/2*s0 + Rational(1,2)*s1
state = TP(two, one)
Print(r'State = $%s$' %fmt(state, which_ltx='sympy'))
rho = OP(state, state)
Print(r'Density matrix is: $\rho=%s$' %latex(rho))
MOP = matrix_as_outer_product(rho)
Print(r'In "outer product form":')
Print(r'$\rho=%s$' %MOP.latex())
MOP.partial_trace(1)
Print('After tracing out the second bit:')
Print(r'$\rho=%s$' %MOP.latex())
Print('And in matrix form:')
Print(r'$\rho=%s$' %latex(MOP.M))
state = 1/sqrt(2)*s00 + half*s01 + half*s11
Print(r'State = $%s$' %fmt(state, which_ltx='mine'))
rho = OP(state, state)
Print(r'Density matrix is: $\rho=%s$' %latex(rho))
MOP = matrix_as_outer_product(rho)
Print(r'In "outer product form":')
Print(r'$\rho=%s$' %MOP.latex())
MOP.partial_trace(1)
Print('After tracing out the second bit:')
Print(r'$\rho=%s$' %MOP.latex())
Print('And in matrix form:')
Print(r'$\rho=%s$' %latex(MOP.M))
Example 1: Take a maximally entangled state
It is not a mixture. So the purity should be maximum and the entropy minimum.
Print(r'$\psi = %s$'%fmt(bell_1))
rho_maxE = OP(bell_1, bell_1)
examine_dm(rho_maxE)
Find the reduced density for the first bit
Since the entanglement was maximal, either individual bit should contain no information about the outcome. This looks just like a maximally mixed state. Both the purity and entropy should reflect this (minimum purity, maximum entropy).
rho = matrix_as_outer_product(rho_maxE)
rho.partial_trace(1)
examine_dm(rho.M)
The reduced density for the second bit should be just like the first
rho = matrix_as_outer_product(rho_maxE)
rho.partial_trace(0)
examine_dm(rho.M)
Example 2: Reproduce Cugini page 13 classical mixture (base e)
This is really just a test that my code works for base $e$.
rho = mat(.5,0,0,0, 0,0,0,0, 0,0,0,0, 0,0,0,.5)
S = Entropy(rho, base=E)
Print(r'$\rho = %s, \;\; S = %s$' %(latex(rho), float(S)))
Example 3: A non-maximally entangled 2 bit state
The "intermediate" values of purity and entropy are a reflection of "some amount" of entanglement.
r3 = 1/sqrt(3)
psi = r3*s00 + r3*s01 + r3*s10
Print(r'$\psi = %s$'%fmt(psi))
rho = OP(psi,psi)
examine_dm(rho)
# Then do the reduced density for both bit. Should be the same.
rho_op = matrix_as_outer_product(rho)
rho_op.partial_trace(1)
examine_dm(rho_op.M)
rho_op = matrix_as_outer_product(rho)
rho_op.partial_trace(0)
examine_dm(rho_op.M)
Example 4: No measurable correlations, yet entangled.
The minimal purity and maximal entropy of one of the bits verifies that this state is maximally entangled.
psi = half*s00 + half*s01 + half*s10 - half*s11
Print(r'$\psi = %s$'%fmt(psi))
rho = OP(psi,psi)
S = Entropy(rho)
Print(r'$\rho = %s, \;\; S = %s$' %(latex(rho), latex(S)))
# Then do the reduced density for bit zero
rho_op = matrix_as_outer_product(rho)
rho_op.partial_trace(1)
examine_dm(rho_op.M)
Now just change the sign on the last term
Giving the product state: $\left(\frac{1}{\sqrt{2}}|0\rangle \frac{1}{\sqrt{2}}|1\rangle\right) \otimes \left(\frac{1}{\sqrt{2}}|0\rangle \frac{1}{\sqrt{2}}|1\rangle\right)$
And we get maximum purity, minimum entropy.
psi = half*s00 + half*s01 + half*s10 + half*s11
Print(r'$\psi = %s$'%fmt(psi))
rho = OP(psi,psi)
S = Entropy(rho)
Print(r'$\rho = %s, \;\; S = %s$' %(latex(rho), latex(S)))
# Then do the reduced density for bit one
rho_op = matrix_as_outer_product(rho)
rho_op.partial_trace(1)
examine_dm(rho_op.M)
P = [] # Purity values
S = [] # Entropy values
M = [] # Percentage of -z mixed in
for percent_mZ in range(51):
frac_mZ = percent_mZ/100.0
rho = (1-frac_mZ)*OP(pZ, pZ) + frac_mZ*OP(mZ,mZ)
M += [ percent_mZ, ]
P += [ Purity(rho), ]
S += [ Entropy(rho), ]
pyplot.plot(M, P, label="Purity", color="green");
pyplot.plot(M, S, label="Entropy", color="red");
pyplot.legend(loc='best');
pyplot.title('One bit purity/entropy vs percent mixed');
psi = alpha*s0 - beta*s1
Print(r'$\psi=%s$' %myltx(psi))
rho = OP(psi,psi)
Print(r'$\rho=%s$' %myltx(rho))
MOP = matrix_as_outer_product(rho)
Print(r'$%s$' %MOP.latex())
Print(r'$\rho=%s$' %latex(MOP.M))
M = mat(1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8)
#M = mat(-1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,-16)
MOP = matrix_as_outer_product(M)
Print('The original matrix:')
Print(r'$%s$' %MOP.latex())
Print('Partial Trace of bit 0:')
MOP.partial_trace(0)
Print(r'$%s$' %MOP.latex())
Print('Another partial trace of bit 0:')
MOP.partial_trace(0)
Print(r'$%s$' %MOP.latex())
#
# How about matrices with symbolic values?
#
psi = TP(1/sqrt(2)*s0 + 1/sqrt(2)*s1, col(alpha, beta))
Print(r'$\psi=%s$'%latex(psi))
Print(r'$\psi=%s$'%fmt(psi))
M = OP(psi,psi)
Print(r'$\rho = %s$' %latex(M))
MOP = matrix_as_outer_product(M)
Print(r'$%s$'%MOP.latex())
MOP.partial_trace(0)
Print('Tracing out bit zero:')
Print(r'$%s$'%MOP.latex())
#
# Now the teleportation state
#
psi = half*TP(s00, col(alpha,beta))\
+ half*TP(s01, col(beta,alpha))\
+ half*TP(s10, col(alpha,-beta))\
+ half*TP(s11, col(-beta,alpha))
Print(r'$\psi=%s$'%fmt(psi)) # THIS IS A PROBLEM !!!
Print(r'$\psi=%s$'%latex(psi))
M = OP(psi,psi)
Print(r'$\rho = %s$' %latex(M))
MOP = matrix_as_outer_product(M)
Print(r'$%s$'%MOP.latex())
MOP.partial_trace(0)
Print('Tracing out Alice\'s first bit:')
Print(r'$%s$'%MOP.latex())
MOP.partial_trace(0)
Print('Tracing out Alice\'s second bit:')
Print(r'$%s$'%MOP.latex())
Print('And the reduced density matrix for Bob is:')
Print(r'$\rho_B = %s$' %latex(MOP.M))