diff options
| author | Denys Vlasenko <vda.linux@googlemail.com> | 2017-04-10 00:24:16 +0200 |
|---|---|---|
| committer | Denys Vlasenko <vda.linux@googlemail.com> | 2017-04-10 00:24:16 +0200 |
| commit | c6476dca5484a9413d20dcbf887cc78055e6965f (patch) | |
| tree | 96c1328084f6658e9d77b1abad056d0bb4feefc6 /coreutils | |
| parent | c804d4ec5cb222c842644bb99d9b077f5c6576f2 (diff) | |
| download | busybox-w32-c6476dca5484a9413d20dcbf887cc78055e6965f.tar.gz busybox-w32-c6476dca5484a9413d20dcbf887cc78055e6965f.tar.bz2 busybox-w32-c6476dca5484a9413d20dcbf887cc78055e6965f.zip | |
factor: simpler isqrt
function old new delta
isqrt_odd 102 87 -15
Signed-off-by: Denys Vlasenko <vda.linux@googlemail.com>
Diffstat (limited to 'coreutils')
| -rw-r--r-- | coreutils/factor.c | 63 |
1 files changed, 39 insertions, 24 deletions
diff --git a/coreutils/factor.c b/coreutils/factor.c index 85284aa27..48833e103 100644 --- a/coreutils/factor.c +++ b/coreutils/factor.c | |||
| @@ -47,33 +47,37 @@ typedef unsigned long half_t; | |||
| 47 | /* Returns such x that x+1 > sqrt(N) */ | 47 | /* Returns such x that x+1 > sqrt(N) */ |
| 48 | static inline half_t isqrt(wide_t N) | 48 | static inline half_t isqrt(wide_t N) |
| 49 | { | 49 | { |
| 50 | wide_t x; | 50 | wide_t mask_2bits; |
| 51 | unsigned c; | 51 | half_t x; |
| 52 | 52 | ||
| 53 | // Never called with N < 1 | 53 | // Never called with N < 1 |
| 54 | // if (N == 0) | 54 | // if (N == 0) |
| 55 | // return 0; | 55 | // return 0; |
| 56 | // | 56 | |
| 57 | /* Count leading zeros */ | 57 | /* First approximation x > sqrt(N) - half as many bits: |
| 58 | c = 0; | 58 | * 1xxxxx -> 111 (six bits to three) |
| 59 | while (!(N & TOPMOST_WIDE_BIT)) { | 59 | * 01xxxx -> 111 |
| 60 | c++; | 60 | * 001xxx -> 011 |
| 61 | N <<= 1; | 61 | * 0001xx -> 011 and so on. |
| 62 | */ | ||
| 63 | x = HALF_MAX; | ||
| 64 | mask_2bits = TOPMOST_WIDE_BIT | (TOPMOST_WIDE_BIT >> 1); | ||
| 65 | while (mask_2bits && !(N & mask_2bits)) { | ||
| 66 | x >>= 1; | ||
| 67 | mask_2bits >>= 2; | ||
| 62 | } | 68 | } |
| 63 | N >>= c; | 69 | dbg("x:%"HALF_FMT"x", x); |
| 64 | 70 | ||
| 65 | /* Make x > sqrt(n) */ | ||
| 66 | x = (wide_t)1 << ((WIDE_BITS + 1 - c) / 2); | ||
| 67 | dbg("x:%llx", x); | ||
| 68 | for (;;) { | 71 | for (;;) { |
| 69 | wide_t y = (x + N/x) / 2; | 72 | half_t y = (x + N/x) / 2; |
| 70 | dbg("y:%llx y^2:%llx (y+1)^2:%llx]", y, y*y, (y+1)*(y+1)); | 73 | dbg("y:%x y^2:%llx", y, (wide_t)y * y); |
| 71 | if (y >= x) { | 74 | /* |
| 72 | /* Handle degenerate case N = 0xffffffffff...fffffff */ | 75 | * "real" y may be one bit wider: 0x100000000 and get truncated to 0. |
| 73 | if (y == (wide_t)HALF_MAX + 1) | 76 | * In this case, "real" y is > x. The first check below is for this case: |
| 74 | y--; | 77 | */ |
| 75 | dbg("isqrt(%llx)=%llx"HALF_FMT, N, y); | 78 | if (y == 0 || y >= x) { |
| 76 | return y; | 79 | dbg("isqrt(%llx)=%"HALF_FMT"x", N, x); |
| 80 | return x; | ||
| 77 | } | 81 | } |
| 78 | x = y; | 82 | x = y; |
| 79 | } | 83 | } |
| @@ -92,6 +96,8 @@ static NOINLINE void factorize(wide_t N) | |||
| 92 | half_t factor; | 96 | half_t factor; |
| 93 | half_t max_factor; | 97 | half_t max_factor; |
| 94 | unsigned count3; | 98 | unsigned count3; |
| 99 | unsigned count5; | ||
| 100 | unsigned count7; | ||
| 95 | 101 | ||
| 96 | if (N < 4) | 102 | if (N < 4) |
| 97 | goto end; | 103 | goto end; |
| @@ -103,6 +109,8 @@ static NOINLINE void factorize(wide_t N) | |||
| 103 | 109 | ||
| 104 | max_factor = isqrt_odd(N); | 110 | max_factor = isqrt_odd(N); |
| 105 | count3 = 3; | 111 | count3 = 3; |
| 112 | count5 = 6; | ||
| 113 | count7 = 9; | ||
| 106 | factor = 3; | 114 | factor = 3; |
| 107 | for (;;) { | 115 | for (;;) { |
| 108 | /* The division is the most costly part of the loop. | 116 | /* The division is the most costly part of the loop. |
| @@ -123,11 +131,18 @@ static NOINLINE void factorize(wide_t N) | |||
| 123 | * ^ ^ ^ ^ ^ ^ ^ _ ^ ^ _ ^ ^ ^ ^ | 131 | * ^ ^ ^ ^ ^ ^ ^ _ ^ ^ _ ^ ^ ^ ^ |
| 124 | * (^ = primes, _ = would-be-primes-if-not-divisible-by-5) | 132 | * (^ = primes, _ = would-be-primes-if-not-divisible-by-5) |
| 125 | */ | 133 | */ |
| 126 | --count3; | 134 | count7--; |
| 127 | if (count3 == 0) { | 135 | count5--; |
| 136 | count3--; | ||
| 137 | if (count3 && count5 && count7) | ||
| 138 | continue; | ||
| 139 | if (count3 == 0) | ||
| 128 | count3 = 3; | 140 | count3 = 3; |
| 129 | goto next_factor; | 141 | if (count5 == 0) |
| 130 | } | 142 | count5 = 5; |
| 143 | if (count7 == 0) | ||
| 144 | count7 = 7; | ||
| 145 | goto next_factor; | ||
| 131 | } | 146 | } |
| 132 | end: | 147 | end: |
| 133 | if (N > 1) | 148 | if (N > 1) |