From 8fc0378f78af2e1ca568bd5cba74778d554d2b6d Mon Sep 17 00:00:00 2001 From: "Kevin J. Sung" Date: Wed, 20 Mar 2024 23:57:28 -0400 Subject: [PATCH] fix bug where givens decomposition modified the original matrix (#104) --- python/ffsim/linalg/givens.py | 2 +- tests/gates/orbital_rotation_test.py | 44 ++++++++++++++++++++++--- tests/linalg/givens_test.py | 33 +++++++++++++++++++ tests/test_data/orbital_rotation-0.npy | Bin 0 -> 1152 bytes 4 files changed, 74 insertions(+), 5 deletions(-) create mode 100644 tests/test_data/orbital_rotation-0.npy diff --git a/python/ffsim/linalg/givens.py b/python/ffsim/linalg/givens.py index fc717a4bf..2c588bfce 100644 --- a/python/ffsim/linalg/givens.py +++ b/python/ffsim/linalg/givens.py @@ -72,7 +72,7 @@ def givens_decomposition( ) -> tuple[list[tuple[tuple[float, complex], tuple[int, int]]], np.ndarray]: """Givens rotation decomposition of a unitary matrix.""" n, _ = mat.shape - current_matrix = mat.astype(complex, copy=False) + current_matrix = mat.astype(complex, copy=True) left_rotations = [] right_rotations: list[tuple[tuple[float, complex], tuple[int, int]]] = [] diff --git a/tests/gates/orbital_rotation_test.py b/tests/gates/orbital_rotation_test.py index b20d79871..b4ac71368 100644 --- a/tests/gates/orbital_rotation_test.py +++ b/tests/gates/orbital_rotation_test.py @@ -13,6 +13,7 @@ from __future__ import annotations import itertools +from pathlib import Path import numpy as np import pyscf @@ -89,7 +90,7 @@ def test_apply_orbital_rotation_spin(dtype: type, atol: float): np.testing.assert_allclose(result, expected, atol=atol) -def test_apply_orbital_rotation_no_side_effects(): +def test_apply_orbital_rotation_no_side_effects_vec(): """Test applying orbital basis change doesn't modify the original vector.""" norb = 5 rng = np.random.default_rng() @@ -262,9 +263,6 @@ def test_apply_orbital_rotation_compose(): np.testing.assert_allclose(result, expected_state) -@pytest.mark.skip( - reason="Flaky test. See https://github.com/qiskit-community/ffsim/issues/97" -) def test_apply_orbital_rotation_nitrogen(): """Test a special case that was found to cause issues.""" mol = pyscf.gto.Mole() @@ -300,3 +298,41 @@ def test_apply_orbital_rotation_nitrogen(): ) expected = scipy.sparse.linalg.expm_multiply(op, original_vec, traceA=1) np.testing.assert_allclose(result, expected, atol=1e-12) + + +def test_apply_orbital_rotation_no_side_effects_mat(): + """Test applying orbital rotation doesn't modify the original matrix.""" + norb = 8 + nelec = 5, 5 + vec = ffsim.hartree_fock_state(norb, nelec) + + rng = np.random.default_rng() + for _ in range(5): + mat = ffsim.random.random_unitary(norb, seed=rng) + original_mat = mat.copy() + _ = ffsim.apply_orbital_rotation(vec, mat, norb, nelec) + + assert ffsim.linalg.is_unitary(original_mat) + assert ffsim.linalg.is_unitary(mat) + np.testing.assert_allclose(mat, original_mat, atol=1e-12) + + +def test_apply_orbital_rotation_no_side_effects_special_case(): + """Test applying orbital rotation doesn't modify the original matrix.""" + datadir = Path(__file__).parent.parent / "test_data" + filepath = datadir / "orbital_rotation-0.npy" + + with open(filepath, "rb") as f: + mat = np.load(f) + assert ffsim.linalg.is_unitary(mat, atol=1e-12) + + norb = 8 + nelec = 5, 5 + vec = ffsim.hartree_fock_state(norb, nelec) + + original_mat = mat.copy() + _ = ffsim.apply_orbital_rotation(vec, mat, norb, nelec) + + assert ffsim.linalg.is_unitary(original_mat) + assert ffsim.linalg.is_unitary(mat) + np.testing.assert_allclose(mat, original_mat, atol=1e-12) diff --git a/tests/linalg/givens_test.py b/tests/linalg/givens_test.py index 14b2f6a66..00468ac94 100644 --- a/tests/linalg/givens_test.py +++ b/tests/linalg/givens_test.py @@ -12,6 +12,8 @@ from __future__ import annotations +from pathlib import Path + import numpy as np from scipy.linalg.lapack import zrot @@ -36,3 +38,34 @@ def test_givens_decomposition(): reconstructed = reconstructed.T np.testing.assert_allclose(reconstructed, mat) + + +def test_givens_decomposition_no_side_effects(): + """Test that the Givens decomposition doesn't modify the original matrix.""" + norb = 8 + rng = np.random.default_rng() + for _ in range(5): + mat = ffsim.random.random_unitary(norb, seed=rng) + original_mat = mat.copy() + _ = givens_decomposition(mat) + + assert ffsim.linalg.is_unitary(original_mat) + assert ffsim.linalg.is_unitary(mat) + np.testing.assert_allclose(mat, original_mat, atol=1e-12) + + +def test_givens_decomposition_no_side_effects_special_case(): + """Test that the Givens decomposition doesn't modify the original matrix.""" + datadir = Path(__file__).parent.parent / "test_data" + filepath = datadir / "orbital_rotation-0.npy" + + with open(filepath, "rb") as f: + mat = np.load(f) + assert ffsim.linalg.is_unitary(mat, atol=1e-12) + + original_mat = mat.copy() + _ = givens_decomposition(mat) + + assert ffsim.linalg.is_unitary(original_mat) + assert ffsim.linalg.is_unitary(mat) + np.testing.assert_allclose(mat, original_mat, atol=1e-12) diff --git a/tests/test_data/orbital_rotation-0.npy b/tests/test_data/orbital_rotation-0.npy new file mode 100644 index 0000000000000000000000000000000000000000..4a570d497f8282fd85c5829499265577c6a7fa0f GIT binary patch literal 1152 zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlWb_FuA`uymS0p-l$aNvUzCyxl5k7RDNY57 z7iT0EqyqUG7CH(RnmP)#3giMV2C~2hd$P1)(R9MPK5gEDXK3avuz@H|*aM*tK=~C= z8r7hNKj$FoS{fmA!x1!f#~|Vi5qr);?AwoO!Y&&pzq5B+tVwqL(M_xUjJ}eSLM2u_78su6n~j~-G-HCM)zO-SNmfn zo?c!0==h!qp3C=0B<~`0dPihLpT@8oP%sFgxy!!LPYDWdjy2PGro4#B7PrRJH zYLe+hn=FI8>;Hu>?0?|lk-g!<3^aWV3=h}lZ+WS-$^P?F4@Y<7J2rEk&-asOJhRM&%NYrpZ&J* zFb|uX({1$PGp$OGF1FuxQ#vZKei5oZU@%BdE?>g*@V$M{`Y&Y%7fINPpFcX&PBkH>5vAr2Eq zggCvCQr`IbH5zqNlggtq_p%`59;-YxqN z5E?DrFfh!IvZ~)ul5Q^>4W{?I28&*Ql73@p!S`B z`U~NUH-X!P%Z=~se+!}Q*XOUdE4h1n|5^wQ(T`Bh5H8_rz4VX4{!I|te%iB(BDFX7 V*e`(4`^#8Zwf4@q0;4hX0|3|-Kso>b literal 0 HcmV?d00001