There is a bug in the block-diagonal schur matrix computation inline 473. It should betmodal = tmodal @ rperm.T. Otherwisetmodal @ amodal @ np.linalg.inv(tmodal)$\neq$a wherea is the input to the function.