H is the Hamiltonian. Here we set µ=0. For grand-canonical calculations add -6µ to the main diagonal.

U	-(Sqrt[2]*t)	-t	0
-(Sqrt[2]*t)	U	0	-2*t
-t	0	2*U	-(Sqrt[2]*t)
0	-2*t	-(Sqrt[2]*t)	2*U