Matlab fft2とMKL DftiComputeForwardを使用した結果が異なる

Question

2次元データのForward Fourie Transformを計算しようとすると、結果が異なります。 Matlabの、リスティング：

fft2([25.6798, 26.0815, 29.0069; 33.5761 37.123 38.4696; 38.6358 38.0078 37.649]) 

Matlabの結果：以下の3x3の行列のため

簡易試験例

ans = 
    1.0e+02 * 
    3.0423 + 0.0000i -0.0528 + 0.0339i -0.0528 - 0.0339i 
    -0.3096 + 0.0444i 0.0112 + 0.0646i -0.0144 + 0.0225i 
    -0.3096 - 0.0444i -0.0144 - 0.0225i 0.0112 - 0.0646i 

MKL、リスト：

DFTI_DESCRIPTOR_HANDLE descriptor1; 
double test[3][3] = {{25.6798, 26.0815, 29.0069}, 
         {33.5761, 37.123, 38.4696}, 
         {38.6358, 38.0078, 37.649}}; 
MKL_LONG status1, l1[2]; l1[0] = 3; l1[1] = 3; 
MKL_Complex16 fftu1[3][3]; 

status1 = DftiCreateDescriptor(&descriptor1, DFTI_DOUBLE, DFTI_REAL, 2, l1); 
status = DftiCommitDescriptor(descriptor1); 
status = DftiComputeForward(descriptor1, test, fftu1); 

MKL結果：

4.02248e-315+2.35325e-310i 6.42285e-323+6.95254e-310i 2.35325e-310+2.35325e-310i 
6.95254e-310+6.95254e-310i 2.35308e-310+2.35325e-310i 0+2.35325e-310i 
2.35325e-310+2.35325e-310i 2.35325e-310+2.35325e-310i 7.41098e-323+1.03754e-322i 

この問題は、MKLの場合の出力ストレージ構成であるディスクリプタによって発生する可能性があることがわかりました。しかし、私はこの記述子を設定する正しい方法を見つけることができません。

私は間違っていますか？私にいくつかのヒントを与えてください。

Answer 1

最後に私はそれを持っています。解決策が誰かに役立つかもしれません。

正確なC++ MKLのリストは、明確化のコメントと共に以下に提供されています。問題のコードに加えて、

#include <complex> 
#define MKL_Complex16 std::complex<double> //This definition should be done before including MKL files. You will need it to use STL functions, for example conj, with MKL_Complex16 data 
#include "mkl.h" 

第二に、DFTI_NOT_INPLACEとDFTI_STORAGE、この場合のために定義する必要があります： まず、次の定義を含む区間で使用されるべきである

DFTI_DESCRIPTOR_HANDLE descriptor1; 
double test[3][3] = {{25.6798, 26.0815, 29.0069}, 
         {33.5761, 37.123, 38.4696}, 
         {38.6358, 38.0078, 37.649}}; 
MKL_LONG status1, l1[2]; l1[0] = 3; l1[1] = 3; 
MKL_Complex16 fftu1[3][3]; 
status1 = DftiCreateDescriptor(&descriptor1, DFTI_DOUBLE, DFTI_REAL, 2, l1); 
status = DftiSetValue(descriptor1, DFTI_PLACEMENT, DFTI_NOT_INPLACE); 
status = DftiSetValue(descriptor1, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX); 
MKL_LONG rs[2]; rs[0] = 0; rs[1] = 3; rs[2] = 1; //describing 2D 3x3 matrix 
status = DftiSetValue(descriptor1, DFTI_INPUT_STRIDES, rs); 
MKL_LONG cs[2]; cs[0] = 0; cs[1] = 3; cs[2] = 1; /*Describing output matrix. Warning! Only N1x(N2/2)+1 half part will contain correct answer. Rest part should be restored!!! According to the manual, cs[1]=N2/2+1, so cs[1] should be 2 in our case. But if cs[1]=2, this leads to results shift in answer.. I hope, this is (making cs[1]=N2} is a correct way to deal with shift, bu i'm not sure there. Need to be checked*/ 
status = DftiSetValue(descriptor1, DFTI_OUTPUT_STRIDES, cs); 
status = DftiCommitDescriptor(descriptor1); 
status = DftiComputeForward(descriptor1, test, fftu1); // Forward Fourie done 
/*Now the complex-valued data sequences in the conjugate-even domain can be reconstructed as described in https://software.intel.com/en-us/mkl-developer-reference-c-dfti-packed-format*/ 
for (size_t ii=0; ii< 3; ii++){ 
    for (size_t jj=3/2+1; jj< 3; jj++){ 
    fftu1[ii][jj] = std::conj((MKL_Complex16)fftu1 [(3-ii)%3] [(3-jj)%3]); 
    } 
} 

今fftu1の結果は、Matlab計算の場合と同じです（小数点第4位まで）。