Hur kan MLK göra så snabba matrismultiplikationer?

Jag behöver tydligen inte normalisera när jag kör FFTPack 5.1, men kör jag FFTPack 5.0 (Det första optimerade) https://github.com/marton78/pffft/blob/master/fftpack.c så måste jag normalisera efter invers FFT.

Vi kan ta ett arbetsexempel.

#define row 3
#define column 3

int main() {
	clock_t start, end;
	float cpu_time_used;
	start = clock();

	/* Create X */
	float X[row * column] = { -0.6485,   0.2851,   1.3475,
					       -0.6743, -0.1499, -1.5549,
						1.4951,   0.4504, -0.4982 };

	/* Do FFT2 */
	fft2(X, row, column);

	/* Print */
	print(X, row, column);

	end = clock();
	cpu_time_used = ((float)(end - start)) / CLOCKS_PER_SEC;
	printf("\nTotal speed  was %f\n", cpu_time_used);

	return EXIT_SUCCESS;
}

Utskriften blir

0.0058111       0.0258111       -0.1242458
0.1611111       -0.2950723      0.3444953
0.3681955       -0.2190056      0.0864390

Kör jag

ifft2(X, row, column);

direkt efter

fft2(X, row, column);

så får jag samma matris X tillbaka.

Men om det inte spelar någon roll, så låter jag det vara.
Det viktigaste är att jag kan återskapa datat.

i7-8700k:
		1 tråd		8 trådar
mkl		0.36 ms		0.14 ms
mkl-fftw	0,36 ms		0.15 ms
fftw		0,28 ms		0.094 ms

i7-12850HX:
		1 tråd		8 trådar
mkl-fftw	0,19 ms		0.11 ms
fftw		0,17 ms		0.09 ms

Jag laddar upp exemplet här: https://easyupload.io/bg2zbv

Du laddar ned FFTPack_5_1.zip och packar upp den. Då finns det två intressanta filer: fft2.c och ifft2.c

Kör med detta exempel

void print(float A[], size_t row, size_t column) {
	size_t i, j;
	for (i = 0; i < row; i++) {
		for (j = 0; j < column; j++) {
			printf("%0.7f\t", *(A++));
		}
		printf("\n");
	}
	printf("\n");
}

#define row 3
#define column 3

int main() {

	/* Create X */
	float X[row * column] = { -0.6485,   0.2851,   1.3475,
							   -0.6743, -0.1499, -1.5549,
								1.4951,   0.4504, -0.4982 };

	/* Do inverse FFT2 */
	fft2(X, row, column);
	ifft2(X, row, column);

	/* Print */
	print(X, row, column);

	return 0;
}

Jag har löst detta problem nu. Jag bara använde cfft2f/b och multiplicerade varje tal med n*n. Då fick jag samma resultat som Matlab. Om det är värt det, är tveksamt. Det är ju bara en skalning.

Hur som helst så kanske det är värt att testa MLK för FFT, nu när FFTpack 5.1 är något långsamare än FFTPack 5.1.

Hur ser deras komplexa rutin ut? Kräver den C99?

Jag kollade upp FFT för MKL och den kräver C99 för FFT complex. Jag försöker alltid att vara halvt besvärlig igenom att strikt hålla mig till ANSI C. Detta har med att ANSI C är liksom oändlig och kod som är skrivet i ANSI C blir för evig.

Men FFTW innehåller så mycket filer. Jag gillar när det är bara några .c filer som man kan bara klistra in. Då slipper man allt görm med att länka hit och dit.

FFTPack är bara en mapp med massa .f filer som man exporterar till .c via f2c. Inge mer.

Jag är nyfiken på vilken mapp du pratar om då i FFTW3.

Jag misstänker att FFTW3 använder sig av hårdvaran?

Nu verkar jag ha problem igen med cfft2f hos FFTPack 5.1.
Jag trodde jag kunde multiplicera alla element med antalet celler i matrisen t.ex. m*n. Men det verkar inte gälla för alla scenarion.

Är det någon som vet hur FFTPack är normaliserad? Enligt manualen så är den normaliserad. Men jag har ett behov utav att få den mer lik matlabs.

Jadu.
Det fungerar bra vid 1D, men vid 2D så fungerar det inte alls.

