The Spectrum of the Transverse Field Ising Model and the Jordan-Wigner Transformation

This example shows how to code up the transverse-field Ising chain and the Jordan-Wigner-equivalent fermion p-wave superconductor:

\[\begin{split}H&=&\sum_{j=0}^{L-1}-J\sigma^z_{j+1}\sigma^z_j - h\sigma^x_j, \nonumber\\ H&=&\sum_{j=0}^{L-1}J\left(-c^\dagger_jc_{j+1} + c_jc^\dagger_{j+1} \right) +J\left( -c^\dagger_jc^\dagger_{j+1} + c_jc_{j+1}\right) + 2h\left(n_j-\frac{1}{2}\right).\end{split}\]

Details about the code below can be found in SciPost Phys. 7, 020 (2019).

Script

download

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
from __future__ import print_function, division
import sys,os
# line 4 and line 5 below are for development purposes and can be removed
qspin_path = os.path.join(os.getcwd(),"../../")
sys.path.insert(0,qspin_path)
#####################################################################
#                            example 4                              #
#    In this script we demonstrate how to construct fermionic       #
#    Hamiltonians, and check the Jordan-Wigner transformation.      #
#####################################################################
from quspin.operators import hamiltonian # Hamiltonians and operators
from quspin.basis import spin_basis_1d, spinless_fermion_basis_1d # Hilbert space spin basis
import numpy as np # generic math functions
import matplotlib.pyplot as plt # plotting library
#
##### define model parameters #####
L=8 # system size
J=1.0 # spin zz interaction
h=np.sqrt(2) # z magnetic field strength
#
# loop over spin inversion symmetry block variable and boundary conditions
for zblock,PBC in zip([-1,1],[1,-1]):
	#
	##### define spin model
	# site-coupling lists (PBC for both spin inversion sectors)
	h_field=[[-h,i] for i in range(L)]
	J_zz=[[-J,i,(i+1)%L] for i in range(L)] # PBC
	# define spin static and dynamic lists
	static_spin =[["zz",J_zz],["x",h_field]] # static part of H
	dynamic_spin=[] # time-dependent part of H
	# construct spin basis in pos/neg spin inversion sector depending on APBC/PBC
	basis_spin = spin_basis_1d(L=L,zblock=zblock) 
	# build spin Hamiltonians
	H_spin=hamiltonian(static_spin,dynamic_spin,basis=basis_spin,dtype=np.float64)
	# calculate spin energy levels
	E_spin=H_spin.eigvalsh()
	#
	##### define fermion model
	# define site-coupling lists for external field
	h_pot=[[2.0*h,i] for i in range(L)]
	if PBC==1: # periodic BC: odd particle number subspace only
		# define site-coupling lists (including boudary couplings)
		J_pm=[[-J,i,(i+1)%L] for i in range(L)] # PBC
		J_mp=[[+J,i,(i+1)%L] for i in range(L)] # PBC
		J_pp=[[-J,i,(i+1)%L] for i in range(L)] # PBC
		J_mm=[[+J,i,(i+1)%L] for i in range(L)] # PBC
		# construct fermion basis in the odd particle number subsector
		basis_fermion = spinless_fermion_basis_1d(L=L,Nf=range(1,L+1,2))
	elif PBC==-1: # anti-periodic BC: even particle number subspace only
		# define bulk site coupling lists
		J_pm=[[-J,i,i+1] for i in range(L-1)]
		J_mp=[[+J,i,i+1] for i in range(L-1)]
		J_pp=[[-J,i,i+1] for i in range(L-1)]
		J_mm=[[+J,i,i+1] for i in range(L-1)]
		# add boundary coupling between sites (L-1,0)
		J_pm.append([+J,L-1,0]) # APBC
		J_mp.append([-J,L-1,0]) # APBC
		J_pp.append([+J,L-1,0]) # APBC
		J_mm.append([-J,L-1,0]) # APBC
		# construct fermion basis in the even particle number subsector
		basis_fermion = spinless_fermion_basis_1d(L=L,Nf=range(0,L+1,2))
	# define fermionic static and dynamic lists
	static_fermion =[["+-",J_pm],["-+",J_mp],["++",J_pp],["--",J_mm],['z',h_pot]]
	dynamic_fermion=[]
	# build fermionic Hamiltonian
	H_fermion=hamiltonian(static_fermion,dynamic_fermion,basis=basis_fermion,
							dtype=np.float64,check_pcon=False,check_symm=False)
	# calculate fermionic energy levels
	E_fermion=H_fermion.eigvalsh()
	#
	##### plot spectra
	plt.plot(np.arange(H_fermion.Ns),E_fermion/L,marker='o'
									,color='b',label='fermion')
	plt.plot(np.arange(H_spin.Ns),E_spin/L,marker='x'
									,color='r',markersize=2,label='spin')
	plt.xlabel('state number',fontsize=16)
	plt.ylabel('energy',fontsize=16)
	plt.xticks(fontsize=16)
	plt.yticks(fontsize=16)
	plt.legend(fontsize=16)
	plt.grid()
	plt.tight_layout()
	plt.savefig('example4.pdf', bbox_inches='tight')
	#plt.show()
	plt.close()