[ci skip] Trade out algorithm from 'The Joy of Factoring' to Wikipedia's version which reduces the number of required temporaries. In fact, the speedup is not large, but the code is more compact, and for larger types, the difference becomes more noticeable.

include/boost/integer/extended_euclidean.hpp

@@ -13,9 +13,8 @@
 
 namespace boost { namespace integer {
 
-// From "The Joy of Factoring", Algorithm 2.7.
+// From "The Joy of Factoring", Algorithm 2.7, with a small optimization to remove tmps from Wikipedia.
 // Solves mx + ny = gcd(m,n). Returns tuple with (gcd(m,n), x, y).
-// Is this the natural ordering?, or must people simply have to read the docs?
 
 template<class Z>
 struct euclidean_result_t {
@@ -38,45 +37,27 @@ euclidean_result_t<Z> extended_euclidean(Z m, Z n)
     {
         BOOST_THROW_EXCEPTION(std::domain_error("Arguments must be strictly positive.\n"));
     }
-    bool swapped = false;
-    if (m < n)
-    {
-        swapped = true;
-        std::swap(m, n);
-    }
-    Z u0 = m;
-    Z u1 = 1;
-    Z u2 = 0;
-    Z v0 = n;
-    Z v1 = 0;
-    Z v2 = 1;
-    Z w0;
-    Z w1;
-    Z w2;
-    while(v0 > 0)
-    {
-        Z q = u0/v0;
-        w0 = u0 - q*v0;
-        w1 = u1 - q*v1;
-        w2 = u2 - q*v2;
-        u0 = v0;
-        u1 = v1;
-        u2 = v2;
-        v0 = w0;
-        v1 = w1;
-        v2 = w2;
-    }
-    if (swapped)
-    {
-        std::swap(u1, u2);
-        BOOST_ASSERT(u2*m+u1*n==u0);
-    }
-    else
-    {
-        BOOST_ASSERT(u1*m+u2*n==u0);
+
+    Z s = 0;
+    Z old_s = 1;
+    Z r = n;
+    Z old_r = m;
+
+    while (r != 0) {
+      Z q = old_r/r;
+      Z tmp = r;
+      r = old_r - q*tmp;
+      old_r = tmp;
+
+      tmp = s;
+      s = old_s - q*tmp;
+      old_s = tmp;
     }
 
-    return {u0, u1, u2};
+    Z y = (old_r - old_s*m)/n;
+
+    BOOST_ASSERT(old_s*m+y*n==old_r);
+    return {old_r, old_s, y};
 }
 
 }}