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

2*h + U	(-I/2)*(-3*I + Sqrt[3])*t	0
(I/2)*(3*I + Sqrt[3])*t	2*h + U	((-3 - I*Sqrt[3])*t)/Sqrt[2]
0	(I*(3*I + Sqrt[3])*t)/Sqrt[2]	2*h + U