complex-real FFT adjusment for MKL

FFTlib/cr1fft.c

@@ -113,8 +113,6 @@ void cr1fft(complex *cdata, REAL *rdata, int n, int sign)
         }
 
         Status = DftiSetValue(handle, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX);
-		//This options is what we would like, but is depreciated in the future
-        //Status = DftiSetValue(handle, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_REAL);
         if (! DftiErrorClass(Status, DFTI_NO_ERROR)) {
             dfti_status_print(Status);
             printf(" DftiSetValue FAIL\n");
@@ -133,9 +131,12 @@ void cr1fft(complex *cdata, REAL *rdata, int n, int sign)
         printf(" DftiComputeBackward FAIL\n");
     }
     rdata[0] = tmp[0];
-    for (i=1; i<n; i++) {
-        rdata[i] = -sign*tmp[n-i];
-    }
+	if (sign < 0) {
+    	for (i=1; i<n; i++) rdata[i] = -sign*tmp[n-i];
+	}
+	else {
+    	for (i=1; i<n; i++) rdata[i] = tmp[i];
+	}
 	free(tmp);
 #else
 	cr1_fft(cdata, rdata, n, sign);

FFTlib/crmfft.c

@@ -170,9 +170,12 @@ void crmfft(complex *cdata, REAL *rdata, int n1, int n2, int ldc, int ldr, int s
     for (j=0; j<n2; j++) {
     	Status = DftiComputeBackward(handle, (MKL_Complex8 *)&cdata[j*ldc], tmp);
     	rdata[j*ldr] = tmp[0];
-    	for (i=1; i<n1; i++) {
-        	rdata[j*ldr+i] = -sign*tmp[n1-i];
-    	}
+		if (sign < 0) {
+    		for (i=1; i<n1; i++) rdata[j*ldr+i] = -sign*tmp[n1-i];
+		}
+		else {
+    		for (i=1; i<n1; i++) rdata[j*ldr+i] = tmp[i];
+		}
 	}
     free(tmp);
     if(! DftiErrorClass(Status, DFTI_NO_ERROR)){

FFTlib/rc1fft.c

@@ -115,6 +115,13 @@ void rc1fft(REAL *rdata, complex *cdata, int n, int sign)
 			dfti_status_print(Status);
 			printf(" DftiSetValue FAIL\n");
 		}
+/*
+		Status = DftiSetValue(handle, DFTI_FORWARD_DOMAIN, DFTI_REAL);
+		if(! DftiErrorClass(Status, DFTI_NO_ERROR)){
+			dfti_status_print(Status);
+			printf(" DftiSetValue FAIL\n");
+		}
+*/
 
 		Status = DftiSetValue(handle, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX);
 		if (! DftiErrorClass(Status, DFTI_NO_ERROR)) {

FFTlib/test/rc1test.c

@@ -24,14 +24,14 @@ int main () {
 
 	N = 16;
 	k = 5.0;
-	sign = 1;
-	isign = -1;
+	sign = -1;
+	isign = 1;
 	while (N <= Nmax) {
 
 		/* Initialize the data */
 
 		for (i=0;i<N;i++) {
-			c_data[i]   = (float)-0.1+0.5*(N/2-i);
+			c_data[i]   = (float)-0.1+0.5*(N/3-i)*sin(i*M_PI/N);
 //			c_data[i]   = 0.0;
 		}
 //		c_data[0]   = 1.0;