R help - May 2005 - eigenvalues of a circulant matrix

Hi,

It is my understanding that the eigenvectors of a circulant matrix are given as
follows:

1,omega,omega^2,....,omega^{p-1}

where the matrix has dimension given by p x p and omega is one of p complex
roots of unity. (See Bellman for an excellent discussion on this). 

The matrix created  by the attached row and obtained using the following
commands
indicates no imaginary parts for the eigenvectors. It appears that the real
values are close, but not exactly so, and there is no imaginary part
whatsoever. 

x<-scan("kinv.dat")       #length(x) = 216
y<-x[c(109:216,1:108)]
X<-toeplitz(y)
eigen(X)$vectors

Note that the eigenvectors are correct, and they are indeed real, because X is
symmetric.

Is this a bug in R? Any insight if not, please!

Many thanks and best wishes!


This is unrelated, but can the R-help archive maintainers please not put e-mail
addresses in the archive? This would really help people like me who would like
to post using their professional e-mail addresses. Just stripping the e-mail
address from everything else would be great, or make it non-spammable by adding
some random number or something which would be obvious to anyone reading it
without the help of a machine. After all, why give spider programs more fodder?

Best wishes!

__________________________________________________




__________________________________________________




__________________________________________________

On Sun, 1 May 2005, someone who didn't give his name wrote:
> It is my understanding that the eigenvectors of a circulant matrix are 
> given as follows:
>
> 1,omega,omega^2,....,omega^{p-1}
>
> where the matrix has dimension given by p x p and omega is one of p complex
> roots of unity. (See Bellman for an excellent discussion on this).
What is the relevance of this?  Also, your reference is useless to us, 
which is important as this all hinges on your definitions.
> The matrix created by the attached row and obtained using the following 
> commands indicates no imaginary parts for the eigenvectors. It appears 
> that the real values are close, but not exactly so, and there is no 
> imaginary part whatsoever.
>
> x<-scan("kinv.dat")       #length(x) = 216
> y<-x[c(109:216,1:108)]
> X<-toeplitz(y)
> eigen(X)$vectors
We don't have "kinv.dat", but X is not circulant as usually
defined.
> Note that the eigenvectors are correct, and they are indeed real, 
> because X is symmetric.
>
> Is this a bug in R? Any insight if not, please!
Well, first R calls LAPACK or EISPACK, so it would be a bug in one of 
those.  But in so far as I understand you, X is a real symmetric matrix, 
and those have real eigenvalues and eigenvectors.

I think you are confused about the meaning of Toeplitz and circulant.
Compare

http://mathworld.wolfram.com/CirculantMatrix.html
http://mathworld.wolfram.com/ToeplitzMatrix.html

and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix.

There is a very regretable tendency here for people to assume their 
lack of understanding is `a bug in R'.

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
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

Dear Professor Ripley:

Lets do this professionally, shall we?

--- Prof Brian Ripley <ripley at stats.ox.ac.uk>
wrote:
> On Sun, 1 May 2005, someone who didn't give his name wrote:
> 
> > It is my understanding that the eigenvectors of a circulant matrix are
> > given as follows:
> >
> > 1,omega,omega^2,....,omega^{p-1}
> >
> > where the matrix has dimension given by p x p and omega is one of p
complex
> > roots of unity. (See Bellman for an excellent discussion on this).
> 
> What is the relevance of this?  Also, your reference is useless to us, 
> which is important as this all hinges on your definitions.
Bellman is an excellent book on the topic and that is what I was alluding to,
that you can calculate the eigendecomposition by hand without any costly
computations, actually.
> 
> > The matrix created by the attached row and obtained using the
following
> > commands indicates no imaginary parts for the eigenvectors. It appears
> > that the real values are close, but not exactly so, and there is no 
> > imaginary part whatsoever.
> >
> > x<-scan("kinv.dat")       #length(x) = 216
> > y<-x[c(109:216,1:108)]
> > X<-toeplitz(y)
> > eigen(X)$vectors
> 
> We don't have "kinv.dat", but X is not circulant as usually
defined.
Sorry about kinv.dat -- in the e-mail that came back to me, it read
"kinv",
btw.

I am unclear why you say that "X is not circulant as usually defined"
-- do you
think you could clarify? It is true I use a Toeplitz matrix to set this up, but
how does that matter? The end result in this case is still a circulant matrix
that is symmetric, is it not? I would like to know why my result is not
circulant here.
> > Note that the eigenvectors are correct, and they are indeed real, 
> > because X is symmetric.
> >
> > Is this a bug in R? Any insight if not, please!
> 
> Well, first R calls LAPACK or EISPACK, so it would be a bug in one of 
> those.  But in so far as I understand you, X is a real symmetric matrix, 
> and those have real eigenvalues and eigenvectors.
Yes, I know that R calls LAPACK (which now contains EISPACK, btw). But I also
know that LAPACK contains complex eigendecomposition routines in addition to
double precision ones and it would need to be used if there is reason to
believe that the result is complex valued. (In particular ZGESDD would do it.)

The eigendecomposition of a matrix is unique. Whatever you think of Bellman,
the book does show how the eigenvectors of a circulant matrix are given by the
complex roots of unity as given above. We have therefore exhibited an
eigendecomposition without actually going through major computations (which is
good, because statistical computing is best when you use it sparingly). 
Why then does the result differ from that in R, and why by so much? (After all,
the eigendecomposition is unique, or is that only fpr real matrices?)
> I think you are confused about the meaning of Toeplitz and circulant.
Unclear, but would like to hear about your views on the actual differences in
this specific example.
> Compare
> 
> http://mathworld.wolfram.com/CirculantMatrix.html
> http://mathworld.wolfram.com/ToeplitzMatrix.html
> 
> and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix.
In my case, my matrix is symmetric and the result is a circulant matrix.
> There is a very regretable tendency here for people to assume their 
> lack of understanding is `a bug in R'.
True, but bugs in software are not exactly rare. Though R does have very few
bugs and which is why I recommend the software to every Tom.

Besides I asked a question here because I was confused....

Many thanks and best wishes!
> -- 
> Brian D. Ripley,                  ripley at stats.ox.ac.uk
> 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
>

OK, lets redo this again, and ensure that we start with a row that will indeed
lead to a symmetric matrix for the circulant matrix:

