lib/gcd.c

   1 #include <linux/kernel.h>
   2 #include <linux/gcd.h>
   3 #include <linux/export.h>
   4
   5 /*
   6  * This implements the binary GCD algorithm. (Often attributed to Stein,
   7  * but as Knuth has noted, appears in a first-century Chinese math text.)
   8  *
   9  * This is faster than the division-based algorithm even on x86, which
  10  * has decent hardware division.
  11  */
  12
  13 #if !defined(CONFIG_CPU_NO_EFFICIENT_FFS) && !defined(CPU_NO_EFFICIENT_FFS)
  14
  15 /* If __ffs is available, the even/odd algorithm benchmarks slower. */
  16
  17 /**
  18  * gcd - calculate and return the greatest common divisor of 2 unsigned longs
  19  * @a: first value
  20  * @b: second value
  21  */
  22 unsigned long gcd(unsigned long a, unsigned long b)
  23 {
  24         unsigned long r = a | b;
  25
  26         if (!a || !b)
  27                 return r;
  28
  29         b >>= __ffs(b);
  30         if (b == 1)
  31                 return r & -r;
  32
  33         for (;;) {
  34                 a >>= __ffs(a);
  35                 if (a == 1)
  36                         return r & -r;
  37                 if (a == b)
  38                         return a << __ffs(r);
  39
  40                 if (a < b)
  41                         swap(a, b);
  42                 a -= b;
  43         }
  44 }
  45
  46 #else
  47
  48 /* If normalization is done by loops, the even/odd algorithm is a win. */
  49 unsigned long gcd(unsigned long a, unsigned long b)
  50 {
  51         unsigned long r = a | b;
  52
  53         if (!a || !b)
  54                 return r;
  55
  56         /* Isolate lsbit of r */
  57         r &= -r;
  58
  59         while (!(b & r))
  60                 b >>= 1;
  61         if (b == r)
  62                 return r;
  63
  64         for (;;) {
  65                 while (!(a & r))
  66                         a >>= 1;
  67                 if (a == r)
  68                         return r;
  69                 if (a == b)
  70                         return a;
  71
  72                 if (a < b)
  73                         swap(a, b);
  74                 a -= b;
  75                 a >>= 1;
  76                 if (a & r)
  77                         a += b;
  78                 a >>= 1;
  79         }
  80 }
  81
  82 #endif
  83
  84 EXPORT_SYMBOL_GPL(gcd);