From b9bc56d7a3f184f61447c38ce0ae4c4081b467f9 Mon Sep 17 00:00:00 2001 From: Eric Giguere Date: Tue, 23 Apr 2024 14:07:02 -0400 Subject: [PATCH] Fix raise on too large min_step in krylov --- qutip/solver/integrator/krylov.py | 14 ++++++++------ qutip/tests/solver/test_sesolve.py | 16 +++++++++++++--- 2 files changed, 21 insertions(+), 9 deletions(-) diff --git a/qutip/solver/integrator/krylov.py b/qutip/solver/integrator/krylov.py index 5559020982..7e2d9e1d79 100644 --- a/qutip/solver/integrator/krylov.py +++ b/qutip/solver/integrator/krylov.py @@ -51,7 +51,7 @@ def _prepare(self): krylov_tridiag, krylov_basis = \ self._lanczos_algorithm(rand_ket(N).data) if ( - krylov_tridiag.shape[0] < self.options["krylov_dim"] + krylov_tridiag.shape[0] < krylov_dim or krylov_tridiag.shape[0] == N ): self._max_step = np.inf @@ -138,20 +138,22 @@ def krylov_error(t): self._compute_psi(t, *reduced_state) ) / self.options["atol"]) - dt = self.options["min_step"] + # Under 0 will cause an infinite loop in the while loop bellow. + dt = max(self.options["min_step"], 1e-14) + max_step = max(self.options["max_step"], dt) err = krylov_error(dt) if err > 0: - ValueError( + raise ValueError( f"With the krylov dim of {self.options['krylov_dim']}, the " f"error with the minimum step {dt} is {err}, higher than the " f"desired tolerance of {self.options['atol']}." ) - while krylov_error(dt * 10) < 0 and dt < self.options["max_step"]: + while krylov_error(dt * 10) < 0 and dt < max_step: dt *= 10 - if dt > self.options["max_step"]: - return self.options["max_step"] + if dt > max_step: + return max_step sol = root_scalar(f=krylov_error, bracket=[dt, dt * 10], method="brentq", xtol=self.options['atol']) diff --git a/qutip/tests/solver/test_sesolve.py b/qutip/tests/solver/test_sesolve.py index cfa5ee86ea..af247c0aec 100644 --- a/qutip/tests/solver/test_sesolve.py +++ b/qutip/tests/solver/test_sesolve.py @@ -301,9 +301,19 @@ def test_krylovsolve(always_compute_step): e_op.dims = H.dims tlist = np.linspace(0, 1, 11) ref = sesolve(H, psi0, tlist, e_ops=[e_op]).expect[0] - options = {"always_compute_step", always_compute_step} - krylov_sol = krylovsolve(H, psi0, tlist, 20, e_ops=[e_op]).expect[0] - np.testing.assert_allclose(ref, krylov_sol) + options = {"always_compute_step": always_compute_step} + krylov_sol = krylovsolve(H, psi0, tlist, 20, e_ops=[e_op], options=options) + np.testing.assert_allclose(ref, krylov_sol.expect[0]) + + +def test_krylovsolve_error(): + H = qutip.rand_herm(256, density=0.2) + psi0 = qutip.basis([256], [255]) + tlist = np.linspace(0, 1, 11) + options = {"min_step": 1e10} + with pytest.raises(ValueError) as err: + krylovsolve(H, psi0, tlist, 20, options=options) + assert "error with the minimum step" in str(err.value) def test_feedback():