Hi. sorry for not replying for so long. Yes, it seems to be the case that the sign arises from anticommutation rules (knowing this, I found a bug in my code and now it works fine).
For now, the current approach - "if a name starts with C, it is fermionic" - works for me. What I think is more needed is a little more documentation. For example, I got a little confused about the order in which the operators are applied. They are sorted according to site index - but I am not sure if in decreasing or increasing order.
Just to explain why I need my own site class - I have a system in which I have two orbitals per site, but there is no tunneling within the site (there is just an exchange interaction). Apparently in such a situation writing the site as two Hubbard sites is not a good idea, because at some point You have to divide Your real physical site in half.
Anyway, thanks for Your reply, it was really helpful.