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

-2*h + U	(-I/12)*(-11*I + Sqrt[3])*t	0	((3*I + Sqrt[3])*t)/2	0	-((9*I + 5*Sqrt[3])*t)/(6*Sqrt[2])	((-3*I + 5*Sqrt[3])*t)/(6*Sqrt[2])	0	0	0	0	((1 - I*Sqrt[3])*t)/(3*Sqrt[2])	0	0	0	((-3*I + Sqrt[3])*t)/4
(I/12)*(11*I + Sqrt[3])*t	-2*h + U	0	0	((9*I - 5*Sqrt[3])*t)/6	0	0	0	0	0	(-I/3)*Sqrt[2]*(-I + Sqrt[3])*t	0	-(((7 + 3*I)*Sqrt[18 - 6*Sqrt[3]] - (3 - 3*I)*Sqrt[6 - 2*Sqrt[3]])*t)/(6*Sqrt[2]*(-3 + Sqrt[3]))	((1/6 + I/6)*(-3 - (2 - 5*I)*Sqrt[3])*t)/Sqrt[3 + Sqrt[3]]	-((3*I + Sqrt[3])*t)/12	0
0	0	-2*(h - U)	((-I/2)*(-5*I + Sqrt[3])*t)/Sqrt[2]	0	0	0	0	((3*I + Sqrt[3])*t)/(2*Sqrt[2])	0	0	0	0	0	0	Sqrt[2]*(1 + I*Sqrt[3])*t
((-3*I + Sqrt[3])*t)/2	0	((I/2)*(5*I + Sqrt[3])*t)/Sqrt[2]	-2*(h - U)	((-1 - (3*I)*Sqrt[3])*t)/4	0	0	((-3*I + Sqrt[3])*t)/(2*Sqrt[2])	0	-((-3*I + Sqrt[3])*t)/4	0	0	0	0	((3 - I*Sqrt[3])*t)/2	0
0	-((9*I + 5*Sqrt[3])*t)/6	0	((-1 + (3*I)*Sqrt[3])*t)/4	-2*(h - U)	-(t/Sqrt[2])	t/Sqrt[2]	0	-((3*I + Sqrt[3])*t)/4	0	0	-2*Sqrt[2/3]*t	0	0	0	(t - I*Sqrt[3]*t)/2
((9*I - 5*Sqrt[3])*t)/(6*Sqrt[2])	0	0	0	-(t/Sqrt[2])	-2*h + U	0	0	0	Sqrt[3/2]*t	(-2*t)/Sqrt[3]	0	0	0	((-I/2)*(-I + Sqrt[3])*t)/Sqrt[2]	0
((3*I + 5*Sqrt[3])*t)/(6*Sqrt[2])	0	0	0	t/Sqrt[2]	0	-2*h + U	0	0	-(Sqrt[3/2]*t)	-((-3*I + Sqrt[3])*t)/3	0	0	0	((I/2)*(3*I + Sqrt[3])*t)/Sqrt[2]	0
0	0	0	((3*I + Sqrt[3])*t)/(2*Sqrt[2])	0	0	0	-2*(h - U)	((I/2)*(3*I + Sqrt[3])*t)/Sqrt[2]	0	0	0	0	0	0	0
0	0	((-3*I + Sqrt[3])*t)/(2*Sqrt[2])	0	-((-3*I + Sqrt[3])*t)/4	0	0	((-I/2)*(-3*I + Sqrt[3])*t)/Sqrt[2]	-2*(h - U)	(-I/4)*(-3*I + Sqrt[3])*t	0	0	0	0	-(Sqrt[3]*t)	0
0	0	0	-((3*I + Sqrt[3])*t)/4	0	Sqrt[3/2]*t	-(Sqrt[3/2]*t)	0	(I/4)*(3*I + Sqrt[3])*t	-2*(h - U)	0	0	0	0	0	0
0	(I/3)*Sqrt[2]*(I + Sqrt[3])*t	0	0	0	(-2*t)/Sqrt[3]	-((3*I + Sqrt[3])*t)/3	0	0	0	-2*h + U	(-I/6)*(-5*I + Sqrt[3])*t	0	0	0	0
((1 + I*Sqrt[3])*t)/(3*Sqrt[2])	0	0	0	-2*Sqrt[2/3]*t	0	0	0	0	0	(I/6)*(5*I + Sqrt[3])*t	-2*h + U	((1 + I)*(-6 - 9*I + (5 + 4*I)*Sqrt[3])*t)/(Sqrt[6]*(3 - Sqrt[3])^(3/2))	-((-3 + 6*I + (2 + 3*I)*Sqrt[3])*t)/(3*Sqrt[2*(3 + Sqrt[3])])	-((3*I + Sqrt[3])*t)/(3*Sqrt[2])	0
0	(((-7 + 3*I)*Sqrt[18 - 6*Sqrt[3]] + (3 + 3*I)*Sqrt[6 - 2*Sqrt[3]])*t)/(6*Sqrt[2]*(-3 + Sqrt[3]))	0	0	0	0	0	0	0	0	0	((1 + I)*(-9 - 6*I + (4 + 5*I)*Sqrt[3])*t)/(Sqrt[18 - 6*Sqrt[3]]*(-3 + Sqrt[3]))	-2*h	0	0	((1 + 3*I)/Sqrt[12 - 4*Sqrt[3]] + ((1 + I)*Sqrt[9 - 3*Sqrt[3]])/(6 - 2*Sqrt[3]))*t
0	((-1/6 + I/6)*(3 + (2 + 5*I)*Sqrt[3])*t)/Sqrt[3 + Sqrt[3]]	0	0	0	0	0	0	0	0	0	((3 + 6*I - (2 - 3*I)*Sqrt[3])*t)/(3*Sqrt[2*(3 + Sqrt[3])])	0	-2*h	0	((-1 + 3*I + (1 - I)*Sqrt[3])*t)/(2*Sqrt[3 + Sqrt[3]])
0	-((-3*I + Sqrt[3])*t)/12	0	((3 + I*Sqrt[3])*t)/2	0	((I/2)*(I + Sqrt[3])*t)/Sqrt[2]	((-I/2)*(-3*I + Sqrt[3])*t)/Sqrt[2]	0	-(Sqrt[3]*t)	0	0	-((-3*I + Sqrt[3])*t)/(3*Sqrt[2])	0	0	-2*h + U	(-I/4)*(-3*I + Sqrt[3])*t
((3*I + Sqrt[3])*t)/4	0	Sqrt[2]*(1 - I*Sqrt[3])*t	0	(t + I*Sqrt[3]*t)/2	0	0	0	0	0	0	0	((-1/2 + I/2)*(Sqrt[18 - 6*Sqrt[3]] + (2 - I)*Sqrt[6 - 2*Sqrt[3]])*t)/(Sqrt[2]*(-3 + Sqrt[3]))	((1/2 + I/2)*(-2 - I + Sqrt[3])*t)/Sqrt[3 + Sqrt[3]]	(I/4)*(3*I + Sqrt[3])*t	-2*h + U