Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
Announcements
Welcome to the Intel Community. If you get an answer you like, please mark it as an Accepted Solution to help others. Thank you!
6433 Discussions

strange values for the FFT of an odd function

vincent-torri
Beginner
152 Views
I want to compute the Fourier coefficients of the 2-periodic function that is equal to x |---> x on [-1,1]. When I display them, they have the correct imaginary value, but a real value which is not equal to 0 (like it should be), but which is equal to -1.

Here is the code:

[cpp]#include 

#include 
#include 

float fct_der(float x)
{
  return x;
}


int main()
{
  float _Complex x_in[32];
  float _Complex x_out[32];
  DFTI_DESCRIPTOR_HANDLE my_desc1_handle;
  MKL_LONG status;
  int i;

  for (i = 0; i < 32; i++)
    {
      float x;

      x = -1.0f + (float)(2 * i) / 32.0f;
      x_in = fct_der(x);
    }

  status = DftiCreateDescriptor( &my_desc1_handle, DFTI_SINGLE, DFTI_COMPLEX, 1, 32);
  status = DftiSetValue( my_desc1_handle, DFTI_PLACEMENT, DFTI_NOT_INPLACE);
  status = DftiCommitDescriptor( my_desc1_handle );
  status = DftiComputeForward( my_desc1_handle, x_in, x_out);
  status = DftiFreeDescriptor(&my_desc1_handle);

  for (i = 0; i < 32; i++)
    {
      printf ("%f + i %f\n", creal(x_out), cimag(x_out));
    }

  return 0;
}[/cpp]

And here is the output:

[plain]-1.000000 + i 0.000000
-1.000000 + i 10.153171
-1.000000 + i 5.027339
-1.000000 + i 3.296558
-1.000000 + i 2.414214
-1.000000 + i 1.870868
-1.000000 + i 1.496606
-1.000000 + i 1.218504
-1.000000 + i 1.000000
-1.000000 + i 0.820679
-1.000000 + i 0.668179
-1.000000 + i 0.534511
-1.000000 + i 0.414214
-1.000000 + i 0.303347
-1.000000 + i 0.198912
-1.000000 + i 0.098491
-1.000000 + i 0.000000
-1.000000 + i -0.098491
-1.000000 + i -0.198912
-1.000000 + i -0.303347
etc...[/plain]

Does someone see why the real part is not 0 ?

thank you
0 Kudos
11 Replies
Todd_R_Intel
Employee
152 Views
I'm trying to read your program and understand what you're trying to do. You mention a real function and in your program you initialize your variable of type float with real data, but then you assign it to your float _Complex data variable and perform a complex FFT.
Is this intended? Is there a reason you don't just use a real data type for x_in and perform a real FFT?
Todd
vincent-torri
Beginner
152 Views
it is intended, indeed. It is the first time that i use the Intel compiler and the MKL. So I am doing tests. Especially tests that can theoretically be computed. It works for an even function (I tried x^2-1) and the imaginary part is 0. But not with that simple odd function.

In the future, i will use the FFT on a complex valued function (more precisely, i want to solve a non linear Schrdinger equation). So I prefer getting familiar with the complex API.
Dmitry_B_Intel
Employee
152 Views
Hi vtorri,

The function you take FFT of is not odd. An odd function must have x_in[0]==0.
Instead of x_in = [-1 + (0:31)*2/32] try [0, -1 + (0:30)*2/30].

Thanks
Dima
vincent-torri
Beginner
152 Views
Of course the function that i took is odd. I discretized x |--> x on [-1,1] which is obviously odd. On the other hand, your remark raises a good point: must I "shift" the vector values, with something like:

int shift = N/2 + 1;
int i;
for (i = 0; i < N; i++)
out2 = out[i+shift)%N];

so that the discretized value of x=0 is the first value of my array ?
Dmitry_B_Intel
Employee
152 Views
Right, the input x[0:N-1] to FFT is expecpected to bedirectly mapped frominterval[0,P] rather than from [-P/2,P/2], so that f(0) is the first value of the input array.

