Matrix H in Subspace 137 with Dimension 16

H is the Hamiltonian. Here we set µ=0. For grand-canonical calculations add -6µ to the main diagonal.
2*h + U (-I/12)*(-13*I + Sqrt[3])*t ((9*I - 5*Sqrt[3])*t)/(6*Sqrt[2]) 0 ((-9*I + Sqrt[3])*t)/(6*Sqrt[2]) 0 0 0 0 0 0 0 0 (I/3)*Sqrt[2]*(I + Sqrt[3])*t 0 ((3*I + Sqrt[3])*t)/12
(I/12)*(13*I + Sqrt[3])*t 2*h + U 0 Sqrt[3/2]*t 0 ((3*I + Sqrt[3])*t)/(2*Sqrt[2]) 0 0 0 0 -((9*I + Sqrt[3])*t)/(6*Sqrt[2]) -((3*I + 7*Sqrt[3])*t)/(6*Sqrt[2]) ((I/3)*(I + Sqrt[3])*t)/Sqrt[2] 0 -((-3*I + Sqrt[3])*t)/4 0
-((9*I + 5*Sqrt[3])*t)/(6*Sqrt[2]) 0 2*(h + U) ((3 + I*Sqrt[3])*t)/2 0 0 0 0 0 0 t -t (-2*t)/Sqrt[3] 0 ((1 - I*Sqrt[3])*t)/(2*Sqrt[2]) 0
0 Sqrt[3/2]*t ((3 - I*Sqrt[3])*t)/2 2*(h + U) 0 0 0 0 0 0 0 0 0 0 0 (-I)*Sqrt[3/2]*t
((9*I + Sqrt[3])*t)/(6*Sqrt[2]) 0 0 0 2*(h + U) (t - I*Sqrt[3]*t)/2 0 t 0 0 0 0 (-2*t)/Sqrt[3] 0 ((3 + I*Sqrt[3])*t)/(2*Sqrt[2]) 0
0 ((-3*I + Sqrt[3])*t)/(2*Sqrt[2]) 0 0 (t + I*Sqrt[3]*t)/2 2*(h + U) -(Sqrt[2]*t) 0 (t - I*Sqrt[3]*t)/2 Sqrt[2]*t 0 0 0 0 0 ((1 - I*Sqrt[3])*t)/(2*Sqrt[2])
0 0 0 0 0 -(Sqrt[2]*t) 2*(h + U) -(Sqrt[2]*t) 0 0 0 0 0 0 (-1 - I*Sqrt[3])*t 0
0 0 0 0 t 0 -(Sqrt[2]*t) 2*(h + U) (t - I*Sqrt[3]*t)/2 ((1 - I*Sqrt[3])*t)/Sqrt[2] 0 0 0 0 0 ((-1 - I*Sqrt[3])*t)/Sqrt[2]
0 0 0 0 0 (t + I*Sqrt[3]*t)/2 0 (t + I*Sqrt[3]*t)/2 2*(h + U) 0 t -t 0 0 0 0
0 0 0 0 0 Sqrt[2]*t 0 ((1 + I*Sqrt[3])*t)/Sqrt[2] 0 2*(h + U) 0 0 0 0 0 0
0 -((-9*I + Sqrt[3])*t)/(6*Sqrt[2]) t 0 0 0 0 0 t 0 2*h + U 0 0 ((3*I + Sqrt[3])*t)/3 0 ((3 - I*Sqrt[3])*t)/(2*Sqrt[2])
0 ((3*I - 7*Sqrt[3])*t)/(6*Sqrt[2]) -t 0 0 0 0 0 -t 0 0 2*h + U 0 (-2*t)/Sqrt[3] 0 ((1 + I*Sqrt[3])*t)/(2*Sqrt[2])
0 ((-I/3)*(-I + Sqrt[3])*t)/Sqrt[2] (-2*t)/Sqrt[3] 0 (-2*t)/Sqrt[3] 0 0 0 0 0 0 0 2*h + U ((-1 - (2*I)*Sqrt[3])*t)/3 0 ((3*I + Sqrt[3])*t)/(3*Sqrt[2])
(-I/3)*Sqrt[2]*(-I + Sqrt[3])*t 0 0 0 0 0 0 0 0 0 ((-3*I + Sqrt[3])*t)/3 (-2*t)/Sqrt[3] ((-1 + (2*I)*Sqrt[3])*t)/3 2*h + U 0 0
0 -((3*I + Sqrt[3])*t)/4 ((1 + I*Sqrt[3])*t)/(2*Sqrt[2]) 0 ((3 - I*Sqrt[3])*t)/(2*Sqrt[2]) 0 I*(I + Sqrt[3])*t 0 0 0 0 0 0 0 2*h + U (-I/4)*(-5*I + Sqrt[3])*t
((-3*I + Sqrt[3])*t)/12 0 0 I*Sqrt[3/2]*t 0 ((1 + I*Sqrt[3])*t)/(2*Sqrt[2]) 0 (I*(I + Sqrt[3])*t)/Sqrt[2] 0 0 ((3 + I*Sqrt[3])*t)/(2*Sqrt[2]) ((1 - I*Sqrt[3])*t)/(2*Sqrt[2]) ((-3*I + Sqrt[3])*t)/(3*Sqrt[2]) 0 (I/4)*(5*I + Sqrt[3])*t 2*h + U