x<-scan("kinv")
y<-x[c(109:1,2:108)]

X=toeplitz(y)
Z=y
for (i in 2:216) Z=rbind(Z,y[c((216-i+2):216,1:(216-i+1))])

range(X-Z)
[1] 0 0

eigen(X) is the same as eigen(Z), but we know that Z is a circulant matrix so
the eigenvectors are complex....

Any thoughts/screams?


--- Prof Brian Ripley <ripley at stats.ox.ac.uk>
wrote:
> On Sun, 1 May 2005, someone who didn't give his name wrote:
> 
> > It is my understanding that the eigenvectors of a circulant matrix are
> > given as follows:
> >
> > 1,omega,omega^2,....,omega^{p-1}
> >
> > where the matrix has dimension given by p x p and omega is one of p
complex
> > roots of unity. (See Bellman for an excellent discussion on this).
> 
> What is the relevance of this?  Also, your reference is useless to us, 
> which is important as this all hinges on your definitions.
> 
> > The matrix created by the attached row and obtained using the
following
> > commands indicates no imaginary parts for the eigenvectors. It appears
> > that the real values are close, but not exactly so, and there is no 
> > imaginary part whatsoever.
> >
> > x<-scan("kinv.dat")       #length(x) = 216
> > y<-x[c(109:216,1:108)]
> > X<-toeplitz(y)
> > eigen(X)$vectors
> 
> We don't have "kinv.dat", but X is not circulant as usually
defined.
> 
> > Note that the eigenvectors are correct, and they are indeed real, 
> > because X is symmetric.
> >
> > Is this a bug in R? Any insight if not, please!
> 
> Well, first R calls LAPACK or EISPACK, so it would be a bug in one of 
> those.  But in so far as I understand you, X is a real symmetric matrix, 
> and those have real eigenvalues and eigenvectors.
> 
> I think you are confused about the meaning of Toeplitz and circulant.
> Compare
> 
> http://mathworld.wolfram.com/CirculantMatrix.html
> http://mathworld.wolfram.com/ToeplitzMatrix.html
> 
> and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix.
> 
> There is a very regretable tendency here for people to assume their 
> lack of understanding is `a bug in R'.
> 
> -- 
> Brian D. Ripley,                  ripley at stats.ox.ac.uk
> 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
>

I just Googled around a bit and found definitions of Toeplitz and
circulant matrices as follows:

A Toeplitz matrix is any n x n matrix with values constant along each
(top-left to lower-right) diagonal.  matrix has the form

	a_0 a_1 .   .  .   .  ... a_{n-1}
	a_{-1} a_0 a_1        ... a_{n-2}
	a_{-2} a_{-1} a_0 a_1 ...    .
	   .      .    .   .   .     .
	   .      .    .   .   .     .
	   .      .    .   .   .     .
	a_{-(n-1)} a_{-(n-2)} ... a_1 a_0

(A Toeplitz matrix ***may*** be symmetric.)

A circulant matrix is an n x n matrix whose rows are composed of
cyclically shifted versions of a length-n vector. For example, the
circulant matrix on the vector (1, 2, 3, 4)  is

	4 1 2 3
	3 4 1 2
	2 3 4 1
	1 2 3 4

So circulant matrices are a special case of Toeplitz matrices.
However a circulant matrix cannot be symmetric.

The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i,
2 - 2i, and 2 --- certainly not roots of unity.  Bellman may have
been talking about the particular (important) case of a circulant
matrix where the vector from which it is constructed is a canonical
basis vector e_i with a 1 in the i-th slot and zeroes elsewhere.

Such a matrix is in fact a unitary matrix (operator), whence its
spectrum is contained in the unit circle; its eigenvalues are indeed
n-th roots of unity.

Such matrices are related to the unilateral shift operator on
Hilbert space (which is the ``primordial'' Toeplitz operator).
It arises as multiplication by z on H^2 --- the ``analytic''
elements of L^2 of the unit circle.

On (infinite dimensional) Hilbert space the unilateral shift
looks like

	0 0 0 0 0 ...
	1 0 0 0 0 ...
	0 1 0 0 0 ...
	0 0 1 0 0 ...
        . . . . . ...
        . . . . . ...

which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ...  on and on
forever.  On (say) 4 dimensional space we can have a unilateral
shift operator/matrix

	0 0 0 0
	1 0 0 0
	0 1 0 0
	0 0 1 0

but its range is a 3 dimensional subspace (e_4 gets ``killed'').

The ``corresponding'' circulant matrix is

	0 0 0 1
	1 0 0 0
	0 1 0 0
	0 0 1 0

which is an onto mapping --- e_4 gets sent back to e_1.

I hope this clears up some of the confusion.

				cheers,

					Rolf Turner
					rolf at math.unb.ca

--- Rolf Turner <rolf at math.unb.ca> wrote:
> I just Googled around a bit and found definitions of Toeplitz and
> circulant matrices as follows:
> 
> A Toeplitz matrix is any n x n matrix with values constant along each
> (top-left to lower-right) diagonal.  matrix has the form
> 
> 	a_0 a_1 .   .  .   .  ... a_{n-1}
> 	a_{-1} a_0 a_1        ... a_{n-2}
> 	a_{-2} a_{-1} a_0 a_1 ...    .
> 	   .      .    .   .   .     .
> 	   .      .    .   .   .     .
> 	   .      .    .   .   .     .
> 	a_{-(n-1)} a_{-(n-2)} ... a_1 a_0
> 
> (A Toeplitz matrix ***may*** be symmetric.)
Agreed. As may a circulant matrix if a_i = a_{p-i+2}
> 
> A circulant matrix is an n x n matrix whose rows are composed of
> cyclically shifted versions of a length-n vector. For example, the
> circulant matrix on the vector (1, 2, 3, 4)  is
> 
> 	4 1 2 3
> 	3 4 1 2
> 	2 3 4 1
> 	1 2 3 4
> 
> So circulant matrices are a special case of Toeplitz matrices.
> However a circulant matrix cannot be symmetric.
> 
> The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i,
> 2 - 2i, and 2 --- certainly not roots of unity. 
The eigenvalues are 4+1*omega+2*omega^2+3*omega^3.
omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above
holds.

 Bellman may have
> been talking about the particular (important) case of a circulant
> matrix where the vector from which it is constructed is a canonical
> basis vector e_i with a 1 in the i-th slot and zeroes elsewhere.
No, that is not true: his result can be verified for any circulant matrix,
directly.
> Such a matrix is in fact a unitary matrix (operator), whence its
> spectrum is contained in the unit circle; its eigenvalues are indeed
> n-th roots of unity.
> 
> Such matrices are related to the unilateral shift operator on
> Hilbert space (which is the ``primordial'' Toeplitz operator).
> It arises as multiplication by z on H^2 --- the ``analytic''
> elements of L^2 of the unit circle.
> 
> On (infinite dimensional) Hilbert space the unilateral shift
> looks like
> 
> 	0 0 0 0 0 ...
> 	1 0 0 0 0 ...
> 	0 1 0 0 0 ...
> 	0 0 1 0 0 ...
>         . . . . . ...
>         . . . . . ...
> 
> which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ...  on and on
> forever.  On (say) 4 dimensional space we can have a unilateral
> shift operator/matrix
> 
> 	0 0 0 0
> 	1 0 0 0
> 	0 1 0 0
> 	0 0 1 0
> 
> but its range is a 3 dimensional subspace (e_4 gets ``killed'').
> 
> The ``corresponding'' circulant matrix is
> 
> 	0 0 0 1
> 	1 0 0 0
> 	0 1 0 0
> 	0 0 1 0
> 
> which is an onto mapping --- e_4 gets sent back to e_1.
> 
> I hope this clears up some of the confusion.
> 
> 				cheers,
> 
> 					Rolf Turner
> 					rolf at math.unb.ca
Many thanks and best wishes!

It's hard to argue against the fact that a real symmetric matrix has real
eigenvalues. The eigenvalues of the circulant matrix with first row v are
*polynomials* (not the roots of 1 themselves, unless as Rolf suggested you
start with a vector with all zeros except one 1) in the roots of 1, with
coefficients equal to the entries in v. This is the finite Fourier transform
of v, by the way, and takes real values when the coefficients are real and
symmetric, ie when the matrix is symmetric.

Reid Huntsinger

-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe Trotter
Sent: Monday, May 02, 2005 10:23 AM
To: Rolf Turner
Cc: r-help at stat.math.ethz.ch
Subject: Re: [R] eigenvalues of a circulant matrix



--- Rolf Turner <rolf at math.unb.ca> wrote:
> I just Googled around a bit and found definitions of Toeplitz and
> circulant matrices as follows:
> 
> A Toeplitz matrix is any n x n matrix with values constant along each
> (top-left to lower-right) diagonal.  matrix has the form
> 
> 	a_0 a_1 .   .  .   .  ... a_{n-1}
> 	a_{-1} a_0 a_1        ... a_{n-2}
> 	a_{-2} a_{-1} a_0 a_1 ...    .
> 	   .      .    .   .   .     .
> 	   .      .    .   .   .     .
> 	   .      .    .   .   .     .
> 	a_{-(n-1)} a_{-(n-2)} ... a_1 a_0
> 
> (A Toeplitz matrix ***may*** be symmetric.)
Agreed. As may a circulant matrix if a_i = a_{p-i+2}
> 
> A circulant matrix is an n x n matrix whose rows are composed of
> cyclically shifted versions of a length-n vector. For example, the
> circulant matrix on the vector (1, 2, 3, 4)  is
> 
> 	4 1 2 3
> 	3 4 1 2
> 	2 3 4 1
> 	1 2 3 4
> 
> So circulant matrices are a special case of Toeplitz matrices.
> However a circulant matrix cannot be symmetric.
> 
> The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i,
> 2 - 2i, and 2 --- certainly not roots of unity. 
The eigenvalues are 4+1*omega+2*omega^2+3*omega^3.
omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above
holds.

 Bellman may have
> been talking about the particular (important) case of a circulant
> matrix where the vector from which it is constructed is a canonical
> basis vector e_i with a 1 in the i-th slot and zeroes elsewhere.
No, that is not true: his result can be verified for any circulant matrix,
directly.
> Such a matrix is in fact a unitary matrix (operator), whence its
> spectrum is contained in the unit circle; its eigenvalues are indeed
> n-th roots of unity.
> 
> Such matrices are related to the unilateral shift operator on
> Hilbert space (which is the ``primordial'' Toeplitz operator).
> It arises as multiplication by z on H^2 --- the ``analytic''
> elements of L^2 of the unit circle.
> 
> On (infinite dimensional) Hilbert space the unilateral shift
> looks like
> 
> 	0 0 0 0 0 ...
> 	1 0 0 0 0 ...
> 	0 1 0 0 0 ...
> 	0 0 1 0 0 ...
>         . . . . . ...
>         . . . . . ...
> 
> which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ...  on and on
> forever.  On (say) 4 dimensional space we can have a unilateral
> shift operator/matrix
> 
> 	0 0 0 0
> 	1 0 0 0
> 	0 1 0 0
> 	0 0 1 0
> 
> but its range is a 3 dimensional subspace (e_4 gets ``killed'').
> 
> The ``corresponding'' circulant matrix is
> 
> 	0 0 0 1
> 	1 0 0 0
> 	0 1 0 0
> 	0 0 1 0
> 
> which is an onto mapping --- e_4 gets sent back to e_1.
> 
> I hope this clears up some of the confusion.
> 
> 				cheers,
> 
> 					Rolf Turner
> 					rolf at math.unb.ca
Many thanks and best wishes!

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

When the matrix is symmetric and omega is not real, omega and its conjugate
(= inverse) give the same eigenvalue, so you have a 2-dimensional
eigenspace. R chooses a real basis of this, which is perfectly fine since
it's not looking for circulant structure.

For example,
> m
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    2    3    3    2
[2,]    2    1    2    3    3
[3,]    3    2    1    2    3
[4,]    3    3    2    1    2
[5,]    2    3    3    2    1
> eigen(m)
$values
[1] 11.000000 -0.381966 -0.381966 -2.618034 -2.618034

$vectors
          [,1]      [,2]       [,3]       [,4]      [,5]
[1,] 0.4472136  0.000000 -0.6324555  0.6324555  0.000000
[2,] 0.4472136  0.371748  0.5116673  0.1954395  0.601501
[3,] 0.4472136 -0.601501 -0.1954395 -0.5116673  0.371748
[4,] 0.4472136  0.601501 -0.1954395 -0.5116673 -0.371748
[5,] 0.4472136 -0.371748  0.5116673  0.1954395 -0.601501

and you can match these columns up with the "canonical" eigenvectors
exp(2*pi*1i*(0:4)*j/5) for j = 0,1,2,3,4. E.g.,
> Im(exp(2*pi*1i*(0:4)*3/5))
[1]  0.0000000 -0.5877853  0.9510565 -0.9510565  0.5877853

which can be seen to be a scalar multiple of column 2. 

Reid Huntsinger

Reid Huntsinger

-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Huntsinger, Reid
Sent: Monday, May 02, 2005 10:43 AM
To: 'Globe Trotter'; Rolf Turner
Cc: r-help at stat.math.ethz.ch
Subject: RE: [R] eigenvalues of a circulant matrix


It's hard to argue against the fact that a real symmetric matrix has real
eigenvalues. The eigenvalues of the circulant matrix with first row v are
*polynomials* (not the roots of 1 themselves, unless as Rolf suggested you
start with a vector with all zeros except one 1) in the roots of 1, with
coefficients equal to the entries in v. This is the finite Fourier transform
of v, by the way, and takes real values when the coefficients are real and
symmetric, ie when the matrix is symmetric.

Reid Huntsinger

-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe Trotter
Sent: Monday, May 02, 2005 10:23 AM
To: Rolf Turner
Cc: r-help at stat.math.ethz.ch
Subject: Re: [R] eigenvalues of a circulant matrix



--- Rolf Turner <rolf at math.unb.ca> wrote:
> I just Googled around a bit and found definitions of Toeplitz and
> circulant matrices as follows:
> 
> A Toeplitz matrix is any n x n matrix with values constant along each
> (top-left to lower-right) diagonal.  matrix has the form
> 
> 	a_0 a_1 .   .  .   .  ... a_{n-1}
> 	a_{-1} a_0 a_1        ... a_{n-2}
> 	a_{-2} a_{-1} a_0 a_1 ...    .
> 	   .      .    .   .   .     .
> 	   .      .    .   .   .     .
> 	   .      .    .   .   .     .
> 	a_{-(n-1)} a_{-(n-2)} ... a_1 a_0
> 
> (A Toeplitz matrix ***may*** be symmetric.)
Agreed. As may a circulant matrix if a_i = a_{p-i+2}
> 
> A circulant matrix is an n x n matrix whose rows are composed of
> cyclically shifted versions of a length-n vector. For example, the
> circulant matrix on the vector (1, 2, 3, 4)  is
> 
> 	4 1 2 3
> 	3 4 1 2
> 	2 3 4 1
> 	1 2 3 4
> 
> So circulant matrices are a special case of Toeplitz matrices.
> However a circulant matrix cannot be symmetric.
> 
> The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i,
> 2 - 2i, and 2 --- certainly not roots of unity. 
The eigenvalues are 4+1*omega+2*omega^2+3*omega^3.
omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above
holds.

 Bellman may have
> been talking about the particular (important) case of a circulant
> matrix where the vector from which it is constructed is a canonical
> basis vector e_i with a 1 in the i-th slot and zeroes elsewhere.
No, that is not true: his result can be verified for any circulant matrix,
directly.
> Such a matrix is in fact a unitary matrix (operator), whence its
> spectrum is contained in the unit circle; its eigenvalues are indeed
> n-th roots of unity.
> 
> Such matrices are related to the unilateral shift operator on
> Hilbert space (which is the ``primordial'' Toeplitz operator).
> It arises as multiplication by z on H^2 --- the ``analytic''
> elements of L^2 of the unit circle.
> 
> On (infinite dimensional) Hilbert space the unilateral shift
> looks like
> 
> 	0 0 0 0 0 ...
> 	1 0 0 0 0 ...
> 	0 1 0 0 0 ...
> 	0 0 1 0 0 ...
>         . . . . . ...
>         . . . . . ...
> 
> which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ...  on and on
> forever.  On (say) 4 dimensional space we can have a unilateral
> shift operator/matrix
> 
> 	0 0 0 0
> 	1 0 0 0
> 	0 1 0 0
> 	0 0 1 0
> 
> but its range is a 3 dimensional subspace (e_4 gets ``killed'').
> 
> The ``corresponding'' circulant matrix is
> 
> 	0 0 0 1
> 	1 0 0 0
> 	0 1 0 0
> 	0 0 1 0
> 
> which is an onto mapping --- e_4 gets sent back to e_1.
> 
> I hope this clears up some of the confusion.
> 
> 				cheers,
> 
> 					Rolf Turner
> 					rolf at math.unb.ca
Many thanks and best wishes!

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

----------------------------------------------------------------------------
--
Notice:  This e-mail message, together with any attachments,...{{dropped}}

The example that I submitted earlier in the day. Would you like me to send
again?

Thanks!


--- "Huntsinger, Reid" <reid_huntsinger at merck.com> wrote:
> For which X?
> 
> Reid Huntsinger
> 
> -----Original Message-----
> From: Globe Trotter [mailto:itsme_410 at yahoo.com] 
> Sent: Monday, May 02, 2005 2:34 PM
> To: Huntsinger, Reid; Rolf Turner
> Cc: r-help at stat.math.ethz.ch
> Subject: RE: [R] eigenvalues of a circulant matrix
> 
> 
> By the way, I just noticed that eigen(X) returns eigenvectors, at least two
> of
> which are NaN's. 
> 
> Best wishes!
> 
> --- "Huntsinger, Reid" <reid_huntsinger at merck.com>
wrote:
> 
> > When the matrix is symmetric and omega is not real, omega and its
> conjugate
> > (= inverse) give the same eigenvalue, so you have a 2-dimensional
> > eigenspace. R chooses a real basis of this, which is perfectly fine
since
> > it's not looking for circulant structure.
> > 
> > For example,
> > 
> > > m
> >      [,1] [,2] [,3] [,4] [,5]
> > [1,]    1    2    3    3    2
> > [2,]    2    1    2    3    3
> > [3,]    3    2    1    2    3
> > [4,]    3    3    2    1    2
> > [5,]    2    3    3    2    1
> > 
> > > eigen(m)
> > $values
> > [1] 11.000000 -0.381966 -0.381966 -2.618034 -2.618034
> > 
> > $vectors
> >           [,1]      [,2]       [,3]       [,4]      [,5]
> > [1,] 0.4472136  0.000000 -0.6324555  0.6324555  0.000000
> > [2,] 0.4472136  0.371748  0.5116673  0.1954395  0.601501
> > [3,] 0.4472136 -0.601501 -0.1954395 -0.5116673  0.371748
> > [4,] 0.4472136  0.601501 -0.1954395 -0.5116673 -0.371748
> > [5,] 0.4472136 -0.371748  0.5116673  0.1954395 -0.601501
> > 
> > and you can match these columns up with the "canonical"
eigenvectors
> > exp(2*pi*1i*(0:4)*j/5) for j = 0,1,2,3,4. E.g.,
> > 
> > > Im(exp(2*pi*1i*(0:4)*3/5))
> > [1]  0.0000000 -0.5877853  0.9510565 -0.9510565  0.5877853
> > 
> > which can be seen to be a scalar multiple of column 2. 
> > 
> > Reid Huntsinger
> > 
> > Reid Huntsinger
> > 
> > -----Original Message-----
> > From: r-help-bounces at stat.math.ethz.ch
> > [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Huntsinger,
Reid
> > Sent: Monday, May 02, 2005 10:43 AM
> > To: 'Globe Trotter'; Rolf Turner
> > Cc: r-help at stat.math.ethz.ch
> > Subject: RE: [R] eigenvalues of a circulant matrix
> > 
> > 
> > It's hard to argue against the fact that a real symmetric matrix
has real
> > eigenvalues. The eigenvalues of the circulant matrix with first row v
are
> > *polynomials* (not the roots of 1 themselves, unless as Rolf suggested
you
> > start with a vector with all zeros except one 1) in the roots of 1,
with
> > coefficients equal to the entries in v. This is the finite Fourier
> transform
> > of v, by the way, and takes real values when the coefficients are real
and
> > symmetric, ie when the matrix is symmetric.
> > 
> > Reid Huntsinger
> > 
> > -----Original Message-----
> > From: r-help-bounces at stat.math.ethz.ch
> > [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe
Trotter
> > Sent: Monday, May 02, 2005 10:23 AM
> > To: Rolf Turner
> > Cc: r-help at stat.math.ethz.ch
> > Subject: Re: [R] eigenvalues of a circulant matrix
> > 
> > 
> > 
> > --- Rolf Turner <rolf at math.unb.ca> wrote:
> > > I just Googled around a bit and found definitions of Toeplitz and
> > > circulant matrices as follows:
> > > 
> > > A Toeplitz matrix is any n x n matrix with values constant along
each
> > > (top-left to lower-right) diagonal.  matrix has the form
> > > 
> > > 	a_0 a_1 .   .  .   .  ... a_{n-1}
> > > 	a_{-1} a_0 a_1        ... a_{n-2}
> > > 	a_{-2} a_{-1} a_0 a_1 ...    .
> > > 	   .      .    .   .   .     .
> > > 	   .      .    .   .   .     .
> > > 	   .      .    .   .   .     .
> > > 	a_{-(n-1)} a_{-(n-2)} ... a_1 a_0
> > > 
> > > (A Toeplitz matrix ***may*** be symmetric.)
> > 
> > Agreed. As may a circulant matrix if a_i = a_{p-i+2}
> > 
> > > 
> > > A circulant matrix is an n x n matrix whose rows are composed of
> > > cyclically shifted versions of a length-n vector. For example,
the
> > > circulant matrix on the vector (1, 2, 3, 4)  is
> > > 
> > > 	4 1 2 3
> > > 	3 4 1 2
> > > 	2 3 4 1
> > > 	1 2 3 4
> > > 
> > > So circulant matrices are a special case of Toeplitz matrices.
> > > However a circulant matrix cannot be symmetric.
> > > 
> > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i,
> > > 2 - 2i, and 2 --- certainly not roots of unity. 
> > 
> > The eigenvalues are 4+1*omega+2*omega^2+3*omega^3.
> > omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the
> above
> > holds.
> > 
> >  Bellman may have
> > > been talking about the particular (important) case of a circulant
> > > matrix where the vector from which it is constructed is a
canonical
> > > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere.
> > 
> > No, that is not true: his result can be verified for any circulant
matrix,
> > directly.
> > 
> > > Such a matrix is in fact a unitary matrix (operator), whence its
> > > spectrum is contained in the unit circle; its eigenvalues are
indeed
> > > n-th roots of unity.
> > > 
> > > Such matrices are related to the unilateral shift operator on
> > > Hilbert space (which is the ``primordial'' Toeplitz
operator).
> > > It arises as multiplication by z on H^2 --- the
``analytic''
> > > elements of L^2 of the unit circle.
> > > 
> > > On (infinite dimensional) Hilbert space the unilateral shift
> > > looks like
> > > 
> > > 	0 0 0 0 0 ...
> > > 	1 0 0 0 0 ...
> > > 	0 1 0 0 0 ...
> > > 	0 0 1 0 0 ...
> > >         . . . . . ...
> > >         . . . . . ...
> > > 
> > > which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ...  on and on
> > > forever.  On (say) 4 dimensional space we can have a unilateral
> > > shift operator/matrix
> > > 
> > > 	0 0 0 0
> > > 	1 0 0 0
> > > 	0 1 0 0
> > > 	0 0 1 0
> > > 
> > > but its range is a 3 dimensional subspace (e_4 gets
``killed'').
> > > 
> > > The ``corresponding'' circulant matrix is
> > > 
> > > 	0 0 0 1
> > > 	1 0 0 0
> > > 	0 1 0 0
> > > 	0 0 1 0
> > > 
> > > which is an onto mapping --- e_4 gets sent back to e_1.
> > > 
> > > I hope this clears up some of the confusion.
> > > 
> > > 				cheers,
> > > 
> > > 					Rolf Turner
> > > 					rolf at math.unb.ca
> > 
> > Many thanks and best wishes!
> > 
> > ______________________________________________
> > R-help at stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide!
> > http://www.R-project.org/posting-guide.html
> > 
> > ______________________________________________
> > R-help at stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide!
> > http://www.R-project.org/posting-guide.html
> > 
> >
>
----------------------------------------------------------------------------
> > --
> > Notice:  This e-mail message, together with any attachments, contains
> > information of Merck & Co., Inc. (One Merck Drive, Whitehouse
Station, New
> > Jersey, USA 08889), and/or its affiliates (which may be known outside
the
> > United States as Merck Frosst, Merck Sharp & Dohme or MSD and in
Japan, as
> > Banyu) that may be confidential, proprietary copyrighted and/or
legally
> > privileged. It is intended solely for the use of the individual or
entity
> > named on this message.  If you are not the intended recipient, and
have
> > received this message in error, please notify us immediately by reply
> e-mail
> > and then delete it from your system.
> >
>
----------------------------------------------------------------------------
> > --
> > 
> > 
> > 
> > 
> >
>
----------------------------------------------------------------------------
> --
> > Notice:  This e-mail message, together with any attachments, contains
> > information of Merck & Co., Inc. (One Merck Drive, Whitehouse
Station, New
> > Jersey, USA 08889), and/or its affiliates (which may be known outside
the
> > United States as Merck Frosst, Merck Sharp & Dohme or MSD and in
Japan, as
> > Banyu) that may be confidential, proprietary copyrighted and/or
legally
> > privileged. It is intended solely for the use of the individual or
entity
> > named on this message.  If you are not the intended recipient, and
have
> > received this message in error, please notify us immediately by reply
> e-mail
> > and then delete it from your system.
> >
>
----------------------------------------------------------------------------
> --
> > 
> 
> 
> __________________________________________________
> Do You Yahoo!?
> http://mail.yahoo.com 
> 
> 
> 
> 
> 
>
------------------------------------------------------------------------------
> Notice:  This e-mail message, together with any attachment...{{dropped}}

The construction was

y<-x[c(109:1,2:108)]

so y is symmetric in the sense of the usual way of writing a function on
integers mod n as a vector with 1-based indexing. I.e., y[i+1] = y[n-(i+1)]
for i=0,1,...,n-1. So the assignment

Z <- toeplitz(y)

*does* create a symmetric circulant matrix. It is diagonalizable but does
not have distinct eigenvalues, hence the eigenspaces may be more than
one-dimensional, so you can't just pick a unit vector and call it
"the"
eigenvector for that eigenvalue. You choose a basis for each eigenspace. R
detects the symmetry:

...
symmetric: if `TRUE', the matrix is assumed to be symmetric (or
          Hermitian if complex) and only its lower triangle is used. If
          `symmetric' is not specified, the matrix is inspected for
          symmetry.

(from help(eigen))
and knows that computations can be done with real arithmetic.

As for why you get NaN, you should submit --along with your example--
details of your platform (machine, R version, how R was built and installed,
etc).

Reid Huntsinger

-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of
Ted.Harding at nessie.mcc.ac.uk
Sent: Monday, May 02, 2005 5:28 PM
To: r-help at stat.math.ethz.ch
Subject: Re: [R] eigenvalues of a circulant matrix


On 02-May-05 Ted Harding wrote:
> On 02-May-05 Rolf Turner wrote:
>> I just Googled around a bit and found definitions of Toeplitz and
>> circulant matrices as follows:
>> [...]
>> A circulant matrix is an n x n matrix whose rows are composed of
>> cyclically shifted versions of a length-n vector. For example, the
>> circulant matrix on the vector (1, 2, 3, 4)  is
>> 
>>       4 1 2 3
>>       3 4 1 2
>>       2 3 4 1
>>       1 2 3 4
>> 
>> So circulant matrices are a special case of Toeplitz matrices.
>> However a circulant matrix cannot be symmetric.
> 
> I suspect the confusion may lie in what's meant by "cyclically
> shifted". In Rolf's example above, each row is shifted right by 1
> and the one that falls off the end is put at the beginning. This
> cannot be symmetric for general values in the fist row.
> 
> However, if you shift left instead, then you get
> 
>         4 1 2 3
>         1 2 3 4
>         2 3 4 1
>         3 4 1 2
> 
> and this *is* symmetric (and indeed will always be so, for
> general values in the first row).
I just had a look at ?toeplitz

(We should have done that earlier!)

toeplitz                package:stats                R Documentation
Form Symmetric Toeplitz Matrix
     *********

Description:
     Forms a symmetric Toeplitz matrix given its first row.
             *********
[...]
Examples:

     x <- 1:5
     toeplitz (x)
> x <- 1:5
>      toeplitz (x)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    2    3    4    5
[2,]    2    1    2    3    4
[3,]    3    2    1    2    3
[4,]    4    3    2    1    2
[5,]    5    4    3    2    1

Since "Globe Trotter's" construction was

  Y<-toeplitz(x)

it's not surprising what he got (and it *certainly* wasn't
a circulant!!!).

Everybody barking up the wring tree here!

Best wishes to all,
Ted.


--------------------------------------------------------------------
E-Mail: (Ted Harding) <Ted.Harding at nessie.mcc.ac.uk>
Fax-to-email: +44 (0)870 094 0861
Date: 02-May-05                                       Time: 22:27:32
------------------------------ XFMail ------------------------------

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

OK, here we go:

I am submitting two attachments. The first is the datafile called kinv used to
create my circulant matrix, using the following commands:


x<-scan("kinv")
y<-x[c(109:1,0:108)]
X=toeplitz(y)
eigen(X)
write(X,ncol=216,file="test.dat")

reports the following columns full of NaN's: 18, 58, 194, 200. (Note that
eigen(X,symmetric=T) makes no difference and I get the same as above).

The second attachment contains only the eigenvectors obtained on calling a
LAPACK routine directly (from C). The eigenvalues are essentially the same as
that obtained using R. Here, I use the LAPACK-recommended double precision
routine dspevd() routine for symmetric matrices in packed storage format. Note
the absence of the NaN's....I would be happy to send my C programs to
whoever
is interested.

I am using 

:~> uname -a
Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 i686 i386 GNU/Linux

and R.2.0.1.

Many thanks and best wishes!

Well since I know nothing about this topic I have lurked so far, but here's
my two bob's worth.

Firstly I tried to make sense of Brian's initial reply. I have got no idea
who Bellman is and you have not referenced (his/her) work in a way I can access
the issues you refer to. So I assumed that's exactly what Brian was talking
about.

Secondly. 

toeplitz(1:4)
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    2    1    2    3
[3,]    3    2    1    2
[4,]    4    3    2    1

require(magic)
 circulant(4)
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    4    1    2    3
[3,]    3    4    1    2
[4,]    2    3    4    1

So they are obviously two different things. Although I think you may have
implied (not stated) that the particular combination you were using resulted in
both being exactly the same.

It does appear as if in this case the (X) matrix is circulant. But then I'm
no expert in even such simple things.

Then I had no idea where I was going. So I tried the variations in eigen.

I ran you code
x<-scan("h:/t.txt")
y<-x[c(109:216,1:108)]
X<-toeplitz(y)
 and then 
> X[is.na(X)]
numeric(0)

So I didn't get any NAs

t1 <- eigen(X)$vectors
t2 <- eigen(X,symmetric = TRUE)$vectors
> identical(t1,t2)
[1] TRUE
> 
Then

t2 <- eigen(X,symmetric = TRUE,EISPACK =
TRUE)$vectors
> identical(t1,t2)
[1] FALSE
> 
So there'e obviously more than one way of getting the vectors. Does the
second one make more sense to you?

I also noticed in the eigen help that there are references to issues such as
"IEEE 754 arithmetic","(They may also differ between methods and
between platforms.)" and "or Hermitian if complex". All of these
are out of my competence but they do signal to me that there are issues which
may relate to hardware, digital arithmetic and other things of that ilk.

I added the comment about complex because I have a vague idea that they are
related to imaginary parts that you refer to.

So not coming to any conclusion that makes sense to me, and given that there are
often threads about supposed inaccuracies that have answers such as the digits
you see are not always what are held by the machine I set my options(digits =
22) and noticed that some of the numbers are still going at the 22 decimal place
suggesting that the machine might be incapable of producing perfectly accurate
results using digital arithmetic.

My other big sphere of ignorance is complex numbers.

So I tried
X<-toeplitz(complex(real = y))
t1 <- eigen(X)$vectors
> t1[1:20]
 [1]  0.068041577278880341+0i -0.068041577140546913+0i  0.068041576864811659+0i
-0.068041576452430155+0i
 [5]  0.068041575907139579+0i -0.068041575231135451+0i  0.068041574435267163+0i
-0.068041573525828514+0i
 [9]  0.068041572538722991+0i -0.068041571498323253+0i  0.068041570619888622+0i
-0.068041570256170081+0i
[13]  0.068041568759931989+0i -0.068041566476633147+0i  0.068041563560502477+0i
-0.068041560000305007+0i
[17]  0.068041555538765813+0i -0.068041549792984865+0i  0.068041544123969511+0i
-0.068041537810956801+0i
> t2[1:20]
 [1]  0.068041381743976906 -0.068041381743976850  0.068041381743976781
-0.068041381743976753  0.068041381743976587
 [6] -0.068041381743976725  0.068041381743976920 -0.068041381743976836 
0.068041381743976892 -0.068041381743976781
[11]  0.068041381743976781 -0.068041381743977392  0.068041381743976725
-0.068041381743976753  0.068041381743976753
[16] -0.068041381743976698  0.068041381743976587 -0.068041381743976642 
0.068041381743976698 -0.068041381743976490
> 

Which is again different. I have no idea what I'm doing but you do seem to
get slightly different answers depending upon which method you use. I do not
know if one is superior to the others or where one draws the line in terms of
accuracy.

Tom
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch
> [mailto:r-help-bounces at stat.math.ethz.ch]On Behalf Of Globe Trotter
> Sent: Tuesday, 3 May 2005 10:51 AM
> To: r-help at stat.math.ethz.ch
> Subject: Re: [R] eigenvalues of a circulant matrix
> 
> 
> OK, here we go:
> 
> I am submitting two attachments. The first is the datafile 
> called kinv used to
> create my circulant matrix, using the following commands:
> 
> 
> x<-scan("kinv")
> y<-x[c(109:1,0:108)]
> X=toeplitz(y)
> eigen(X)
> write(X,ncol=216,file="test.dat")
> 
> reports the following columns full of NaN's: 18, 58, 194, 
> 200. (Note that
> eigen(X,symmetric=T) makes no difference and I get the same as above).
> 
> The second attachment contains only the eigenvectors obtained 
> on calling a
> LAPACK routine directly (from C). The eigenvalues are 
> essentially the same as
> that obtained using R. Here, I use the LAPACK-recommended 
> double precision
> routine dspevd() routine for symmetric matrices in packed 
> storage format. Note
> the absence of the NaN's....I would be happy to send my C 
> programs to whoever
> is interested.
> 
> I am using 
> 
> :~> uname -a
> Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 
> i686 i386 GNU/Linux
> 
> and R.2.0.1.
> 
> Many thanks and best wishes!
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! 
http://www.R-project.org/posting-guide.html

How was your R 2.0.1 built? Which Lapack did it link to, and which one does
it load? Which BLAS? Is the BLAS threaded? Does it link to the optimized
pthreads library? Etc. (I'm not a Fedora Core 3 user so I'm not sure
what
the default setup is, and of course I don't know if that's what you
have.)

As I pointed out in my previous post, in eigen() the matrix is "inspected
for symmetry" so symmetric=TRUE is the same as not specifying this at all.
You could try symmetric=FALSE...

Reid Huntsinger

-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe Trotter
Sent: Monday, May 02, 2005 10:51 PM
To: r-help at stat.math.ethz.ch
Subject: Re: [R] eigenvalues of a circulant matrix


OK, here we go:

I am submitting two attachments. The first is the datafile called kinv used
to
create my circulant matrix, using the following commands:


x<-scan("kinv")
y<-x[c(109:1,0:108)]
X=toeplitz(y)
eigen(X)
write(X,ncol=216,file="test.dat")

reports the following columns full of NaN's: 18, 58, 194, 200. (Note that
eigen(X,symmetric=T) makes no difference and I get the same as above).

The second attachment contains only the eigenvectors obtained on calling a
LAPACK routine directly (from C). The eigenvalues are essentially the same
as
that obtained using R. Here, I use the LAPACK-recommended double precision
routine dspevd() routine for symmetric matrices in packed storage format.
Note
the absence of the NaN's....I would be happy to send my C programs to
whoever
is interested.

I am using 

:~> uname -a
Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 i686 i386
GNU/Linux

and R.2.0.1.

Many thanks and best wishes!

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

Hi,

The R was downloaded in binary form (Fedora Core 3 RPM) from a CRAN mirror. I
do not know which LAPACK it links to, but the only LAPACK on my machine is
lapack-3.0-28 (RPM installation). How does one figure out which BLAS or whether
it is threaded?

Many thanks and best wishes!

--- "Huntsinger, Reid" <reid_huntsinger at merck.com>
wrote:
> How was your R 2.0.1 built? Which Lapack did it link to, and which one does
> it load? Which BLAS? Is the BLAS threaded? Does it link to the optimized
> pthreads library? Etc. (I'm not a Fedora Core 3 user so I'm not
sure what
> the default setup is, and of course I don't know if that's what you
have.)
> 
> As I pointed out in my previous post, in eigen() the matrix is
"inspected
> for symmetry" so symmetric=TRUE is the same as not specifying this at
all.
> You could try symmetric=FALSE...
> 
> Reid Huntsinger
> 
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch
> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe Trotter
> Sent: Monday, May 02, 2005 10:51 PM
> To: r-help at stat.math.ethz.ch
> Subject: Re: [R] eigenvalues of a circulant matrix
> 
> 
> OK, here we go:
> 
> I am submitting two attachments. The first is the datafile called kinv used
> to
> create my circulant matrix, using the following commands:
> 
> 
> x<-scan("kinv")
> y<-x[c(109:1,0:108)]
> X=toeplitz(y)
> eigen(X)
> write(X,ncol=216,file="test.dat")
> 
> reports the following columns full of NaN's: 18, 58, 194, 200. (Note
that
> eigen(X,symmetric=T) makes no difference and I get the same as above).
> 
> The second attachment contains only the eigenvectors obtained on calling a
> LAPACK routine directly (from C). The eigenvalues are essentially the same
> as
> that obtained using R. Here, I use the LAPACK-recommended double precision
> routine dspevd() routine for symmetric matrices in packed storage format.
> Note
> the absence of the NaN's....I would be happy to send my C programs to
> whoever
> is interested.
> 
> I am using 
> 
> :~> uname -a
> Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 i686 i386
> GNU/Linux
> 
> and R.2.0.1.
> 
> Many thanks and best wishes!
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html
> 
> 
> 
> 
> 
>
------------------------------------------------------------------------------
> Notice:  This e-mail message, together with any attachment...{{dropped}}

Under Linux run "ldd" on the binary to see what shared libraries the
binary
has been linked against and will attempt to load. The command you run is a
shell script which sets and exports LD_LIBRARY_PATH and then runs R, so to
be sure first start R and execute
> system("echo $LD_LIBRARY_PATH")
and then using this value for LD_LIBRARY_PATH do

$ export LD_LIBRARY_PATH=<R's value for this>
$ ldd `R RHOME`/bin/exec/R

which will list the shared libraries R will try to load to resolve links to
shared libraries. 

Having said that, I have tried your example on several platforms with
several configurations, and I get NaNs only when R is calling the Lapack
routine dsyevr (R's builtin Lapack) with an external (optimized) blas
(either Goto or ATLAS's blas) but not R's builtin blas. Moreover even
with
the optimized blas, if I use dysev instead of dsyevr, (eigen() executes a
.Call("La_rs",x,only.values,"dsyevr",PACKAGE="base"),
I execute this from R
with "dsyev" in place of "dsyevr") I get no NaNs. 

So it looks like a problem with a blas routine used by dsyevr but not dsyev,
but I have yet to confirm.

I note that your C program does not use either of these lapack routines.

Reid Huntsinger


-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe Trotter
Sent: Tuesday, May 03, 2005 12:21 PM
To: Huntsinger, Reid; r-help at stat.math.ethz.ch
Subject: RE: [R] eigenvalues of a circulant matrix


Hi,

The R was downloaded in binary form (Fedora Core 3 RPM) from a CRAN mirror.
I
do not know which LAPACK it links to, but the only LAPACK on my machine is
lapack-3.0-28 (RPM installation). How does one figure out which BLAS or
whether
it is threaded?

Many thanks and best wishes!

--- "Huntsinger, Reid" <reid_huntsinger at merck.com>
wrote:
> How was your R 2.0.1 built? Which Lapack did it link to, and which one
does
> it load? Which BLAS? Is the BLAS threaded? Does it link to the optimized
> pthreads library? Etc. (I'm not a Fedora Core 3 user so I'm not
sure what
> the default setup is, and of course I don't know if that's what you
have.)
> 
> As I pointed out in my previous post, in eigen() the matrix is
"inspected
> for symmetry" so symmetric=TRUE is the same as not specifying this at
all.
> You could try symmetric=FALSE...
> 
> Reid Huntsinger
> 
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch
> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Globe Trotter
> Sent: Monday, May 02, 2005 10:51 PM
> To: r-help at stat.math.ethz.ch
> Subject: Re: [R] eigenvalues of a circulant matrix
> 
> 
> OK, here we go:
> 
> I am submitting two attachments. The first is the datafile called kinv
used
> to
> create my circulant matrix, using the following commands:
> 
> 
> x<-scan("kinv")
> y<-x[c(109:1,0:108)]
> X=toeplitz(y)
> eigen(X)
> write(X,ncol=216,file="test.dat")
> 
> reports the following columns full of NaN's: 18, 58, 194, 200. (Note
that
> eigen(X,symmetric=T) makes no difference and I get the same as above).
> 
> The second attachment contains only the eigenvectors obtained on calling a
> LAPACK routine directly (from C). The eigenvalues are essentially the same
> as
> that obtained using R. Here, I use the LAPACK-recommended double precision
> routine dspevd() routine for symmetric matrices in packed storage format.
> Note
> the absence of the NaN's....I would be happy to send my C programs to
> whoever
> is interested.
> 
> I am using 
> 
> :~> uname -a
> Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 i686 i386
> GNU/Linux
> 
> and R.2.0.1.
> 
> Many thanks and best wishes!
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html
> 
> 
> 
> 
> 
>
----------------------------------------------------------------------------
--
> Notice:  This e-mail message, together with any attachment...{{dropped}}
______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html