On Mar 14, 2009, at 10:47 PM, Jason Evans wrote:
> Robert Bradshaw wrote:
>> http://trac.cython.org/cython_trac/ticket/127 I only special-cased
>> 0-3, the rest is generic (repeated squaring) code.
>
> Here's some code I wrote a while back that uses binary
> decomposition to
> improve performance. I don't know how much it really matters...
>
> Jason
>
> // Use binary decomposition to reduce the number of multiplication
> operations
> // when calculating x^n. See "Hackers Delight", pp 212-213 for
> further
> // information.
> LkmpInline int64_t
> LkpIntExp(int64_t x, int64_t n) {
> int64_t ret, p;
> bool neg = (n < 0);
>
> // The main loop only works for non-negative exponents, so
> transform
> // negative exponents outside the loop, and clean up at the
> end with a
> // division.
> if (neg) {
> // This can overflow if n is -2^63, but the results are
> explicitly
> // undefined on underflow/overflow, so there is no need to
> check.
> n = -n;
> }
>
> for (ret = 1, p = x;
> n != 0;
> n >>= 1, p *= p) {
> if (n & 1) {
> ret *= p;
> }
> }
>
> if (neg) {
> ret = 1 / ret;
> }
>
> return ret;
Yep, this is basically what I implemented. The switch at the top is
handy for small exponents, and is primarily for constant exponents.
The messy expression at the top of the loop is to avoid branching
(may or may not be worth it--depends on how randomly the bits are set
in the exponent).
----------
static INLINE %(type)s %(func_name)s(%(type)s b, %(type)s e) {
%(type)s t = b;
switch (e) {
case 3:
t *= b;
case 2:
t *= b;
case 1:
return t;
case 0:
return 1;
}
if (unlikely(e<0)) return 0;
t = 1;
while (likely(e)) {
t *= (b * (e&1)) | ((~e)&1); /* 1 or b */
b *= b;
e >>= 1;
}
return t;
_______________________________________________
Cython-dev mailing list
[email protected]
http://codespeak.net/mailman/listinfo/cython-dev
