From 17d39896b07979ecba40eecfc252f8f5e56f4a61 Mon Sep 17 00:00:00 2001 From: Paul Menczel Date: Thu, 25 Sep 2025 16:36:40 +0900 Subject: [PATCH 1/6] Hopefully fix loky MPI deadlock --- qutip/solver/parallel.py | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/qutip/solver/parallel.py b/qutip/solver/parallel.py index 425a317f21..dd9e84d41d 100644 --- a/qutip/solver/parallel.py +++ b/qutip/solver/parallel.py @@ -483,6 +483,26 @@ def mpi_pmap(task, values, task_args=None, task_kwargs=None, from mpi4py.futures import MPIPoolExecutor + # If loky is used (see loky_pmap), it starts two "resource tracker" child + # processes that are normally never shut down. However, using mpi4py + # when these child processes are still alive may lead to a deadlock. + # -- WE RECOMMEND NOT TO MIX loky_pmap AND mpi_pmap IN THE SAME PROGRAM -- + # However, it happens in our test runs. For this reason, we here try to + # shut down the resource tracker processes manually + try: + import loky + loky_resource_tracker = loky.backend.resource_tracker._resource_tracker + if loky_resource_tracker: + loky_resource_tracker._stop() + + if hasattr(multiprocessing, 'resource_tracker'): + mp_resource_tracker =\ + multiprocessing.resource_tracker._resource_tracker + if mp_resource_tracker: + mp_resource_tracker._stop() + except ImportError: + pass + # If the provided num_cpus is None, we use the default value instead. # We thus intentionally make it impossible to call # MPIPoolExecutor(max_workers=None, ...) From 3a936e4bbab51c729fe6004887fa82895cfb1d98 Mon Sep 17 00:00:00 2001 From: Paul Menczel Date: Thu, 25 Sep 2025 16:37:05 +0900 Subject: [PATCH 2/6] Correct calculation of root mean square error --- qutip/utilities.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/qutip/utilities.py b/qutip/utilities.py index 52ceaa55fd..79f18a9afb 100644 --- a/qutip/utilities.py +++ b/qutip/utilities.py @@ -490,7 +490,7 @@ def _rmse(fun, xdata, ydata, params): if (yhat == ydata).all(): return 0 return ( - np.sqrt(np.mean((yhat - ydata) ** 2) / len(ydata)) + np.sqrt(np.mean(np.abs(yhat - ydata) ** 2)) / (np.max(ydata) - np.min(ydata)) ) From 9c7318453c71212067373abd9ed0f576dff2ca8e Mon Sep 17 00:00:00 2001 From: Paul Menczel Date: Thu, 25 Sep 2025 16:37:43 +0900 Subject: [PATCH 3/6] Avoid calling np.random.seed in tests, update fitting tolerances for fixed rmse calculation --- qutip/tests/core/test_brtools.py | 11 +++++----- qutip/tests/core/test_environment.py | 32 ++++++++++++++-------------- qutip/tests/solver/test_results.py | 1 + qutip/tests/test_mkl.py | 4 ++-- qutip/tests/test_utilities.py | 23 ++++++++++++-------- 5 files changed, 38 insertions(+), 33 deletions(-) diff --git a/qutip/tests/core/test_brtools.py b/qutip/tests/core/test_brtools.py index 095e40644a..77822366b2 100644 --- a/qutip/tests/core/test_brtools.py +++ b/qutip/tests/core/test_brtools.py @@ -10,9 +10,8 @@ ) -def _make_rand_data(shape): - np.random.seed(11) - array = np.random.rand(*shape) + 1j*np.random.rand(*shape) +def _make_rand_data(shape, rng): + array = rng.random(shape) + 1j*rng.random(shape) return qutip.data.Dense(array) @@ -31,9 +30,9 @@ def _make_rand_data(shape): ids=['', 'transpose', 'conj', 'dag']) def test_matmul_var(datatype, transleft, transright): shape = (5, 5) - np.random.seed(11) - left = qutip.data.to(datatype, _make_rand_data(shape)) - right = qutip.data.to(datatype, _make_rand_data(shape)) + rng = np.random.default_rng(seed=11) + left = qutip.data.to(datatype, _make_rand_data(shape, rng)) + right = qutip.data.to(datatype, _make_rand_data(shape, rng)) expected = qutip.data.matmul( transform[transleft](left), diff --git a/qutip/tests/core/test_environment.py b/qutip/tests/core/test_environment.py index 70f8253fa7..22399a9a43 100644 --- a/qutip/tests/core/test_environment.py +++ b/qutip/tests/core/test_environment.py @@ -532,8 +532,8 @@ def test_fixed_cf_fit(self, reference, tMax, full_ansatz, tol): assert info["Nr"] == 2 assert info["Ni"] == 2 - assert info["rmse_real"] < 5e-3 - assert info["rmse_imag"] < 5e-3 + assert info["rmse_real"] < 5e-2 + assert info["rmse_imag"] < 5e-2 for key in ["fit_time_real", "fit_time_imag", "params_real", "params_imag", "summary"]: assert key in info @@ -551,7 +551,7 @@ def test_dynamic_cf_fit(self, reference, tMax, tol, full_ansatz): ) tlist = np.linspace(0, tMax, 100)[1:] # exclude t=0 fit, info = env.approximate( - "cf", tlist, target_rmse=0.01, Nr_max=3, Ni_max=3, + "cf", tlist, target_rmse=0.1, Nr_max=3, Ni_max=3, full_ansatz=full_ansatz ) @@ -564,8 +564,8 @@ def test_dynamic_cf_fit(self, reference, tMax, tol, full_ansatz): assert info["Nr"] == 1 assert info["Ni"] == 1 - assert info["rmse_real"] <= 0.01 - assert info["rmse_imag"] <= 0.01 + assert info["rmse_real"] <= 0.1 + assert info["rmse_imag"] <= 0.1 for key in ["fit_time_real", "fit_time_imag", "params_real", "params_imag", "summary"]: assert key in info @@ -594,7 +594,7 @@ def test_fixed_sd_fit(self, reference, wMax, tol): assert info["N"] == 4 assert info["Nk"] == 1 - assert info["rmse"] < 5e-3 + assert info["rmse"] < 5e-2 for key in ["fit_time", "params", "summary"]: assert key in info @@ -611,7 +611,7 @@ def test_dynamic_sd_fit(self, reference, wMax, tol, params): ) wlist = np.linspace(0, wMax, 100) fit, info = env.approximate( - "sd", wlist, Nk=1, target_rmse=0.01, Nmax=5, **params + "sd", wlist, Nk=1, target_rmse=0.1, Nmax=5, **params ) assert isinstance(fit, ExponentialBosonicEnvironment) @@ -623,14 +623,14 @@ def test_dynamic_sd_fit(self, reference, wMax, tol, params): assert info["N"] < 5 assert info["Nk"] == 1 - assert info["rmse"] < 0.01 + assert info["rmse"] < 0.1 for key in ["fit_time", "params", "summary"]: assert key in info @pytest.mark.parametrize(["reference", "tMax", "N", "tol"], [ pytest.param(OhmicReference(3, .75, 10, 1), - 15, 8, 1e-3, id="Ohmic Example"), - pytest.param(UDReference(1, .5, .1, 1), 2, 4, 1e-3, id='UD Example'), + 15, 8, 1e-2, id="Ohmic Example"), + pytest.param(UDReference(1, .5, .1, 1), 2, 4, 1e-2, id='UD Example'), ]) @pytest.mark.parametrize("separate", [True, False]) def test_fixed_prony_fit(self, reference, tMax, N, tol, separate): @@ -665,8 +665,8 @@ def test_fixed_prony_fit(self, reference, tMax, N, tol, separate): @pytest.mark.parametrize(["reference", "tMax", "N", "tol"], [ pytest.param(OhmicReference(3, .75, 10, 1), - 15, 8, 1e-3, id="Ohmic Example"), - pytest.param(UDReference(1, .5, .1, 1), 2, 4, 1e-3, id='UD Example'), + 15, 8, 1e-2, id="Ohmic Example"), + pytest.param(UDReference(1, .5, .1, 1), 2, 4, 1e-2, id='UD Example'), ]) @pytest.mark.parametrize("separate", [True, False]) def test_fixed_esprit_fit(self, reference, tMax, N, tol, separate): @@ -701,8 +701,8 @@ def test_fixed_esprit_fit(self, reference, tMax, N, tol, separate): @pytest.mark.parametrize(["reference", "tMax", "N", "tol"], [ pytest.param(OhmicReference(3, .75, 10, 1), - 15, 8, 1e-3, id="Ohmic Example"), - pytest.param(UDReference(1, .5, .1, 1), 2, 2, 1e-3, id='UD Example'), + 15, 8, 1e-2, id="Ohmic Example"), + pytest.param(UDReference(1, .5, .1, 1), 2, 2, 1e-2, id='UD Example'), ]) @pytest.mark.parametrize("separate", [True, False]) def test_fixed_espira1_fit(self, reference, tMax, N, tol, separate): @@ -737,8 +737,8 @@ def test_fixed_espira1_fit(self, reference, tMax, N, tol, separate): @pytest.mark.parametrize(["reference", "tMax", "N", "tol"], [ pytest.param(OhmicReference(3, .75, 10, 1), - 15, 8, 1e-3, id="Ohmic Example"), - pytest.param(UDReference(1, .5, .1, 1), 2, 2, 1e-3, id='UD Example'), + 15, 8, 1e-2, id="Ohmic Example"), + pytest.param(UDReference(1, .5, .1, 1), 2, 2, 1e-2, id='UD Example'), ]) @pytest.mark.parametrize("separate", [True, False]) def test_fixed_espira2_fit(self, reference, tMax, N, tol, separate): diff --git a/qutip/tests/solver/test_results.py b/qutip/tests/solver/test_results.py index fd3752dce1..2835857d6e 100644 --- a/qutip/tests/solver/test_results.py +++ b/qutip/tests/solver/test_results.py @@ -400,6 +400,7 @@ def test_repr(self, keep_runs_results): if keep_runs_results: assert "Trajectories saved." in repr + @pytest.mark.flaky(reruns=2) @pytest.mark.parametrize('keep_runs_results1', [True, False]) @pytest.mark.parametrize('keep_runs_results2', [True, False]) def test_merge_result(self, keep_runs_results1, keep_runs_results2): diff --git a/qutip/tests/test_mkl.py b/qutip/tests/test_mkl.py index 734f33e25e..647d3225c5 100644 --- a/qutip/tests/test_mkl.py +++ b/qutip/tests/test_mkl.py @@ -19,8 +19,8 @@ def test_single_rhs_vector_real(self): [1, 0, 1], [0, 0, 1]]) As = scipy.sparse.csr_matrix(Adense) - np.random.seed(1234) - x = np.random.randn(3) + rng = np.random.default_rng(seed=1234) + x = rng.standard_normal(3) b = As * x x2 = mkl_spsolve(As, b, verbose=True) np.testing.assert_allclose(x, x2) diff --git a/qutip/tests/test_utilities.py b/qutip/tests/test_utilities.py index dab481ef30..8fd8755dc9 100644 --- a/qutip/tests/test_utilities.py +++ b/qutip/tests/test_utilities.py @@ -124,6 +124,8 @@ def test_cpu_count(monkeypatch): class TestFitting: + rng = np.random.default_rng(seed=42) + def model(self, x, a, b, c): return np.real(a * np.exp(-(b + 1j * c) * x)) @@ -145,8 +147,7 @@ def generate_data(self, request): fparams2 = [3, 2, .5] y = self.model(x, *fparams1) + self.model(x, *fparams2) if noisy: - np.random.seed(42) - noise = np.random.normal(0, 0.01, len(x)) + noise = self.rng.normal(0, 0.01, len(x)) y += noise return x, y, fparams1, fparams2, noisy @@ -156,11 +157,15 @@ def test_fit(self, generate_data): self.model, num_params=3, xdata=x, ydata=y, lower=[-np.inf, -np.inf, 0], target_rmse=1e-8, Nmax=2 ) + fit_y = self.model(x, *params[0]) + self.model(x, *params[1]) + assert rmse == pytest.approx( + np.sqrt(np.mean((y - fit_y)**2)) / (np.max(y) - np.min(y)) + ) + if noisy: - assert rmse < 1e-3 - # The atol is the std of noise times the maximum of the signal - assert (np.all(np.isclose(params, [fparams1, fparams2], atol=1e-2*np.max(y))) or - np.all(np.isclose(params, [fparams2, fparams1], atol=1e-2*np.max(y)))) + assert rmse < 1e-2 + assert (np.all(np.isclose(params, [fparams1, fparams2], atol=.2)) or + np.all(np.isclose(params, [fparams2, fparams1], atol=.2))) else: assert rmse < 1e-8 assert (np.all(np.isclose(params, [fparams1, fparams2], atol=1e-3)) or @@ -181,7 +186,7 @@ def test_espira_I(self, generate_data): x, y, _, _, noisy = generate_data rmse, params = utils.espira1(y, 4, tol=1e-16) if noisy: - assert rmse < 1e-3 + assert rmse < 1e-2 np.testing.assert_allclose(self.eval_prony(len(x), params), y, atol=1e-2*np.max(y)) else: assert rmse < 1e-8 @@ -191,7 +196,7 @@ def test_espira_II(self, generate_data): x, y, _, _, noisy = generate_data rmse, params = utils.espira2(y, 4, tol=1e-16) if noisy: - assert rmse < 1e-3 + assert rmse < 1e-2 np.testing.assert_allclose(self.eval_prony(len(x), params), y, atol=1e-2*np.max(y)) else: assert rmse < 1e-8 @@ -202,7 +207,7 @@ def test_prony_methods(self, generate_data, method): x, y, _, _, noisy = generate_data rmse, params = utils.prony_methods(method, y, 4) if noisy: - assert rmse < 1e-3 + assert rmse < 1e-2 np.testing.assert_allclose(self.eval_prony(len(x), params), y, atol=2e-2*np.max(y)) else: assert rmse < 1e-8 From 4d1d03bfc7faba5b8334659f6ed92f3dd3fdb8f9 Mon Sep 17 00:00:00 2001 From: Paul Menczel Date: Fri, 26 Sep 2025 11:01:17 +0900 Subject: [PATCH 4/6] Add another flaky mark --- qutip/tests/solver/test_results.py | 1 + 1 file changed, 1 insertion(+) diff --git a/qutip/tests/solver/test_results.py b/qutip/tests/solver/test_results.py index 2835857d6e..7a33316919 100644 --- a/qutip/tests/solver/test_results.py +++ b/qutip/tests/solver/test_results.py @@ -283,6 +283,7 @@ def test_NmmcResult(self, include_no_jump, martingale, result_trace): for i, (s1, s2) in enumerate(zip(m_res.average_states, result_trace)): assert s1 == s2 * qutip.fock_dm(10, i) + @pytest.mark.flaky(reruns=2) @pytest.mark.parametrize('keep_runs_results', [True, False]) @pytest.mark.parametrize('include_no_jump', [True, False]) @pytest.mark.parametrize(["e_ops", "results"], [ From d201a03e12a652ce4cf67fa5384394ce4290f863 Mon Sep 17 00:00:00 2001 From: Paul Menczel Date: Fri, 26 Sep 2025 13:11:06 +0900 Subject: [PATCH 5/6] Update tolerances --- qutip/tests/core/test_environment.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/qutip/tests/core/test_environment.py b/qutip/tests/core/test_environment.py index 22399a9a43..2c60a3e45f 100644 --- a/qutip/tests/core/test_environment.py +++ b/qutip/tests/core/test_environment.py @@ -775,7 +775,7 @@ def test_fixed_espira2_fit(self, reference, tMax, N, tol, separate): @pytest.mark.parametrize(["reference", "wMax", "tol"], [ pytest.param(OhmicReference(3, .75, 10, 1), 15, .2, id="Ohmic Example"), - pytest.param(UDReference(1, .5, .1, 1), 2, 1e-4, id='UD Example'), + pytest.param(UDReference(1, .5, .1, 1), 2, 1e-3, id='UD Example'), ]) def test_fixed_aaa_fit(self, reference, wMax, tol): env = BosonicEnvironment.from_spectral_density( @@ -797,7 +797,7 @@ def test_fixed_aaa_fit(self, reference, wMax, tol): @pytest.mark.parametrize(["reference", "wMax", "tol"], [ pytest.param(OhmicReference(3, .75, 10, 1), 15, .2, id="DL Example"), - pytest.param(UDReference(1, .5, .1, 1), 2, 1e-4, id='UD Example'), + pytest.param(UDReference(1, .5, .1, 1), 2, 1e-3, id='UD Example'), ]) def test_fixed_ps_fit(self, reference, wMax, tol): env = BosonicEnvironment.from_spectral_density( From 2295f5a4657f46e4ca70dc42cc14323bbdbab412 Mon Sep 17 00:00:00 2001 From: Paul Menczel Date: Fri, 26 Sep 2025 15:08:37 +0900 Subject: [PATCH 6/6] Update power spectrum fitting test again since it failed on MacOS --- qutip/tests/core/test_environment.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/qutip/tests/core/test_environment.py b/qutip/tests/core/test_environment.py index 2c60a3e45f..34b3d16959 100644 --- a/qutip/tests/core/test_environment.py +++ b/qutip/tests/core/test_environment.py @@ -804,7 +804,7 @@ def test_fixed_ps_fit(self, reference, wMax, tol): reference.spectral_density, T=reference.T, tag="test" ) wlist = np.linspace(-wMax, wMax, 200) - fit, info = env.approximate("ps", wlist,Nmax=6) + fit, info = env.approximate("ps", wlist, Nmax=6, maxfev=1e5) assert isinstance(fit, ExponentialBosonicEnvironment) assert fit.T == env.T