Particle conserving Hamiltonians dont need this term - just specify the sector you want by the initial state as you said. If they, for example, only conserve parity due to a pinning term (like equation 17 in your arxiv link), then you’ll need that term to control filling.
Regardless, with OpSum
that term is easy to add, e.g. a 1d Hubbard model
os = OpSum()
for b in 1:(N - 1)
os += -t, "Cdagup", b, "Cup", b + 1
os += -t, "Cdagup", b + 1, "Cup", b
os += -t, "Cdagdn", b, "Cdn", b + 1
os += -t, "Cdagdn", b + 1, "Cdn", b
end
for i in 1:N
os += U, "Nupdn", i
os += -mu, "Ntot", i # The term you're asking about
end
H = MPO(os, sites)