Jag har en liten idé!
FFTpack som jag använde först, den där optimerade från Pretty Fast FFT på GitHub. Den kan bara göra 1D FFT och den gav exakt samma resultat som Matlab gjorde. Dessutom var den riktigt optimerad, trots bara en .c och .h fil. Orsaken har med att skaparen har själv tagit bort lite onödiga saker så som massa gotos osv.

Men orsaken varför jag gick över till FFTPack 5.1 har med att den har 2D beräkningar för komplexa tal, trots ANSI C.

Matlab's FFT2 fungerar så att den tar FFT på alla kolumner, sedan tar den transponatet av resultatet och sedan tar FFT på alla rader. Så här alltså.

fft(fft(X)')'

Pretty Fast FFT https://github.com/marton78/pffft/blob/master/fftpack.h
Denna har cffti, cfftf, cfftb och med dessa så kanske jag kan göra en FFT 1D på en 2D matris.
Men då är frågan hur jag skapar min komplexa array?

Hur ska arrayen c vara? Jag vet att den ska vara komplex, men om den ska vara n lång, så måste den väll vara 2*n lång om den ska vara komplex?

Jag har svårt att tolka hur arrayen c ska skapas. Det är av datatypen float. Så jag misstänker att antingen ska arrayen vara 2*n, eller så så blir det bara halva frekvensen?

Tackar. Jag ska testa detta.
Jag hoppas att det ska gå bättre att göra 2D FFT av 1D FFT. Jag gillar inte att det är icke-samma resultat från FFTPack 5.1, trots att man kan återskapa signalen med invers FFT.

Återkommer.

Nu har jag testat ytterligare.

Jag testade göra FFT 2D av två stycken rffti -> rfftf -> cffti -> cfftf
r står för tids-data och c står för komplex data.

Med den där optimerade fftpack.h och fftpack.c från Pretty Fast FFT så tog det 3.15 sekunder att göra FFT 2D på en 5000*5000 matris.
Med FFTpack 5.1 så tog det 5.3 sekunder att göra FFT 2D på en 5000*5000 matris.

Dessutom gav fftpack.h och fftpack.c filen från Pretty Fast FFT 100% exakt samma resultat som Matlab.
Notera att detta är utan någon optimeringsflagga eller liknande. Bara ren ANSI C kod som är plattformsoberoende.

Det var trixigt att få till rffti -> rfftf -> cffti -> cfftf då det inte finns någon manual eller ens folk som verkar använda det, trots att FFTpack är snabbt om man väljer den från Pretty Fast FFT.

Jo, fast FFTW använder säkert trådar. fftpack.h och fftpack.c kommer från PFFFT och det är en modifierad FFTpack version. Troligtvis så har dom själv manuellt tagit bort massa saker för att det ska gå snabbare.

FFTW är säkert jättebra, men det verkar vara mycket filer. fftpack.c/h från PFFFT är bara två filer och är riktigt snabbt. Så jag känner mig nöjd.

Tror jag inte kan få det snabbare. FFTW3 sägs inte vara den snabbaste heller. Tror mycket handlar om vilken hårdvara och kompilator man använder. Jag kör MSVC.

Testa och kör dessa exempel: https://easyupload.io/b7kr10

Kolla vad du får för tid vid 5000*5000 matris.
Det tog 3.15 sekunder för mig.

Det går självklart snabbare I Octave

>> tic; fft2(randn(5000,5000)); toc
Elapsed time is 0.671597 seconds.
>>

Så 1.7 sekunder är inte illa, med tanke på att det är bara en enkel C fil Jag har försökt optimera detta så bra som möjligt.

Så:
Din dator FFT2 med min kod: 1.7 sekunder
Min dator FFT2 med min kod: 3.15 sekunder
Din dator FFT2 med octave: 0.13 sekunder
Min dator FFT2 med min kod: 0.67 sekunder

Detta betyder att din dator är ungefär 2 gånger snabbare än min.

Ja, det gick mycket snabbare.

>> m=single(randn(5000,5000)); tic; fft2(m); toc
Elapsed time is 0.191178 seconds.

Problemet är att jag vet inte vilka filer jag ska ha. FFTW är inte pedagogiskt skrivet direkt.
Sedan så kräver FFTW3 att man ska vara licensierad mot MIT om det ska användas i kommersiella produkter.