home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD1.mdf / magazine / drdobbs / 1991 / 02 / divisor.asc < prev    next >
Text File  |  1990-12-26  |  4KB  |  126 lines
  1. _OPTIMIZING INTEGER DIVISION BY A CONSTANT DIVISOR_
  2. by Robert Grappel
  3.  
  4.  
  5.  
  6. [LISTING ONE]
  7.  
  8. #include stdio.h
  9.  
  10. /* Program to generate "star-chain-sequence" for division of an unsigned 
  11. 16-bit integer numerator by an unsigned 16-bit integer constant denominator.
  12. It assumes the existence of two 32-bit integer "registers", add, subtract, and
  13. shift instructions. It uses the "star-chain" multiplication routine from DDJ
  14. #125, March 1987 page 35 */
  15.  
  16. static unsigned int mult; /* global variable for trim_trailing() & binmul() */
  17.  
  18. /* Support subroutine to trim trailing 0s or 1s from global variable "mult". */
  19. int trim_trailing(one_zero) int one_zero;
  20. { /* if one_zero    == 0, trim trailing zeros in "mult", return "count"
  21.             == 1, ones                                       */
  22.     int c; /* bit counter */
  23.     for (c = 0; ((mult & 1) == one_zero); c++, mult >>= 1) ;
  24.     return c;
  25. }
  26.  
  27. /* Slightly modified version of multiplication routine */
  28. binmul(m) long m;
  29. {
  30.     int    last_shift,    /* final shift count */
  31.         last_cnt,    /* count of low-order zeros */
  32.         stkptr,     /* pointer to "stack[]" */
  33.         cnt,        /* bit counter */
  34.         ts,        /* top-of-stack element */
  35.         flag,        /* flag for special-case */
  36.         stack[16];    /* stack for shift-add/subs */
  37.     mult = m;
  38.     stkptr = last_cnt = 0;    /* init. stack ptr. and count */
  39.     last_shift = trim_trailing(0);    /* trim trailing 0s */
  40.     while (1)
  41.     { /* loop to decompose "mult", building stack */
  42.         cnt = trim_trailing(1); /* count low-order 1s */
  43.         if (cnt > 1)
  44.         { /* more than 1 bit, shift-subtract */
  45.             flag = 0;
  46.             if (last_cnt == 1) 
  47.      /* shift "k",sub,shift 1,add --> shift "k+1", sub */ 
  48.                 stack[stkptr-1] = -(cnt+1); /* overwrite */
  49.             else
  50.                 stack[stkptr++] = -cnt;
  51.         }
  52.         else 
  53.             flag = 1; /* need another shift-add */
  54.         if (mult == 0) break; /* decomp. "mult", now output */
  55.         last_cnt = trim_trailing(0) + flag; /* low-order 0s */
  56.         stack[stkptr++] = last_cnt; /* shift-add */
  57.     }
  58.     while (stkptr > 0)
  59.     { /* output code from stack */
  60.         ts = stack[--stkptr]; /* get top stack element */
  61.         if (ts < 0)
  62.         { /* generate shift-subtract instructions */
  63.             printf("\nRw <<= %d",-ts);
  64.             printf("\nRw -= R1");
  65.         }
  66.         else
  67.         { /* generate shift-add instructions */
  68.             printf("\nRw <<= %d",ts);
  69.             printf("\nRw += R1");
  70.         }
  71.     }
  72.     if (last_shift != 0) printf("\nRw <<= %d",last_shift);
  73. }
  74.  
  75. main()
  76. { /* generate pseudo-instructions for star-chain division */
  77.     int    a,b,r, /* computed multiplier, addend, and remainder */
  78.         i2,    /* number of bits to scale divisor */
  79.         denom,    /* intended divisor */
  80.         denom2; /* divisor scaled by powers-of-2 */
  81.     int z = 65536;    /* 2^16 */
  82.     printf("\nEnter positive integer denominator: ");
  83.     scanf("%d",&denom);
  84.     if (denom != 0)
  85.     { /* scale denominator by powers-of-2 */
  86.         mult = denom;
  87.         i2 = trim_trailing(0); /* how many powers-of-2? */
  88.         denom2 = mult;
  89.         if (denom2 == 1)
  90.         { /* divisor was power-of-2, simply scale it */
  91.             printf("\nRw = R1");
  92.             if (i2 > 0) printf("\nRw >>= %d",i2);
  93.         }
  94.         else
  95.         { /* divisor not power-of-2, scale and more */
  96.             if (i2 > 0) 
  97.                 printf("\nR1 >>= %d",i2); /* handle scaling */
  98.             printf("\nRw = R1");    /* load work register */
  99.             a = z / denom2;     /* scaled reciprocal */
  100.              r = z % denom2;    /* remainder of recip. */
  101.             b = a + r - 1;
  102.             binmul(a);
  103.             printf("\nRw += %d",b);
  104.             printf("\nRw >>= 16");
  105.         }
  106.     }
  107.     else
  108.         printf("\nCannot divide by zero\n"); /* special case */
  109. }
  110.  
  111.  
  112.  
  113. Example 1: Division by 102
  114.  
  115. R1 >>= 1        scale 102 down to 51 (odd)
  116. Rw = R1             mult. by a = (65536/51) = 1285
  117. Rw <<= 2
  118. Rw += R1
  119. Rw <<= 6
  120. Rw += R1
  121. Rw <<= 2
  122. Rw += R1        
  123. Rw += 1285            add b = (a + r - 1) = 1285
  124. Rw >>= 16        divide by z = 65536
  125.  
  126.