Using fftw with basic square matrix columns (armadillo library)

Question

I find the Armadillo C ++ library very handy for matrix computing. How to perform two-dimensional FFT on armadillo matrices using the FFTW library?

I understand that the armadillo matrix class stores data in the main column order. How to transfer this to FFTW? The fftw 3.3.3 documentation says:

If you have an array stored in column order and you want to convert it using FFTW, this is pretty easy to do. When creating a plan, simply pass the dimensions of the array to the scheduler in the reverse order. For example, if your array is a three-dimensional matrix N x M x L ​​in the main column order, you should pass the dimensions of the array as if it were a matrix L x M x N (which, from the point of view of FFTW)

I cannot fully understand what this means, given the syntax for creating a plan as follows.

fftw_plan fftw_plan_dft_2d(int n0, int n1,
                            fftw_complex *in, fftw_complex *out,
                            int sign, unsigned flags);

Can anyone explain this?

Answer 1

It just means you can try both

  fftw_plan plan=fftw_plan_dft_2d(4, 2,(double(*)[2])&AAA(0,0), (double(*)[2])&BBB(0,0), FFTW_FORWARD, FFTW_ESTIMATE);

and

 fftw_plan plan=fftw_plan_dft_2d(2, 4,(double(*)[2])&AAA(0,0), (double(*)[2])&BBB(0,0), FFTW_FORWARD, FFTW_ESTIMATE);

And keep the correct order.

Typically, inand outstand fftw_malloc()for the preservation aligned memory. Small tests show that this is already the case with the matrix cx_matfrom Armadillo (without a terrible segmentation error or incorrect values ​​...).

Here is the test code (and the correct order ...):

    #include <iostream>
#include <fftw3.h>
#include "armadillo"

using namespace arma;
using namespace std;


int main(int argc, char** argv)
{
    cout << "Armadillo version: " << arma_version::as_string() << endl;


    cx_mat AAA = eye<cx_mat>(2,4);
    AAA(0,0)=0;
    AAA(0,1)=1;
    AAA(0,2)=2;
    AAA(0,3)=3;

    AAA(1,0)=0;
    AAA(1,1)=1;
    AAA(1,2)=2;
    AAA(1,3)=3;


    cx_mat BBB = eye<cx_mat>(2,4);

    fftw_plan plan=fftw_plan_dft_2d(4, 2,(double(*)[2])&AAA(0,0), (double(*)[2])&BBB(0,0), FFTW_FORWARD, FFTW_ESTIMATE);

    fftw_execute(plan);

    BBB.print("BBB:");
    return 0;
}

Compile it with g++ -O2 -o example1 example1.cpp -larmadillo -llapack -lblas -lfftw3!

Till,

Francis