BSD 4 release

[unix-history] / usr / src / cmd / lisp / divbig.c

static char *sccsid = "@(#)divbig.c     34.1 10/3/80";

#include "global.h"

#define b 0x40000000
#define toint(p) ((int) (p))

divbig(dividend, divisor, quotient, remainder)
lispval dividend, divisor, *quotient, *remainder;
{
       register *ujp, *vip;
       int *sp(), *alloca(), d, negflag = 0, m, n, carry, rem, qhat, j;
       int borrow, negrem = 0;
       int *utop = sp(), *ubot, *vbot, *qbot;
       register lispval work; lispval export();
       snpand(3);

       /* copy dividend */
       for(work = dividend; work; work = work ->s.CDR)
               stack(work->s.I);
       ubot = sp();
       if(*ubot < 0) {         /* knuth's division alg works only for pos
                                       bignums                         */
               negflag ^= 1;
               negrem = 1;
               dsmult(utop-1,ubot,-1);
       }
       stack(0);
       ubot = sp();

       
       /*copy divisor */
       for(work = divisor; work; work = work->s.CDR)
               stack(work->s.I);

       vbot = sp();
       stack(0);
       if(*vbot < 0) {
               negflag ^= 1;
               dsmult(ubot-1,vbot,-1);
       }

       /* check validity of data */
       n = ubot - vbot;
       m = utop - ubot - n - 1;
       if (n == 1) {
               /* do destructive division by  a single. */
               rem = dsdiv(utop-1,ubot,*vbot);
               if(negrem)
                       rem = -rem;
               if(negflag)
                       dsmult(utop-1,ubot,-1);
               if(remainder)
                       *remainder = inewint(rem);
               if(quotient)
                       *quotient = export(utop,ubot);
               return;
       }
       if (m < 0) {
               if (remainder)
                       *remainder = dividend;
               if(quotient)
                       *quotient = inewint(0);
               return;
       }
       qbot = alloca(toint(utop) + toint(vbot) - 2 * toint(ubot));
d1:
       d = b /(*vbot +1);
       dsmult(utop-1,ubot,d);
       dsmult(ubot-1,vbot,d);

d2:     for(j=0,ujp=ubot; j <= m; j++,ujp++) {

       d3:     
               qhat = calqhat(ujp,vbot);
       d4:
               if((borrow = mlsb(ujp + n, ujp, ubot, -qhat)) < 0) {
                       adback(ujp + n, ujp, ubot);
                       qhat--;
               }
               qbot[j] = qhat;
       }
d8:     if(remainder) {
               dsdiv(utop, utop - n, d);
               if(negrem) dsmult(utop-1,utop-n,-1);
               *remainder = export(utop,utop-n);
       }
       if(quotient) {
               if(negflag)
                       dsmult(qbot+m,qbot,-1);
               *quotient = export(qbot + m + 1, qbot);
       }
}
/*
* asm code commented out due to optimizer bug
calqhat(ujp,v1p)
register int *ujp, *v1p;
{
asm("   movl    $0x3fffffff,r0");
asm("   cmpl    (r10),(r11)");
asm("   beql    on1");
asm("   emul    (r11),$0x40000000,4(r11),r1");
asm("   ediv    (r10),r1,r0,r5");
asm("on1:");
asm("   emul    r0,4(r10),$0,r1");
asm("   emul    r5,$0x40000000,8(r11),r3");
asm("   subl2   r3,r1");
asm("   sbwc    r4,r2");
asm("   bleq    out1");
asm("   decl    r0");
asm("out1:");
}
mlsb(utop,ubot,vtop,nqhat)
register int *utop, *ubot, *vtop;
register int nqhat;
{
asm("   clrl    r0");
asm("loop2:     addl2   (r11),r0");
asm("   emul    r8,-(r9),r0,r2");
asm("   extzv   $0,$30,r2,(r11)");
asm("   extv    $30,$32,r2,r0");
asm("   acbl    r10,$-4,r11,loop2");
}
adback(utop,ubot,vtop)
register int *utop, *ubot, *vtop;
{
asm("   clrl    r0");
asm("loop3:     addl2   -(r9),r0");
asm("   addl2   (r11),r0");
asm("   extzv   $0,$30,r0,(r11)");
asm("   extv    $30,$2,r0,r0");
asm("   acbl    r10,$-4,r11,loop3");
}
dsdiv(top,bot,div)
register int* bot;
{
asm("   clrl    r0");
asm("loop4:     emul    r0,$0x40000000,(r11),r1");
asm("   ediv    12(ap),r1,(r11),r0");
asm("   acbl    4(ap),$4,r11,loop4");
}
dsmult(top,bot,mult)
register int* top;
{
asm("   clrl    r0");
asm("loop5:     emul    12(ap),(r11),r0,r1");
asm("   extzv   $0,$30,r1,(r11)");
asm("   extv    $30,$32,r1,r0");
asm("   acbl    8(ap),$-4,r11,loop5");
asm("   movl    r1,4(r11)");
}
lispval export(top,bot)
register lispval bot;
{
       register r10, r9, r8, r7, r6;
asm("   movl    4(ap),r10");
asm("   movl    $0xC0000000,r4");
asm("   jmp     Bexport");
}
*/

