of the same order of magnitude ≈ O(N ). As a reference, a straightforward approach requires two steps: (I) evaluating the determinant of the canonical orbital overlap matrix for every S ab ij matrix (Eq. ( )), which takes O(N 5 o N 2 v ); (II) computing the summation in Eq. ( ) for