eigen()

classic Classic list List threaded Threaded
9 messages Options
Reply | Threaded
Open this post in threaded view
|

eigen()

Robin Hankin
Hi

I am having difficulty with eigen() on   R-devel_2006-01-05.tar.gz

Specifically,  in R-2.2.0 I get expected behaviour:


 > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
[1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
[3] -4.805412e-15+0.000000e+00i  1.347691e-15+4.487511e-15i
[5]  1.347691e-15-4.487511e-15i -4.269863e-16+0.000000e+00i
[7]  1.364748e-16+0.000000e+00i -1.269735e-16+0.000000e+00i
[9] -1.878758e-18+5.031259e-17i -1.878758e-18-5.031259e-17i
 >


The same command gives different results in the development version:


 > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
[1]  3.903094e-118 -3.903094e-118 -2.610848e-312 -2.995687e-313  
-2.748516e-313
[6] -1.073138e-314 -1.061000e-314 -1.060998e-314  4.940656e-324    
0.000000e+00
 > R.version()
Error: attempt to apply non-function
 > R.version
                _
platform       powerpc-apple-darwin8.3.0
arch           powerpc
os             darwin8.3.0
system         powerpc, darwin8.3.0
status         Under development (unstable)
major          2
minor          3.0
year           2006
month          01
day            04
svn rev        36984
language       R
version.string Version 2.3.0 Under development (unstable) (2006-01-04  
r36984)
 >


Note the strange magnitude of the output.

[
I need this to work because one of my packages fails under R-devel
]


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

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Peter Dalgaard
Robin Hankin <[hidden email]> writes:

> Hi
>
> I am having difficulty with eigen() on   R-devel_2006-01-05.tar.gz
>
> Specifically,  in R-2.2.0 I get expected behaviour:
>
>
>  > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
> [1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
> [3] -4.805412e-15+0.000000e+00i  1.347691e-15+4.487511e-15i
> [5]  1.347691e-15-4.487511e-15i -4.269863e-16+0.000000e+00i
> [7]  1.364748e-16+0.000000e+00i -1.269735e-16+0.000000e+00i
> [9] -1.878758e-18+5.031259e-17i -1.878758e-18-5.031259e-17i
>  >
>
>
> The same command gives different results in the development version:
>
>
>  > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
> [1]  3.903094e-118 -3.903094e-118 -2.610848e-312 -2.995687e-313  
> -2.748516e-313
> [6] -1.073138e-314 -1.061000e-314 -1.060998e-314  4.940656e-324    
> 0.000000e+00
>  > R.version()
> Error: attempt to apply non-function
>  > R.version

Strange and semi-random results on SuSE 9.3 as well:

>  eigen(matrix(1:100,10,10))$values
 [1] -5.393552e+194   3.512001e-68   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
>  eigen(matrix(1:100,10,10))$values
 [1]  1.526259e-311 -1.041529e-311  1.181720e-313   0.000000e+00   0.000000e+00
 [6]   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00
>  eigen(matrix(1:100,10,10))$values
 [1] -9.338774e+93  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
>  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



--
   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark          Ph:  (+45) 35327918
~~~~~~~~~~ - ([hidden email])                  FAX: (+45) 35327907

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Peter Dalgaard
Peter Dalgaard <[hidden email]> writes:

> Robin Hankin <[hidden email]> writes:
>
> > Hi
> >
> > I am having difficulty with eigen() on   R-devel_2006-01-05.tar.gz
> >
> > Specifically,  in R-2.2.0 I get expected behaviour:
> >
> >
> >  > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
> > [1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
> > [3] -4.805412e-15+0.000000e+00i  1.347691e-15+4.487511e-15i
> > [5]  1.347691e-15-4.487511e-15i -4.269863e-16+0.000000e+00i
> > [7]  1.364748e-16+0.000000e+00i -1.269735e-16+0.000000e+00i
> > [9] -1.878758e-18+5.031259e-17i -1.878758e-18-5.031259e-17i
> >  >
> >
> >
> > The same command gives different results in the development version:
> >
> >
> >  > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
> > [1]  3.903094e-118 -3.903094e-118 -2.610848e-312 -2.995687e-313  
> > -2.748516e-313
> > [6] -1.073138e-314 -1.061000e-314 -1.060998e-314  4.940656e-324    
> > 0.000000e+00
> >  > R.version()
> > Error: attempt to apply non-function
> >  > R.version
>
> Strange and semi-random results on SuSE 9.3 as well:
>
> >  eigen(matrix(1:100,10,10))$values
>  [1] -5.393552e+194   3.512001e-68   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
> >  eigen(matrix(1:100,10,10))$values
>  [1]  1.526259e-311 -1.041529e-311  1.181720e-313   0.000000e+00   0.000000e+00
>  [6]   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00
> >  eigen(matrix(1:100,10,10))$values
>  [1] -9.338774e+93  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
> >  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


On closer inspection:

> eigen(matrix(as.numeric(1:100),10),FALSE,TRUE)
$values
 [1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
 [3]  9.895951e-15+3.370683e-15i  9.895951e-15-3.370683e-15i
 [5] -7.018984e-15+2.881924e-15i -7.018984e-15-2.881924e-15i
 [7] -7.978136e-16+2.629350e-15i -7.978136e-16-2.629350e-15i
 [9]  1.818143e-16+6.007106e-16i  1.818143e-16-6.007106e-16i

$vectors
NULL


I.e. there's a bug in that eigen() doesn't check the storage mode, but
there's also a simple workaround.

--
   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark          Ph:  (+45) 35327918
~~~~~~~~~~ - ([hidden email])                  FAX: (+45) 35327907

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Hin-Tak Leung-2
In reply to this post by Peter Dalgaard
Peter Dalgaard wrote:

> Robin Hankin <[hidden email]> writes:
>
>
>>Hi
>>
>>I am having difficulty with eigen() on   R-devel_2006-01-05.tar.gz
>>
>>Specifically,  in R-2.2.0 I get expected behaviour:
>>
>>
>> > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
>>[1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
>>[3] -4.805412e-15+0.000000e+00i  1.347691e-15+4.487511e-15i
>>[5]  1.347691e-15-4.487511e-15i -4.269863e-16+0.000000e+00i
>>[7]  1.364748e-16+0.000000e+00i -1.269735e-16+0.000000e+00i
>>[9] -1.878758e-18+5.031259e-17i -1.878758e-18-5.031259e-17i
>> >
>>
>>
>>The same command gives different results in the development version:
>>
>>
>> > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
>>[1]  3.903094e-118 -3.903094e-118 -2.610848e-312 -2.995687e-313  
>>-2.748516e-313
>>[6] -1.073138e-314 -1.061000e-314 -1.060998e-314  4.940656e-324    
>>0.000000e+00
>> > R.version()
>>Error: attempt to apply non-function
>> > R.version
>
>
> Strange and semi-random results on SuSE 9.3 as well:
>
>
>> eigen(matrix(1:100,10,10))$values
>
>  [1] -5.393552e+194   3.512001e-68   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
>
>> eigen(matrix(1:100,10,10))$values
>
>  [1]  1.526259e-311 -1.041529e-311  1.181720e-313   0.000000e+00   0.000000e+00
>  [6]   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00
>
>> eigen(matrix(1:100,10,10))$values
>
>  [1] -9.338774e+93  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
>
>> 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

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Peter Dalgaard
Hin-Tak Leung <[hidden email]> writes:

> Peter Dalgaard wrote:
> > Robin Hankin <[hidden email]> writes:
> >
> >>Hi
> >>
> >>I am having difficulty with eigen() on   R-devel_2006-01-05.tar.gz
> >>
> >>Specifically,  in R-2.2.0 I get expected behaviour:
> >>
> >>
> >> > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
> >>[1]  5.208398e+02+0.000000e+00i -1.583980e+01+0.000000e+00i
> >>[3] -4.805412e-15+0.000000e+00i  1.347691e-15+4.487511e-15i
> >>[5]  1.347691e-15-4.487511e-15i -4.269863e-16+0.000000e+00i
> >>[7]  1.364748e-16+0.000000e+00i -1.269735e-16+0.000000e+00i
> >>[9] -1.878758e-18+5.031259e-17i -1.878758e-18-5.031259e-17i
> >> >
> >>
> >>
> >>The same command gives different results in the development version:
> >>
> >>
> >> > eigen(matrix(1:100,10,10),FALSE,TRUE)$values
> >> [1]  3.903094e-118 -3.903094e-118 -2.610848e-312 -2.995687e-313
> >> -2.748516e-313
> >> [6] -1.073138e-314 -1.061000e-314 -1.060998e-314  4.940656e-324
> >> 0.000000e+00
> >> > R.version()
> >>Error: attempt to apply non-function
> >> > R.version
> > Strange and semi-random results on SuSE 9.3 as well:
> >
> >> eigen(matrix(1:100,10,10))$values
> >  [1] -5.393552e+194   3.512001e-68   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
> >
> >> eigen(matrix(1:100,10,10))$values
> >  [1]  1.526259e-311 -1.041529e-311  1.181720e-313   0.000000e+00
> > 0.000000e+00
> >  [6]   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00
> >
> >> eigen(matrix(1:100,10,10))$values
> >  [1] -9.338774e+93  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
> >
> >> 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.

--
   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark          Ph:  (+45) 35327918
~~~~~~~~~~ - ([hidden email])                  FAX: (+45) 35327907

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: 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

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Prof Brian Ripley
I haven't seen most of this thread, but this is a classic case of passing
integers instead of doubles.  And indeed

     else if(is.numeric(x)) {
  storage.mode(x) <- "double"

has been removed from eigen.R in R-devel in r36952.  So that's the
culprit.

[BTW, x86 is not 64-bit, which is x86_64.  I'd take x86 to be ix86, what
Intel calls ia32 (and AMD does call x86, given what 'i' stands for here).]

On Tue, 10 Jan 2006, Robin Hankin wrote:

>
>
> 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
>
> ______________________________________________
> [hidden email] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
>

--
Brian D. Ripley,                  [hidden email]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Peter Dalgaard
Prof Brian Ripley <[hidden email]> writes:

> I haven't seen most of this thread, but this is a classic case of
> passing integers instead of doubles.  And indeed
>
>      else if(is.numeric(x)) {
>   storage.mode(x) <- "double"
>
> has been removed from eigen.R in R-devel in r36952.  So that's the
> culprit.

Thought as much...
 
> [BTW, x86 is not 64-bit, which is x86_64.  I'd take x86 to be ix86,
> what Intel calls ia32 (and AMD does call x86, given what 'i' stands
> for here).]

Actually, I have

> version
               _
platform       x86_64-unknown-linux-gnu
arch           x86_64
os             linux-gnu
system         x86_64, linux-gnu
status         Under development (unstable)
major          2
minor          3.0
year           2006
month          01
day            10
svn rev        37041
language       R
version.string Version 2.3.0 Under development (unstable) (2006-01-10 r37041)

Which is in fact Pentium 4 EMT (not to be confused with ia64 which is
Itanium, AFAIR).

I don't see how anyone could figure that out from the information
given, though...
   
        -p

> On Tue, 10 Jan 2006, Robin Hankin wrote:
>
> >
> >
> > 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
> >
> > ______________________________________________
> > [hidden email] mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-devel
> >
> >
>
> --
> Brian D. Ripley,                  [hidden email]
> Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
> University of Oxford,             Tel:  +44 1865 272861 (self)
> 1 South Parks Road,                     +44 1865 272866 (PA)
> Oxford OX1 3TG, UK                Fax:  +44 1865 272595
>

--
   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark          Ph:  (+45) 35327918
~~~~~~~~~~ - ([hidden email])                  FAX: (+45) 35327907

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Reply | Threaded
Open this post in threaded view
|

Re: eigen()

Martin Maechler
In reply to this post by Prof Brian Ripley
>>>>> "BDR" == Prof Brian Ripley <[hidden email]>
>>>>>     on Tue, 10 Jan 2006 15:01:00 +0000 (GMT) writes:

    BDR> I haven't seen most of this thread, but this is a classic case of passing
    BDR> integers instead of doubles.  And indeed

    BDR> else if(is.numeric(x)) {
    BDR> storage.mode(x) <- "double"

    BDR> has been removed from eigen.R in R-devel in r36952.  So that's the
    BDR> culprit.

and I am the culprit of that revision.  I'll fix this ASAP.
Martin

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel