ffsim.linalg.double_factorized

ffsim.linalg.double_factorized(two_body_tensor, *, tol=1e-08, max_vecs=None, optimize=False, method='L-BFGS-B', callback=None, options=None, diag_coulomb_indices=None, return_optimize_result=False, cholesky=True)[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=0}^{L - 1} \sum_{k\ell} Z^{(t)}_{k\ell} U^{(t)}_{pk} U^{(t)}_{qk} U^{(t)}_{r\ell} U^{(t)}_{s\ell},\]

where each \(Z^{(t)}\) is a real symmetric matrix representing a diagonal Coulomb operator and each \(U^{(t)}\) is a unitary matrix representing an orbital rotation.

The number of terms \(L\) 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 \(L\). 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 default behavior of this routine is to perform a straightforward “exact” factorization of the two-body tensor based on a nested eigenvalue decomposition, and then truncate the terms based on the values of tol and max_vecs. If optimize is set to True, then the entries of the resulting tensors (the diagonal Coulomb matrices and orbital rotations) are further optimized with scipy.optimize.minimize to reduce the error in the factorization. The diagonal Coulomb matrices returned by the optimization can be optionally constrained to have only certain elements allowed to be nonzero. This is achieved by passing the diag_coulomb_indices parameter, which is a list of matrix entry indices (integer pairs) specifying where the diagonal Coulomb matrices are allowed to be nonzero. Since the diagonal Coulomb matrices are symmetric, only upper triangular indices should be given, i.e., pairs \((i, j)\) where \(i \leq j\).

References

Note: Currently, only real-valued two-body tensors are supported.

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

  • tol (float) – Tolerance for 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. This argument overrides tol.

  • optimize (bool) – Whether to optimize the tensors returned by the decomposition to to minimize the error in the factorization.

  • method (str) – The optimization method. See the documentation of scipy.optimize.minimize for possible values. This argument is ignored if optimize is set to False.

  • callback – Callback function for the optimization. See the documentation of scipy.optimize.minimize for usage. This argument is ignored if optimize is set to False.

  • options (dict | None) – Options for the optimization. See the documentation of scipy.optimize.minimize for usage. This argument is ignored if optimize is set to False.

  • diag_coulomb_indices (list[tuple[int, int]] | None) – Allowed indices for nonzero values of the diagonal Coulomb matrices. Matrix entries corresponding to indices not in this list will be set to zero. This list should contain only upper triangular indices, i.e., pairs \((i, j)\) where \(i \leq j\). This argument is ignored if optimize is set to False.

  • return_optimize_result (bool) – Whether to also return the OptimizeResult returned by scipy.optimize.minimize. This argument is ignored if optimize is set to False.

  • cholesky (bool) – Whether to perform the factorization using a modified Cholesky decomposition. If False, a full eigenvalue decomposition is used instead, which can be much more expensive. This argument is ignored if optimize is set to True.

Returns:

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 (L, n, n) where L is the number of terms in the decomposition and n is the number of orbitals. If optimize and return_optimize_result are both set to True, the OptimizeResult returned by scipy.optimize.minimize is also returned.

Raises:

ValueError – diag_coulomb_indices contains lower triangular indices.