SVD issue calculating null basis

I am tyring to reproduce the image results for V. I have written the following code. However, the V produced by svd command is not matching what I should see. I have double checked my G matrix and the S output is correct. Any advice appreciated.

%create the matrix, run svd
G=[1 0 0 1 0 0 1 0 0; 0 1 0 0 1 0 0 1 0; 0 0 1 0 0 1 0 0 1; 1 1 1 0 0 0 0 0 0; 
    0 0 0 1 1 1 0 0 0; 0 0 0 0 0 0 1 1 1; sqrt(2) 0 0 0 sqrt(2) 0 0 0 sqrt(2); 0 0 0 0 0 0 0 0 sqrt(2)];
[U,S,V] = svd(G,0);
s_d=diag(S);

%The model null space, N(G), is spanned by the two orthonormal vectors
%that form the 8th and 9th columns of V. An orthonormal basis for the null space is
s = diag(S);
column_basis = U(:,logical(s));
rank_G= nnz(s);
null_basis=V(:,~s);

[Walter Roberson]

rank(G) shows 7.

You are using a proxy for the rank when you look at the number of nonzeros of the diagonal. As they discuss, there is numeric noise in one of the elements.

[L]

Hi Walter,

Thanks, the rank makes sense; I'm just not sure why my V matrix is not correct.

[L]

Is it because the null space basis is nonunique, so my results are not incorrect, just different? Or should I be seeing the same thing?

Answer 1

[John D'Errico]

The problem is that while the last element of s is small, it is not technically zero.

s
s =
          3.17982071419273
                         2
          1.73205080756888
          1.73205080756888
          1.73205080756888
          1.60697051494866
         0.553521444640095
      8.80988299920657e-17
>> nnz(s)
ans =
     8

So what does nnz tell you? Hey it is not zero after all. 8.8e-17 is tiny, but is it zero?

s == 0
ans =
  8×1 logical array
   0
   0
   0
   0
   0
   0
   0
   0

In fact, that least element is close enough to zero for our purposes that here, we consider it to be so. But it is not in fact zero, and nnz is rather strict. As far as nnz cares, it simply is not zero.

You might do something like this:

sum(s >= eps*max(s))
ans =
     7

Here we see that what matters is if the smallest number really is tiny, when compared to the maximum element in the vector s. Note that I am not worried about whether elements of s are negative, because they are singular values! Otherwise, I might have needed a call to abs in there.

[L]

Thank you, that makes sense. In terms of the V vector, is it because the null space basis is nonunique, so my results are not incorrect, just different? Or should I be seeing the same thing and the treatment of the last element of s as nonzero is causing my question?

[John D'Errico]

size(G)
ans =
     8     9
>> null(G)
ans =
         0.335734498771089         0.232269268310709
         0.232269268310709        -0.335734498771089
        -0.568003767081798          0.10346523046038
        -0.232269268310709         0.335734498771089
        -0.335734498771089        -0.232269268310709
         0.568003767081798         -0.10346523046038
         -0.10346523046038        -0.568003767081798
          0.10346523046038         0.568003767081798
     -1.80411241501588e-16      2.77555756156289e-17

So, G is 8x9, with numerical rank 7. So the null space of G is 2-dimensional. Thus the set of vectors in 9 dimensions such that when you right multiply them times G gives you zero, those vectors live in a 2-dimensional subspace.

It is indeed true that we could rotate those two vectors arbitrarily, choosing different sets of two vectors, but still spanning that 2-dimensional suspace. (See my comments below.) HOWEVER, what you did is not valid:

rank_G= nnz(s);
null_basis=V(:,~s);

As I said, nnz does not count 7 non-zeros. It counts 8 non-zeros. Similarly, when you used ~s in there, that was again wrong.

~s
ans =
  8×1 logical array
   0
   0
   0
   0
   0
   0
   0
   0

We can use null(G) to compute that orthogonal basis. Or the last two columns of V, which are in fact the same as returned by null(G).

V(:,sum(s >= eps*max(s))+1:end)
ans =
         0.335734498771089         0.232269268310709
         0.232269268310709        -0.335734498771089
        -0.568003767081798          0.10346523046038
        -0.232269268310709         0.335734498771089
        -0.335734498771089        -0.232269268310709
         0.568003767081798         -0.10346523046038
         -0.10346523046038        -0.568003767081798
          0.10346523046038         0.568003767081798
     -1.80411241501588e-16      2.77555756156289e-17

I cannot copy and paste the array V0 in, since you give us only a picture. But let me see how my typing does.

V0 = [-.0620 -.4035;-.4035 .0620;.4655 .3415; .4035 -.0620; .0620 .4035;-.4655 -.3415;-.3415 .4655;.3415 -.4655; 0 0];
V0
V0 =
   -0.0620   -0.4035
   -0.4035    0.0620
    0.4655    0.3415
    0.4035   -0.0620
    0.0620    0.4035
   -0.4655   -0.3415
   -0.3415    0.4655
    0.3415   -0.4655
         0         0
>> G*V0
ans =
     0     0
     0     0
     0     0
     0     0
     0     0
     0     0
     0     0
     0     0

However, if you look at null(G), it was clearly different, but null(G) still kills off G, as it must.

G*null(G)
ans =
   5.5511e-17            0
   1.1102e-16  -3.3307e-16
  -4.0246e-16    1.249e-16
            0   3.0531e-16
  -2.2204e-16  -2.3592e-16
  -9.7145e-17  -1.9429e-16
  -3.3095e-17   2.0579e-16
  -2.5514e-16   3.9252e-17

The latter result is zero as far as this discussion is concerned. Therefore, I could have chosen a different basis for the nullspace of G, such as V0. I don't know how they built V0 in that text fragment, but V0 is somewhat 2-dimensionally arbitrary.

[L]

This was a really good explanation, thank you