Implementing Basic Arithmetic for Large Integers: Division: Implementation

[This page is auto-generated; please prefer to read the article’s PDF Form.]

Implementation

To summarize, if the divisor and dividend have been normalized, then the process of computing d is the following. Find the ﬁrst estimate d₁ using (7). Since we use Δ = 2, the desired d must be one among d₁ − 2, d₁ − 1, d₁. Speciﬁcally, d is the maximum among these three values with the property that (d₁ − i)b ≤ x, i = 0,1,2.

Note that when b is a single digit divisor, IDD x is simply the two-digits (yz), and so the ﬁrst estimate d₁ = ⌊(yz)∕e⌋ must be the correct digit d always. So, no normalization is needed in this case.

Below we implement the long-division for bigints, while using the above approach for computing a quotient-digit. Please refer to [1] and [2] for implementation of methods referred by these methods.

/* a: dividend, b: divisor, qt: quotient, rm: remainder */
int divide(bigint *a, bigint *b, bigint *qt, bigint *rm)
{
  int    na, nb, m;
  ulong  yz, d;
  bigint x, t1, t2, af, bf;
  uint   e, f = 1;

  if(BINT_ISZERO(b))
  {
    printf("divide: divisor is 0\n");
    return -1;
  }

  if(BINT_LEN(a) < BINT_LEN(b))
  {
    BINT_INIT(qt, 0);
    copy(rm, a);
    return 0;
  }

  e = b->digits[b->msd];

  if((BINT_LEN(b) > 1) && (e < B/2))
  {
    /* normalization */
    f = B/(e + 1);
    multiply_digit(a, f, &af);
    multiply_digit(b, f, &bf);
    a = ⁡
    b = &bf;
    e = b->digits[b->msd];
  }

  na = BINT_LEN(a);
  nb = BINT_LEN(b);
  /* na >= nb holds. */

  /* quotient can have maximum (na-nb+1) digits  */
  qt->msd = WIDTH - (na-nb+1);

  /* loop-invariant P: first m digits of ’a’ have been brought-down
     and processed. */

  m = 0;

  while(m < na)
  {
    if(m == 0)
    {
      m = nb;
      prefix(&x, a, m);
    }
    else
    {
      m++;
      shift_left(&x, 1);
      x.digits[WIDTH-1] = a->digits[a->msd + m-1];
    }

    yz = yz_from_x(&x, nb+1);
    d = yz/e;

    if(nb > 1)
      correct_d_and_subtract(&x, b, &d);
    else
    {
      /* we don’t need correction in d if nb=1 */
      yz = yz - e*d;
      /* above remainder is less than e, so must be single digit */
      BINT_INIT(&x, yz);
    }

    qt->digits[qt->msd + m-nb] = d;
  }

  /* (loop-invariant P) AND (m=na) holds. */
  /* Now x contains the remainder. */

  rm_leading_0s(qt);

  if(f > 1)
  {
    BINT_INIT(&t1, f);
    divide(&x, &t1, rm, &t2);
  }
  else
    copy(rm, &x);

  return 0;
}

/* treat x as having width w, by prefixing some 0s if necessary */
static ulong yz_from_x(bigint *x, int w)
{
  int  iy = WIDTH-w;  /* index of y */
  uint y, z;

  if(iy >= x->msd)
    y = x->digits[iy];
  else
    y = 0;

  if(iy+1 >= x->msd)
    z = x->digits[iy+1];
  else
    z = 0;

  return ((ulong)y)*B + z;
}

static void correct_d_and_subtract(bigint *x, bigint *b, ulong *d)
{
  bigint t;

  if(*d > B-1)
    *d = B-1;

  multiply_digit(b, *d, &t);
  /* Now (t = bd) holds. */

  if(compare(&t, x) == 1)
  {
    *d = *d - 1;
    subtract(&t, b, &t);
    /* Now (t = bd) holds again. */

    if(compare(&t, x) == 1)
    {
      *d = *d - 1;
      subtract(&t, b, &t);
      /* Now (t = bd) holds again. */
      /* We must now have t <= x. */
    }
  }

  /* ((t = bd) AND (t <= x)) holds (implies bd <= x). */
  subtract(x, &t, x);
}

References

[1] Nitin Verma. Implementing Basic Arithmetic for Large Integers: Introduction. https://mathsanew.com/articles/
implementing_large_integers_introduction.pdf (2021).

[2] Nitin Verma. Implementing Basic Arithmetic for Large Integers: Multiplication. https://mathsanew.com/articles/
implementing_large_integers_multiplication.pdf (2021).

[3] D. A. Pope, M. L. Stein. Multiple Precision Arithmetic. C. ACM, Vol 3 (12) (1960), 652–654.

[4] D. E. Knuth. The Art of Computer Programming, Vol 2, Third Edition. Addison-Wesley (1997).

[5] P. B. Hansen. Multiple-Length Division Revisited: A Tour of the Mineﬁeld. Syracuse University EECS - Technical Reports, 166 (1992). https://surface.syr.edu/eecs_techreports/166.

[6] E. V. Krishnamurthy, S. K. Nandi. On the Normalization Requirement of Divisor in Divide-and-Correct Methods. C. ACM, Vol 10 (12) (1967), 809–813.