Re: Frequency of integer division in source code?

chase@Think.COM (David Chase)
Tue, 15 Feb 1994 18:57:55 GMT

          From comp.compilers

Related articles
Frequency of integer division in source code? Peter-Lawrence.Montgomery@cwi.nl (1994-02-07)
Re: Frequency of integer division in source code? hbaker@netcom.com (1994-02-08)
Re: Frequency of integer division in source code? meissner@osf.org (1994-02-09)
Re: Frequency of integer division in source code? tfj@apusapus.demon.co.uk (1994-02-10)
Re: Frequency of integer division in source code? glew@ichips.intel.com (1994-02-11)
Re: Frequency of integer division in source code? bazyar@netcom.com (1994-02-10)
Re: Frequency of integer division in source code? cliffc@dawn.cs.rice.edu (1994-02-14)
Re: Frequency of integer division in source code? chase@Think.COM (1994-02-15)
Re: Frequency of integer division in source code? pardo@cs.washington.edu (1994-02-17)
| List of all articles for this month |
Newsgroups: comp.compilers
From: chase@Think.COM (David Chase)
Keywords: architecture, optimize
Organization: Thinking Machines Corporation, Cambridge MA, USA
References: 94-02-058 94-02-088
Date: Tue, 15 Feb 1994 18:57:55 GMT

cliffc@dawn.cs.rice.edu (Cliff Click) writes:
|> I would be more worried about code like:
|>
|> int offset;
|> struct foo *a, *b;
|> ...
|> offset = a - b;
|>
|> where sizeof(struct foo) is not a power of two.


No, if this is a problem, then you employ my favorite stupid compiler
trick, replacing a constant (odd) division with multiplication by the ring
reciprocal (this only works when you know that the remainder will be zero,
but that is the case for address arithmetic, or else the result is either
un- or implementation-defined).


For example, if you have a structure that is 12 bytes long, you might
compute "a-b" as:


    sub a,b,t ! get byte difference
    sra t,2,u ! u = t / 4;
    add t,u,b ! Multiplication by 2863311531 starts here:
    sll b,4,t
    add t,b,b
    sll b,8,t
    add b,t,b
    sll b,16,t
    add b,t,b
    add b,b,t
    sub b,t,t


(11 instructions -- not necessarily a win. There are other ways
  to expand the constant multiplication that produce more parallelism,
  if it is needed by a particular chip.)


Note that there is no need to monkey with absolute values -- it works
for negative numbers as well as positive numbers.


Or, in C:


    int t = (int) a - (int) b;
    int u = t >> 2;
    Rb = t + u; /* was Ra + (Ra << 2); */
    Rb = (Rb << 4) + Rb;
    Rb = (Rb << 8) + Rb;
    Rb = (Rb << 16) + Rb;
    t = Rb - (Rb << 1);




This can also be used to replace a rest for (unsigned)
remainder-equal-to-zero with a bounds test on the "quotient" (if the
"quotient" is greater than floor((2**wordsize)/divisor), then the remainder
cannot be zero).


The code to compute inverses in a ring follows. You need a compiler that
supports 64 bit integers if you want to generate results for 32 bits, and
you need Ansi C in either case. The relevant formatting and type
definitions/macros appear at the beginning of the file. The magic number
above was obtained thusly:


% ./invert
Enter number to invert and modulus of ring to invert in: 3 4294967296
                invert(4294967296, 3) returns 2863311531
Inverse of 3 (3) in 4294967296 (100000000) is 2863311531 (aaaaaaab),
product is 1
Enter number to invert and modulus of ring to invert in:^D


===================================================================
typedef unsigned long long int cell;
const cell word = 0x100000000LL;
#define SIGNED_CELL_FMT "%lld"
#define UNSIGNED_CELL_FMT "%llu"
#define HEX_CELL_FMT "%llx"


/* Find inverse of a in ring of integers modulo b. */
cell invert(cell a, cell b) {
    cell q = b/a;
    cell r = b - q*a;
    cell result;


    if (r == 1) {
        /* That means that q*a is the square root of
              the inverse. */
        result = (((q*q)%b)*a) % b;


    } else if (a-r == 1) {
        /* otherwise we were one away from the other end. */
        result = q+1;


    } else if (2*r > a) { /* Always cut problem in half. */
        q = q + 1;
        r = a - r;
        {
            /* Obtain inverse of r in a. */
            cell nq = invert(r, a);
            /* Overshoot, and then step back in the larger
ring. */
            result = (q*nq - (nq*r)/a) %b;
        }


    } else {
        cell nq = invert(r, a);
        /* Undershoot, then step forward to negative one,
              and then correct for square of negative one. */
        result = (q*nq + (nq*r)/a) %b;
        result = (((result*result) % b) * a) %b;
    }


    printf("\tinvert("UNSIGNED_CELL_FMT", "UNSIGNED_CELL_FMT") returns "UNSIGNED_CELL_FMT"\n", b, a, result);
    return result;
}


main() {
    cell a,b,c;
    unsigned int x,y,z;
    int rc;


    while (1) {
        printf("Enter number to invert and modulus of ring to invert in: ");
        rc = scanf(SIGNED_CELL_FMT SIGNED_CELL_FMT, &a, &c);
        if (rc != 2)
            break;


        b = invert(a, c);


        printf("Inverse of "UNSIGNED_CELL_FMT" ("HEX_CELL_FMT
                      ") in "UNSIGNED_CELL_FMT" ("HEX_CELL_FMT") is "
                      UNSIGNED_CELL_FMT" ("HEX_CELL_FMT
                      "),\n\tproduct is "UNSIGNED_CELL_FMT"\n",
a, a, c, c, b, b, a*b % c);
    }
}
===================================================================


David Chase, Thinking Machines Corp
--


Post a followup to this message

Return to the comp.compilers page.
Search the comp.compilers archives again.