FFT functions

[sebpiq]

I don't know if that stays in the scope of mathjs.
If you're interested I can implement them.

[josdejong]

Somewhere there will be a border between "basic" functions included in the math.js library, and "advanced" functions located in separate math.js modules. I'm not yet sure where this line will be exactly, but I think this will become clear over time. Let's just start with adding these functions (we definitely need them!), and when some section grows too big/advanced, move the functionality to a separate module.

[sebpiq]

Yes, definitely. Once again taking example on numpy, there is only very basic stuff there, and advanced functionality is provided in scipy which is just an extension for numpy.

But FFT is extremely basic stuff. In fact it is in numpy.

I can take care of that.

[josdejong]

Nice to have the road paved already by applications like numpy, so we don't have to invent everything ourselves :)

[drom]

Nice Library. I desperately need FFT in it. What is the progress with implementation? Can I help?

[josdejong]

FFT is still not implemented, sorry. Help would be welcome here. For now you could use the FFT of for example http://numericjs.com/, you can import functions from this library into math.js if your like.

[drom]

FFT/iFFT would need to operate on the real/complex vector.
We probably can start with input:Float64Array -> output:Float64Array function.
Do you have recommendation how to store Float64Array array of complex numbers, interleaved / half-half?

[josdejong]

I would start with a basic fft implementation using regular numbers and a regular Array (and with support for the Matrix type of math.js, but that will probably be a wrapper). For ifft I would start with an implementation using the Complex type of math.js, and using a regular Array and with support for Matrix.

For a second iteration we could optimize performance by using typed arrays. However typed arrays are not (yet) supported by all browsers so the library must be able to fall back on regular Arrays. Support for typed arrays is a needed for all functions of math.js, not just ftt. The changes currently being made for v2 will pave the way for library-wide support for typed arrays, resulting in way better performance.

One future matrix type could be a ComplexMatrix. This matrix could store a series of complex numbers in two typed arrays (one for im, one for re), which I expect to best performance.

[drom]

I have started small repo here: https://github.com/drom/fourier
Testing with Travis encountered float differences between node 0.12 / iojs and node 0.10
How are you dealing with floating point differences in testing of Math operations?

[josdejong]

To prevent round-off errors from breaking on round-off errors, we have a little utility approx.js, which compares whether two numbers are approximately equal.

https://github.com/josdejong/mathjs/blob/master/tools/approx.js

[drom]

I have done KISS implementation of dft / idft. Will do the radix 2 / 4 / 8 / 16 next. For the accuracy will follow the following methodology http://www.fftw.org/accuracy/method.html in order to compare with other FFTs.

[josdejong]

sounds good :)

[drom]

I have implemented power-of-2 radix-2 decimation in time version with normal arrays and typed array.
Performance and accuracy different between node 0.10, 0.12 and iojs.
https://travis-ci.org/drom/fourier
Surprised to see only 10-15% performance gain with typed arrays.

[josdejong]

Interesting to see what typed arrays do for the performance of such an algorithm. If it's just 10-15% it may not even be worth all extra work. Could there be more asm.js related techniques allowing the JS engine to better optimize this? (or is the current implementation already close to what you could get from a bare-bone C implementation?...)

[drom]

I have opened relevant issue here: drom/fourier#1
Any feedback is welcome!

[josdejong]

Thanks. So if I understand it correctly:

• The speed improvements for using typed arrays are actually significant, but the differences in performance differs a lot on different JavaScript engines.
• The JavaScript implementation performs way slower than other implementations like in C

[drom]

I have done some more experiments with FFT code. Results are here: drom/fourier#1

[echo66]

@josdejong take a look at the emscripten compilation of FFTW that I upload to this repository: https://github.com/echo66/FFTW3-emscripten/

[drom]

It does not look complete. How can one use it? I have implemented asm.js version of my FFT. You can see the results.

[echo66]

@drom , look at this file: https://github.com/echo66/FFTW3-emscripten/blob/master/test_code.js

[drom]

@echo66 As i mentioned above, it does not look like FFT implementation. How do you use this code?
BTW. It looks like use asm is good for Firefox and bad for V8.
drom/fourier#1
https://jsperf.com/fft/2

[echo66]

Oh, I just noticed that I forgot to upload the compiled files.

Done.

Just open the simple_1d_wrapper.html and see the results.

[echo66]

The initial delay is due to the loading of the "simple_1d_wrapper.js" file (1.4 MB).

[drom]

@echo66 it is better now. But, how dose one simply use it? Can you benchmark it for 64K with jsperf, or benchmark.js ?

[echo66]

http://jsperf.com/benchmark-for-fftw-js

[echo66]

A comparison between several FFT solutions: http://jsperf.com/comparing-two-fft-functions/3

I will add @drom FFT code to it. @drom, can you upload a "browserified" version of your library?

