Re: FFT proposal for gsoc 2021

From: Janek Kozicki (janek.listy.mailowe_at_[hidden])
Date: 2021-04-07 13:57:24

Eduardo Quintana via Boost said:     (by the date of Fri, 02 Apr 2021 08:18:07 +0000)

> the Gsoc proposal. That you can find here:
> https://github.com/Lagrang3/gsoc-2021/releases/download/v2.0/gsoc.pdf

Hi, I only have a couple of remarks after reading the proposal.

> Also dft can be computed for a limited number of real number representation types: float and double.

Not true. fftw supports also long double and float128, but not boost::mpfr or boost::cpp_bin_float.
See e.g. debian packages: libfftw3-long3 , libfftw3-quad3. And a super short example:

#include <iostream>
#include <boost/cstdfloat.hpp>
#include <fftw3.h>

int main() {
        using Real = boost::float128_t;
        using Complex = std::complex<Real>;

        std::vector<Complex> in {{1,-1}, { 0 , 0 }, { 1  , 0.5} , { 0  , 0.5 } , {0.3Q, 1.1Q} };
        std::vector<Complex> out{{0, 0}, { 0 , 0 }, { 0  , 0  } , { 0  , 0   } , { 0  ,  0  } };
        size_t N=in.size();
        fftwq_plan p;
        p = fftwq_plan_dft_1d(N , (Real (*)[2])(in.data()) , (Real (*)[2])(out.data()) , FFTW_FORWARD, FFTW_ESTIMATE);
        fftwq_execute(p);
        fftwq_destroy_plan(p);

        int digs = std::numeric_limits<Real>::digits10;
        std::cout << "precision digits=" << digs << "\n";
        std::cout << std::setprecision(digs+1);
        for(size_t i = 0 ; i < N ; ++i) { std::cout << in [i] << "\n"; }
        std::cout << "---\n";
        for(size_t i = 0 ; i < N ; ++i) { std::cout << (out[i]/std::sqrt(Real(N))) << "\n"; }
        std::cout << "---\n";
        std::cout << std::sqrt(Real(3.Q)) << "\n";

}

You need to link -lfftw3q. Same with -lfftw3l for long double. Just
replace s/fftwq/fftwl/gc in function names. And using Real = long double.

> double[2], which lacks of the algebraic semantics that
> std::complex<double> has;

not true. The C++ standard requires that any std::complex<T> number
can be reinterpret_cast down to an array[2] of that type, see:
https://en.cppreference.com/w/cpp/numeric/complex
and https://isocpp.org/files/papers/n4296.pdf paragraph 26.4 page 926

> In astrophysics fft is used to solve the differential field equations

It would be good to make sure that your API design is also compatible with the
native boost library for solving ordinary differential equations:

  https://www.boost.org/doc/libs/1_75_0/libs/numeric/odeint/doc/html/index.html

> 3.2	Future outlook
> The D-dimensional mpi dft is an issue we look forward to tackle
> using a scalable (D − 1)-dimensional domain decomposition

Since the main focus will be on API (a good choice), please make sure
that your multi-dimensional API is compatible with
boost::ublas::tensor. Regular users will want to use other
ublas::tensor capabilities. Maybe even some of them, such as
dimensional slicing (you discuss it in details in the proposal), will
be useful for you.

> https://github.com/Lagrang3/P2.16_Parallel_FFT/tree/master/solution

this link does not work.

> 4 Timeline

It is very nicely detailed. It should be easy to track progress with
it. I am still worried that there is a little bit too much in it. But
we will see. I'm sure it can be adjusted during run-time :)

best regards
Janek

--
Janek Kozicki, PhD. DSc. Arch. Assoc. Prof.
Gdańsk University of Technology
Faculty of Applied Physics and Mathematics
Department of Theoretical Physics and Quantum Information
--
http://yade-dem.org/
http://pg.edu.pl/jkozicki (click English flag on top right)