How to create local identity operator at site j? Indirectly, one can use the identity "S^zjS^zj" or "projUp+projDn", which have built-in functions like
op("projDn",sites,j). Is there a more standard way?
Moreover, how about high-spin sites, fermion/electron sites?