Implement `lcm` and `gcd` from library fundamentals V2

[mclow.lists]

Summary

These are probably going to end up in C++17 as well.

This is a straight addition; no changes to existing code.

[howard.hinnant]

I created a top-level gcd which did nothing but make everything unsigned and do the abs once, and then called a __gcd specialized for unsigned and only one kind of unsigned, and got measurably faster results (about 10%).

template<class _Tp>
constexpr
std::enable_if_t<is_unsigned<_Tp>{}, _Tp> _LIBCPP_INLINE_VISIBILITY
__gcd(_Tp __m, _Tp __n)
{
	return __n == 0 ? __m : __gcd(__n, __m % __n);
}

template<class _Tp, class _Up>
constexpr
common_type_t<_Tp,_Up> _LIBCPP_INLINE_VISIBILITY
gcd(_Tp __m, _Up __n)
{
	static_assert((is_integral<_Tp>::value && is_integral<_Up>::value), "Arguments to gcd must be integer types");
	using _Rp = common_type_t<_Tp,_Up>;
	using _Wp = make_unsigned_t<_Rp>;
	return static_cast<_Rp>(__gcd(static_cast<_Wp>(__abs<_Tp>()(__m)),
	                              static_cast<_Wp>(__abs<_Up>()(__n))));
}

Here is the test driver I used:

#include <chrono>
#include <iostream>

int
main()
{
    std::uint64_t x = 0;
    auto t0 = std::chrono::steady_clock::now();
    for (auto i = 1'999'990'000; i <= 2'000'000'000; ++i)
    {
        for (auto j = 1'999'990'000; j <= 2'000'000'000; ++j)
        {
            x += gcd(i, j);
        }
    }
    auto t1 = std::chrono::steady_clock::now();
    std::cout << std::chrono::duration<double>(t1-t0).count() << '\n';
    return x;
}

