Фурье

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "omp.h"
#include <complex.h>

#define N 8
#define PI 3.141592653589793

void swap(double complex* a, int i, int j)
{
    double complex temp;
    temp = a[i];
    a[i] = a[j];
    a[j] = temp;
}

void shuffle(double complex* a)
{
    int j = N / 2, t;
    for(int i = 1; i < N-1; i++)
    {
        if (i < j) swap(a, i, j);
        t = N / 2;
        while(j / t == 1)
        {
            j = j - t;
            t /= 2;
        }
        j += t;
    }
}



void get_omega(double complex* W)
{
    for(int k = 0; k < N / 2; k++)
        W[k] = cexp(2*PI*I*k/N);
}

void Furhe(double complex* x)
{
    double complex W[N/2];
    double complex omega;
    double complex temp;
    int i;
    int l = N / 2;
    int t;
    get_omega(W);
    //for(int i = 0; i <= N/2; i++)
        //printf("%f + %f", creal(W[i]), cimag(W[i]));
    shuffle(x);
    int k = 1;
    l = N / 2;
    while(k < N)
    {
        t = 0;
        omega = W[0];
        for(int j = 0; j < k; j++)
        {
            i = j;
            while(i < N)
            {
                temp = x[i];
                x[i] = temp + omega*x[i+k];
                x[i + k] = temp - omega*x[i+k];
                i = i +2*k;
            }
            t += l;
            omega = W[t];
        }
        k = 2*k;
        l = l /2;
    }


}

int main()
{
    double complex x[N] = {0, 1, 0, 1};

    Furhe(x);
    for(int i = 0; i < N; i++)
        printf("%f + %f\n", creal(x[i]), cimag(x[i]));
    return 0;
}