[Rcpp-devel] Is declaring an argument as a reference the same as declaring it a pointer?

Whit Armstrong armstrong.whit at gmail.com
Sun Mar 28 00:56:53 CET 2010


All the newer c++ matrix template libraries are very good.  I like
Armadillo in particular with MTL as a close second.

As much of a pain as it is to learn a new library, you code will be
infinitely more readable if you use a matrix library.

here are the operator examples from MTL4 and Armadillo for reference:

http://osl.iu.edu/research/mtl/mtl4/doc/matrix__vector__expr.html
http://arma.sourceforge.net/docs.html#operators

My secret desire is to implement R using Armadillo as the core data
structures...

just my 2c.

-Whit


On Sat, Mar 27, 2010 at 7:22 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
> On Sat, Mar 27, 2010 at 6:12 PM, Dirk Eddelbuettel <edd at debian.org> wrote:
>>
>> That's a long with a number of questions, I'll just cherry-pick:
>>
>> On 27 March 2010 at 17:10, Douglas Bates wrote:
>> | I am writing code to interface from Rcpp to the Lapack linear algebra
>> | code which is written in Fortran.
>> |
>> | Because Lapack is written in Fortran all arguments must be passed from
>> | C as pointers.  There are many examples of the C versions of the
>> | declarations in the file R_HOME/include/R_ext/Lapack.h where you will
>> | see that a common idiom is to pass a pointer to the contents of a
>> | matrix, a pointer to an integer giving the number of columns, the same
>> | for the number of rows, and a third pointer which is the "leading
>> | dimension of the array as declared in the calling program".  (If all
>> | this seems bizarre, be grateful that you never needed to learn to
>> | program in Fortran.)
>> |
>> | For example, the C declaration of dgemm, which computes C := alpha * A
>> | %*% B + beta * C, is
>>
>> (Everybody's favourite Blas-3 function :)
>>
>> | /* DGEMM - perform one of the matrix-matrix operations    */
>> | /* C := alpha*op( A )*op( B ) + beta*C */
>> | void
>> | F77_NAME(dgemm)(const char *transa, const char *transb, const int *m,
>> |               const int *n, const int *k, const double *alpha,
>> |               const double *a, const int *lda,
>> |               const double *b, const int *ldb,
>> |               const double *beta, double *c, const int *ldc);
>> |
>> | The arguments m, n, k, alpha, lda, ldb, and beta are all int's or
>> | double's that must be passed as pointers.  If A, B and C are
>> | NumericMatrix objects then I need to write code like
>> |
>> |      const char *transa = "N", transb = "N";
>> |      const int m = A.nrow(), n = A.ncol(), ...
>> |      const double one = 1.0, zero = 0.0;
>> |
>> |      F77_CALL(dgemm)(transa, transb, &m, &n, &one, A.begin(), &m, ...)
>> |
>> | which is somewhat roundabout because I need to declare int's and
>> | double's and initialize them to fixed values so I can pass those fixed
>> | values as pointers for the Fortran routine.
>> |
>> | Can I declare a C++ interface replacing the const int *m, const double
>> | *alpha, etc. by
>> | const int &m, const double &alpha, etc. so that I can call the Fortran
>> | subroutine as
>> |
>> |     F77_CALL(dgemm)(transa, transb, A.nrow(), A.ncol(), 1.0,
>> | A.begin(), A.nrow(), ...
>>
>> I played with something related the other day when I (re-)worked the two MPI
>> examples contributed to RInside. I didn't end up getting this sorted out. It
>> seems one really needs the "C approach" here with 'const int m = N.nrow()'
>> and passing &m.
>>
>> | In other words, do the formal argument declarations const int &m and
>> | const int *m both end up passing an address of an int, with the only
>> | difference being in what the form of the actual argument is?
>> |
>> | If so, is there any way to avoid the dance with the character strings?
>> |  It turns out that only the first character is ever examined in the
>> | Lapack code by I still need to pass an address to that character, not
>> | the character itself.
>>
>> I fear it may be the same issue, so you may need the char[].
>>
>> Would the C++ interface to Lapack offer help?  Otherwise, you could play
>> "just don't do it" and use a modern shim like Armadillo, Eigen, ... or even
>> uBlas from Boost which all end up calling Lapack functions for you.  Or you
>> sweat the distateful bits, write'em once and hide them behind other functions
>> you call.
>
> I think a C++ interface would help in that it makes it easier to read.
>  Because everything in Fortran is passed as an address there is no
> real distinction between scalars and vectors (or pointers in C).  I
> wrote the R_HOME/include/R_ext/Lapack.h file and went to some trouble
> to differentiate between a declaration of int *m and int iwork[].  The
> first would be a scalar (passed as a pointer) and the second would be
> a vector, even though the two declarations are equivalent in C.
> _______________________________________________
> 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
>


More information about the Rcpp-devel mailing list