[Rd] eigen()

[Robin Hankin]

On 10 Jan 2006, at 14:14, Peter Dalgaard wrote:
>


>>> Strange and semi-random results on SuSE 9.3 as well:
>>>
>>>
>>>> eigen(matrix(1:100,10,10))$values
>>>  [1]  5.4e-311+ 0.0e+00i -2.5e-311+3.7e-311i -2.5e-311-3.7e-311i
>>>  [4]  2.5e-312+ 0.0e+00i -2.4e-312+ 0.0e+00i  3.2e-317+ 0.0e+00i
>>>  [7]   0.0e+00+ 0.0e+00i   0.0e+00+ 0.0e+00i   0.0e+00+ 0.0e+00i
>>> [10]   0.0e+00+ 0.0e+00i
>>>
>>
>> Mine is closer to Robin's, but not the same (EL4 x86).
>>
>>> eigen(matrix(1:100,10,10))$values
>>   [1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
>>   [3]  6.292457e-16+2.785369e-15i  6.292457e-16-2.785369e-15i
>>   [5] -1.055022e-15+0.000000e+00i  3.629676e-16+0.000000e+00i
>>   [7]  1.356222e-16+2.682405e-16i  1.356222e-16-2.682405e-16i
>>   [9]  1.029077e-16+0.000000e+00i -1.269181e-17+0.000000e+00i
>>>
>>
>> But surely, my matrix algebra is a bit rusty, I think this matrix is
>> solveable analytically? Most of the eigenvalues shown are almost
>> exactly zero, except the first two, actually, which is about 521
>> and -16 to the closest integer.
>>
>> I think the difference between mine and Robin's are rounding errors
>> (the matrix is simple enough I expect the solution to be simple  
>> integers
>> or easily expressible analystical expressions, so 8 e-values being  
>> zero
>> is fine). Peter's number seems to be all 10 e-values are zero or one
>> being a huge number! So Peter's is odd... and Peter's machine also
>> seems
>> to be of a different archtecture (64-bit machine)?
>>
>> HTL
>
> Notice that Robin got something completely different in _R-devel_
> which is where I did my check too.  In R 2.2.1 I get the expected two
> non-zero eigenvalues.
>
> I'm not sure whether (and how) you can work out the eigenvalues
> analytically, but since all columns are linear progressions, it is
> at least obvious that the matrix must have column rank two.
>




For everyone's entertainment, here's an example where the analytic  
solution
is known.

fact 1:  the first eigenvalue of a magic square is equal to its constant
fact 2: the sum of the other eigenvalues of a magic square is zero
fact 3: the constant of a magic square of order 10 is 505.

R-2.2.0:

 > library(magic)
 > round(Re(eigen(magic(10),F,T)$values))
[1]  505  170 -170 -105  105   -3    3    0    0    0
 >

answers as expected.


R-devel:



 > a <- structure(c(68, 66, 92, 90, 16, 14, 37, 38, 41, 43, 65, 67, 89,
91, 13, 15, 40, 39, 44, 42, 96, 94, 20, 18, 24, 22, 45, 46, 69,
71, 93, 95, 17, 19, 21, 23, 48, 47, 72, 70, 4, 2, 28, 26, 49,
50, 76, 74, 97, 99, 1, 3, 25, 27, 52, 51, 73, 75, 100, 98, 32,
30, 56, 54, 80, 78, 81, 82, 5, 7, 29, 31, 53, 55, 77, 79, 84,
83, 8, 6, 60, 58, 64, 62, 88, 86, 9, 10, 33, 35, 57, 59, 61,
63, 85, 87, 12, 11, 36, 34), .Dim = c(10, 10))

[no magic package!  it fails R CMD check !]

 > round(Re(eigen(magic(10),F,T)$values))
[1] 7.544456e+165  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e 
+00
[6]  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e 
+00
 >


not as expected.



--
Robin Hankin
Uncertainty Analyst
National Oceanography Centre, Southampton
European Way, Southampton SO14 3ZH, UK
  tel  023-8059-7743