# How to solve a linear system with a close to zero eigenvalue?

The following is a small part of my research:

I got a 4 by 4 matrix $$\mathbf{M}$$ of rank 3:

M={{1, 1, 1, 1}, {2, -2, 2.0006219879401477, -2.0006219879401477}, {4.193503236170054, 0.07680669095242888, 4.195764689057942, 0.0767525731461281}, {19.781537502364305, -0.36231149755725733, 19.793967330980973, -0.36208844823571257}}


I need to solve the system $$\mathbf{MC=0}$$ where $$\mathbf{C}=[C_1,C_2,C_3,C_4]^\top$$. I know that $$\mathbf{M}$$ has rank equal to 3, but due to computation error, Mathematica gives its rank equal to 4. (One can check using the command "Eigenvalues[M]" to find that one of its eigenvalues is very close to zero.)

I set $$C_1=1$$, and want to solve for $$C_2,C_3,C_4$$. I'm wondering if there is a way I can get a good approximation of $$C_2,C_3,C_4$$. I thought of using Gaussian Elimination, but I have no idea how to proceed, as I cannot change the fact that the rank is 4 in Mathemtaica.

Any help or hint will be much appreciated!

• The Tolerance option ofNullSpace is good for this situation. In[2220]:= mat = {{1, 1, 1, 1}, {2, -2, 2.0006219879401477, -2.0006219879401477}, {4.193503236170054, 0.07680669095242888, 4.195764689057942, 0.0767525731461281}, {19.781537502364305, -0.36231149755725733, 19.793967330980973, -0.36208844823571257}}; NullSpace[mat, Tolerance -> 10^(-6)] Out[2221]= {{0.23992, -0.665256, -0.23978, 0.665115}} Dec 20, 2023 at 22:02
• @DanielLichtblau NullSpace gives exactly the same result as a SingularValueDecomposition  (see update in my solution) and seems much simpler! Dec 21, 2023 at 7:32
• @Roman As you may have guessed, indeed NullSpace for approximate values does use the SVD. Dec 21, 2023 at 14:57

M = {{1, 1, 1, 1},
{2, -2, 2.0006219879401477, -2.0006219879401477},
{4.193503236170054, 0.07680669095242888, 4.195764689057942, 0.0767525731461281},
{19.781537502364305, -0.36231149755725733, 19.793967330980973, -0.36208844823571257}};


Look for the eigenvector associated with the smallest eigenvalue (by absolute value):

c = MinimalBy[Transpose[Eigensystem[M]], Abs@*First][[1, 2]]
(*    {0.240408, -0.66508, -0.240268, 0.664939}    *)


The result is pretty close to $$M\cdot c=0$$:

M . c
(*    {-3.98559*10^-7, 1.1026*10^-6, 3.98325*10^-7, -1.10236*10^-6}    *)


Normalizing by the first element $$c_1=1$$:

cn = Normalize[c, First]
(*    {1., -2.76646, -0.999414, 2.76587}    *)


If your matrix is large, you could speed things up by trying to compute only the smallest eigenvalue/eigenvector pair:

Eigensystem[M, -1]
(*    {{-1.65784*10^-6}, {{0.240408, -0.66508, -0.240268, 0.664939}}}    *)


and therefore

c = Eigenvectors[M, -1][[1]]
(*    {0.240408, -0.66508, -0.240268, 0.664939}    *)


# update: singular-value decomposition

On second thought, an eigenvalue decomposition may not be the best way forward and a more general singular-value decomposition could give a slightly more accurate result.

The singular value decomposition is

{u, σ, v} = SingularValueDecomposition[M];


Look for the smallest singular value by absolute magnitude:

s = OrderingBy[Diagonal[σ], Abs, 1][[1]]
(*    4    *)


The minimizing vector is the corresponding vector in $$v$$:

d = Transpose[v][[s]]
(*    {0.23992, -0.665256, -0.23978, 0.665115}    *)


See that it fits better:

{Norm[M . c], Norm[M . d]}
(*    {1.65784*10^-6, 5.37312*10^-7}    *)


Faster by computing only the smallest singular value by Arnoldi iteration:

d = SingularValueDecomposition[M, -1,
Method -> {"Arnoldi",
"Criteria" -> "Magnitude",
"MaxIterations" -> 10^6}][[3, All, 1]]
(*    {-0.23992, 0.665256, 0.23978, -0.665115}    *)


Even faster and fully equivalent: compute the smallest eigenvalue and its associated eigenvector of $$M^T\cdot M$$: this is the vector that minimizes $$\|M\cdot\vec{d}\|^2=\vec{d}\cdot (M^T\cdot M)\cdot \vec{d}$$ under the constraint that $$\|\vec{d}\|^2=\vec{d}\cdot\vec{d}=1$$:

d = Eigenvectors[Transpose[M] . M, -1,
Method -> {"Arnoldi",
"Criteria" -> "Magnitude",
"MaxIterations" -> 10^6}][[1]]
(*    {-0.23992, 0.665256, 0.23978, -0.665115}    *)