From 1fb4e70a86966ab8b6d6980cd91184572a4f9df7 Mon Sep 17 00:00:00 2001
From: Ammon Fischer <ammon.fischer@rwth-aachen.de>
Date: Fri, 24 Jun 2022 09:53:36 -0400
Subject: [PATCH] Added Bandstructure with correct Sz polarization and file to
 compute spin texture from skretch

---
 .DS_Store                       | Bin 0 -> 6148 bytes
 __pycache__/tSnS.cpython-38.pyc | Bin 5153 -> 4786 bytes
 spin_texture.py                 |  59 ++++++++++++++++++++++++++
 tSnS.py                         |  72 +++++++++++++-------------------
 4 files changed, 87 insertions(+), 44 deletions(-)
 create mode 100644 .DS_Store
 create mode 100644 spin_texture.py

diff --git a/.DS_Store b/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..5008ddfcf53c02e82d7eee2e57c38e5672ef89f6
GIT binary patch
literal 6148
zcmeH~Jr2S!425mzP>H1@V-^m;4Wg<&0T*E43hX&L&p$$qDprKhvt+--jT7}7np#A3
zem<@ulZcFPQ@L2!n>{z**<q8>++&mCkOWA81W14cNZ<zv;LbK1Poaz?KmsK2CSc!(
z0ynLxE!0092;Krf2c+FF_Fe*7ECH>lEfg7;MkzE(HCqgga^y>{tEnwC%0;vJ&^%eQ
zLs35+`xjp>T0<F0fCPF1$Cyrb|F7^5{eNG?83~ZUUlGt@xh*qZDeu<Z%US-OSsOPv
j)R!Z4KLME7ReXlK;d!wEw5GODWMKRea10D2@KpjYNUI8I

literal 0
HcmV?d00001

diff --git a/__pycache__/tSnS.cpython-38.pyc b/__pycache__/tSnS.cpython-38.pyc
index 45462ff99d850960856b55c316c47e6591ea6ea7..18e5191016c906076dfecd7be4d8892488f30e6a 100644
GIT binary patch
literal 4786
zcmd5AZEPG@aldYF_wMZYD{<_EWRtigm)3SZiXx!leCLKVRJc|1fmju`n|tf$wfA;6
zyK5)5?gTWcREk0YQJ{aK_y?6JNL7WPf<G4|3P==DT8hdF0a7c_6#-FyK?DsmZ|^Rb
zbIcEvqQ296Z+2$h%)HroZ)fggPfvz}mb>!3(_h=jac|&5XA|LJ3VOH`2Eq~E;q+wY
zb>S#S1R`GGh<IMmMZzEF@=^o%xS2AJovkfrFJLWB=Lx3^gojm3f-Vt>q%LrJieyNd
zC>MBLCMw9w#!__LJ-%`xqaL@SP+Ro1XwK58<yEXGu*;esc)k@xwiRmLtQMWM!^k2U
zEYVDt0H13GW@Kq4$1Bfi_n8&ji9FXf-3cwDmiZ17C4g1v_d^di0VU!t^MptHsJO~k
zIYVF=u(OR38Ac=`k()dcVyfJ{0KQ5;CX8?Zi!GQ6t9%Wl%%?8PRfKVUx|&E=7+wWT
z2G`@dRp0@=S%S3O5X#DKLA#J6(XA4dWKj32fOMm0EiYOGTunZd-Q?{BxCeUgFqg%(
z8*Z|IQJNO?ku1>r2?uW^H$%Vi#_MD9!}m_ne`#-S-+O2ORDh18&)xCn`Af(8u2B6m
zzb`-e;Pe#bU*WIJrVdXvu0v}G$ON6WKPS`GLIWM2ABY9lk5gt4n9H$H+z<VLE));O
zvSYiZLu1i)mtr*vOgHqs(2DsJvE+I|B~FzCFAQ~g$#WJfRxJ2-EF6XPM`8WZgSwjF
zjglVoGkGbNLd%(r<%P3mKw(wUuwCPvE}z4>_$@9N25YWOPo4~|Ae=NS70;cu!}6RJ
zOx|mS3z6qfMmC*`#!uTWv0XY|@rdP&H{)(Rbi+xE*9m_)Rt%77M~2Y@BE!wNUvYqk
zpSY9H@i}pnzYBhS{HVCn!`8I;%jl6-O97CEJ`Fv@!&Al3;NWG8F!Mv)Fjtcz;W8XQ
z7^Ny-mGJeh;$ak4@)=cqaLx=ZZPqMDUce5CH9KpSBYVlxPLwoeKXCmML}U9b8@a3X
zd@8`>5a9g~parpLo%3UncoC4`uHa&&0;kZ87ArCD#C%21I%X8vWy@Hy%FL0;fLo7c
zf6g|X${4^Q+n+qduVfQnF-m~gGTsxdhB|tov#n>qm@3>D3v!r?xJUVjhcgCe6Mh1f
zaAYSHq2z>0O{z>KDn%kmUj+9M1y91N2q%x9hr1@h*CdcEEC|8lBujFnhxC#Sq>uEI
z0WwI2$VRdWxMecDg`+*SlsiT?lMymXwveslePkQiPBgNE>?FGw*I2?;g<E2TtI06q
zEs!GW2q$p0bWH*H4u(+}My+K4^L~ce%`n+o4lsKdW{hEaYQ2ETgS-te?`4!b8Ktk*
z516}vGEkMsI2nHoebSs~r4$FeauvKZ0wHsV6S;f1HlNVs6V~BVt9U!H8*o*`fJ>y*
z2HO(4H4S!Y=%PpbWu-c}%FnA^eBG9G%ju%8kpuT3vlOD>J#H!7s!g9VOFMkowa4mU
z%isUqS1n(3akuoyiAOs0biI>x^>p#AqX)hr>+Gj%{cX~-c0ae*ZyVZ?1h>R{I_BYv
ze8cl?cw$u?W*>)j&uZ(muXo18<Thx%@Ag^EwLJro+lJpYeI!-w8{)_$#5LTX0PO<^
z?$9A_GxS@;^Z&T7_4P16w(b=fqv5t{?Qibs7Q1V|P4WMwwPUXWxy^fpND`oDz|#Ym
zm_YUn->6Q^Wo4D40?Uqvs%f&n9bI4vvux;g`o3+Cw8ow$<|ccT&bQHq-xbTM5@aX`
zpe~RglgL&jcm~hRtK=ZL`=a<IHIy;)sY;OAePHDoNZOzz!is23aif_LBz8$|rum-j
zM%uze1{RFB8o<SXLc;7u;fJ&ZZCpD@v?6BT+F3iA)5@L~Kz<9Q6H6Iy76yw|RnxT7
zP-2B(!XD3$Y!A}(Jjf|Ek;f*F9RtR^(7@t=Wp}|2gxKRKK+DiO+SNjU*<YTIQzTl3
zqE3b<wk|Iv%P~W|0nyk0`ohK^q+Y30efH^6W6CekYp!+a?h(5{r@u4#i|4G-I(_KK
z?|PpI{W^W|fw#|G+41o@ed5C(`}~h8C+hV4GvB2vzr4RdUs5(j&wdjyKRfmOuBWyi
zt<#@=aPP|(N~Jn|)q3kZ{r=&4V;iWZSKt2I=+vXXuG7H{hhF;jZyqVozdSJi+QOMn
z7pU{Y=|5ch$>st*`P%c@ul{POK)<{@_vEX8Di-KhUe?u5zHp>KH~##@sfTwTtkbP8
zeqzUOcLDJ1KQ5m7g8fn8dFtBVJ!;nJwtqkV#cPjzs@_mQVFUfqzzU#&zGz?;-oT{1
zfk}BIi$e~EnZ`R%Ur~Ud{juHx3IQ(73yC990_}fm{Buu#{mfx-vSdb{lh4Jy#Sm&8
ze1AuRzzbNaAE5R)UG^%zW1WlBffde~z7?ldEGV{uexL>#vZ!UxonNF7CIM#Z!4^a+
zrsI?yDAMv9lRAN4h||qj(Nn&;=-9eAXI6B{hrrS0$S4+#eQYkZ=lxAQ)Mkvntfyfv
z8I=~`FE=Tb=Ik?$qHzE@9zx!IEdYdJ`B5x^UeDKsg-TMHV6hjl@}#$PpV9g601(0-
z!@#LIVS*Q~%e(;TyKr3+dj&;Mk|$1R&7AU%EZ+(%r4=ms*bgB^*~9OF(of^}U>Dj@
zb^tIF!kY00ubix4V_IXYC%VV<1kWUb&0~y~5w;2|gRIQDd51H3J~>-yEL2^yVi`u9
zF$|~&7afFihH+-mbefd3VGyru7{LhmguVOV5JJp@HwN@&fOS4bUWieoh*C-!Rx%6=
k^|u%<hjGeXtoX}Gs3G^Q#Kg)yi7g+&3xm~+3a^O&0faYGdjJ3c

literal 5153
zcmd5AZHya7b-%p!+P=@+M=lNR=9<7ct?zPoX_7`EK3~!#C+M~0t|Bt!mW^j`ym$8Q
zW_RxsTUYW?A=C<~?GGxXC{{%TMF<I02?+^D;zOXq4@y;`Z7G2UMJ_D_1VTln;mvsM
zIKG%4=uXO5GxPS%oA<SE-p=zAEiGXN+CyKyc)HxfFz=9u`o_b<7<8u}2E-sX%cLqZ
zo8pc$h(r8&2Jz>(6pz>`CMJ}Dj+so8vA)&#>;|mHq*%nHIK;v(${{5n0R_%8sQ?P2
zAPSvlQzDXpU&IsKDPwANIxJ0Tt|Kp4vYXLlteH8@we_@Y*`}r0uC6(<IV-zMy5nj{
zhApycAfPieTXi*gJZq*i@(DGkXI<0KRbx;NOKG-FMFC(5`lp~fJ+R_3msrFS{V-o-
z^Nhk#4Cq<I@D#%%9*I7UKrne`jssr>FK~Sfz<d?PT#+pRm$|?tF;8H!KbWs@4^g}X
zm@wH-_7#B!<VFx#i*^xSk$Z?<grXu_LJ|rS*+q`5D=q73iAF%>(-Gc2Z7aYn&|7<%
z2-(|qn+9~n*PtCmV7&t|@J6DO_-pU{cBlB{(J}m?eC@v74|I&##F6;ByRMzPIMIFu
zr`Fz1fA#SrW0-x7y)qm4<XHJ8v@%DCpwo7}GTlg&iQ{9P9%ooyK(%dk#p9C0&__~S
za?}&Ex}j#V$Lq$TC%Lw2IF{*X9y{#`hH2-#K-xAPCnYYL*@c|uahC3J$6^0**nfO9
zB~|dsc*<iZV}d6*T6WeG=a*C)!>**F8_IG@Tqbk!Ta4Wc+Wh;-&=ZbkJ40$NXBxA*
zlg?=N&{55qcTH=^)p5p+pVkee8#tabk(Q17aTj+CXNbh>pta(K6yT}5iqZlsovp-P
zzLSL?-_J(bD8G$;5Pt3KHh#5*?y1tJiASm}06-A>5$Fy%Jb4lt47_X}X10szWeS4J
zU4r8WqmXCw0(refau~VQSXh!C&!~<j&#GzHwCN$yW@ojut1oKu^tep*2dWQ%HEB=L
zk-1Td1#EI0Y?2>rqJhV2%a+F@(}ksOf(Yf{6dLYA&SSG4n@dHqs_W`$O<B~^)RC2e
znC^*IMpv@Aod7%3KCz2kjZ}O^S8Re!lRQyvkdK%)=+yMRphq4u#=Pi6UFIz7vT(-W
zY{HMj0vYK^B~ayr3jV5y1uVEc3O)<&K_PMyiaeY=b`Ek)1z!+=GdIuKPooHmq88MO
z+E6>{K%Hm{>O$S92lb+@s1LUl0!BaDhVDVz(Z|re=sqN)9jG7Ok9MM6=;P=UXn@e<
zAw%>r1rcVF0Rl+sV>qS|EQA1#QA~(pq(T@lyD8=Yiis4WfO(K&;uO<TXa&q5#SBqQ
zTcI5=dnjg@Vmb<)fEl5fQHt48=*kOdFWP&7_{g7Ug#ZIIVjjHK2O+ed;h95BjgIr_
zxF&Q#^iZCG%*6snvB>w*AgpjIbl134rB%6fH<7|uq_`|(!0sYD7xL*sz6{X22q8}8
zahoVls>uL(<XAPe$k#aEXi?STW)@Z7L!*Bi`82Gr^9*Et-}~HMu4!A%r&QCr=xDGp
z-(B^0lPK|{sOpOb>Z%+OHOgstCr#xv&^3|M-qfCk{Tn1_y*(T2w+`PPiJY;!I_5iT
zF%Mb$uGSQ8Nc414pVr{b%ad?F^K&D_ehuEYX?Z#7%Mk1Gvbc$P8E&VyvRw_{H2Y8>
zzom;o`ys=?HyfY>0Ks>=i`fdjE`wE1{QvoVYc7X=PpdXV+-$JaF5jj*SXYZ4nhiG8
z)+SQxY^C4!jkGl{z&-C_sIG*(0GHr{tC!(=3&kU;9oLnCvg*PWxN?h?EmJpKd44bq
z8{*XlaMH$5S{SZ#P@b3L@+gv%q)w5SbT=cXP1A--0E$d{xxSqlY?dThmQO=5>VOIZ
zrse7;6f`m5GwvfNh9)L}FvgWh8Ab2G_6`Ek$87+uK(A{z5;nOf$5=0b+!ZK)MR=n7
zVnSs*xl7-I@b_Y`zVQ<OS_z+Q?L7IJgq6UD;*Z`R|20nF?Z0}gbNil)Kf=$P{_`_y
zFIpvhb@eaj-fiENz(;pXF8uCid<nlhb7pPnACH#sPu@EA)|(SkC7l0p$BCA;QziV1
ztJ$AF_Woo7zx@^K@)v*8Tfzg<Q?GP<|5yo+e(BUdFT4^@;OAcWLI3+)GJy-v_k8X<
zXJ0SjSD#$SmA;`SaOW>SKlAObQwjXSf381U?CVS5m%skk>wkECF@dlA@x8U5es;Ko
z%YQF?dite<34Hm@%wg?~L;~Y;4-f5<jwLYn-htm24vd!Y^pUg4zwHO$Wz{?r9{JN@
zyklnRU;91+nh$^Y^x1Dd53-*70v=IjW)fxMyE3V(%EWJFQrnb6WFYlUnbbSw2pOVe
zAjMm`6^2QhZutl*SqzWRu9m7mBF0tT0bv-Eyw;=x^&@%HkJz?p)AGqCBD`SQ%vo7&
z*$dj5lTj_r3#@8Tz}hldmC~j$w}2s*Y$7t&Q>nAqc`xWkP%2=l3t2tIXVhFuuppdL
zqN^m6${sq8*XCg#549~_Ur7aF9#?Y9K44XWTy@D>@o97ZY*EsrGEyN0fRLUfZ1z-v
zVNk`i8mqn5@e@k@pDkdFa}5SYigJT2cT;3JsED|m0^iDoIDtIfxql!={@1Dy4Uvla
z)-N=`4)FJ~GCRNxkP)hE80m8jigW^IkXtw2iB_!aT{i}0%5!kt80J>D(5m)!4@}0`
z$^i|M{{;*+rzwgTRum||7qSG7D$1D!HS4bg6$P1TMX~$9XEZkKeFP%^PSg#gvbM?J
zDkL<W5DXz+2!wh=VTy&nJ9uZs@dCy|&RVHNA)&70KpD~kJ#?s|%VCmrXo((TL;QaM
DBW$93

diff --git a/spin_texture.py b/spin_texture.py
new file mode 100644
index 0000000..a9f5635
--- /dev/null
+++ b/spin_texture.py
@@ -0,0 +1,59 @@
+import numpy as np
+import matplotlib.pyplot as plt
+
+from momentum import Momentum
+from tSnS import tSnS
+
+
+def spin_texture(Nk=30, bandset=2): 
+
+    """
+    Compute Spin texture in BZ
+    """
+
+    #Pauli matrices and orbital identity
+    sigma_x_tau_0 = np.array([[0,0,1,0],[0,0,0,1],[1,0,0,0],[0,1,0,0]])
+    sigma_y_tau_0 = 1j*np.array([[0,0,-1,0],[0,0,0,-1],[1,0,0,0],[0,1,0,0]])
+    sigma_z_tau_0 = np.array([[1,0,0,0],[0,1,0,0],[0,0,-1,0],[0,0,0,-1]])
+
+
+    #Create Mesh and compute spin polarization
+    mesh        = BZ.Bravais_mesh_simple(Nk)
+    ekf         = np.zeros((Nk**2, model.N))
+    ukf         = np.zeros((Nk**2, model.N, model.N), dtype=complex)
+    spin_pol    = np.zeros((Nk**2, 3, model.N), dtype=complex)
+
+    for i in range(len(mesh)):
+        ekf, ukf        = np.linalg.eigh(model.BH_Hamiltonian(mesh[i], bandset=bandset))
+        spin_pol[i,0,:] = (np.diagonal(ukf.T.conjugate()@sigma_x_tau_0@ukf)).real
+        spin_pol[i,1,:] = (np.diagonal(ukf.T.conjugate()@sigma_y_tau_0@ukf)).real
+        spin_pol[i,2,:] = (np.diagonal(ukf.T.conjugate()@sigma_z_tau_0@ukf)).real
+
+
+    fig, axes = plt.subplots(3, model.N, sharex=True, sharey=True, figsize=(4,3))
+    for i in range(3): 
+        for j in range(model.N): 
+            axes[i,j].scatter(mesh[:,0], mesh[:,1], c=spin_pol[:,i,j], cmap=plt.cm.RdBu, vmin=-1, vmax=1, s=1)
+            axes[i,j].axis('scaled')
+            axes[i,j].set_xticks([])
+            axes[i,j].set_yticks([])
+
+    axes[0,0].set_ylabel(r'$\langle S_x \rangle$')
+    axes[1,0].set_ylabel(r'$\langle S_y \rangle$')
+    axes[2,0].set_ylabel(r'$\langle S_z \rangle$')
+
+    axes[0,0].set_xlabel(r'$n_b = 1$', labelpad=-70)
+    axes[0,1].set_xlabel(r'$n_b = 2$', labelpad=-70)
+    axes[0,2].set_xlabel(r'$n_b = 3$', labelpad=-70)
+    axes[0,3].set_xlabel(r'$n_b = 4$', labelpad=-70)
+
+
+    plt.subplots_adjust(hspace=0, wspace=0)
+    plt.show()
+
+
+if __name__ == '__main__':
+
+    model = tSnS()
+    BZ    = Momentum(model.G1, model.G2)
+    spin_texture()
diff --git a/tSnS.py b/tSnS.py
index 7bfe6be..9f1ed2f 100644
--- a/tSnS.py
+++ b/tSnS.py
@@ -61,25 +61,15 @@ class tSnS:
 
 		if bandset == 1: 
 
-			popt = [-0.28147216773927675, 0.003966634424731179, 
-					0.0029073602294481787, -4.4983518626626615e-07, 
-					-0.004052637332926278, 0.00016072142349528914, 
-					-3.7721284264216035e-05, -0.0009141524777857822, 
-					-0.0011772772860306072, -0.0014951352789262786, 
-					0.0009438170971867384, -9.088937950755711e-05, 
-					-0.000692733427968305, 0.0002266622565665427, 
-					0.0005191476715863555]
+			popt = [-0.28159930344528794, 0.003968154118035716, 0.0030637767514870547, -0.00010069737918143397, 
+					-0.004016436852695664, -0.0008890332864148689, -0.001246870994836976, 0.000848707931621465, 
+					-0.0012602009377081112, -0.0007056692333075815, -0.0005572460653001037, -0.00045337702969447147, np.pi/30]
 
 		else: 
 
-			popt = [-0.32047590155917033, 0.0021542300052420408, 
-					0.001469771721045236, 9.574936436961549e-05, 
-					0.005574114270507838, -0.00026463100214217314, 
-					0.0004064637000665132, 0.0004858014275756446, 
-					0.0006724769177247459, -0.00029328415914492434, 
-					-0.0002932514194152098, 1.0106788563322311e-05, 
-					0.0007016914741439542, 0.0021456042488447995, 
-					0.0018535546929481545]
+			popt = [-0.32044985418524075,   0.002055488150036208,   0.0018580964148542333,  9.451241186629052e-05, 
+					0.0055709631952728005, 0.0005463154794551496, 0.0006627259014835454, -0.000292263072270747, 
+					-0.00029205370932945825, 0.0006839685785136236, -0.0021245180394823893, -0.0017931410640271215, np.pi/8]
 
 
 
@@ -87,11 +77,11 @@ class tSnS:
 		t_AA_0				= popt[0]
 		t_AB_0, t_AB_1 		= popt[1], popt[2]
 		t_AA_x, t_AA_y 		= popt[3], popt[4]
-		t_AB_mx, t_AB_px 	= popt[5], popt[6]
-		t_AB_my, t_AB_py 	= popt[7], popt[8]
-		t_AA_m, t_AA_p 		= popt[9], popt[10]
-		t_AA_2x, t_AA_2y 	= popt[11], popt[12]
-		tR_AB_0, tR_AB_1 	= popt[13], popt[14] 
+		t_AB_my, t_AB_py 	= popt[5], popt[6]
+		t_AA_m, t_AA_p 		= popt[7], popt[8]
+		t_AA_2y 			= popt[9]
+		tR_AB_0, tR_AB_1 	= popt[10], popt[11]
+		phi_AB 				= popt[12] 
 
 
 		##################################################
@@ -102,15 +92,6 @@ class tSnS:
 		ham[:,0,:,0] += t_AA_0*pauli[0]
 		ham[:,1,:,1] += t_AA_0*pauli[0]
 
-		#Nearest-Neighbor AB/BA coupling (1. nn)
-		hop = t_AB_0*(1.+self.phi_lm(k,0,-1))*pauli[0]
-		ham[:,0,:,1] += hop
-		ham[:,1,:,0] += np.conjugate(hop).T
-
-		hop = t_AB_1*(self.phi_lm(k,1,0)+self.phi_lm(k,1,-1))*pauli[0]
-		ham[:,0,:,1] += hop
-		ham[:,1,:,0] += np.conjugate(hop).T
-
 		#Nearest-Neighbor AA/BB coupling (2. nn)
 		hop = t_AA_x*(self.phi_lm(k,1,0))*pauli[0]
 		ham[:,0,:,0] += hop
@@ -150,12 +131,6 @@ class tSnS:
 		ham[:,1,:,1] += np.conjugate(hop)*pauli[0]
 
 		#Second shell (5. nn)
-		hop = t_AA_2x*self.phi_lm(k,2,0)*pauli[0]
-		ham[:,0,:,0] += hop
-		ham[:,0,:,0] += np.conjugate(hop).T
-		ham[:,1,:,1] += hop
-		ham[:,1,:,1] += np.conjugate(hop).T
-
 		hop = t_AA_2y*self.phi_lm(k,0,2)*pauli[0]
 		ham[:,0,:,0] += hop
 		ham[:,0,:,0] += np.conjugate(hop).T
@@ -166,18 +141,10 @@ class tSnS:
 		#AB coupling
 		################
 		#Edge states of second shell (4. nn)
-		hop = t_AB_mx*(self.phi_lm(k,2,0)+self.phi_lm(k,2,-1))*pauli[0]
-		ham[:,0,:,1] += hop
-		ham[:,1,:,0] += np.conjugate(hop).T
-
 		hop = t_AB_my*(self.phi_lm(k,1,1)+self.phi_lm(k,1,-2))*pauli[0]
 		ham[:,0,:,1] += hop
 		ham[:,1,:,0] += np.conjugate(hop).T
 
-		hop = t_AB_px*(self.phi_lm(k,-1,0)+self.phi_lm(k,-1,-1))*pauli[0]
-		ham[:,0,:,1] += hop
-		ham[:,1,:,0] += np.conjugate(hop).T
-
 		hop = t_AB_py*(self.phi_lm(k,0,1)+self.phi_lm(k,0,-2))*pauli[0]
 		ham[:,0,:,1] += hop
 		ham[:,1,:,0] += np.conjugate(hop).T			
@@ -203,6 +170,23 @@ class tSnS:
 		ham[:,0,:,1] += hop
 		ham[:,1,:,0] += np.conjugate(hop).T
 
+
+		######################################
+		#SOC Kane-Mele Coupling
+		######################################
+		#Nearest-neighbors (1. nn)
+		km = np.array([[np.exp(1j*phi_AB),0],[0,np.exp(1j*phi_AB).conjugate()]])
+		hop = t_AB_0*(1.+self.phi_lm(k,0,-1))*km
+		ham[:,0,:,1] += hop
+		ham[:,1,:,0] += np.conjugate(hop).T
+
+		km = np.array([[np.exp(-1j*phi_AB),0],[0,np.exp(-1j*phi_AB).conjugate()]])
+		hop = t_AB_1*(self.phi_lm(k,1,0)+self.phi_lm(k,1,-1))*km
+		ham[:,0,:,1] += hop
+		ham[:,1,:,0] += np.conjugate(hop).T
+
 		ham = ham.reshape((4,4))
 
+		assert np.allclose(ham, ham.conjugate().T, atol=1e-12)
+
 		return ham
-- 
GitLab