-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 I'm trying to understand some C code in an R package I'm using. I'm address this question here as it's matrix algebra...and I'm no pro at that! the C command reads: double alpha = 1.0, beta = 0.0; dsyrk_("L", "N", nGenes, nGenes, & alpha, mat1, nGenes, & beta, mat2, nGenes); - From google, I've found out that dsyrk is for performing one of the symmetric rank k operations - whatever that means!? From here: http://linux.die.net/man/l/dsyrk I've found that the calculation being performed is: alpha*A*A' + beta*C However, since alpha is 1 and beta is 0, this reduces to: => 1*A*A' + 0*C => A*A' Which is simply the cross product....am I correct? Cheers, Nath - -- - -------------------------------------------------------- Dr. Nathan S. Watson-Haigh OCE Post Doctoral Fellow CSIRO Livestock Industries Queensland Bioscience Precinct St Lucia, QLD 4067 Australia Tel: +61 (0)7 3214 2922 Fax: +61 (0)7 3214 2900 Web: http://www.csiro.au/people/Nathan.Watson-Haigh.html - -------------------------------------------------------- -----BEGIN PGP SIGNATURE----- Version: GnuPG v1.4.9 (MingW32) Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org iEYEARECAAYFAklkGVYACgkQ9gTv6QYzVL73HgCgvx4OCxcuczv8nd0n6gOEPFYa w3UAnAnDIkvPDen9p7ahz+BdG47V/D/S =gSGC -----END PGP SIGNATURE-----
On Wed, 7 Jan 2009, Nathan S. Watson-Haigh wrote:> -----BEGIN PGP SIGNED MESSAGE----- > Hash: SHA1 > > I'm trying to understand some C code in an R package I'm using. I'm address this question here as > it's matrix algebra...and I'm no pro at that! > > the C command reads: > > double alpha = 1.0, beta = 0.0; > dsyrk_("L", "N", nGenes, nGenes, & alpha, mat1, nGenes, > & beta, mat2, nGenes);That's a Fortran subroutine. The code is in blas.f (as 'DSYRK') and is commented. So, you should be able to work thru it. If you want to see how it is used, grep the sources for 'dsyrk' (lower case). I think array.c uses it for symmetric crossproducts along the lines of the alpha=1.0 and beta = 0.0 example you cite. It is often a good idea to search the R sources for routines whose function is puzzling and find examples of their use to help one understand what they do. HTH, Chuck> > - From google, I've found out that dsyrk is for performing one of the symmetric rank k operations - > whatever that means!? From here: > http://linux.die.net/man/l/dsyrk > > I've found that the calculation being performed is: > alpha*A*A' + beta*C > > However, since alpha is 1 and beta is 0, this reduces to: > => 1*A*A' + 0*C > => A*A' > > Which is simply the cross product....am I correct? > > Cheers, > Nath > > - -- > - -------------------------------------------------------- > Dr. Nathan S. Watson-Haigh > OCE Post Doctoral Fellow > CSIRO Livestock Industries > Queensland Bioscience Precinct > St Lucia, QLD 4067 > Australia > > Tel: +61 (0)7 3214 2922 > Fax: +61 (0)7 3214 2900 > Web: http://www.csiro.au/people/Nathan.Watson-Haigh.html > - -------------------------------------------------------- > > -----BEGIN PGP SIGNATURE----- > Version: GnuPG v1.4.9 (MingW32) > Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org > > iEYEARECAAYFAklkGVYACgkQ9gTv6QYzVL73HgCgvx4OCxcuczv8nd0n6gOEPFYa > w3UAnAnDIkvPDen9p7ahz+BdG47V/D/S > =gSGC > -----END PGP SIGNATURE----- > > ______________________________________________ > R-help at r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. >Charles C. Berry (858) 534-2098 Dept of Family/Preventive Medicine E mailto:cberry at tajo.ucsd.edu UC San Diego http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
On Tue, Jan 6, 2009 at 8:54 PM, Nathan S. Watson-Haigh <nathan.watson-haigh at csiro.au> wrote:> -----BEGIN PGP SIGNED MESSAGE----- > Hash: SHA1 > > I'm trying to understand some C code in an R package I'm using. I'm address this question here as > it's matrix algebra...and I'm no pro at that! > > the C command reads: > > double alpha = 1.0, beta = 0.0; > dsyrk_("L", "N", nGenes, nGenes, & alpha, mat1, nGenes, > & beta, mat2, nGenes); > > - From google, I've found out that dsyrk is for performing one of the symmetric rank k operations - > whatever that means!? From here: > http://linux.die.net/man/l/dsyrk > > I've found that the calculation being performed is: > alpha*A*A' + beta*C > > However, since alpha is 1 and beta is 0, this reduces to: > => 1*A*A' + 0*C > => A*A' > > Which is simply the cross product....am I correct?Not quite. The crossprod function in R computes A'A. This particular operation is what is computed by the tcrossprod function in R. Another difference is that this call modifies only the diagonal and lower triangle of mat2. The tcrossprod function calls the same underlying code then copies the lower triangle to the upper triangle so as to produce a symmetric matrix. When you are calling Lapack or BLAS routines from C it is helpful to look at the declarations in include/R_ext/BLAS.h and include/R_ext/Lapack.h /* DSYRK - perform one of the symmetric rank k operations */ /* C := alpha*A*A' + beta*C or C := alpha*A'*A + beta*C */ BLAS_extern void F77_NAME(dsyrk)(const char *uplo, const char *trans, const int *n, const int *k, const double *alpha, const double *a, const int *lda, const double *beta, double *c, const int *ldc); This helps to explain the mysterious "L" and "N" as the first two arguments in the call. By the way, the use of the literal name dsyrk_ in that code is risky. It is preferable to use the macro F77_SUB(dsyrk) which is defined in the include file include/R_ext/RS.h If you are wondering what the difference between F77_SUB and F77_NAME might be, there isn't a difference in R. Long, long ago there was a difference between the two macros in S and S-PLUS related to a peculiar Fortran compiler, if I recall correctly. Most of the time these macros simply append an underscore to the name but there may still be cases where something else occurs and it is best to allow for those cases.