diff options
-rw-r--r-- | coreutils/factor.c | 44 |
1 files changed, 10 insertions, 34 deletions
diff --git a/coreutils/factor.c b/coreutils/factor.c index 11097c12d..f910fdb44 100644 --- a/coreutils/factor.c +++ b/coreutils/factor.c | |||
@@ -48,41 +48,17 @@ typedef unsigned long half_t; | |||
48 | static inline half_t isqrt(wide_t N) | 48 | static inline half_t isqrt(wide_t N) |
49 | { | 49 | { |
50 | half_t x; | 50 | half_t x; |
51 | unsigned shift; | ||
51 | 52 | ||
52 | // Never called with N < 1 | 53 | shift = WIDE_BITS - 2; |
53 | //if (N == 0) | 54 | x = 0; |
54 | // return 0; | 55 | do { |
55 | 56 | x = (x << 1) + 1; | |
56 | x = HALF_MAX; | 57 | if ((wide_t)x * x > (N >> shift)) |
57 | /* First approximation of x+1 > sqrt(N) - all-ones, half as many bits: | 58 | x--; /* whoops, that +1 was too much */ |
58 | * 1xxxxx -> 111 (six bits to three) | 59 | shift -= 2; |
59 | * 01xxxx -> 111 | 60 | } while ((int)shift >= 0); |
60 | * 001xxx -> 011 | 61 | return x; |
61 | * 0001xx -> 011 and so on. | ||
62 | */ | ||
63 | // It is actually not performance-critical at all. | ||
64 | // Can simply start Newton loop with very conservative x=0xffffffff. | ||
65 | //wide_t mask_2bits; | ||
66 | //mask_2bits = TOPMOST_WIDE_BIT | (TOPMOST_WIDE_BIT >> 1); | ||
67 | //while (!(N & mask_2bits)) { | ||
68 | // x >>= 1; | ||
69 | // mask_2bits >>= 2; | ||
70 | //} | ||
71 | dbg("x:%"HALF_FMT"x", x); | ||
72 | |||
73 | for (;;) { | ||
74 | half_t y = (x + N/x) / 2; | ||
75 | dbg("y:%x y^2:%llx", y, (wide_t)y * y); | ||
76 | /* | ||
77 | * "real" y may be one bit wider: 0x100000000 and get truncated to 0. | ||
78 | * In this case, "real" y is > x. The first check below is for this case: | ||
79 | */ | ||
80 | if (y == 0 || y >= x) { | ||
81 | dbg("isqrt(%llx)=%"HALF_FMT"x", N, x); | ||
82 | return x; | ||
83 | } | ||
84 | x = y; | ||
85 | } | ||
86 | } | 62 | } |
87 | 63 | ||
88 | static NOINLINE half_t isqrt_odd(wide_t N) | 64 | static NOINLINE half_t isqrt_odd(wide_t N) |