#define MAXINT 0x8000000L

Ihau(fix)
register int fix;
{
       register count;
       if(fix==MAXINT)
               return(32);
       if(fix < 0)
               fix = -fix;
       for(count = 0; fix; count++)
               fix /= 2;
       return(count);
}
lispval
Lhau()
{
       register count;
       register lispval handy;
       register dum1,dum2;
       register struct argent *lbot, *np;
       lispval Labsval();

       handy = lbot->val;
top:
       switch(TYPE(handy)) {
       case INT:
               count = Ihau(handy->i);
               break;
       case SDOT:
               lbot->val = Labsval();
               for(count = 0; handy->s.CDR!=((lispval) 0); handy = handy->s.CDR)
                       count += 30;
               count += Ihau(handy->s.I);
               break;
       default:
               handy = errorh(Vermisc,"Haulong: bad argument",nil,
                              TRUE,997,handy);
               goto top;
       }
       return(inewint(count));
}
lispval
Lhaipar()
{
       int *sp();
       register lispval work;
       register n;
       register int *top = sp() - 1;
       register int *bot;
       register struct argent *lbot, *np;
       int mylen;

       /*chkarg(2);*/
       work = lbot->val;
                                       /* copy data onto stack */
on1:
       switch(TYPE(work)) {
       case INT:
               stack(work->i);
               break;
       case SDOT:
               for(; work!=((lispval) 0); work = work->s.CDR)
                       stack(work->s.I);
               break;
       default:
               work = errorh(Vermisc,"Haipart: bad first argument",nil,
                               TRUE,996,work);
               goto on1;
       }
       bot = sp();
       if(*bot < 0) {
               stack(0);
               dsmult(top,bot,-1);
               bot--;
       }
       for(; *bot==0 && bot < top; bot++);
                               /* recalculate haulong internally */
       mylen = (top - bot) * 30 + Ihau(*bot);
                               /* get second argument            */
       work = lbot[1].val;
       while(TYPE(work)!=INT)
               work = errorh(Vermisc,"Haipart: 2nd arg not int",nil,
                               TRUE,995,work);
       n = work->i;
       if(n >= mylen || -n >= mylen)
               goto done;
       if(n==0) return(inewint(0));
       if(n > 0) {
                               /* Here we want n most significant bits
                                  so chop off mylen - n bits */
               stack(0);
               n = mylen - n;
               for(n; n >= 30; n -= 30)
                       top--;
               if(top < bot)
                       error("Internal error in haipart #1",FALSE);
               dsdiv(top,bot,1<<n);

       } else {
                               /* here we want abs(n) low order bits */
               stack(0);
               bot = top + 1;
               for(; n <= 0; n += 30)
                       bot--;
               n = 30 - n;
               *bot &= ~ (-1<<n);
       }
done:
       return(export(top + 1,bot));
}