The minimal example is a chain of qubits of 4, and then followed by a 2-site gate. The strange thing is that: the answer is incorrect when the 2 site gate is applied to qubits (1,2) using MPS.apply_local_op", but it is correct when applying to qubits (0,1). What is even stranger is that, the correctness depends on the specific applied "2-site gate". Below I will provide the example to reproduce it.
I am using the version of
tenpy 0.10.0 (compiled without HAVE_MKL),
git revision 0ada42d89dcf186a22ff820159dec378f2cde077 using
python 3.11.4 (main, Jul 5 2023, 14:15:25) [GCC 11.2.0]
numpy 1.24.3, scipy 1.11.4
Python: Select all
import tenpy
import numpy as np
from tenpy.linalg import np_conserved as npc
# construct a simple |0000> state
qubit_site = tenpy.networks.site.SpinHalfSite(conserve=None,)
legs=[f'p{i}' for i in range(4)]
vec=np.zeros((2**4,),dtype=complex)
vec[0]=1
wavefunction_tenpy = npc.Array.from_ndarray(vec.reshape((2,2,2,2)), [qubit_site.leg] * 4,labels=legs,)
# convert mps back to state vector
def toarray(psi_mps):
wf=psi_mps.get_B(0)
for mps in psi_mps._B[1:]:
wf=npc.tensordot(wf,mps,axes=([-1],[0]))
wf_contracted = npc.trace(wf, leg1=0, leg2=-1)
return wf_contracted.to_ndarray().flatten()
wf_mps=tenpy.networks.mps.MPS.from_full( [qubit_site]*4, wavefunction_tenpy,form='B',normalize=True)
toarray(wf_mps)
# The output is correct.
array([1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j,
0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j])
Python: Select all
U1=(1+1j)*np.arange(1,17).reshape((2,2,2,2))
U_npc1=npc.Array.from_ndarray(U1, [qubit_site.leg] * 4,labels=['p0','p1','p0*','p1*'],)
wf_mps=tenpy.networks.mps.MPS.from_full( [qubit_site]*4, wavefunction_tenpy,form='B',normalize=True)
wf_mps.apply_local_op(1,U_npc1)
toarray(wf_mps)
# Output:
array([0.04256283+0.04256283j, 0. +0.j ,
0.21281413+0.21281413j, 0. +0.j ,
0.38306544+0.38306544j, 0. +0.j ,
0.55331674+0.55331674j, 0. +0.j ,
0. +0.j , 0. +0.j ,
0. +0.j , 0. +0.j ,
0. +0.j , 0. +0.j ,
0. +0.j , 0. +0.j ])
wf_mps=tenpy.networks.mps.MPS.from_full( [qubit_site]*4, wavefunction_tenpy,form='B',normalize=True)
wf_mps.apply_local_op(0,U_npc1)
toarray(wf_mps)
# Output:
array([0.04256283+0.04256283j, 0. +0.j ,
0. +0.j , 0. +0.j ,
0.21281413+0.21281413j, 0. +0.j ,
0. +0.j , 0. +0.j ,
0.38306544+0.38306544j, 0. +0.j ,
0. +0.j , 0. +0.j ,
0.55331674+0.55331674j, 0. +0.j ,
0. +0.j , 0. +0.j ])
Now comes to the strange part: If the applied gate is the following (which I attach it in the pickle), then the result is wrong.
Python: Select all
import pickle
with open('UU.pkl','rb') as f:
U2=pickle.load(f)
U2
# Output:
array([[[[ 0.04578135-0.19817758j, -0.0670271 -0.11887706j],
[ 0.21649947+0.26349682j, -0.17076864+0.89137018j]],
[[-0.1950499 -0.04680251j, 0.236489 +0.71568491j],
[ 0.56152975-0.10577503j, -0.25196233-0.04054729j]]],
[[[-0.25624666+0.32897489j, -0.62812977-0.04047919j],
[ 0.02579847-0.56901467j, -0.29334038+0.13941218j]],
[[-0.84659876-0.16666881j, 0.08502642+0.09913931j],
[-0.44095103+0.19573532j, -0.02838949+0.06984177j]]]])
Python: Select all
U_npc2=npc.Array.from_ndarray(U2, [qubit_site.leg] * 4,labels=['p0','p1','p0*','p1*'],)
wf_mps=tenpy.networks.mps.MPS.from_full( [qubit_site]*4, wavefunction_tenpy,form='B',normalize=True)
wf_mps.apply_local_op(1,U_npc2)
toarray(wf_mps)
# Output:
array([ 0.29762525-0.81477149j, 0. +0.j ,
-0.49749284-0.00821637j, 0. +0.j ,
-0.26306051+0.42233375j, 0. +0.j ,
-0.84588209-0.19213763j, 0. +0.j ,
0. +0.j , 0. +0.j ,
0. +0.j , 0. +0.j ,
0. +0.j , 0. +0.j ,
0. +0.j , 0. +0.j ])
toarray(wf_mps)
This is already wrong. The vector components do not even exist in the original U2 matrix.
However, if I just apply the gate to (0,1), then it is correct again:
Python: Select all
wf_mps=tenpy.networks.mps.MPS.from_full( [qubit_site]*4, wavefunction_tenpy,form='B',normalize=True)
wf_mps.apply_local_op(0,U_npc2)
toarray(wf_mps)
# Output:
array([ 0.04578135-0.19817758j, 0. +0.j ,
0. +0.j , 0. +0.j ,
-0.1950499 -0.04680251j, 0. +0.j ,
0. +0.j , 0. +0.j ,
-0.25624666+0.32897489j, 0. +0.j ,
0. +0.j , 0. +0.j ,
-0.84659876-0.16666881j, 0. +0.j ,
0. +0.j , 0. +0.j ])
I have also provide the ipynb which can directly run to reproduce it.