A FermionOperator that conserves particle number and spin Z can be converted to a LinearOperator, and the current implementation is here:

Currently, the implementation uses the functions des_a, des_b, cre_a, cre_b from pyscf.fci.addons. These functions don't preserve symmetries, and indeed they change the shape of the vector. The efficiency can be vastly improved by implementing a function that can directly apply a number- and spin-Z-conserving term to a vector. This should be implemented in Rust, but it can take as an argument a linkstr_index produced by pyscf.fci.cistring.gen_linkstr_index (or rather, the cached wrapper in ffsim).

Recommended approach:

1. Rewrite _apply_fermion_term in Python to apply the fermion term in a more efficient manner.
2. Port the function to Rust.

After doing this, we can increase the maximum norb covered by some tests.