[Rcpp-devel] symmatu does not give symmetric matrix

[Douglas Bates]

First, this is not a good way of constructing a matrix to test the Cholesky
decomposition because the result is rank-deficient. It is a 2 by 2 matrix
of rank 1.  In fact it is the very definition of a rank-1 matrix - a column
vector multiplied by its transpose.

I am surprised at the result because I can't think of a way in computer
arithmetic that the off-diagonals would be different.  They should both be

prod(x)

The dput function will print the values out to full precision so try

dput(prod(x))
dput(x %*% t(x))

to see the full precision results.

Thirdly, R has two functions, crossprod and tcrossprod, to compute X'X and
XX', respectively.  They compute one triangle of the result and copy it
into the other, so the result is guaranteed to be symmetric. Try

dput(tcrossprod(x))

Finally, what are you going to do with the Cholesky decomposition?  If you
are doing this simply to get the decomposition, compiling C++ code, even
with all the help provided by Rcpp, is doing things the hard way. Just use
the chol function in R or the cholesky decomposition in the Matrix package
if you want more control.

On Wed, Jul 27, 2016 at 7:40 AM Dirk Eddelbuettel <edd at debian.org> wrote:

>
> Sorry for the delay in responding.  Mail delivery was stuck at R-Forge for
> a
> few days it seems.
>
> On 21 July 2016 at 10:46, Jendoubi, Takoua wrote:
> | Dear all,
> |
> | I am using RcppArmadillo to deal with some matrix computations.
> Specifically I
> | need to cholesky factorization of some symmetric matrices.
> |
> | I am generating random vectors using Rcpp and using them to construct
> symmetric
> | matrices.
> |
> | I have an error stating that my matrix is not symmetric although it
> definitely
> | should be. Here is the example I am working with (in R):
> |
> | >x
> | [1] -1.6683320 -0.8597148
> | >x%*%t(x)
> |          [,1]      [,2]
> | [1,] 2.783332 1.4342896
> | [2,] 1.434290 0.7391095
> |
> | Apparently, it is a rounding-off error. Is there any way to ensure that
> x%*%t
> | (x) gives an exactly symmetric matrix to use for cholesky factorization?
>
> I don't know of an automatic way.  R surely has no native type for
> symmetric
> matrices (at least not at the SEXP level as it really only has vectors with
> dim attributes).  I would probably start with some helper functions at the
> R
> or C++ level.  Some may have better tricks...
>
> Dirk
>
>
> --
> http://dirk.eddelbuettel.com | @eddelbuettel | edd at debian.org
> _______________________________________________
> Rcpp-devel mailing list
> Rcpp-devel at lists.r-forge.r-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20160727/443e97cc/attachment.html>