Repository: pymor/pymor

Issue:
Better default tolerances for method of snapshots?
See discussion in #1012.

@lbalicki, as I understand https://www.netlib.org/lapack/lug/node89.html, then every EV smaller than the largest eigenvalue times machine precision should be considered as potentially zero. Everything else should be guaranteed to be a non-zero EV. So maybe it would make sense to throw everything away below that threshold (or a fixed constant times this threshold) and trust that the corresponding eigenvectors are at least not linearly dependent? 

Discussion:

Patch:
diff --git a/src/pymor/algorithms/svd_va.py b/src/pymor/algorithms/svd_va.py
index d8e364358..d636dada0 100644
--- a/src/pymor/algorithms/svd_va.py
+++ b/src/pymor/algorithms/svd_va.py
@@ -96,6 +96,10 @@ def method_of_snapshots(A, product=None, modes=None, rtol=1e-7, atol=0., l2_err=
         if modes is not None:
             selected_modes = min(selected_modes, modes)
 
+        if selected_modes > A.dim:
+            logger.warning('Number of computed singular vectors larger than array dimension! Truncating ...')
+            selected_modes = A.dim
+
         s = np.sqrt(evals[:selected_modes])
         V = V[:selected_modes]
         Vh = V.conj()


Test patch:
diff --git a/src/pymortests/algorithms/svd_va.py b/src/pymortests/algorithms/svd_va.py
index 43b0d326c..175258c81 100644
--- a/src/pymortests/algorithms/svd_va.py
+++ b/src/pymortests/algorithms/svd_va.py
@@ -11,6 +11,7 @@ from pymor.algorithms.basic import almost_equal
 from pymor.algorithms.svd_va import method_of_snapshots, qr_svd
 from pymor.algorithms.basic import contains_zero_vector
 from pymor.core.logger import log_levels
+from pymor.vectorarrays.numpy import NumpyVectorSpace
 from pymortests.base import runmodule
 from pymortests.fixtures.operator import operator_with_arrays_and_products
 from pymortests.strategies import given_vector_arrays
@@ -57,5 +58,12 @@ def test_method_of_snapshots_with_product(operator_with_arrays_and_products, met
     assert np.all(almost_equal(A, UsVh, rtol=4e-8))
 
 
+@pytest.mark.parametrize('method', methods)
+def test_not_too_many_modes(method):
+    vec_array = NumpyVectorSpace.from_numpy(np.logspace(-5, 0, 10).reshape((-1, 1)))
+    U, s, V = method(vec_array, atol=0, rtol=0)
+    assert len(U) == len(s) == len(V) == 1
+
+
 if __name__ == "__main__":
     runmodule(filename=__file__)