double_factorized(two_body_tensor, *, error_threshold=1e-08, max_vecs=None)[source]#

Double-factorized decomposition of a two-body tensor.

The double-factorized decomposition is a representation of a two-body tensor \(h_{pqrs}\) as

\[h_{pqrs} = \sum_{t=1}^N \sum_{k\ell} U^{t}_{pk} U^{t}_{qk} Z^{t}_{k\ell} U^{t}_{r\ell} U^{t}_{s\ell}\]

Here each \(Z^{(t)}\) is a real symmetric matrix, referred to as a “diagonal Coulomb matrix,” and each \(U^{t}\) is a unitary matrix, referred to as an “orbital rotation.”

The number of terms \(N\) in the decomposition depends on the allowed error threshold. A larger error threshold may yield a smaller number of terms. Furthermore, the max_vecs parameter specifies an optional upper bound on \(N\). The max_vecs parameter is always respected, so if it is too small, then the error of the decomposition may exceed the specified error threshold.

The input tensor is assumed to be positive definite when reshaped into a matrix.


No checking is performed to verify whether the input tensor is positive definite. If the input tensor is not positive definite, then the decomposition returned will be invalid.


  • two_body_tensor (np.ndarray) – The two-body tensor to decompose.

  • error_threshold (float) – Threshold for allowed error in the decomposition. The error is defined as the maximum absolute difference between an element of the original tensor and the corresponding element of the reconstructed tensor.

  • max_vecs (int | None) – An optional limit on the number of terms to keep in the decomposition of the two-body tensor.


The diagonal Coulomb matrices and the orbital rotations. Each list of matrices is collected into a numpy array, so this method returns a tuple of two numpy arrays, the first containing the diagonal Coulomb matrices and the second containing the orbital rotations. Each numpy array will have shape (t, n, n) where t is the rank of the decomposition and n is the number of orbitals.

Return type:

tuple[np.ndarray, np.ndarray]