Thanks
Dima
vincent-torri
Beginner
152 Views
I shifted the values:

[bash]  for (i = 0; i < 32; i++)
    {
      float x;

      x = -1.0f + (float)(2 * i) / 32.0f;
      tmp = fct_der(x);
    }

  shift = 32 / 2;
  for (i = 0; i < 32; i++)
    {
      x_in = tmp[(i+shift)%32];
    }[/bash]

but the result is still wrong (1st column : i, 2nd : x_in, 3rd and 4th x_out) :

[bash]0 ** 0.000000 ** -1.000000 + i 0.000000
1 ** 0.062500 ** 1.000000 + i -10.153171
2 ** 0.125000 ** -1.000000 + i 5.027339
3 ** 0.187500 ** 1.000000 + i -3.296558
4 ** 0.250000 ** -1.000000 + i 2.414214
5 ** 0.312500 ** 1.000000 + i -1.870868
6 ** 0.375000 ** -1.000000 + i 1.496606
7 ** 0.437500 ** 1.000000 + i -1.218504
8 ** 0.500000 ** -1.000000 + i 1.000000
9 ** 0.562500 ** 1.000000 + i -0.820679
10 ** 0.625000 ** -1.000000 + i 0.668179
11 ** 0.687500 ** 1.000000 + i -0.534511
12 ** 0.750000 ** -1.000000 + i 0.414214
13 ** 0.812500 ** 1.000000 + i -0.303347
14 ** 0.875000 ** -1.000000 + i 0.198912
15 ** 0.937500 ** 1.000000 + i -0.098491
16 ** -1.000000 ** -1.000000 + i 0.000000
17 ** -0.937500 ** 1.000000 + i 0.098491
18 ** -0.875000 ** -1.000000 + i -0.198912
19 ** -0.812500 ** 1.000000 + i 0.303347
20 ** -0.750000 ** -1.000000 + i -0.414214
21 ** -0.687500 ** 1.000000 + i 0.534511
22 ** -0.625000 ** -1.000000 + i -0.668179
23 ** -0.562500 ** 1.000000 + i 0.820679
24 ** -0.500000 ** -1.000000 + i -1.000000
25 ** -0.437500 ** 1.000000 + i 1.218504
26 ** -0.375000 ** -1.000000 + i -1.496606
27 ** -0.312500 ** 1.000000 + i 1.870868
28 ** -0.250000 ** -1.000000 + i -2.414214
29 ** -0.187500 ** 1.000000 + i 3.296558
30 ** -0.125000 ** -1.000000 + i -5.027339
31 ** -0.062500 ** 1.000000 + i 10.153171
[/bash]

Any idea ?
Dmitry_B_Intel
Employee
152 Views
Hi vtorri,

To be odd, the discretization should obey the property:
x_in[0] = 0
x_in = -x_in[N-k], k=1...N-1
In your case this property doesn't hold for k=16.

Thanks
Dima
vincent-torri
Beginner
152 Views
Ok. So the shift is not what should be done, actually. It's more a symetry. And it's certainly why it works without problem with en even function. I'll try that.
vincent-torri
Beginner
152 Views
it works well for an even function because the number of coefficients are even. The problem is the number of coefficient for an odd function if I have to satisfy : x_in[0]=0 and x_in=x_in[N-i] i=1...N-1 (the number of coef is then odd)
Dmitry_B_Intel
Employee
152 Views

Consider even N = 2L.

For even function this should hold:
x_in[0] = any
x_in = x_in[N-k]
In particular, x_in = x_in which means x_in is arbitrary.

For odd function this should hold:
x_in[0] = any
x_in = -x_in[N-k]
In particular, x_in = -x_in which means x_in must be zero.

Thanks
Dima

vincent-torri
Beginner
152 Views
I have found my error. The way I discretized x |----> x is wrong.The value ax x=-1 should be 0 and not -1. Thank you for your answers !
Reply