From: Nuralanur@...
Date: 2006-08-07T19:56:40+09:00
Subject: Re: Setting precision in Ruby-GSL calculations ?

-------------------------------1154948191
Content-Type: text/plain; charset="US-ASCII"
Content-Transfer-Encoding: 7bit

Dear Yoshiki, dear Ed,
 
thank you for your replies, over which I have been meditating for a  weekend.
 
I have tried to implement a singular value decomposition with U and V  always
square - and maybe of different sizes. I was following a calculated  example 
from a
textbook to do this.
The idea is the following.
 
Assume we want to calculate the SVD of matrix a (not necessarily  square).
I was considering a=[[10 5 -10],[2 -11 10]] as an example, so no  isssues with
respect to 64 bit precision can matter here.
I'll describe  the procedure and draw a line how far I could follow
it (to Float's  precision).
It was suggested to calculate b=a*a.tranpose, to get the  eigenvalues
of b : lambda_1=360 and lambda_2=90. 
The singular values of  matrix a are then the square roots of these numbers, 
i.e.
sigma_1=6*sqrt(10)  and sigma_2=3*sqrt(10).
One may then calculate normed eigenvectors to  lambda_1: [-1/sqrt(2)  
1/sqrt(2) ]
and to lambda_2 := [ 1/sqrt(2)  1/sqrt(2)], these two vectors form the matrix 
V.
 
Now, to calculate U, one uses the fact that 
 
 
a.tranpose*v.column(0)=sigma_1*u_column(0)
and  a.tranpose*v.column(1)=sigma_2*u_column(1).


This gives two vectors, but as U is 3x3, I need a third vector, which is  
constructed using
Gram-Schmidt orthogonalization.
In the textbook, this  third vector is determined to be  [ 1/3 2/3 
2/3].transpose , giving
U as  [ [-2/sqrt(45) 2/sqrt(5) 1/3],[-4/sqrt(45) -1/sqrt(5) 2/3],[5/sqrt(45) 
0  2/3]].
 
I tried to construct U from a random vector, but I got something like  [0.304 
0.657 0.657] instead
and determinants for U between 0.96 and 1.03 ,  where they should be 1 or 
quite close to 1 
all the time.
Up to  calculating U, all the calculations I got were in Float-perfect 
harmony with  
the textbook example.
So, I've checked my Gram-Schmidt implementation,  but I can't find any errors 
in there, and
I can reproduce some worked-out  examples from textbooks ..
So, my guess was that  maybe the sigmas I  calculated from GSL eigenvalue 
calculations
might be to blame, but you tell  me that that's not the case ....
 
Best regards,
 
Axel 


-------------------------------1154948191--