Closed martijnende closed 5 years ago
I implemented option 1 in the matlab wrapper, only for MESHDIM=1. For MESHDIM=2, this needs to be done in the fortran code. I also modified the user's manual (the google doc) accordingly.
I implemented option 1 also in the fortran code. Option 2 does not work with our periodic kernels expressed analytically in spectral domain (wavenumber). Hence we can close this issue.
Thanks Pablo. I've had a look at the error check you implemented in fftsg.f90
, and I have a few remarks:
my_rdft
is called. Although this is by far the safest way to go, I propose to take a little more risk in favour of performance and do the sanity check only once at the initiation stage (i.e. during/after reading the input file).(n .AND. (n-1)) == 0
(an equivalent expression works in C++ and Python, never tried it in Fortran). This expression only requires 2 clock cycles or so, as opposed to several tens of cycles required to calculate the next power of 2.Have you tested your code changes before you pushed them to master
?
n>1 .and. popcnt(n)==1
. I implemented that too, but commented it out in case some users don't have a Fortran 2008 compiler (I have seen that happen in another software project).I didn't test the code. I am editing in a browser. No compiler at hand.
I would avoid using bit manipulation, the 2^n check is called only once each model run thus not a big burden for computation
All right, I think Pablo's implementation is satisfactory, so I'll close this issue now.
As @gpercy and @jpampuero mentioned to me, the current FFT routine requires a vector of size
2^N
, though there are no internal checks for this requirement, possibly resulting in corrupt/inaccurate data. Possible workarounds:The latter is more user-friendly, while the former is computationally more efficient.