[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: complex arrays in octfiles
From: |
John W. Eaton |
Subject: |
Re: complex arrays in octfiles |
Date: |
Wed, 29 Jun 2011 10:03:50 -0400 |
On 29-Jun-2011, c. wrote:
| Here is an example of how to build A ComplexNDArray:
|
| #include <octave/oct.h>
|
| DEFUN_DLD (pippo, args, , "test return statement")
| {
|
| dim_vector dims;
| dims.resize (3);
| dims(2) = 4;
| dims(1) = 3;
| dims(0) = 2;
|
|
| ComplexNDArray a (dims, std::complex<double> (0.0, 0.0));
| // or
| //Array<std::complex<double> > a (dims, std::complex<double> (0.0, 0.0));
| for (octave_idx_type i = 0; i < dims(0); i++)
| for (octave_idx_type j = 0; j < dims(1); j++)
| for (octave_idx_type k = 0; k < dims(2); k++)
| a (i, j, k) = std::complex<double> (i, j + k);
|
| return octave_value (a);
|
| }
This should work, but using the indexing operator for assignment is
relatively slow because of Octave's copy-on-write semantics. So if
you want better performance, you can write it this way:
ComplexNDArray a (dims, Complex (0.0, 0.0));
Complex *cptr = a.fortran_vec ();
for (octave_idx_type k = 0; k < dims(2); k++)
for (octave_idx_type j = 0; j < dims(1); j++)
for (octave_idx_type i = 0; i < dims(0); i++)
*cptr++ = Complex (i, j + k);
Note that the data is stored in column-major order, so I reversed the
loop order.
Complex is a typedef for std::complex<double> in the Octave sources.
jwe