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 unsigned long gcd(unsigned long a, unsigned long b)
  17 {
  18         unsigned long r = a | b;
  19
  20         if (!a || !b)
  21                 return r;
  22
  23         b >>= __ffs(b);
  24         if (b == 1)
  25                 return r & -r;
  26
  27         for (;;) {
  28                 a >>= __ffs(a);
  29                 if (a == 1)
  30                         return r & -r;
  31                 if (a == b)
  32                         return a << __ffs(r);
  33
  34                 if (a < b)
  35                         swap(a, b);
  36                 a -= b;
  37         }
  38 }
  39
  40 #else
  41
  42 /* If normalization is done by loops, the even/odd algorithm is a win. */
  43 unsigned long gcd(unsigned long a, unsigned long b)
  44 {
  45         unsigned long r = a | b;
  46
  47         if (!a || !b)
  48                 return r;
  49
  50         /* Isolate lsbit of r */
  51         r &= -r;
  52
  53         while (!(b & r))
  54                 b >>= 1;
  55         if (b == r)
  56                 return r;
  57
  58         for (;;) {
  59                 while (!(a & r))
  60                         a >>= 1;
  61                 if (a == r)
  62                         return r;
  63                 if (a == b)
  64                         return a;
  65
  66                 if (a < b)
  67                         swap(a, b);
  68                 a -= b;
  69                 a >>= 1;
  70                 if (a & r)
  71                         a += b;
  72                 a >>= 1;
  73         }
  74 }
  75
  76 #endif
  77
  78 EXPORT_SYMBOL_GPL(gcd);