diff --git a/.gitignore b/.gitignore index 7dea02f3..967aa5dd 100644 --- a/.gitignore +++ b/.gitignore @@ -44,3 +44,5 @@ dist .DS_Store .eslintcache + +CPU.* \ No newline at end of file diff --git a/benchmark/reimFFT.ts b/benchmark/reimFFT.ts new file mode 100644 index 00000000..ee4d11d9 --- /dev/null +++ b/benchmark/reimFFT.ts @@ -0,0 +1,28 @@ +import { reimFFT } from '../src/reim/reimFFT.ts'; + +const size = 2 ** 16; +const re = new Float64Array(size); +const im = new Float64Array(size); +for (let i = 0; i < size; i++) { + re[i] = Math.random(); + im[i] = Math.random(); +} +// Warmup +reimFFT({ re, im }); + +// Benchmark: repeat until ~5s +const targetMs = 5000; +let iterations = 0; + +console.time('reimFFT'); +const start = performance.now(); +while (performance.now() - start < targetMs) { + reimFFT({ re, im }); + iterations++; +} +console.timeEnd('reimFFT'); + +console.log(`${iterations} iterations, ${size} points each`); +console.log( + `${((performance.now() - start) / iterations).toFixed(2)} ms per iteration`, +); diff --git a/eslint.config.js b/eslint.config.js index 00cc43d2..f7512d55 100644 --- a/eslint.config.js +++ b/eslint.config.js @@ -1,4 +1,7 @@ import { defineConfig, globalIgnores } from 'eslint/config'; import cheminfo from 'eslint-config-cheminfo-typescript'; -export default defineConfig(globalIgnores(['coverage', 'lib']), cheminfo); +export default defineConfig( + globalIgnores(['coverage', 'lib', 'benchmark']), + cheminfo, +); diff --git a/src/reim/__tests__/reimFFT.test.ts b/src/reim/__tests__/reimFFT.test.ts index 0af28fd4..346247ad 100644 --- a/src/reim/__tests__/reimFFT.test.ts +++ b/src/reim/__tests__/reimFFT.test.ts @@ -12,4 +12,17 @@ test('reimFFT', () => { }); expect(inverse.re).toStrictEqual(re); + expect(inverse.im).toStrictEqual(im); + // check pointer are different + expect(inverse.re).not.toBe(re); + expect(inverse.im).not.toBe(im); + + const transformed2 = reimFFT({ re, im }, { inPlace: true }); + const inverse2 = reimFFT(transformed2, { inverse: true, inPlace: true }); + + expect(inverse2.re).toStrictEqual(re); + expect(inverse2.im).toStrictEqual(im); + // check pointer are the same + expect(inverse2.re).toBe(re); + expect(inverse2.im).toBe(im); }); diff --git a/src/reim/reimFFT.ts b/src/reim/reimFFT.ts index f1cd6c7b..2bcedbee 100644 --- a/src/reim/reimFFT.ts +++ b/src/reim/reimFFT.ts @@ -3,9 +3,13 @@ import FFT from 'fft.js'; import type { DataReIm } from '../types/index.ts'; import { xRotate } from '../x/index.ts'; +let output: Float64Array | undefined; +let complexArray: Float64Array | undefined; +let fft: FFT | undefined; export interface ReimFFTOptions { inverse?: boolean; applyZeroShift?: boolean; + inPlace?: boolean; } /** @@ -18,20 +22,21 @@ export function reimFFT( data: DataReIm, options: ReimFFTOptions = {}, ): DataReIm { - const { inverse = false, applyZeroShift = false } = options; - const { re, im } = data; const size = re.length; const csize = size << 1; - let complexArray = new Float64Array(csize); + const { inPlace = false, inverse = false, applyZeroShift = false } = options; + + if (fft?.size !== size) fft = new FFT(size); + + if (complexArray?.length !== csize) complexArray = new Float64Array(csize); for (let i = 0; i < csize; i += 2) { complexArray[i] = re[i >>> 1]; complexArray[i + 1] = im[i >>> 1]; } - const fft = new FFT(size); - let output = new Float64Array(csize); + if (output?.length !== csize) output = new Float64Array(csize); if (inverse) { if (applyZeroShift) complexArray = zeroShift(complexArray, true); fft.inverseTransform(output, complexArray); @@ -40,14 +45,21 @@ export function reimFFT( if (applyZeroShift) output = zeroShift(output); } - const newRe = new Float64Array(size); - const newIm = new Float64Array(size); - for (let i = 0; i < csize; i += 2) { - newRe[i >>> 1] = output[i]; - newIm[i >>> 1] = output[i + 1]; + if (inPlace) { + for (let i = 0; i < csize; i += 2) { + re[i >>> 1] = output[i]; + im[i >>> 1] = output[i + 1]; + } + return data as DataReIm; + } else { + const newRe = new Float64Array(size); + const newIm = new Float64Array(size); + for (let i = 0; i < csize; i += 2) { + newRe[i >>> 1] = output[i]; + newIm[i >>> 1] = output[i + 1]; + } + return { re: newRe, im: newIm }; } - - return { re: newRe, im: newIm }; } function zeroShift(