Btw, It is kind of ironic to see that chrome is faster than firefox when dealing with ASM code. x)

[echo66]

It is kind of ironic to see that chrome is faster than firefox when dealing with ASM code. x)

[drom]

@echo66 I have added browserified version of my FFT. And updated your benchmark to use it.
http://jsperf.com/comparing-two-fft-functions/5
p.s. Comparing bigger FFTs (like 64K), make more sense, IMHO.

[echo66]

An exaustive benchmark for FFTWJS and drom's FFT:
http://jsperf.com/benchmark-for-fftw-js-and-drom-s-fft/2

Can you check if I'm using fourier.js in the right way?

@drom, can you describe the details of your implementation, regarding optimization and code templates? I'm really interested in knowing more about those.

[drom]

@echo66 your benchmark look fine to me.
My current implementation uses:

• classic Radix-2;
• Decimation in Time (DIT) algorithm;
• in-place, (performs data permutations before first stage);
• one typed data array, complex data, all real point, then all imaginary points.
• generator: https://github.com/drom/fourier/blob/master/bin/fft-gen.js generates individual code for each FFT size using lodash template: https://github.com/drom/fourier/blob/master/src/fft-template.js (can be done dynamically if needed)
• I use asm.js conventions on function interface.
• all twiddle factors (1.25 of SIN period) calculated upfront ( .init() function ) and stored in the twiddle array.
• in both asm and raw versions I perform i32 casting with (x | 0) hint, u32 casting with (x >>> 0)
• in asm version I perform f32 casting using fround(x) hint, and f64 casting using +(x) hint.
• asm version has 'use asm'; pragma and passes Firefox AOT validation (as it can be seeing in the log).

More details are here: drom/fourier#1

[echo66]

If you create a multidimensional version of fourier.js, I see no reason to use FFTWJS.

[drom]

@echo66 I have never done Nd FFT in my life. Can you suggest good description? reference code? or pull request? ;)

[echo66]

Sorry, I never implemented a FFT before. :/

Still, you have this implementation: https://github.com/scijs/ndarray-fft

[echo66]

Note 1: for each execution of forward and inverse in the FFTW wrapper, I'm creating a new FFTW plan. And it is always the same type of plan. I probably should take a look at it.

Note 2: for all browsers in the benchmark, the FFTWJS has a stable, but sub-optimal, performance.

Note 3: I didn't include the other FFT implementations because they presented some memory problems in previous tests, when dealing with "big" vectors (32768, 65536).

[mrquincle]

Hi Jos,

Created a simple FFT implementation with radix-2. This is a decimation-in-time version. It is out-of-place, memory/space complexity is not optimized. Implementation details are common knowledge, for example on https://en.wikipedia.org/wiki/Cooley%E2%80%93Tukey_FFT_algorithm

https://gist.github.com/mrquincle/b11fff96209c9d1396b0

I haven't visualized it, but this should be nice to combine with vis.js.

Cheers,

Anne

[josdejong]

Thanks Anne, that looks neat!

Should be quite straight forward to add this code to math.js, though we have to add docs and unit tests of course. Would you like to pick that up too?

[mrquincle]

Sure, remind me if I forget.

[josdejong]

great, thanks!

[arkajitmandal]

Hi is this implemented yet?

[josdejong]

No this hasn't been implemented. Anyone interested in picking this up?

[arkajitmandal]

I can do this once I am done the eig functions

[josdejong]

Great, thanks for your offer Arkajit!

[mklilley]

I would like to add a +1 for adding this feature. I'm a physicist trying to create an ObservableHQ notebook to solve the Schrödinger equation. Currently I'm using numjs because it's got FFT, but it doesn't have complex arithmetic so I'm having the create that functionality myself....and probably not doing a good job at it.

[cshaa]

A cool comparison of different implementations of FFT in JavaScript: https://github.com/scijs/fourier-transform/blob/master/benchmark.md

[HanchaiN]

Which category should the function be implemented in? (Just in case someone will actually implement it.)

Also, the multidimensional Fourier transform seems to be a good generalization as the module normally works with multidimensional arrays/matrices. (It may or may not affect the performance to implement this.)

[josdejong]

Which category should the function be implemented in?

I think we can add fft to the matrix category. If you have an other suggestion just let me know.

Also, the multidimensional Fourier transform seems to be a good generalization as the module normally works with multidimensional arrays/matrices. (It may or may not affect the performance to implement this.)

That sounds good. If there is performance impact I suppose we could have two implementations under the hood (1 dimensional and multi dimensional), and automatically pick the right one.

[HanchaiN]

As for now, the function can only operate on a power-of-two--size array. It still needs to be fixed to support the input of arbitrary sizes. (Either open a new issue or reopen this one would be fine.)

[josdejong]

Yes that is true. Can you open a new issue about fft support for arbitrary size matrices?