2N/A/* Copyright (c) 1984, 1986, 1987, 1988, 1989 AT&T */
2N/A/* All Rights Reserved */
2N/A
2N/A
2N/A/*
2N/A * Copyright (c) 1980 Regents of the University of California.
2N/A * All rights reserved. The Berkeley software License Agreement
2N/A * specifies the terms and conditions for redistribution.
2N/A */
2N/A/* Portions Copyright(c) 1988, Sun Microsystems Inc. */
2N/A/* All Rights Reserved */
2N/A
2N/A/*
2N/A * Copyright (c) 1997, by Sun Microsystems, Inc.
2N/A * All rights reserved.
2N/A */
2N/A
2N/A#ident "%Z%%M% %I% %E% SMI" /* SVr4.0 1.1 */
2N/A
2N/A/* LINTLIBRARY */
2N/A
2N/A#include <mp.h>
2N/A#include <stdio.h>
2N/A#include <stdlib.h>
2N/A#include <sys/types.h>
2N/A#include "libmp.h"
2N/A
2N/Astatic void m_div(MINT *, MINT *, MINT *, MINT *);
2N/A
2N/Avoid
2N/Amp_mdiv(MINT *a, MINT *b, MINT *q, MINT *r)
2N/A{
2N/A MINT x, y;
2N/A int sign;
2N/A
2N/A sign = 1;
2N/A x.len = y.len = 0;
2N/A _mp_move(a, &x);
2N/A _mp_move(b, &y);
2N/A if (x.len < 0) {
2N/A sign = -1;
2N/A x.len = -x.len;
2N/A }
2N/A if (y.len < 0) {
2N/A sign = -sign;
2N/A y.len = -y.len;
2N/A }
2N/A _mp_xfree(q);
2N/A _mp_xfree(r);
2N/A m_div(&x, &y, q, r);
2N/A if (sign == -1) {
2N/A q->len = -q->len;
2N/A r->len = -r->len;
2N/A }
2N/A _mp_xfree(&x);
2N/A _mp_xfree(&y);
2N/A}
2N/A
2N/Astatic int
2N/Am_dsb(int qx, int n, short *a, short *b)
2N/A{
2N/A int borrow;
2N/A int s3b2shit;
2N/A int j;
2N/A short fifteen = 15;
2N/A short *aptr, *bptr;
2N/A#ifdef DEBUGDSB
2N/A (void) printf("m_dsb %d %d %d %d\n", qx, n, *a, *b);
2N/A#endif
2N/A
2N/A borrow = 0;
2N/A aptr = a;
2N/A bptr = b;
2N/A for (j = n; j > 0; j--) {
2N/A#ifdef DEBUGDSB
2N/A (void) printf("1 borrow=%x %d %d %d\n", borrow, (*aptr * qx),
2N/A *bptr, *aptr);
2N/A#endif
2N/A borrow -= (*aptr++) * qx - *bptr;
2N/A#ifdef DEBUGDSB
2N/A (void) printf("2 borrow=%x %d %d %d\n", borrow, (*aptr * qx),
2N/A *bptr, *aptr);
2N/A#endif
2N/A *bptr++ = (short)(borrow & 077777);
2N/A#ifdef DEBUGDSB
2N/A (void) printf("3 borrow=%x %d %d %d\n", borrow, (*aptr * qx),
2N/A *bptr, *aptr);
2N/A#endif
2N/A if (borrow >= 0) borrow >>= fifteen; /* 3b2 */
2N/A else borrow = 0xfffe0000 | (borrow >> fifteen);
2N/A#ifdef DEBUGDSB
2N/A (void) printf("4 borrow=%x %d %d %d\n", borrow, (*aptr * qx),
2N/A *bptr, *aptr);
2N/A#endif
2N/A }
2N/A borrow += *bptr;
2N/A *bptr = (short)(borrow & 077777);
2N/A if (borrow >= 0) s3b2shit = borrow >> fifteen; /* 3b2 */
2N/A else s3b2shit = 0xfffe0000 | (borrow >> fifteen);
2N/A if (s3b2shit == 0) {
2N/A#ifdef DEBUGDSB
2N/A (void) printf("mdsb 0\n");
2N/A#endif
2N/A return (0);
2N/A }
2N/A borrow = 0;
2N/A aptr = a;
2N/A bptr = b;
2N/A for (j = n; j > 0; j--) {
2N/A borrow += *aptr++ + *bptr;
2N/A *bptr++ = (short)(borrow & 077777);
2N/A if (borrow >= 0) borrow >>= fifteen; /* 3b2 */
2N/A else borrow = 0xfffe0000 | (borrow >>fifteen);
2N/A }
2N/A#ifdef DEBUGDSB
2N/A (void) printf("mdsb 1\n");
2N/A#endif
2N/A return (1);
2N/A}
2N/A
2N/Astatic int
2N/Am_trq(short v1, short v2, short u1, short u2, short u3)
2N/A{
2N/A short d;
2N/A int x1;
2N/A int c1;
2N/A
2N/A c1 = u1 * 0100000 + u2;
2N/A if (u1 == v1) {
2N/A d = 077777;
2N/A } else {
2N/A d = (short)(c1 / v1);
2N/A }
2N/A do {
2N/A x1 = c1 - v1 * d;
2N/A x1 = x1 * 0100000 + u3 - v2 * d;
2N/A --d;
2N/A } while (x1 < 0);
2N/A#ifdef DEBUGMTRQ
2N/A (void) printf("mtrq %d %d %d %d %d %d\n", v1, v2, u1, u2, u3, (d+1));
2N/A#endif
2N/A return ((int)d + 1);
2N/A}
2N/A
2N/Astatic void
2N/Am_div(MINT *a, MINT *b, MINT *q, MINT *r)
2N/A{
2N/A MINT u, v, x, w;
2N/A short d;
2N/A short *qval;
2N/A short *uval;
2N/A int j;
2N/A int qq;
2N/A int n;
2N/A short v1;
2N/A short v2;
2N/A
2N/A u.len = v.len = x.len = w.len = 0;
2N/A if (b->len == 0) {
2N/A _mp_fatal("mdiv divide by zero");
2N/A return;
2N/A }
2N/A if (b->len == 1) {
2N/A r->val = _mp_xalloc(1, "m_div1");
2N/A mp_sdiv(a, b->val[0], q, r->val);
2N/A if (r->val[0] == 0) {
2N/A free(r->val);
2N/A r->len = 0;
2N/A } else {
2N/A r->len = 1;
2N/A }
2N/A return;
2N/A }
2N/A if (a -> len < b -> len) {
2N/A q->len = 0;
2N/A r->len = a->len;
2N/A r->val = _mp_xalloc(r->len, "m_div2");
2N/A for (qq = 0; qq < r->len; qq++) {
2N/A r->val[qq] = a->val[qq];
2N/A }
2N/A return;
2N/A }
2N/A x.len = 1;
2N/A x.val = &d;
2N/A n = b->len;
2N/A d = 0100000 / (b->val[n - 1] + 1);
2N/A mp_mult(a, &x, &u); /* subtle: relies on mult allocing extra space */
2N/A mp_mult(b, &x, &v);
2N/A#ifdef DEBUG_MDIV
2N/A (void) printf(" u=%s\n", mtox(&u));
2N/A (void) printf(" v=%s\n", mtox(&v));
2N/A#endif
2N/A v1 = v.val[n - 1];
2N/A v2 = v.val[n - 2];
2N/A qval = _mp_xalloc(a -> len - n + 1, "m_div3");
2N/A uval = u.val;
2N/A for (j = a->len - n; j >= 0; j--) {
2N/A qq = m_trq(v1, v2, uval[j + n], uval[j + n - 1],
2N/A uval[j + n - 2]);
2N/A if (m_dsb(qq, n, v.val, uval + j))
2N/A qq -= 1;
2N/A qval[j] = (short)qq;
2N/A }
2N/A x.len = n;
2N/A x.val = u.val;
2N/A _mp_mcan(&x);
2N/A#ifdef DEBUG_MDIV
2N/A (void) printf(" x=%s\n", mtox(&x));
2N/A (void) printf(" d(in)=%d\n", (d));
2N/A#endif
2N/A mp_sdiv(&x, d, &w, &d);
2N/A#ifdef DEBUG_MDIV
2N/A (void) printf(" w=%s\n", mtox(&w));
2N/A (void) printf(" d(out)=%d\n", (d));
2N/A#endif
2N/A r->len = w.len;
2N/A r->val = w.val;
2N/A q->val = qval;
2N/A qq = a->len - n + 1;
2N/A if (qq > 0 && qval[qq - 1] == 0)
2N/A qq -= 1;
2N/A q->len = qq;
2N/A if (qq == 0)
2N/A free(qval);
2N/A if (x.len != 0)
2N/A _mp_xfree(&u);
2N/A _mp_xfree(&v);
2N/A}