Error: Using EigenFFT in for loop

[JamieAl]

I am testing EigenFFT with some Eigen tensors and matrices. I have been so far successful with implementing 1D FFTs along different directions and comparing my results to MATLAB. So, for example, my MATLAB code looks like:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

clearvars; clc; close all; %3D FFT test Nx = 8; Ny = 4;%8; Nz= 6;%8; Lx =16;%2*pi;% Ly = 6;% 2*pi;% dx = Lx/Nx; % Need to calculate dx dy = Ly/Ny;

xi_x = (2*pi)/Lx; yi_y = (2*pi)/Ly; xi = ((0:Nx-1)/Nx)*(2*pi); yi = ((0:Ny-1)/Ny)*(2*pi); x = xi/xi_x; y = yi/yi_y; zlow = 0; %a zupp =6; %b Lz = (zupp-zlow); zgl = cos ( pi * ( 0 : Nz ) / Nz )'; zgl = (1/2)*(((zupp-zlow)*zgl) + (zupp+zlow)); [X,Z,Y] = meshgrid(x,zgl,y); %this gives 3d grid with z-by-x-by-y size %ICs A = 2*pi / Lx; B = 2*pi / Ly; uFun = (Z-zlow) .* (Z-zupp) .* sin(A*X).* sin(B*Y); %Take 1D FFTs along each direction uh1 =(fft(u,[],1));%each columns numRows=Nz+1;numCols=Nx; for j= 1:numCols uk1(:,j,:) = fft(u(:,j,:)); end uh2 =(fft(u,[],2));%each row for j= 1:numRows uk2(j,:,:) = fft(u(j,:,:)); end uh3 =(fft(u,[],3));%over 3rd dimension (i.e. Y) for i= 1:numRows for j= 1:numCols uk(i,j,:) = fft(u(i,j,:)); end end

Now, in C/C++ I am able to reproduce both uh1 =(fft(u,[],1)); and uh2 =(fft(u,[],2)); really the third dimension that I am struggling with the syntax of Eigen and getting errors.

My C/C++ code looks like:

static const int nx = 8; 
static const int ny = 4;
static const int nz = 6;

int main(){
Eigen::Tensor<double, 3> uFun((nz+1),nx,ny);//tensor(row,col,matrix)
uFun.setZero();
for(int x = 0; x< nx; x++){
		    for(int y = 0; y< ny; y++){
                for(int z = 0; z< nz+1; z++){
                uFun(z,x,y) = //define function u same as MATLAB
            }
        }
    }
//define matrix and copy the tensor into it
//map tensor to matrix

Eigen::MatrixXd dummy((nz+1),nx); //to match tensor slices
dummy.setZero();
Eigen::Tensor<double, 2> Tensor2D((nz+1),nx);   //tensor(row,col,matrix)
Tensor2D.setZero();

Tensor2D = uFun.slice(std::array{0, 0, 3}, std::array{nz + 1, nx, 1}).reshape(std::array{nz + 1, nx}) ;
dummy = static_cast<MatrixType<double>>(Eigen::Map<const MatrixType<double>>(Tensor2D.data(), (nz+1), nx));
Eigen::MatrixXcd tempOut1((nz+1),nx);
//uh1= fft(u,[],1);
	for (int k=0; k< dummy.cols(); k++){
		tempOut1.col(k) = fft.fwd(dummy.col(k));		
	}
//uh2= fft(u,[],2);
Eigen::MatrixXcd tempOut((nz+1),nx);
Eigen::VectorXcd row(nx);
	for (int k=0; k< (dummy.rows()); k++){
		fft.fwd(row, dummy.row(k));
		tempOut.row(k) = row;
	}

}

The above code works, but now when trying to reproduce this part from MATLAB code:

for i= 1:numRows
   for j= 1:numCols
      uk(i,j,:) = fft(u(i,j,:));
   end
end

I get the error:

36)  (14.790,-7.123)  (10.458,-5.036)
test.out: /mnt/c/Users/Documents/eigen-3.4.0/eigen-3.4.0/Eigen/src/Core/Block.h:120: Eigen::Block<XprType, BlockRows, BlockCols, InnerPanel>::Block(XprType&, Eigen::Index) [with XprType = Eigen::Matrix<double, -1, -1>; int BlockRows = 1; int BlockCols = -1; bool InnerPanel = false; Eigen::Index = long int]: Assertion `(i>=0) && ( ((BlockRows==1) && (BlockCols==XprType::ColsAtCompileTime) && i<xpr.rows()) ||((BlockRows==XprType::RowsAtCompileTime) && (BlockCols==1) && i<xpr.cols()))' failed.
Aborted (core dumped) 

My attempt in C/C++ looks like:

for (int i = 0; i < dummy.rows(); ++i) {
        for (int j = 0; j < dummy.cols(); ++j) {
            int index =  i * (dummy.cols()) + j; // Flattened index
	//int index =  i * (dummy.rows()) + j; // Flattened index
            VectorXcd slice = dummy.row(index);
            tempOut.row(index) = fft.fwd(slice);

        }
    }

std::cout<<tempOut(index)<<endl;

How can I loop over # of rows and columns for each matrix similarly to the MATLAB code?