I also tried the traditional iterative solution (as opposed to recursive) and was surprised that it was not faster. I then tried "binary" gcd (https://en.wikipedia.org/wiki/Binary_GCD_algorithm) and was again surprised it wasn't faster. Both of these alternatives were in the ball park though, and it could be that they might measure faster over other domains. The binary variant used __m >>= std::__ctz(__m); instead of a while loop to remove factors of 2 (in two places).

[EricWF]

include/experimental/numeric

45

This should be _LIBCPP_BEGIN_NAMESPACE_LFTS.

47

_VSTD:: is redundant here.

62

It might be beneficial to calculate the absolute values during the initial call instead of at every level of recursion, since most of those calls will be redundant.

80

Copy/Paste error

test/std/experimental/numeric/numeric.ops.overview/nothing_to_do.pass.cpp

10

Since it includes the header it should have: // UNSUPPORTED: c++98, c++03, c++11.

test/std/experimental/numeric/numeric.ops/nothing_to_do.pass.cpp

10

Since it includes the header it should have: // UNSUPPORTED: c++98, c++03, c++11.

[howard.hinnant]

Can't call std::abs as it isn't constexpr. You'll have to roll your own.

[EricWF]

Gosh darn it Howard. Your like 10 minutes ahead of me.

Here is a constexpr gcd test.
https://gist.github.com/EricWF/46ea4489f405ec9f83325a278fe99535

include/experimental/numeric

51

Both the call operator and abs are not constexpr.

56

Call operator is not constexpr.

[mclow.lists]

Updated based on comments here and elsewhere.

• Do the abs once, rather than at every level of recursion
• Add tests for bool.
• Constexpr

Also, add an assertion to detect overflow in lcm.

I have not yet incorporated Eric's constexpr tests.

[mclow.lists]

Updated with Eric's constexpr tests.
Once this is approved, it will be applied twice - once in <experimental/numeric>, and once in <numeric>

[EricWF]

LGTM.

[Eugene.Zelenko]

Committed in r276750.

include/experimental/numeric

+// -*- C++ -*-
+//===--------------------------- numeric ----------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef _LIBCPP_EXPERIMENTAL_NUMERIC
+#define _LIBCPP_EXPERIMENTAL_NUMERIC
+/*
+    experimental/numeric synopsis
+
+// C++1z
+namespace std {
+namespace experimental {
+inline namespace fundamentals_v2 {
+
+  // 13.1.2, Greatest common divisor
+  template<class M, class N>
+  constexpr common_type_t<M,N> gcd(M m, N n);
+
+  // 13.1.3, Least common multiple
+  template<class M, class N>
+  constexpr common_type_t<M,N> lcm(M m, N n);
+
+} // namespace fundamentals_v2
+} // namespace experimental
+} // namespace std
+
+ */
+
+#include <experimental/__config>
+#include <numeric>
+#include <type_traits>              // is_integral
+#include <limits>                   // numeric_limits
+
+#if !defined(_LIBCPP_HAS_NO_PRAGMA_SYSTEM_HEADER)
+#pragma GCC system_header
+#endif
+
+#if _LIBCPP_STD_VER > 11
+
+_LIBCPP_BEGIN_NAMESPACE_LFTS
+
+template <typename _Tp, bool _IsSigned = is_signed<_Tp>::value> struct __abs;
+
+template <typename _Tp>
+struct __abs<_Tp, true> {
+    _LIBCPP_CONSTEXPR _LIBCPP_INLINE_VISIBILITY
+    _Tp operator()(_Tp __t) const noexcept { return __t >= 0 ? __t : -__t; }
+};
+
+template <typename _Tp>
+struct __abs<_Tp, false> {
+    _LIBCPP_CONSTEXPR _LIBCPP_INLINE_VISIBILITY
+    _Tp operator()(_Tp __t) const noexcept { return __t; }
+};
+
+
+template<class _Tp>
+_LIBCPP_CONSTEXPR _LIBCPP_INLINE_VISIBILITY
+_Tp __gcd(_Tp __m, _Tp __n)
+{
+    static_assert((!is_signed<_Tp>::value), "" );
+    return __n == 0 ? __m : __gcd<_Tp>(__n, __m % __n);
+}
+
+
+template<class _Tp, class _Up>
+_LIBCPP_CONSTEXPR _LIBCPP_INLINE_VISIBILITY
+common_type_t<_Tp,_Up>
+gcd(_Tp __m, _Up __n)
+{
+    static_assert((is_integral<_Tp>::value && is_integral<_Up>::value), "Arguments to gcd must be integer types");
+    using _Rp = common_type_t<_Tp,_Up>;
+    using _Wp = make_unsigned_t<_Rp>;
+    return static_cast<_Rp>(__gcd(static_cast<_Wp>(__abs<_Tp>()(__m)),
+                                  static_cast<_Wp>(__abs<_Up>()(__n))));
+}
+
+template<class _Tp, class _Up>
+_LIBCPP_CONSTEXPR _LIBCPP_INLINE_VISIBILITY
+common_type_t<_Tp,_Up>
+lcm(_Tp __m, _Up __n)
+{
+    static_assert((is_integral<_Tp>::value && is_integral<_Up>::value), "Arguments to lcm must be integer types");
+    if (__m == 0 || __n == 0)
+        return 0;
+
+    using _Rp = common_type_t<_Tp,_Up>;
+    _Rp __val1 = __abs<_Tp>()(__m) / gcd(__m,__n);
+    _Up __val2 = __abs<_Up>()(__n);
+    _LIBCPP_ASSERT((numeric_limits<_Rp>::max() / __val1 > __val2), "Overflow in lcm");
+    return __val1 * __val2;
+}
+
+_LIBCPP_END_NAMESPACE_LFTS
+
+#endif /* _LIBCPP_STD_VER > 11 */
+#endif /* _LIBCPP_EXPERIMENTAL_NUMERIC */

test/std/experimental/numeric/numeric.ops.overview/nothing_to_do.pass.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+//
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+#include <experimental/numeric>
+
+int main () {}

test/std/experimental/numeric/numeric.ops/nothing_to_do.pass.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+//
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+#include <experimental/numeric>
+
+int main () {}

test/std/experimental/numeric/numeric.ops/numeric.ops.gcd/gcd.not_integral1.fail.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+// 
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+// template<class _M, class _N>
+// constexpr common_type_t<_M,_N> gcd(_M __m, _N __n)
+
+// Remarks: If either M or N is not an integer type, the program is ill-formed.
+
+#include <experimental/numeric>
+
+
+int main()
+{
+    std::experimental::gcd(2.0, 4);
+}

test/std/experimental/numeric/numeric.ops/numeric.ops.gcd/gcd.not_integral2.fail.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+// 
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+// template<class _M, class _N>
+// constexpr common_type_t<_M,_N> gcd(_M __m, _N __n)
+
+// Remarks: If either M or N is not an integer type, the program is ill-formed.
+
+#include <experimental/numeric>
+
+
+int main()
+{
+    std::experimental::gcd(4, 6.0);
+}

test/std/experimental/numeric/numeric.ops/numeric.ops.gcd/gcd.pass.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+// 
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+// template<class _M, class _N>
+// constexpr common_type_t<_M,_N> gcd(_M __m, _N __n)
+
+#include <experimental/numeric>
+#include <cassert>
+#include <cstdlib>    // for rand()
+#include <iostream>
+
+constexpr struct {
+  int x;
+  int y;
+  int expect;
+} Cases[] = {
+    {0, 0, 0},
+    {1, 0, 1},
+    {0, 1, 1},
+    {1, 1, 1},
+    {2, 3, 1},
+    {2, 4, 2},
+    {36, 17, 1},
+    {36, 18, 18}
+};
+
+
+template <typename Input1, typename Input2, typename Output>
+constexpr bool test0(Input1 in1, Input2 in2, Output out)
+{
+    static_assert((std::is_same<Output, decltype(std::experimental::gcd(in1, in2))>::value), "" );
+    static_assert((std::is_same<Output, decltype(std::experimental::gcd(in2, in1))>::value), "" );
+    return out == std::experimental::gcd(in1, in2) ? true : (std::abort(), false);
+}
+
+
+template <typename Input1, typename Input2 = Input1>
+constexpr bool do_test(int dummy = 0)
+{
+    using S1 = typename std::make_signed<Input1>::type;
+    using S2 = typename std::make_signed<Input2>::type;
+    using U1 = typename std::make_unsigned<Input1>::type;
+    using U2 = typename std::make_unsigned<Input2>::type;
+    bool accumulate = true;
+    for (auto TC : Cases) {
+        { // Test with two signed types
+            using Output = std::common_type_t<S1, S2>;
+            accumulate &= test0<S1, S2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S1, S2, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<S1, S2, Output>(TC.x, -TC.y, TC.expect);
+            accumulate &= test0<S1, S2, Output>(-TC.x, -TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(TC.x, -TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(-TC.x, -TC.y, TC.expect);
+        }
+        { // test with two unsigned types
+            using Output = std::common_type_t<U1, U2>;
+            accumulate &= test0<U1, U2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<U2, U1, Output>(TC.x, TC.y, TC.expect);
+        }
+        { // Test with mixed signs
+            using Output = std::common_type_t<S1, U2>;
+            accumulate &= test0<S1, U2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<U2, S1, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S1, U2, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<U2, S1, Output>(TC.x, -TC.y, TC.expect);
+        }
+        { // Test with mixed signs
+            using Output = std::common_type_t<S2, U1>;
+            accumulate &= test0<S2, U1, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<U1, S2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S2, U1, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<U1, S2, Output>(TC.x, -TC.y, TC.expect);
+        }
+    }
+    return accumulate;
+}
+
+int main()
+{
+    auto non_cce = std::rand(); // a value that can't possibly be constexpr
+
+    static_assert(do_test<signed char>(), "");
+    static_assert(do_test<short>(), "");
+    static_assert(do_test<int>(), "");
+    static_assert(do_test<long>(), "");
+    static_assert(do_test<long long>(), "");
+
+    assert(do_test<signed char>(non_cce));
+    assert(do_test<short>(non_cce));
+    assert(do_test<int>(non_cce));
+    assert(do_test<long>(non_cce));
+    assert(do_test<long long>(non_cce));
+
+    static_assert(do_test< int8_t>(), "");
+    static_assert(do_test<int16_t>(), "");
+    static_assert(do_test<int32_t>(), "");
+    static_assert(do_test<int64_t>(), "");
+
+    assert(do_test< int8_t>(non_cce));
+    assert(do_test<int16_t>(non_cce));
+    assert(do_test<int32_t>(non_cce));
+    assert(do_test<int64_t>(non_cce));
+
+    static_assert(do_test<signed char, int>(), "");
+    static_assert(do_test<int, signed char>(), "");
+    static_assert(do_test<short, int>(), "");
+    static_assert(do_test<int, short>(), "");
+    static_assert(do_test<int, long>(), "");
+    static_assert(do_test<long, int>(), "");
+    static_assert(do_test<int, long long>(), "");
+    static_assert(do_test<long long, int>(), "");
+
+    assert((do_test<signed char, int>(non_cce)));
+    assert((do_test<int, signed char>(non_cce)));
+    assert((do_test<short, int>(non_cce)));
+    assert((do_test<int, short>(non_cce)));
+    assert((do_test<int, long>(non_cce)));
+    assert((do_test<long, int>(non_cce)));
+    assert((do_test<int, long long>(non_cce)));
+    assert((do_test<long long, int>(non_cce)));
+}

test/std/experimental/numeric/numeric.ops/numeric.ops.lcm/lcm.not_integral1.fail.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+//
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+// template<class _M, class _N>
+// constexpr common_type_t<_M,_N> lcm(_M __m, _N __n)
+
+// Remarks: If either M or N is not an integer type, the program is ill-formed.
+
+#include <experimental/numeric>
+
+
+int main()
+{
+    std::experimental::lcm(2.0, 4);
+}

test/std/experimental/numeric/numeric.ops/numeric.ops.lcm/lcm.not_integral2.fail.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+// 
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+// template<class _M, class _N>
+// constexpr common_type_t<_M,_N> lcm(_M __m, _N __n)
+
+// Remarks: If either M or N is not an integer type, the program is ill-formed.
+
+#include <experimental/numeric>
+
+
+int main()
+{
+    std::experimental::lcm(4, 6.0);
+}

test/std/experimental/numeric/numeric.ops/numeric.ops.lcm/lcm.pass.cpp

+//===----------------------------------------------------------------------===//
+//
+//                     The LLVM Compiler Infrastructure
+//
+// This file is dual licensed under the MIT and the University of Illinois Open
+// Source Licenses. See LICENSE.TXT for details.
+//
+//===----------------------------------------------------------------------===//
+// 
+// UNSUPPORTED: c++98, c++03, c++11
+// <numeric>
+
+// template<class _M, class _N>
+// constexpr common_type_t<_M,_N> gcd(_M __m, _N __n)
+
+#include <experimental/numeric>
+#include <cassert>
+#include <cstdlib>
+#include <iostream>
+
+constexpr struct {
+  int x;
+  int y;
+  int expect;
+} Cases[] = {
+    {0, 0, 0},
+    {1, 0, 0},
+    {0, 1, 0},
+    {1, 1, 1},
+    {2, 3, 6},
+    {2, 4, 4},
+    {3, 17, 51},
+    {36, 18, 36}
+};
+
+template <typename Input1, typename Input2, typename Output>
+constexpr bool test0(Input1 in1, Input2 in2, Output out)
+{
+    static_assert((std::is_same<Output, decltype(std::experimental::lcm(Input1(0), Input2(0)))>::value), "" );
+    static_assert((std::is_same<Output, decltype(std::experimental::lcm(Input2(0), Input1(0)))>::value), "" );
+    return out == std::experimental::lcm(in1, in2) ? true : (std::abort(), false);
+}
+
+
+template <typename Input1, typename Input2 = Input1>
+constexpr bool do_test(int dummy = 0)
+{
+    using S1 = typename std::make_signed<Input1>::type;
+    using S2 = typename std::make_signed<Input2>::type;
+    using U1 = typename std::make_unsigned<Input1>::type;
+    using U2 = typename std::make_unsigned<Input2>::type;
+    bool accumulate = true;
+    for (auto TC : Cases) {
+        { // Test with two signed types
+            using Output = std::common_type_t<S1, S2>;
+            accumulate &= test0<S1, S2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S1, S2, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<S1, S2, Output>(TC.x, -TC.y, TC.expect);
+            accumulate &= test0<S1, S2, Output>(-TC.x, -TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(TC.x, -TC.y, TC.expect);
+            accumulate &= test0<S2, S1, Output>(-TC.x, -TC.y, TC.expect);
+        }
+        { // test with two unsigned types
+            using Output = std::common_type_t<U1, U2>;
+            accumulate &= test0<U1, U2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<U2, U1, Output>(TC.x, TC.y, TC.expect);
+        }
+        { // Test with mixed signs
+            using Output = std::common_type_t<S1, U2>;
+            accumulate &= test0<S1, U2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<U2, S1, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S1, U2, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<U2, S1, Output>(TC.x, -TC.y, TC.expect);
+        }
+        { // Test with mixed signs
+            using Output = std::common_type_t<S2, U1>;
+            accumulate &= test0<S2, U1, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<U1, S2, Output>(TC.x, TC.y, TC.expect);
+            accumulate &= test0<S2, U1, Output>(-TC.x, TC.y, TC.expect);
+            accumulate &= test0<U1, S2, Output>(TC.x, -TC.y, TC.expect);
+        }
+    }
+    return accumulate;
+}
+
+int main()
+{
+    auto non_cce = std::rand(); // a value that can't possibly be constexpr
+
+    static_assert(do_test<signed char>(), "");
+    static_assert(do_test<short>(), "");
+    static_assert(do_test<int>(), "");
+    static_assert(do_test<long>(), "");
+    static_assert(do_test<long long>(), "");
+
+    assert(do_test<signed char>(non_cce));
+    assert(do_test<short>(non_cce));
+    assert(do_test<int>(non_cce));
+    assert(do_test<long>(non_cce));
+    assert(do_test<long long>(non_cce));
+
+    static_assert(do_test< int8_t>(), "");
+    static_assert(do_test<int16_t>(), "");
+    static_assert(do_test<int32_t>(), "");
+    static_assert(do_test<int64_t>(), "");
+
+    assert(do_test< int8_t>(non_cce));
+    assert(do_test<int16_t>(non_cce));
+    assert(do_test<int32_t>(non_cce));
+    assert(do_test<int64_t>(non_cce));
+
+    static_assert(do_test<signed char, int>(), "");
+    static_assert(do_test<int, signed char>(), "");
+    static_assert(do_test<short, int>(), "");
+    static_assert(do_test<int, short>(), "");
+    static_assert(do_test<int, long>(), "");
+    static_assert(do_test<long, int>(), "");
+    static_assert(do_test<int, long long>(), "");
+    static_assert(do_test<long long, int>(), "");
+
+    assert((do_test<signed char, int>(non_cce)));
+    assert((do_test<int, signed char>(non_cce)));
+    assert((do_test<short, int>(non_cce)));
+    assert((do_test<int, short>(non_cce)));
+    assert((do_test<int, long>(non_cce)));
+    assert((do_test<long, int>(non_cce)));
+    assert((do_test<int, long long>(non_cce)));
+    assert((do_test<long long, int>(non_cce)));
+}