# KrylovKit.svdsolve() terminating too early, with non-orthogonal singular vectors

**URL:** <https://discourse.julialang.org/t/krylovkit-svdsolve-terminating-too-early-with-non-orthogonal-singular-vectors/73109>\
**Category:** Performance\
**Created:** [December 15, 2021, 2:36am UTC](https://discourse.julialang.org/t/krylovkit-svdsolve-terminating-too-early-with-non-orthogonal-singular-vectors/73109 "2021-12-15T02:36:28Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![iamsuddhasattwa](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/iamsuddhasattwa/32/7441_2.png) [@iamsuddhasattwa](https://discourse.julialang.org/u/iamsuddhasattwa)\
**Post date:** [December 15, 2021, 2:36am UTC](https://discourse.julialang.org/t/krylovkit-svdsolve-terminating-too-early-with-non-orthogonal-singular-vectors/73109/1 "2021-12-15T02:36:28Z")

</div>

I have for long had troubles with the `KrylovKit.svdsolve()` method. I am using Julia version **1.7.0** , and I am trying to find the SVD of a sparse 1001X1001 matrix. It seems that if I try to find L=150 svd pairs, the U, V matrices are not even orthogonal. Here is the simple code

```julia
using DelimitedFiles;
using Serialization;
using KrylovKit;
using RecursiveArrayTools;
using LinearAlgebra;
K = deserialize( "K.txt" );
L=150;
λ, U, V, info = KrylovKit.svdsolve( K, L, krylovdim=2L+1, tol=1E-15, maxiter=10000, verbosity=2 );
L=size(λ,1);
U = convert(Array,VectorOfArray(U[1:L]));
M = U'*U; M-=I; print("\n", maximum(abs, M));

```

As the output below shows, the left singular vectors are not even orthogonal, the very basic requirement of SVD. The M matrix is meant to test this orthogonality.

```julia
[ Info: GKL svdsolve in iter 1, krylovdim 301: 18 values converged, normres = (0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 4.05e-66, 0.00e+00, 0.00e+00, 6.30e-59, 0.00e+00, 4.03e-62, 6.30e-114, 2.03e-65, 0.00e+00, 0.00e+00, 5.28e-108, 0.00e+00, 0.00e+00, 3.90e-30, 2.07e-02, 3.09e-13, 0.00e+00, 0.00e+00, 0.00e+00, 9.65e-14, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 9.36e-29, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 3.02e-60, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 7.97e-108, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 3.17e-61, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 5.72e-108, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00)
┌ Info: GKL svdsolve finished after 2 iterations:
│ * 263 singular values converged
│ * norm of residuals = (0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.258099826319623e-16, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.2616591976142803e-16, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)
└ * number of operations = 830

0.9442596238752043

```

Also notice that the procedure had just one iteration. In all trials of this experiment, it always ends after one iteration.Also, the output says

```julia
18 values converged ... 

```

So obviously, some unknown factor is causing the program to terminate before at least 150 values have converged. Could you tell me what is going wrong ?

By the way, this same piece of code on the same matrix, runs properly on a larger computer with the same version of Julia. There is no memory overflow in my current machine, so I do not understand why this could be a problem.
