Annotation of embedaddon/php/ext/bcmath/libbcmath/src/div.c, revision 1.1

1.1     ! misho       1: /* div.c: bcmath library file. */
        !             2: /*
        !             3:     Copyright (C) 1991, 1992, 1993, 1994, 1997 Free Software Foundation, Inc.
        !             4:     Copyright (C) 2000 Philip A. Nelson
        !             5: 
        !             6:     This library is free software; you can redistribute it and/or
        !             7:     modify it under the terms of the GNU Lesser General Public
        !             8:     License as published by the Free Software Foundation; either
        !             9:     version 2 of the License, or (at your option) any later version.
        !            10: 
        !            11:     This library is distributed in the hope that it will be useful,
        !            12:     but WITHOUT ANY WARRANTY; without even the implied warranty of
        !            13:     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
        !            14:     Lesser General Public License for more details.  (COPYING.LIB)
        !            15: 
        !            16:     You should have received a copy of the GNU Lesser General Public
        !            17:     License along with this library; if not, write to:
        !            18: 
        !            19:       The Free Software Foundation, Inc.
        !            20:       59 Temple Place, Suite 330
        !            21:       Boston, MA 02111-1307 USA.
        !            22: 
        !            23:     You may contact the author by:
        !            24:        e-mail:  philnelson@acm.org
        !            25:       us-mail:  Philip A. Nelson
        !            26:                 Computer Science Department, 9062
        !            27:                 Western Washington University
        !            28:                 Bellingham, WA 98226-9062
        !            29:        
        !            30: *************************************************************************/
        !            31: 
        !            32: #include <config.h>
        !            33: #include <stdio.h>
        !            34: #include <assert.h>
        !            35: #include <stdlib.h>
        !            36: #include <ctype.h>
        !            37: #include <stdarg.h>
        !            38: #include "bcmath.h"
        !            39: #include "private.h"
        !            40: 
        !            41: 
        !            42: /* Some utility routines for the divide:  First a one digit multiply.
        !            43:    NUM (with SIZE digits) is multiplied by DIGIT and the result is
        !            44:    placed into RESULT.  It is written so that NUM and RESULT can be
        !            45:    the same pointers.  */
        !            46: 
        !            47: static void
        !            48: _one_mult (num, size, digit, result)
        !            49:      unsigned char *num;
        !            50:      int size, digit;
        !            51:      unsigned char *result;
        !            52: {
        !            53:   int carry, value;
        !            54:   unsigned char *nptr, *rptr;
        !            55: 
        !            56:   if (digit == 0)
        !            57:     memset (result, 0, size);
        !            58:   else
        !            59:     {
        !            60:       if (digit == 1)
        !            61:        memcpy (result, num, size);
        !            62:       else
        !            63:        {
        !            64:          /* Initialize */
        !            65:          nptr = (unsigned char *) (num+size-1);
        !            66:          rptr = (unsigned char *) (result+size-1);
        !            67:          carry = 0;
        !            68: 
        !            69:          while (size-- > 0)
        !            70:            {
        !            71:              value = *nptr-- * digit + carry;
        !            72:              *rptr-- = value % BASE;
        !            73:              carry = value / BASE;
        !            74:            }
        !            75: 
        !            76:          if (carry != 0) *rptr = carry;
        !            77:        }
        !            78:     }
        !            79: }
        !            80: 
        !            81: 
        !            82: /* The full division routine. This computes N1 / N2.  It returns
        !            83:    0 if the division is ok and the result is in QUOT.  The number of
        !            84:    digits after the decimal point is SCALE. It returns -1 if division
        !            85:    by zero is tried.  The algorithm is found in Knuth Vol 2. p237. */
        !            86: 
        !            87: int
        !            88: bc_divide (bc_num n1, bc_num n2, bc_num *quot, int scale TSRMLS_DC)
        !            89: {
        !            90:   bc_num qval;
        !            91:   unsigned char *num1, *num2;
        !            92:   unsigned char *ptr1, *ptr2, *n2ptr, *qptr;
        !            93:   int  scale1, val;
        !            94:   unsigned int  len1, len2, scale2, qdigits, extra, count;
        !            95:   unsigned int  qdig, qguess, borrow, carry;
        !            96:   unsigned char *mval;
        !            97:   char zero;
        !            98:   unsigned int  norm;
        !            99: 
        !           100:   /* Test for divide by zero. */
        !           101:   if (bc_is_zero (n2 TSRMLS_CC)) return -1;
        !           102: 
        !           103:   /* Test for divide by 1.  If it is we must truncate. */
        !           104:   if (n2->n_scale == 0)
        !           105:     {
        !           106:       if (n2->n_len == 1 && *n2->n_value == 1)
        !           107:        {
        !           108:          qval = bc_new_num (n1->n_len, scale);
        !           109:          qval->n_sign = (n1->n_sign == n2->n_sign ? PLUS : MINUS);
        !           110:          memset (&qval->n_value[n1->n_len],0,scale);
        !           111:          memcpy (qval->n_value, n1->n_value,
        !           112:                  n1->n_len + MIN(n1->n_scale,scale));
        !           113:          bc_free_num (quot);
        !           114:          *quot = qval;
        !           115:        }
        !           116:     }
        !           117: 
        !           118:   /* Set up the divide.  Move the decimal point on n1 by n2's scale.
        !           119:      Remember, zeros on the end of num2 are wasted effort for dividing. */
        !           120:   scale2 = n2->n_scale;
        !           121:   n2ptr = (unsigned char *) n2->n_value+n2->n_len+scale2-1;
        !           122:   while ((scale2 > 0) && (*n2ptr-- == 0)) scale2--;
        !           123: 
        !           124:   len1 = n1->n_len + scale2;
        !           125:   scale1 = n1->n_scale - scale2;
        !           126:   if (scale1 < scale)
        !           127:     extra = scale - scale1;
        !           128:   else
        !           129:     extra = 0;
        !           130:   num1 = (unsigned char *) safe_emalloc (1, n1->n_len+n1->n_scale, extra+2);
        !           131:   if (num1 == NULL) bc_out_of_memory();
        !           132:   memset (num1, 0, n1->n_len+n1->n_scale+extra+2);
        !           133:   memcpy (num1+1, n1->n_value, n1->n_len+n1->n_scale);
        !           134: 
        !           135:   len2 = n2->n_len + scale2;
        !           136:   num2 = (unsigned char *) safe_emalloc (1, len2, 1);
        !           137:   if (num2 == NULL) bc_out_of_memory();
        !           138:   memcpy (num2, n2->n_value, len2);
        !           139:   *(num2+len2) = 0;
        !           140:   n2ptr = num2;
        !           141:   while (*n2ptr == 0)
        !           142:     {
        !           143:       n2ptr++;
        !           144:       len2--;
        !           145:     }
        !           146: 
        !           147:   /* Calculate the number of quotient digits. */
        !           148:   if (len2 > len1+scale)
        !           149:     {
        !           150:       qdigits = scale+1;
        !           151:       zero = TRUE;
        !           152:     }
        !           153:   else
        !           154:     {
        !           155:       zero = FALSE;
        !           156:       if (len2>len1)
        !           157:        qdigits = scale+1;      /* One for the zero integer part. */
        !           158:       else
        !           159:        qdigits = len1-len2+scale+1;
        !           160:     }
        !           161: 
        !           162:   /* Allocate and zero the storage for the quotient. */
        !           163:   qval = bc_new_num (qdigits-scale,scale);
        !           164:   memset (qval->n_value, 0, qdigits);
        !           165: 
        !           166:   /* Allocate storage for the temporary storage mval. */
        !           167:   mval = (unsigned char *) safe_emalloc (1, len2, 1);
        !           168:   if (mval == NULL) bc_out_of_memory ();
        !           169: 
        !           170:   /* Now for the full divide algorithm. */
        !           171:   if (!zero)
        !           172:     {
        !           173:       /* Normalize */
        !           174:       norm =  10 / ((int)*n2ptr + 1);
        !           175:       if (norm != 1)
        !           176:        {
        !           177:          _one_mult (num1, len1+scale1+extra+1, norm, num1);
        !           178:          _one_mult (n2ptr, len2, norm, n2ptr);
        !           179:        }
        !           180: 
        !           181:       /* Initialize divide loop. */
        !           182:       qdig = 0;
        !           183:       if (len2 > len1)
        !           184:        qptr = (unsigned char *) qval->n_value+len2-len1;
        !           185:       else
        !           186:        qptr = (unsigned char *) qval->n_value;
        !           187: 
        !           188:       /* Loop */
        !           189:       while (qdig <= len1+scale-len2)
        !           190:        {
        !           191:          /* Calculate the quotient digit guess. */
        !           192:          if (*n2ptr == num1[qdig])
        !           193:            qguess = 9;
        !           194:          else
        !           195:            qguess = (num1[qdig]*10 + num1[qdig+1]) / *n2ptr;
        !           196: 
        !           197:          /* Test qguess. */
        !           198:          if (n2ptr[1]*qguess >
        !           199:              (num1[qdig]*10 + num1[qdig+1] - *n2ptr*qguess)*10
        !           200:               + num1[qdig+2])
        !           201:            {
        !           202:              qguess--;
        !           203:              /* And again. */
        !           204:              if (n2ptr[1]*qguess >
        !           205:                  (num1[qdig]*10 + num1[qdig+1] - *n2ptr*qguess)*10
        !           206:                  + num1[qdig+2])
        !           207:                qguess--;
        !           208:            }
        !           209: 
        !           210:          /* Multiply and subtract. */
        !           211:          borrow = 0;
        !           212:          if (qguess != 0)
        !           213:            {
        !           214:              *mval = 0;
        !           215:              _one_mult (n2ptr, len2, qguess, mval+1);
        !           216:              ptr1 = (unsigned char *) num1+qdig+len2;
        !           217:              ptr2 = (unsigned char *) mval+len2;
        !           218:              for (count = 0; count < len2+1; count++)
        !           219:                {
        !           220:                  val = (int) *ptr1 - (int) *ptr2-- - borrow;
        !           221:                  if (val < 0)
        !           222:                    {
        !           223:                      val += 10;
        !           224:                      borrow = 1;
        !           225:                    }
        !           226:                  else
        !           227:                    borrow = 0;
        !           228:                  *ptr1-- = val;
        !           229:                }
        !           230:            }
        !           231: 
        !           232:          /* Test for negative result. */
        !           233:          if (borrow == 1)
        !           234:            {
        !           235:              qguess--;
        !           236:              ptr1 = (unsigned char *) num1+qdig+len2;
        !           237:              ptr2 = (unsigned char *) n2ptr+len2-1;
        !           238:              carry = 0;
        !           239:              for (count = 0; count < len2; count++)
        !           240:                {
        !           241:                  val = (int) *ptr1 + (int) *ptr2-- + carry;
        !           242:                  if (val > 9)
        !           243:                    {
        !           244:                      val -= 10;
        !           245:                      carry = 1;
        !           246:                    }
        !           247:                  else
        !           248:                    carry = 0;
        !           249:                  *ptr1-- = val;
        !           250:                }
        !           251:              if (carry == 1) *ptr1 = (*ptr1 + 1) % 10;
        !           252:            }
        !           253: 
        !           254:          /* We now know the quotient digit. */
        !           255:          *qptr++ =  qguess;
        !           256:          qdig++;
        !           257:        }
        !           258:     }
        !           259: 
        !           260:   /* Clean up and return the number. */
        !           261:   qval->n_sign = ( n1->n_sign == n2->n_sign ? PLUS : MINUS );
        !           262:   if (bc_is_zero (qval TSRMLS_CC)) qval->n_sign = PLUS;
        !           263:   _bc_rm_leading_zeros (qval);
        !           264:   bc_free_num (quot);
        !           265:   *quot = qval;
        !           266: 
        !           267:   /* Clean up temporary storage. */
        !           268:   efree (mval);
        !           269:   efree (num1);
        !           270:   efree (num2);
        !           271: 
        !           272:   return 0;    /* Everything is OK. */
        !           273: }
        !           274: 

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>