How to efficiently extract real/imaginary parts of complex matrix in Eigen3 library?

Question

I have some complex, dense vectors/matrices in the Eigen3 library, and I want to extract the real and imaginary parts into separate arrays. In Matlab, I could do something like

cplxFoo = [1, 1i; -1i -1]
re = real(cplxFoo)
im = imag(cplxFoo)

which expectantly yields

cplxFoo =
   1.0000 + 0.0000i   0.0000 + 1.0000i
   0.0000 - 1.0000i  -1.0000 + 0.0000i
re =
     1     0
     0    -1
im =
     0     1
    -1     0

Is there anything like the real() and imag() Matlab functions in Eigen3?

Right now, the only thing I know will work is something akin to

MatrixXcd cplxFoo = ...;
MatrixXd re(cplxFoo.rows(), cplxFoo.cols());
MatrixXd im(cplxFoo.rows(), cplxFoo.cols());

for(size_t j=0; j<cplxFoo.cols(); ++j) {
    for(size_t i=0; i<cplxFoo.rows(); ++i) {
        re(i, j) = cplxFoo(i,j).real();
        im(i, j) = cplxFoo(i,j).imag();
    }
}

It works, and I can put it in a function even, but then I'm stuck having to do my own loop vectorizing, unrolling, etc., and I have to make an extra copy.

What I would like to be able to do is wrap a couple of Map<MatrixXd> with appropriate strides around cplxFoo to get the real and imaginary parts. But the problem is that the elements of MatrixXcd are std::complex<double>, and I'm not sure what the layout of that is. My guess is that std::complex<T> is essentially laid out like struct {T real; T imag;}; so that real and imaginary parts are tightly packed and interleaved when you make an array of std::complex<T> (and that also seems to be the consensus at this SO question), but is that guaranteed by the C++ standard? AFAICT, a compliant C++ compiler could lay it out like struct {T imag; T real;}; (note the changed order), or something more exotic like

class {
    T radius;
    T angle;

public:
    T real() const { return radius * cos(angle); }
    T imag() const { return radius * sin(angle); }
    /* ... */
};

So, is it ok to wrap a couple of Map<MatrixXd> with appropriate strides around cplxFoo? If so, how do I set up the strides correctly?

Alternatively, is there any way to get Eigen's complex data types to use separate chunks of memory for the real and imaginary parts?

For what it's worth, the reason I need to do this is because I need to interface the Eigen library with MATLAB which can only handle separate arrays for real and imaginary parts, not interleaved in any way.

The layout of `std::complex` is guaranteed to be compatible with an array of two `T` values (real followed by imaginary parts), and for arrays to generalise as you'd expect. (But I don't know enough about Eigen to answer the other questions). — Mike Seymour, Mar 11 '14 at 10:53
If I don't get the docs wrong `cplxFoo.imag()` and `cplxFoo.real()` should work. Have you tried that? — anderas, Mar 11 '14 at 10:55
@anderas Hot damn, it works! I wish I would have seen it in the docs, could have saved me a ton of trouble. — Nicu Stiurca, Mar 11 '14 at 11:37
@anderas The combination of your comment combined with Mike's fully answers my question. — Nicu Stiurca, Mar 12 '14 at 04:44

score 15 · Accepted Answer · edited Feb 25 '20 at 13:03

15

That's easy, just use the .real() and .imag() views:

MatrixXcd M;
MatrixXd r, i;
r = M.real();
i = M.imag();

Note that you can use M.real() into an expression without copying it into a MatrixXd.

edited Feb 25 '20 at 13:03

chtz

17,329
4
26
56

answered Mar 12 '14 at 07:57

ggael

28,425
2
65
71

Is there a way to do the exact same thing but for an Eigen tensor? – Jamie Oct 07 '22 at 23:01

How to efficiently extract real/imaginary parts of complex matrix in Eigen3 library?

1 Answers1