 Solution

Use $\hslash =k_{\mathrm{B}}=1$ and dimensionless spin operators. At one site, the canonical anticommutation relations imply the on-site spin-square normal-ordering identityIndeed the empty and doubly occupied states have spin zero, while either singly occupied state has ${\mathbf{S}}^2=3/4$. Thus the interaction in normal-ordered form is $3U{\sum }_n{c}_{n\uparrow }^{†}{c}_{n\downarrow }^{†}{c}_{n\downarrow }{c}_{n\uparrow }$. There is no residual $3Un/2$ chemical-potential shift: that contraction has cancelled.

Insert fermionic coherent state resolutions between imaginary-time evolution factors. With ${\xi }_k={ϵ}_k-\mu$ and $D_n={\bar{\psi }}_{n\uparrow }{\bar{\psi }}_{n\downarrow }{\psi }_{n\downarrow }{\psi }_{n\uparrow }$, the grand-canonical coherent-state action isThe independent barred and unbarred Grassmann variables satisfy antiperiodic boundary conditions at $0,\beta$. Adjacent time-slice labels in the normal-ordered symbol fix the continuum prescription and thermal normalization. The interaction coefficient is $3U$, not the $U$ of a differently normalized Hubbard model.