@zoziha has done a great job modernizing it by adding fpm support, CI, modern interfaces, etc. So it would be at least he and I maintaining it. And we are hoping others would join also.

Fast Fourir Transform (FFT) is a well-known and practical algorithm. fftpack is the native fft implementation of Fortran stored in netlib. I believe we have the need to maintain and use it. My original intention was to make fftpack easier-to-use when fortran-lang was rebuilding the Fortran community. We have some APIs we can refer to:

@certik I still have 2 tasks waiting to be submitted. Please wait until I finish these 2 tasks today or tomorrow before we migrate certik/fftpack to fortran-lang/fftpack.

I also have a modern Fortran reimplementation of FFTPACK here: hfsolver/fourier.f90 at master · certik/hfsolver · GitHub, but that should go into a separate package. For fortran-lang/fftpack I would like to preserve the original algorithm, for benchmarking and reference.

Yes please! I haven’t actually used the version at certik/fftpack yet, but I’ve looked at the interface (as explained in the new spec files) and I like what I see. I am definitely in favour of bringing this into fortran-lang, and while I’m rather overloaded right now, in principle I’m willing to help maintain it in the future.

You are right, the kind parameters in stdlib are one of the two modules which are currently undocumented in the internal API documentation. It’s a low-hanging fruit for a new contributor to fix or a straight-forward task for the existing developers to take care of, if nobody takes this up, I’ll send a patch by end of the week.

Ah, now I see the issue. Thank you. So if we generate the list of kinds at build time, and then use that list to generate our generic implementations of stdlib functions, how would you propose handling symbolic shorthands like dp, etc? Should we just omit those for the purpose of stdlib, and assume users of the library will define their own based on whatever kinds they need?

And just out of curiosity, the possible five real kinds you’re referring to… are are they the common 32/64/128 bit kinds, plus 16 bit and 80 bit?

Thanks @kargl for the feedback. I think what you are suggesting makes sense. Please keep the messages up for a while until we document your suggestion.

@zoziha, anything else you want to do before we move fftpack under fortran-lang?

2DECOMP&FFT: scalable library for Fourier transforms, derived from 2decomp&fft, by BenMql. Compatible with in-core Chebyshev transforms.

ffte: computes Discrete Fourier Transforms of 1-, 2- and 3- dimensional sequences of length (2^p) (3^q) (5^r).

fftlog-f90: extended Fortran 90 version of the FFTLog code by Andrew Hamilton to convolve a tabulated function with a Bessel functions, by coccoinomane

fft-overlap: efficient implementations of ffts on multiple GPUs and across multiple nodes, by dappelha. Overlapping data transfer on multiple levels.

fftpack: double precision version of original fftpack, from fortran-lang

FFTPack: easily usable package of functions using wrapping the Fortran 77 FFTPack library, by keurfonluu

modern_fftpack: Fortran 2008 library of fast Fourier transforms – an object-oriented modernization of NCAR’s FFTPACK5.1., by jlokimlin

nufft: code for the 1D, 2D and 3D non-uniform FFTs, from biotrump

Thank you both for your work on this! I’ve already downloaded it, built it with FPM, and played around with the modern interfaces. Quick question (probably for @zoziha) - the documentation says "The contents of wsave must not be changed between calls of dfftf or dfftb" - can I re-use wsave in between multiple forward-backward transforms, or do I need to call dffti before every forward transform?

I may not understand the specific details. I just transferred the instructions in the doc file from netlib/fftpack to fftpack.md. According to the instructions in the doc file, it seems best to call again, of course you can try it boldly.

I’ve not studied FYPP but this preprocessor appears to be employed actively in stdlib. With FYPP, is it possible somehow to make use of REAL_KINDS array named constant in ISO_FORTRAN_ENV in the code implementations of algorithms (and containers) in stdlib, rather than the current use of REAL32, REAL64, and REAL128, none of which are guaranteed to be supported by a conforming Fortran processor.

Note with REAL_KINDS constant, a Fortran standard conforming processor can be expected to have its size be at least 2, often 3 (e.g., IFORT), a few others like gfortran and NagFor to be 4, and in the near future, some processors may have it of be of size 5 or more. As you all know, the precise size of this array and the values of its elements would not matter to a truly generic algorithm. And this is something I hope will taken into close consideration during the Generics development with Fortran 202Y.

In the mean time, it will be good if there is a FYPP approach that can be used with this REAL_KINDS constant. For example, consider a (hypothetical and simple-minded) module toward a generic function to compute the statistical mean of rank-1 arrays:

module xx_stats_rank1_mean
! Module for a generic algorithm; should NOT have any dependence on REAL32, REAL64, etc.
use, intrinsic :: iso_fortran_env, only : rk => real_kinds
generic :: rank1_mean => rank1_mean_1, rank1_mean_2, .. rank1_mean_n ! n is the size of real_kinds
contains
function rank1_mean_1( x ) result(mean)
real(rk(1)), intent(in) :: x(:)
real(rk(1)) :: mean
mean = sum(x) / real(size(x), rk(1))
end function rank1_mean_1
function rank1_mean_2( x ) result(mean)
real(rk(2)), intent(in) :: x(:)
real(rk(2)) :: mean
mean = sum(x) / real(size(x), rk(2))
end function rank1_mean_2
..
function rank1_mean_n( x ) result(mean)
real(rk(n)), intent(in) :: x(:)
real(rk(n)) :: mean
mean = sum(x) / real(size(x), rk(n))
end function rank1_mean_n
end module

Can the macro language be setup using FYPP to expand out the generic function and setup the generic interface based simply on the size of REAL_KINDS constant array?

Oh, yes, the current fortran-lang/fftpack is based on netlib/fftpack4.0 (public domain).
We also have a nacr/fftpack5.1 (FFTPACK license, look like MIT and BSD-3, see comment and FFTPACK5 (the “Software”) License), which natively supports multi-dimensional FFT, but gfortran compilation requires --fallow-argument-mismatch (see Build failure with gfortran 10 · Issue #123 · NCAR/ncl · GitHub) (This issue will be solved by that fpm supports native flag settings, see my discourse help).

My initial idea was to develop certik(netlib)/fftpack4.0 and nacr/fftpack5.1, first do fftpack4.0 and then fftpack5.1. I do not involve john/fftpack5.1 by default, it uses the GPL license, which may be less friendly to users.

But in fact we better maintain only one fortran-lang/fftpack? My question is: do we need to develop based on nacr/fftpack5.1? Or let’s leave it to time, let’s finish fftpack4.0 first.

My thoughts:
First, I will try my best to improve fftpack4.0, and when it has modern features, I am willing to see whether it is necessary to maintain fftpack5.1 (may consider the needs of multi-dimensional fft), if necessary, I will hope to achieve the modern features of fftpack5.1 in another branch (I believe there is fftpack4.0 experience, this work is not difficult, it is just a matter of time). In general, the two fftpacks try to keep the same modern interface.