Почему FFTW в Windows работает быстрее, чем в Linux?

16 FFT очень мало. Вы обнаружите, что БПФ меньше, чем, скажем, 64, будут жестко закодированы на ассемблере без циклов, чтобы получить максимально возможную производительность. Это означает, что они могут быть очень восприимчивы к изменениям в наборах инструкций, оптимизации компилятора и даже к 64- или 32-битным словам.

Что происходит, когда вы запускаете тест размеров БПФ от 16 до 1048576 в степени двойки? Я говорю это, потому что конкретная жестко запрограммированная процедура asm в Linux может быть не лучшим образом оптимизирована для вашей машины, тогда как вам может повезти в реализации Windows для этого конкретного размера. Сравнение всех размеров в этом диапазоне даст вам лучшее представление о производительности Linux и Windows.

Вы калибровали FFTW? При первом запуске FFTW угадывает самую быструю реализацию для каждой машины, однако, если у вас есть специальные наборы инструкций, кеш-память определенного размера или другие функции процессора, это может существенно повлиять на скорость выполнения. В результате выполнение калибровки проверит скорость различных процедур БПФ и выберет самую быструю для вашего конкретного оборудования. Калибровка включает в себя повторное вычисление планов и сохранение сгенерированного файла «Мудрость» FFTW. Сохраненные данные калибровки (это длительный процесс) можно затем использовать повторно. Я предлагаю сделать это один раз при запуске вашего программного обеспечения и каждый раз повторно использовать файл. Я заметил 4-10-кратное улучшение производительности для определенных размеров после калибровки!

Ниже приведен фрагмент кода, который я использовал для калибровки FFTW для определенных размеров. Обратите внимание, что этот код дословно вставлен из библиотеки DSP, над которой я работал, поэтому некоторые вызовы функций специфичны для моей библиотеки. Я надеюсь, что конкретные звонки FFTW будут полезны.

