芯耀
7790
一、FFTW簡(jiǎn)介
FFT（快速傅里葉變換）是數(shù)字信號(hào)處理的超級(jí)經(jīng)典算法，在量子物理、光譜分析、音視頻流信號(hào)處理、石油勘探、地震預(yù)報(bào)、天氣預(yù)報(bào)、概率論、編碼理論、醫(yī)學(xué)斷層診斷等領(lǐng)域有著廣泛的應(yīng)用。在數(shù)學(xué)中，傅里葉級(jí)數(shù)（Fourier series）是把類(lèi)似波的函數(shù)表示成簡(jiǎn)單正弦波的方式。
快速傅里葉變換的算法有許多，既有FFTW這種開(kāi)源庫(kù)，一些芯片企業(yè)也針對(duì)自家芯片設(shè)計(jì)FFT優(yōu)化庫(kù)，例如Intel的MKL_FFT，華為的KML_FFT。
打住，這里有一個(gè)問(wèn)題？
既然有FFTW這種通用的、開(kāi)源的還免費(fèi)的庫(kù)，為什么芯片公司還要費(fèi)人費(fèi)錢(qián)去自己搞FFT優(yōu)化庫(kù)呢？
通過(guò)這篇文章，也許你就會(huì)知其一二。
FFTW（Fast Fourier Transform in the West）是一款用于快速傅里葉變換（FFT）的開(kāi)源庫(kù)，它是自由軟件許可證（類(lèi)似于GPL）下發(fā)布的，您可以在其官方網(wǎng)站或代碼倉(cāng)庫(kù)上找到最新版本的FFTW，然后根據(jù)許可證的規(guī)定使用它。它可以用于在多種計(jì)算平臺(tái)上進(jìn)行高效的數(shù)值計(jì)算。FFTW通過(guò)使用預(yù)分解算法和緩存技術(shù)，可以大大提高FFT的計(jì)算速度，與傳統(tǒng)的FFT算法相比，可以節(jié)省數(shù)倍的計(jì)算時(shí)間。FFTW不僅可以用于科學(xué)計(jì)算和信號(hào)處理，還可以用于圖像處理、音頻處理等領(lǐng)域。FFTW的優(yōu)秀性能使其成為許多科學(xué)計(jì)算和工程應(yīng)用中不可或缺的工具。官網(wǎng)地址：http://fftw.org/ github地址：https://github.com/FFTW/fftw3
本實(shí)驗(yàn)僅是將FFTW搬到RISC-V硬件平臺(tái)上編譯和運(yùn)行，文中將SG2042和FT2000進(jìn)行了單核性能的對(duì)比，只是為了先建立一個(gè)優(yōu)化對(duì)標(biāo)，并不代表這兩款芯片性能的差距。通常一些軟件測(cè)試出來(lái)的性能，與以下幾個(gè)因素相關(guān)：
芯片的硬件算力水平
編譯器（RISC-V的GCC屬于初代版本）
編譯選項(xiàng)（尤其是是否啟用SIMD加速）
操作系統(tǒng)和依賴(lài)庫(kù)等其他因素
二、平臺(tái)環(huán)境
[硬件參數(shù)]
處理器: 算能SG2042 x 1
核心數(shù): 64核
L1 Cache: I:64KB and D:64KB
L2 Cache: 1MB/Cluster
L3 Cache: 64MB System Cache
[軟件環(huán)境]
linux版本: 22.10
gcc版本: 10.2.0
三、安裝流程
我們從官網(wǎng)上下載最新的程序包，當(dāng)前最新版本為：3.3.10 下載鏈接：http://fftw.org/download.html
我們下載的文件為：fftw-3.3.10.tar.gz
解壓：
ubuntu@perfxlab:~$ tar zxvf fftw-3.3.10.tar.gz
配置：
ubuntu@perfxlab:~/fftw-3.3.10$ ./configure --prefix=/usr/fftw3
安裝：
ubuntu@perfxlab:~/fftw-3.3.10$ make install
配置環(huán)境變量：
我們安裝后需要配置相關(guān)的環(huán)境變量，打開(kāi)~/.bashrc, 在文件在最后面添加如下代碼：
export LD_LIBRARY_PATH=/usr/fftw3/lib:$LD_LIBRARY_PATH
驗(yàn)證是否安裝成功：
查詢(xún)FFTW版本：
ubuntu@perfxlab:~/fftw-3.3.10$ fftw-wisdom --version
出現(xiàn)以下信息即為安裝成功：
fftw-wisdom tool for FFTW version 3.3.10.
Copyright (c) 2003, 2007-14 Matteo Frigo
Copyright (c) 2003, 2007-14 Massachusetts Institute of Technology
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program; if not, write to the Free Software
Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
我們寫(xiě)一個(gè)簡(jiǎn)單的測(cè)試代碼驗(yàn)證是否安裝成功：
#include <stdio.h>
#include <stdlib.h>
#include <fftw3.h>
int main() {
int N = 8; // 信號(hào)長(zhǎng)度
fftw_complex *in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
fftw_complex *out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
fftw_plan plan = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
// 初始化輸入信號(hào)
for (int i = 0; i < N; i++) {
in[i][0] = i; // 實(shí)部
in[i][1] = 0; // 虛部
}
// 執(zhí)行傅里葉變換
fftw_execute(plan);
// 打印結(jié)果
printf("Input Signal:n");
for (int i = 0; i < N; i++) {
printf("%f + %fin", in[i][0], in[i][1]);
}
printf("Transformed Signal:n");
for (int i = 0; i < N; i++) {
printf("%f + %fin", out[i][0], out[i][1]);
}
fftw_destroy_plan(plan);
fftw_free(in);
fftw_free(out);
return 0;
}
編譯代碼：
ubuntu@perfxlab:~/cheng$ gcc -o fftw_test testfftw.c -lfftw3 -lm
執(zhí)行程序運(yùn)行結(jié)果如下：
ubuntu@perfxlab:~/cheng$ ./fftw_test
Input Signal:
0.000000 + 0.000000i
1.000000 + 0.000000i
2.000000 + 0.000000i
3.000000 + 0.000000i
4.000000 + 0.000000i
5.000000 + 0.000000i
6.000000 + 0.000000i
7.000000 + 0.000000i
Transformed Signal:
28.000000 + 0.000000i
-4.000000 + 9.656854i
-4.000000 + 4.000000i
-4.000000 + 1.656854i
-4.000000 + 0.000000i
-4.000000 + -1.656854i
-4.000000 + -4.000000i
-4.000000 + -9.656854i
好了，結(jié)果正確，目前已經(jīng)安裝測(cè)試完成，可以愉快的使用了。
四、性能測(cè)試和對(duì)比
下面，我們將基于FFTW的開(kāi)源項(xiàng)目，對(duì)SG2042進(jìn)行性能評(píng)估，并和FT2000做簡(jiǎn)單對(duì)比。
測(cè)試代碼如下
我們主要測(cè)試了FFTW（Fastest Fourier Transform in the West）庫(kù)中的核心函數(shù)fftw_execute的性能。我們通過(guò)多次執(zhí)行fftw_execute函數(shù)并取平均值來(lái)評(píng)估其性能，并將結(jié)果以GFlops（每秒十億次浮點(diǎn)運(yùn)算數(shù)）為單位進(jìn)行了測(cè)量。in_place 、out_place代表是否使用了同一內(nèi)存。
#include <stdio.h>
#include <time.h>
#include <fftw3.h>
#include <math.h> // 包含math.h以使用log2函數(shù)
int main() {
int signal_lengths[] = {32, 64, 128, 256, 512, 1024, 2048, 4096, 8192, 16384}; // 不同信號(hào)長(zhǎng)度
int num_lengths = sizeof(signal_lengths) / sizeof(signal_lengths[0]);
for (int i = 0; i < num_lengths; i++) {
int N = signal_lengths[i];
fftw_complex *in = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * N);
fftw_complex *out = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * N);
// 填充輸入數(shù)據(jù) (使用雙精度浮點(diǎn)數(shù))
for (int j = 0; j < N; j++) {
in[j][0] = 1.0; // 實(shí)部 (double)
in[j][1] = 0.0; // 虛部 (double)
}
double total_elapsed_time = 0.0;
int num_iterations = 20; // 每個(gè)長(zhǎng)度的信號(hào)執(zhí)行20次
// 執(zhí)行雙精度浮點(diǎn)數(shù)傅立葉變換
fftw_plan plan = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
for (int j = 0; j < num_iterations; j++) {
clock_t start_time = clock();
fftw_execute(plan);
clock_t end_time = clock();
double elapsed_time = (double)(end_time - start_time) / CLOCKS_PER_SEC;
total_elapsed_time += elapsed_time;
}
fftw_destroy_plan(plan);
double average_time = total_elapsed_time / num_iterations;
// 計(jì)算FLOPS（每秒浮點(diǎn)運(yùn)算數(shù)）
double flops = 5.0 * N * log2(N) / average_time; // 假設(shè)FFT算法的復(fù)雜度為5*N*log2(N)
// 計(jì)算GFLOPS（每秒十億次浮點(diǎn)運(yùn)算數(shù)）
double gflops = flops / 1e9;
printf("平均FFTW執(zhí)行時(shí)間（長(zhǎng)度 %d）: %f秒 GFLOPS : %fn", N, average_time , gflops);
fftw_free(in);
fftw_free(out);
}
return 0;
}
如果要測(cè)試float類(lèi)型：則需要修改以下地方：
將fftw_complex改為fftwf_complex，將fftw_plan_dft_1d改為fftwf_plan_dft_1d。
編譯選項(xiàng)
以下是編譯double類(lèi)型，如果編譯float 需要用-lfftw3f。
SG2042編譯選項(xiàng)（注意，F(xiàn)FTW目前沒(méi)有支持risc-v的simd版本）
-march=rv64gcv0p7_zfh_xtheadc -mabi=lp64d -mtune=c920 -O3
FT2000 編譯選項(xiàng)
-Wall -Wextra -funroll-loops -O3 -DNDEBUG -DNOPEN_HASWELL_OPT
性能圖和對(duì)比
數(shù)據(jù)類(lèi)型為float in_place：
數(shù)據(jù)類(lèi)型為float out_place：
數(shù)據(jù)類(lèi)型為 double in_place：
數(shù)據(jù)類(lèi)型為 double out_place
五、簡(jiǎn)單總結(jié)
將FFTW遷移到RISC-V處理器SG2042上的實(shí)驗(yàn)過(guò)程還是相對(duì)簡(jiǎn)單，這得益于RISC-V軟件生態(tài)的不斷完善。但同時(shí)也依然存在一些問(wèn)題：
1：從上面的測(cè)試結(jié)果來(lái)看，數(shù)據(jù)類(lèi)型為float時(shí)， SG2042的表現(xiàn)相對(duì)FT2000 有比較明顯的差距；數(shù)據(jù)類(lèi)型為double時(shí)，隨著信號(hào)長(zhǎng)度的增加，性能差距逐漸減小，趨于接近。
2：FFTW目前沒(méi)有支持RISC-V的Vector版本，仍有極大性能潛力待挖掘。我們或者等待FFTW發(fā)行版支持RISC-V；或者有人自研RISC-V FFT優(yōu)化算法庫(kù)?；氐搅宋恼伦铋_(kāi)始的問(wèn)題“既然有FFTW這種通用的、開(kāi)源的還免費(fèi)的庫(kù)，為什么芯片公司還要費(fèi)人費(fèi)錢(qián)去自己搞FFT優(yōu)化庫(kù)呢？”，實(shí)際上還有許多類(lèi)似這樣的計(jì)算庫(kù)，存在類(lèi)似的問(wèn)題。
下載ECAD模型
