From d39f4c2d9be3e904b3bd5861e8191e44b6671667 Mon Sep 17 00:00:00 2001 From: Stefano Zampini Date: Wed, 8 Oct 2025 08:47:34 +0300 Subject: [PATCH] LinearEigensolver: default to eps_target = 0 is conceptually wrong the presence of the option, together with the sinvert for the bcs case can lead to confusion Better to be consistent with whatever SLEPc believe is the proper default --- firedrake/eigensolver.py | 9 ++++---- .../firedrake/regression/test_eigensolver.py | 21 ++++++++++++------- 2 files changed, 18 insertions(+), 12 deletions(-) diff --git a/firedrake/eigensolver.py b/firedrake/eigensolver.py index bb8a11a5ab..568536e9a6 100644 --- a/firedrake/eigensolver.py +++ b/firedrake/eigensolver.py @@ -142,9 +142,10 @@ class LinearEigensolver(OptionsManager): "eps_largest_real": None """ - DEFAULT_EPS_PARAMETERS = {"eps_type": "krylovschur", - "eps_tol": 1e-10, - "eps_target": 0.0} + DEFAULT_EPS_PARAMETERS = { + "eps_type": "krylovschur", + "eps_tol": 1e-10, + } def __init__(self, problem, n_evals, *, options_prefix=None, solver_parameters=None, ncv=None, mpd=None): @@ -158,8 +159,6 @@ def __init__(self, problem, n_evals, *, options_prefix=None, for key in self.DEFAULT_EPS_PARAMETERS: value = self.DEFAULT_EPS_PARAMETERS[key] solver_parameters.setdefault(key, value) - if self._problem.bcs: - solver_parameters.setdefault("st_type", "sinvert") super().__init__(solver_parameters, options_prefix) self.set_from_options(self.es) diff --git a/tests/firedrake/regression/test_eigensolver.py b/tests/firedrake/regression/test_eigensolver.py index f4065083bd..1f34a3e72d 100644 --- a/tests/firedrake/regression/test_eigensolver.py +++ b/tests/firedrake/regression/test_eigensolver.py @@ -20,10 +20,13 @@ def evals(n, degree=1, mesh=None, restrict=False): bc = DirichletBC(V, 0.0, "on_boundary") eigenprob = LinearEigenproblem(a, bcs=bc, bc_shift=-6666., restrict=restrict) - # Create corresponding eigensolver, looking for n eigenvalues - eigensolver = LinearEigensolver( - eigenprob, n, solver_parameters={"eps_largest_real": None} - ) + # Create corresponding eigensolver, looking for n eigenvalues close to 0 + # We use shift-and-invert as spectral transform (SLEPc's default is shift) + solver_parameters = { + "eps_target": 0, + "st_type": "sinvert", + } + eigensolver = LinearEigensolver(eigenprob, n, solver_parameters=solver_parameters) ncov = eigensolver.solve() # boffi solns @@ -67,8 +70,12 @@ def poisson_eigenvalue_2d(i): ep = LinearEigenproblem(inner(grad(u), grad(v)) * dx, bcs=bc, bc_shift=666.0) - es = LinearEigensolver(ep, 1, solver_parameters={"eps_gen_hermitian": None, - "eps_largest_real": None}) + solver_parameters = { + "eps_gen_hermitian": None, + "eps_target": 0, + "st_type": "sinvert", + } + es = LinearEigensolver(ep, 1, solver_parameters=solver_parameters) es.solve() return es.eigenvalue(0)-2.0 @@ -77,7 +84,7 @@ def poisson_eigenvalue_2d(i): @pytest.mark.skipslepc def test_evals_2d(): """2D Eigenvalue convergence test. As with Boffi, we observe that the - convergence rate convergest to 2 from above.""" + convergence rate converges to 2 from above.""" errors = np.array([poisson_eigenvalue_2d(i) for i in range(5)]) convergence = np.log(errors[:-1]/errors[1:])/np.log(2.0)