2

R uses, when interfacing with other languages, the header R_ext/Complex.h which includes the type Rcomplex which seems to be an implementation of std::complex<double>. The standard way of using it would be for a complex vector x_ in R: Rcomplex *px = COMPLEX(x_); However, since I need to pass it to armadillo I then do: arma::cx_vec x(px, nrows(x_), false, false);

but armadillo does not accept Rcomplex types. I have tried doing instead: std::complex<double> *px = COMPLEX(x_); but get the following error: cannot convert ‘Rcomplex*’ to ‘std::complex<double>*’ in initialization

Do you have any clue for passing a complex vector in R to std::complex<double> type? I am aware of Rcpp but would like to have a direct solution relying on base R.

EDIT: Following one of the comments, I wanted to clarify that Rcomplex is of C type but that it is compatible with std::complex<double> according to the answer by @Stephen Canon.

EDIT2: Why does Dirk's answer have more votes than the accepted answer if it is not answering the "without dependencies" question. In addition, I have been downvoted apparently because if one wants preferably to use base R with C or C++ somebody does not like it. Anyway, I have better things to do but this is not the first time that I get no answer to my original question when asking something related to base R interfacing with C or C++ and get a Rcpp related answer that I have not asked for.

Community
  • 1
  • 1
nopeva
  • 1,583
  • 5
  • 22
  • 38
  • Is there a specific reason why you don't want to use RcppArmadillo? – Roland Sep 29 '16 at 10:27
  • @Roland For simple things I usually prefer avoiding dependencies and my function is extremely easy but I want to speed it up. – nopeva Sep 29 '16 at 11:20
  • @nopeva: That's called reinventing the wheel. By the same logic you wouldn't use Armadillo. Or R. – Dirk Eddelbuettel Sep 29 '16 at 13:29
  • Your premise is also wrong: _which seems to be an implementation of std::complex_. Certainly not, as R is implemented only in C and not in C++. – Dirk Eddelbuettel Sep 29 '16 at 13:35
  • Dirk Eddelbuettel@ You are right in the sense that `Rcomplex` is of type `double complex`. However, the `std::complex`type has been developed with care to make them compatible. See [link](http://stackoverflow.com/questions/10540228/c-complex-numbers-in-c) – nopeva Sep 29 '16 at 13:58
  • May I encourage you to add some of the justification and explanation given in the comments to your question (e.g., wanting to avoid dependencies) *and* remove the complaint in your question? Being clearer why you want to avoid dependencies is good; if the SO community thinks that's a bad idea and upvotes Dirk's answer anyway as a better way to solve the underlying problem, that (though it may annoy you) isn't *your* problem. Upvote and accept your preferred answer and move on ... – Ben Bolker Sep 30 '16 at 18:59

2 Answers2

4

Complex numbers are less common in statistics so that has not been an initial focus. However there are use cases Baptiste has one or two packages which pushes to add features to the interface given the existing support in Armadillo and R.

So all the work is done for you in templates -- here is the simplest possibly example of passing complex-values matrix and returning its sum:

#include <RcppArmadillo.h>

// [[Rcpp::depends(RcppArmadillo)]]

// [[Rcpp::export]]
arma::cx_mat foo(const arma::cx_mat & C) {
    return C + C;
}

/*** R
C <- matrix( (1+1i) * 1:4, 2, 2)
C
foo(C)
*/

It does what you'd expect:

R> sourceCpp("/tmp/armaComplex.cpp")

R> C <- matrix( (1+1i) * 1:4, 2, 2)

R> C
     [,1] [,2]
[1,] 1+1i 3+3i
[2,] 2+2i 4+4i

R> foo(C)
     [,1] [,2]
[1,] 2+2i 6+6i
[2,] 4+4i 8+8i
R> 

So we start with complex values in R, pass them via Rcpp and RcppArmadillo to Armadillo and get them back to R. Without writing an additional line of code, and no discernible overhead.

Dirk Eddelbuettel
  • 360,940
  • 56
  • 644
  • 725
2

One type can always be forced to another type using reinterpret_cast. Generally this is a bad idea, but if you can guarantee that the two complex types are indeed compatible, you can do:

 Rcomplex* px = COMPLEX(x_);
 arma::cx_vec x( reinterpret_cast<arma::cx_double*>(px), nrows(x_), false, false );

The type arma::cx_double is a shorthand for std::complex<double>

mtall
  • 3,574
  • 15
  • 23