// Calibration FFTW
void DSP::forceCalibration(void)
{
// Try to import FFTw Wisdom for fast plan creation
FILE *fftw_wisdom = fopen("DSPDLL.ftw", "r");

// If wisdom does not exist, ask user to calibrate
if (fftw_wisdom == 0)
{
    int iStatus2 = AfxMessageBox("FFTw not calibrated on this machine."\
        "Would you like to perform a one-time calibration?\n\n"\
        "Note:\tMay take 40 minutes (on P4 3GHz), but speeds all subsequent FFT-based filtering & convolution by up to 100%.\n"\
        "\tResults are saved to disk (DSPDLL.ftw) and need only be performed once per machine.\n\n"\
        "\tMAKE SURE YOU REALLY WANT TO DO THIS, THERE IS NO WAY TO CANCEL CALIBRATION PART-WAY!", 
        MB_YESNO | MB_ICONSTOP, 0);

    if (iStatus2 == IDYES)
    {
        // Perform calibration for all powers of 2 from 8 to 4194304
        // (most heavily used FFTs - for signal processing)
        AfxMessageBox("About to perform calibration.\n"\
            "Close all programs, turn off your screensaver and do not move the mouse in this time!\n"\
            "Note:\tThis program will appear to be unresponsive until the calibration ends.\n\n"
            "\tA MESSAGEBOX WILL BE SHOWN ONCE THE CALIBRATION IS COMPLETE.\n");
        startTimer();

        // Create a whole load of FFTw Plans (wisdom accumulates automatically)
        for (int i = 8; i <= 4194304; i *= 2)
        {
            // Create new buffers and fill
            DSP::cFFTin = new fftw_complex[i];
            DSP::cFFTout = new fftw_complex[i];
            DSP::fconv_FULL_Real_FFT_rdat = new double[i];
            DSP::fconv_FULL_Real_FFT_cdat = new fftw_complex[(i/2)+1];
            for(int j = 0; j < i; j++)
            {
                DSP::fconv_FULL_Real_FFT_rdat[j] = j;
                DSP::cFFTin[j][0] = j;
                DSP::cFFTin[j][1] = j;
                DSP::cFFTout[j][0] = 0.0;
                DSP::cFFTout[j][1] = 0.0;
            }

            // Create a plan for complex FFT.
            // Use the measure flag to get the best possible FFT for this size
            // FFTw "remembers" which FFTs were the fastest during this test. 
            // at the end of the test, the results are saved to disk and re-used
            // upon every initialisation of the DSP Library
            DSP::pCF = fftw_plan_dft_1d
                (i, DSP::cFFTin, DSP::cFFTout, FFTW_FORWARD, FFTW_MEASURE);

            // Destroy the plan
            fftw_destroy_plan(DSP::pCF);

            // Create a plan for real forward FFT
            DSP::pCF = fftw_plan_dft_r2c_1d
                (i, fconv_FULL_Real_FFT_rdat, fconv_FULL_Real_FFT_cdat, FFTW_MEASURE);

            // Destroy the plan
            fftw_destroy_plan(DSP::pCF);

            // Create a plan for real inverse FFT
            DSP::pCF = fftw_plan_dft_c2r_1d
                (i, fconv_FULL_Real_FFT_cdat, fconv_FULL_Real_FFT_rdat, FFTW_MEASURE);

            // Destroy the plan
            fftw_destroy_plan(DSP::pCF);

            // Destroy the buffers. Repeat for each size
            delete [] DSP::cFFTin;
            delete [] DSP::cFFTout;
            delete [] DSP::fconv_FULL_Real_FFT_rdat;
            delete [] DSP::fconv_FULL_Real_FFT_cdat;
        }

        double time = stopTimer();

        char * strOutput;
        strOutput = (char*) malloc (100);
        sprintf(strOutput, "DSP.DLL Calibration complete in %d minutes, %d seconds\n"\
            "Please keep a copy of the DSPDLL.ftw file in the root directory of your application\n"\
            "to avoid re-calibration in the future\n", (int)time/(int)60, (int)time%(int)60);
        AfxMessageBox(strOutput);

        isCalibrated = 1;

        // Save accumulated wisdom
        char * strWisdom = fftw_export_wisdom_to_string();  
        FILE *fftw_wisdomsave = fopen("DSPDLL.ftw", "w");
        fprintf(fftw_wisdomsave, "%s", strWisdom);

        fclose(fftw_wisdomsave);
        DSP::pCF = NULL;
        DSP::cFFTin = NULL;
        DSP::cFFTout = NULL;
        fconv_FULL_Real_FFT_cdat = NULL;
        fconv_FULL_Real_FFT_rdat = NULL;
        free(strOutput);
    }
}
else 
{
    // obtain file size.
    fseek (fftw_wisdom , 0 , SEEK_END);
    long lSize = ftell (fftw_wisdom);
    rewind (fftw_wisdom);

    // allocate memory to contain the whole file.
    char * strWisdom = (char*) malloc (lSize);

    // copy the file into the buffer.
    fread (strWisdom,1,lSize,fftw_wisdom);

    // import the buffer to fftw wisdom
    fftw_import_wisdom_from_string(strWisdom);

    fclose(fftw_wisdom);
    free(strWisdom);

    isCalibrated = 1;

    return;
}
}

Секрет в том, чтобы создать план с использованием флага FFTW_MEASURE, который специально измеряет сотни процедур, чтобы найти самую быструю для вашего конкретного типа БПФ (действительного, комплексного, 1D, 2D) и размера:

DSP::pCF = fftw_plan_dft_1d (i, DSP::cFFTin, DSP::cFFTout, 
   FFTW_FORWARD, FFTW_MEASURE);

Наконец, все эталонные тесты также должны выполняться с одной стадией планирования БПФ вне выполнения, вызываемой из кода, скомпилированного в режиме выпуска с включенными оптимизациями и отсоединенными от отладчика. Бенчмарки должны выполняться в цикле со многими тысячами (или даже миллионами) итераций, а затем брать среднее время выполнения для вычисления результата. Как вы, наверное, знаете, этап планирования занимает значительное количество времени, а выполнение рассчитано на многократное выполнение одного плана.