tyler
849306d65e
Build system - Add a single hand-written Makefile (GNU make + gcc, pure C). Builds raylib from vendored source and links rspektrum with no premake/lua. Targets: all/run/test/bench/clean; release default, DEBUG=1 for debug; ARCH overridable (defaults to -march=x86-64-v3). - Remove premake entirely: rspektrum.make, raylib.make, build/premake5* binaries, build/premake5.lua, build/ecc/*. The generated top-level Makefile was gitignored, so hand-edits to it were silently lost. - Vendor raylib src/ into the repo (was gitignored -> fresh clones could not build). Commit only src/ (~16MB); examples/projects stay local. Verified: a build from the git-tracked tree alone succeeds offline. - Release flags bumped to -O3 -ffast-math with a portable arch baseline (x86-64-v3 = AVX2+FMA on x64, SSE2 on x86). Confirmed FMA/AVX codegen in fft.o. FFT optimization (src/fft.c) - Precompute twiddle factors and the bit-reversal permutation once per size, cached as a small plan table (FFTW's idea, lightweight). Removes the per-butterfly cexpf() and per-element bit-twiddling that dominated. - 3.6x faster on the mlnl_samples.wav STFT workload (2048-pt, -O2 same flags both sides): 81us -> 22us per FFT. With the new -O3/-ffast-math/ AVX2 release flags stacked: ~15us (5.5x vs the old -O2 baseline). - Verified vs a double-precision reference DFT: 1e-6 relative error, round-trip 2.4e-7. Drop-in: same FFT() signature. Tests/bench (bench/) - fft_verify.c: FFT vs reference DFT + round-trip check (make test). - fft_bench.c: times the real STFT workload (make bench). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
53 lines
2.0 KiB
C
53 lines
2.0 KiB
C
// fft_verify.c - correctness check for FFT(): compare vs naive DFT, and check
|
|
// the forward->inverse round-trip recovers the input.
|
|
#include <stdio.h>
|
|
#include <stdlib.h>
|
|
#include <math.h>
|
|
#include <complex.h>
|
|
#include "../src/fft.h"
|
|
|
|
// Reference DFT accumulated in double precision so the *reference* isn't the
|
|
// source of error when comparing a float FFT at large n.
|
|
static void naive_dft(const float complex* x, float complex* X, int n) {
|
|
for (int k = 0; k < n; k++) {
|
|
double complex s = 0;
|
|
for (int t = 0; t < n; t++)
|
|
s += x[t] * cexp(-2.0 * M_PI * I * k * t / n);
|
|
X[k] = (float complex)s;
|
|
}
|
|
}
|
|
|
|
int main(void) {
|
|
int sizes[] = {2, 4, 8, 16, 64, 256, 1024, 2048};
|
|
double worst_dft = 0, worst_rt = 0;
|
|
for (unsigned si = 0; si < sizeof(sizes)/sizeof(sizes[0]); si++) {
|
|
int n = sizes[si];
|
|
float complex* x = malloc(n*sizeof(float complex));
|
|
float complex* X = malloc(n*sizeof(float complex));
|
|
float complex* Xn = malloc(n*sizeof(float complex));
|
|
float complex* xrt = malloc(n*sizeof(float complex));
|
|
srand(1234 + n);
|
|
for (int i = 0; i < n; i++)
|
|
x[i] = (rand()/(float)RAND_MAX - 0.5f) + (rand()/(float)RAND_MAX - 0.5f)*I;
|
|
|
|
FFT(x, X, n, false);
|
|
naive_dft(x, Xn, n);
|
|
FFT(X, xrt, n, true);
|
|
|
|
double norm = 0;
|
|
for (int i = 0; i < n; i++) norm += cabsf(Xn[i]);
|
|
norm /= n; // mean magnitude, for a relative tolerance
|
|
for (int i = 0; i < n; i++) {
|
|
double e1 = cabsf(X[i] - Xn[i]) / norm; // relative to spectrum scale
|
|
double e2 = cabsf(xrt[i] - x[i]);
|
|
if (e1 > worst_dft) worst_dft = e1;
|
|
if (e2 > worst_rt) worst_rt = e2;
|
|
}
|
|
free(x); free(X); free(Xn); free(xrt);
|
|
}
|
|
printf("max relative |FFT - naive DFT| = %.3e\n", worst_dft);
|
|
printf("max round-trip error = %.3e\n", worst_rt);
|
|
printf("%s\n", (worst_dft < 1e-4 && worst_rt < 1e-4) ? "PASS" : "FAIL");
|
|
return (worst_dft < 1e-4 && worst_rt < 1e-4) ? 0 : 1;
|
|
}
|