Avoid overflow in good_size... functions

[crisluengo]

This should fix issue #15 in the manner discussed there.

I tried to follow the existing code style, but feel free to reformat or edit in any way.

In my own version of these functions, I combined the code for _cmplx and _real by making the maximum factor a template parameter as well. That reduces code duplication a bit. I'm happy to do that here too if you like.

This was my test program, for posterity:

#include <iostream>
#include "pocketfft_hdronly.h"

int main() {
    std::size_t res;
    res = pocketfft::detail::util::good_size_real(10);
    std::cout << "good_size_real(10): " << res << ", correct: " << (10 == res) << '\n';
    res = pocketfft::detail::util::good_size_real(11);
    std::cout << "good_size_real(11): " << res << ", correct: " << (12 == res) << '\n';
    res = pocketfft::detail::util::good_size_real(13);
    std::cout << "good_size_real(13): " << res << ", correct: " << (15 == res) << '\n';
    res = pocketfft::detail::util::good_size_real(101);
    std::cout << "good_size_real(101): " << res << ", correct: " << (108 == res) << '\n';
    res = pocketfft::detail::util::good_size_real(2109375001u);
    std::cout << "good_size_real(2109375001): " << res << ", correct: " << (2123366400u == res) << '\n';
    res = pocketfft::detail::util::good_size_real(4294967295u/5/2); // does the work in the 32-bit algorithm
    std::cout << "good_size_real(4294967295/5/2): " << res << ", correct: " << (429981696 == res) << '\n';
    res = pocketfft::detail::util::good_size_real(4294967295u/5/2+1); // should be elevated to the 64-bit algorithm
    std::cout << "good_size_real(4294967295/5/2+1): " << res << ", correct: " << (429981696 == res) << '\n';
    // Throws in the 32-bit case:
    try {
        res = pocketfft::detail::util::good_size_real(4294967295);
        std::cout << "good_size_real(4294967295): " << res << ", correct: " << (4294967296 == res) << '\n';
    } catch (std::runtime_error const& err) {
        std::cout << "good_size_real(4294967295) throws.\n";
    }

    res = pocketfft::detail::util::good_size_cmplx(10);
    std::cout << "good_size_cmplx(10): " << res << ", correct: " << (10 == res) << '\n';
    res = pocketfft::detail::util::good_size_cmplx(11);
    std::cout << "good_size_cmplx(11): " << res << ", correct: " << (11 == res) << '\n';
    res = pocketfft::detail::util::good_size_cmplx(13);
    std::cout << "good_size_cmplx(13): " << res << ", correct: " << (14 == res) << '\n';
    res = pocketfft::detail::util::good_size_cmplx(101);
    std::cout << "good_size_cmplx(101): " << res << ", correct: " << (105 == res) << '\n';
    res = pocketfft::detail::util::good_size_cmplx(2109375001u);
    std::cout << "good_size_cmplx(2109375001): " << res << ", correct: " << (2112000000u == res) << '\n';
    res = pocketfft::detail::util::good_size_cmplx(4294967295u/11/2); // does the work in the 32-bit algorithm
    std::cout << "good_size_cmplx(4294967295/11/2): " << res << ", correct: " << (195230112 == res) << '\n';
    res = pocketfft::detail::util::good_size_cmplx(4294967295u/11/2+1); // should be elevated to the 64-bit algorithm
    std::cout << "good_size_cmplx(4294967295/11/2+1): " << res << ", correct: " << (195230112 == res) << '\n';
    // Throws in the 32-bit case:
    try {
        res = pocketfft::detail::util::good_size_cmplx(4294967295);
        std::cout << "good_size_cmplx(4294967295): " << res << ", correct: " << (4294967296 == res) << '\n';
    } catch (std::runtime_error const& err) {
        std::cout << "good_size_real(4294967295) throws.\n";
    }

    res = pocketfft::detail::util::prev_good_size_real(10);
    std::cout << "prev_good_size_real(10): " << res << ", correct: " << (10 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_real(11);
    std::cout << "prev_good_size_real(11): " << res << ", correct: " << (10 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_real(13);
    std::cout << "prev_good_size_real(13): " << res << ", correct: " << (12 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_real(107);
    std::cout << "prev_good_size_real(107): " << res << ", correct: " << (100 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_real(2123366399u);
    std::cout << "prev_good_size_real(2123366399): " << res << ", correct: " << (2109375000u == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_real(4294967295u/5); // does the work in the 32-bit algorithm
    std::cout << "prev_good_size_real(4294967295/5): " << res << ", correct: " << (854296875 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_real(4294967295u/5+1); // should be elevated to the 64-bit algorithm
    std::cout << "prev_good_size_real(4294967295/5+1): " << res << ", correct: " << (854296875 == res) << '\n';

    res = pocketfft::detail::util::prev_good_size_cmplx(10);
    std::cout << "prev_good_size_cmplx(10): " << res << ", correct: " << (10 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(11);
    std::cout << "prev_good_size_cmplx(11): " << res << ", correct: " << (11 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(13);
    std::cout << "prev_good_size_cmplx(13): " << res << ", correct: " << (12 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(107);
    std::cout << "prev_good_size_cmplx(107): " << res << ", correct: " << (105 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(2123366399u);
    std::cout << "prev_good_size_cmplx(2123366399): " << res << ", correct: " << (2122312500u == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(4294967295u/11); // does the work in the 32-bit algorithm
    std::cout << "prev_good_size_cmplx(4294967295/11): " << res << ", correct: " << (390297600 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(4294967295u/11+1); // should be elevated to the 64-bit algorithm
    std::cout << "prev_good_size_cmplx(4294967295/11+1): " << res << ", correct: " << (390297600 == res) << '\n';
    res = pocketfft::detail::util::prev_good_size_cmplx(4294967295u);
    std::cout << "prev_good_size_cmplx(4294967295): " << res << ", correct: " << (4293273600u == res) << '\n';
}

[crisluengo]

Also, I would have put the templated functions in a namespace, but util is a struct, not a namespace, and all of this is already in the detail namespace, so I just made the function name a bit longer... :)

[mreineck]

What do you think about #17? This is yuor commit plus some tweaks that allow us to go back to the original naming.

In your test code I changed the signed integers to unsigned ones, and then everything works again.

[crisluengo]

I guess it's not meant to be called by the user? If it's only an internal function, then this is fine. I had made the stub functions so that the input is automatically cast to unsigned int, which is more user-friendly.

[mreineck]

It's no problem to suport signed integers as well, but then we'd have (strictly speaking) to check for negative input to be safe, and I thought we could perhaps avoid this kind of silliness already by using an appropriate input type.

In your example, if you do

int good_size = good_size_cmplx((1u<<31)+5)

you'll end up with a nasty surprise as well, I think, so it doesn't save you from being careful with data types.

[crisluengo]

As long as the function is called with something like vector.size() as input, none of this really matters. I don't know how this function is used, but I would expect it's never with some hard-coded constant.

But be aware that, on a 64-bit system, good_size_cmplx(4294967295u) throws (because the result doesn't fit in an unsigned int), whereas good_size_cmplx(4294967295ull) doesn't.

[mreineck]

Sorry for taking so long with this ...

I think your approch is the most reasonable after all; I only made a few pedantic changes (instead of comparing the size of the type to 4, I check whether it is smaller than uint64?t).

I've tried to push to your branch directly as suggested by the PR, but I couldn't figure out how this works, so I updated #17 accordingly.

[crisluengo]

No worries. I'll just close this PR, you commit the fix in #17 however you think is best.

By the way, thanks for all the work you put into this library, it's fantastic!