[BACK]Return to ddN.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / asm

Diff for /OpenXM_contrib2/asir2000/asm/ddN.c between version 1.5 and 1.14

version 1.5, 2006/08/09 02:43:37 version 1.14, 2018/03/29 01:32:50
Line 45 
Line 45 
  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,   * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.   * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
  *   *
  * $OpenXM: OpenXM_contrib2/asir2000/asm/ddN.c,v 1.4 2001/10/09 01:36:03 noro Exp $   * $OpenXM: OpenXM_contrib2/asir2000/asm/ddN.c,v 1.13 2015/08/29 04:15:04 fujimoto Exp $
 */  */
 #ifndef FBASE  #ifndef FBASE
 #define FBASE  #define FBASE
Line 55 
Line 55 
 #include "base.h"  #include "base.h"
 #include "inline.h"  #include "inline.h"
   
   #if defined(__GNUC__)
   unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)  __attribute__ ((noinline));
   void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)  __attribute__ ((noinline));
   #endif
   
 void divn(N n1,N n2,N *nq,N *nr)  void divn(N n1,N n2,N *nq,N *nr)
 {  {
         int tmp,b;    int tmp,b;
         int i,j,d1,d2,dd;    int i,j,d1,d2,dd;
         unsigned int *m1,*m2;    unsigned int *m1,*m2;
         unsigned int uq,ur,msw;    unsigned int uq,ur,msw;
         N q,r,t1,t2;    N q,r,t1,t2;
   
         if ( !n2 )    if ( !n2 )
                 error("divn: divsion by 0");      error("divn: divsion by 0");
         else if ( !n1 ) {    else if ( !n1 ) {
                 *nq = 0; *nr = 0;      *nq = 0; *nr = 0;
         } else if ( UNIN(n2) ) {    } else if ( UNIN(n2) ) {
                 COPY(n1,*nq); *nr = 0;      COPY(n1,*nq); *nr = 0;
         } else if ( (d1 = PL(n1)) < (d2 = PL(n2)) ) {    } else if ( (d1 = PL(n1)) < (d2 = PL(n2)) ) {
                 COPY(n1,*nr); *nq = 0;      COPY(n1,*nr); *nq = 0;
         } else if ( d2 == 1 ) {    } else if ( d2 == 1 ) {
                 if ( d1 == 1 ) {      if ( d1 == 1 ) {
                         DQR(BD(n1)[0],BD(n2)[0],uq,ur)        DQR(BD(n1)[0],BD(n2)[0],uq,ur)
                         STON(uq,*nq); STON(ur,*nr);        STON(uq,*nq); STON(ur,*nr);
                 } else {      } else {
                         ur = divin(n1,BD(n2)[0],nq); STON(ur,*nr);        ur = divin(n1,BD(n2)[0],nq); STON(ur,*nr);
                 }      }
                 return;      return;
         } else if ( (d1 == d2) && ((tmp = cmpn(n1,n2)) <= 0) )    } else if ( (d1 == d2) && ((tmp = cmpn(n1,n2)) <= 0) )
                 if ( !tmp ) {      if ( !tmp ) {
                         *nr = 0; COPY(ONEN,*nq);        *nr = 0; COPY(ONEN,*nq);
                 } else {      } else {
                         COPY(n1,*nr); *nq = 0;        COPY(n1,*nr); *nq = 0;
                 }      }
         else {    else {
                 msw = BD(n2)[d2-1];      msw = BD(n2)[d2-1];
                 for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;      for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
                 b = BSH-i;      b = BSH-i;
                 if ( b ) {      if ( b ) {
                         bshiftn(n1,-b,&t1); bshiftn(n2,-b,&t2);        bshiftn(n1,-b,&t1); bshiftn(n2,-b,&t2);
                 } else {      } else {
                         COPY(n1,t1); COPY(n2,t2);        COPY(n1,t1); COPY(n2,t2);
                 }      }
                 m1 = (unsigned int *)ALLOCA((PL(t1)+1)*sizeof(unsigned int));      m1 = (unsigned int *)ALLOCA((PL(t1)+1)*sizeof(unsigned int));
                 bcopy(BD(t1),m1,PL(t1)*sizeof(unsigned int)); m1[PL(t1)] = 0;      bcopy(BD(t1),m1,PL(t1)*sizeof(unsigned int)); m1[PL(t1)] = 0;
                 m2 = BD(t2); dd = d1-d2;      m2 = BD(t2); dd = d1-d2;
                 *nq = q = NALLOC(dd+1); INITRC(q);      *nq = q = NALLOC(dd+1); INITRC(q);
                 divnmain(d1,d2,m1,m2,BD(q)); FREEN(t1); FREEN(t2);      divnmain(d1,d2,m1,m2,BD(q)); FREEN(t1); FREEN(t2);
                 PL(q) = (BD(q)[dd]?dd+1:dd);      PL(q) = (BD(q)[dd]?dd+1:dd);
                 if ( !PL(q) ) {      if ( !PL(q) ) {
                         FREEN(q);        FREEN(q);
                 }      }
                 for ( j = d2-1; (j >= 0) && (m1[j] == 0); j--);      for ( j = d2-1; (j >= 0) && (m1[j] == 0); j--);
                 if ( j == -1 )      if ( j == -1 )
                         *nr = 0;        *nr = 0;
                 else {      else {
                         r = NALLOC(j+1); INITRC(r); PL(r) = j+1;        r = NALLOC(j+1); INITRC(r); PL(r) = j+1;
                         bcopy(m1,BD(r),(j+1)*sizeof(unsigned int));        bcopy(m1,BD(r),(j+1)*sizeof(unsigned int));
                         if ( b ) {        if ( b ) {
                                 bshiftn(r,b,nr); FREEN(r);          bshiftn(r,b,nr); FREEN(r);
                         } else        } else
                                 *nr = r;          *nr = r;
                 }      }
         }    }
 }  }
   
 void divsn(N n1,N n2,N *nq)  void divsn(N n1,N n2,N *nq)
 {  {
         int d1,d2,dd,i,b;    int d1,d2,dd,i,b;
         unsigned int *m1,*m2;    unsigned int *m1,*m2;
         unsigned int uq,msw;    unsigned int uq,msw;
         N q,t1,t2;    N q,t1,t2;
   
         if ( !n1 )    if ( !n1 )
                 *nq = 0;      *nq = 0;
         else if ( UNIN(n2) ) {    else if ( UNIN(n2) ) {
                 COPY(n1,*nq);      COPY(n1,*nq);
         } else if ( (d1 = PL(n1)) < (d2 = PL(n2)) )    } else if ( (d1 = PL(n1)) < (d2 = PL(n2)) )
                 error("divsn: cannot happen");      error("divsn: cannot happen");
         else if ( d2 == 1 ) {    else if ( d2 == 1 ) {
                 if ( d1 == 1 ) {      if ( d1 == 1 ) {
                         uq = BD(n1)[0] / BD(n2)[0]; STON(uq,*nq);        uq = BD(n1)[0] / BD(n2)[0]; STON(uq,*nq);
                 } else      } else
                         divin(n1,BD(n2)[0],nq);        divin(n1,BD(n2)[0],nq);
         } else {    } else {
                 msw = BD(n2)[d2-1];      msw = BD(n2)[d2-1];
                 for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;      for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
                 b = BSH-i;      b = BSH-i;
                 if ( b ) {      if ( b ) {
                         bshiftn(n1,-b,&t1); bshiftn(n2,-b,&t2);        bshiftn(n1,-b,&t1); bshiftn(n2,-b,&t2);
                 } else {      } else {
                         COPY(n1,t1); COPY(n2,t2);        COPY(n1,t1); COPY(n2,t2);
                 }      }
                 m1 = (unsigned int *)ALLOCA((PL(t1)+1)*sizeof(unsigned int));      m1 = (unsigned int *)ALLOCA((PL(t1)+1)*sizeof(unsigned int));
                 bcopy(BD(t1),m1,PL(t1)*sizeof(unsigned int)); m1[PL(t1)] = 0;      bcopy(BD(t1),m1,PL(t1)*sizeof(unsigned int)); m1[PL(t1)] = 0;
                 m2 = BD(t2); dd = d1-d2;      m2 = BD(t2); dd = d1-d2;
                 *nq = q = NALLOC(dd+1); INITRC(q);      *nq = q = NALLOC(dd+1); INITRC(q);
                 divnmain(d1,d2,m1,m2,BD(q));      divnmain(d1,d2,m1,m2,BD(q));
                 FREEN(t1); FREEN(t2);      FREEN(t1); FREEN(t2);
                 PL(q) = (BD(q)[dd]?dd+1:dd);      PL(q) = (BD(q)[dd]?dd+1:dd);
                 if ( !PL(q) )      if ( !PL(q) )
                         FREEN(q);        FREEN(q);
         }    }
 }  }
   
 void remn(N n1,N n2,N *nr)  void remn(N n1,N n2,N *nr)
 {  {
         int d1,d2,tmp;    int d1,d2,tmp;
         unsigned int uq,ur;    unsigned int uq,ur;
         N q;    N q;
   
         if ( !n2 )    if ( !n2 )
                 error("remn: divsion by 0");      error("remn: divsion by 0");
         else if ( !n1 || UNIN(n2) )    else if ( !n1 || UNIN(n2) )
                 *nr = 0;      *nr = 0;
         else if ( (d1 = PL(n1)) < (d2 = PL(n2)) ) {    else if ( (d1 = PL(n1)) < (d2 = PL(n2)) ) {
                 COPY(n1,*nr);      COPY(n1,*nr);
         } else if ( d2 == 1 ) {    } else if ( d2 == 1 ) {
                 if ( d1 == 1 ) {      if ( d1 == 1 ) {
                         DQR(BD(n1)[0],BD(n2)[0],uq,ur)        DQR(BD(n1)[0],BD(n2)[0],uq,ur)
                         STON(ur,*nr);        STON(ur,*nr);
                 } else {      } else {
                         ur = rem(n1,BD(n2)[0]); UTON(ur,*nr);        ur = rem(n1,BD(n2)[0]); UTON(ur,*nr);
                 }      }
                 return;      return;
         } else if ( (d1 == d2) && ((tmp = cmpn(n1,n2)) <= 0) )    } else if ( (d1 == d2) && ((tmp = cmpn(n1,n2)) <= 0) )
                 if ( !tmp ) {      if ( !tmp ) {
                         *nr = 0;        *nr = 0;
                 } else {      } else {
                         COPY(n1,*nr);        COPY(n1,*nr);
                 }      }
         else    else
                 divn(n1,n2,&q,nr);      divn(n1,n2,&q,nr);
 }  }
   
 /* d = 2^(32*words)-lower */  /* d = 2^(32*words)-lower */
   
 void remn_special(N a,N d,int bits,unsigned int lower,N *b)  void remn_special(N a,N d,int bits,unsigned int lower,N *b)
 {  {
         int words;    int words;
         N r;    N r;
         int rl,i;    int rl,i;
         unsigned int c,tmp;    unsigned int c,tmp;
         unsigned int *rb,*src,*dst;    unsigned int *rb,*src,*dst;
   
         if ( cmpn(a,d) < 0 ) {    if ( cmpn(a,d) < 0 ) {
                 *b = a;      *b = a;
                 return;      return;
         }    }
         words = bits>>5;    words = bits>>5;
         r = NALLOC(PL(a)); dupn(a,r);    r = NALLOC(PL(a)); dupn(a,r);
         rl = PL(r); rb = BD(r);    rl = PL(r); rb = BD(r);
         while ( rl > words ) {    while ( rl > words ) {
                 src = rb+words;      src = rb+words;
                 dst = rb;      dst = rb;
                 i = rl-words;      i = rl-words;
                 for ( c = 0; i > 0; i--, src++, dst++ ) {      for ( c = 0; i > 0; i--, src++, dst++ ) {
                         DMA2(*src,lower,c,*dst,c,*dst)        DMA2(*src,lower,c,*dst,c,*dst)
                         *src = 0;        *src = 0;
                 }      }
                 for ( ; c; dst++ ) {      for ( ; c; dst++ ) {
                         tmp = *dst + c;        tmp = *dst + c;
                         c = tmp < c ? 1 : 0;        c = tmp < c ? 1 : 0;
                         *dst = tmp;        *dst = tmp;
                 }      }
                 rl = dst-rb;      rl = dst-rb;
         }    }
         for ( i = words-1; i >= 0 && !rb[i]; i-- );    for ( i = words-1; i >= 0 && !rb[i]; i-- );
         PL(r) = i + 1;    PL(r) = i + 1;
         if ( cmpn(r,d) >= 0 ) {    if ( cmpn(r,d) >= 0 ) {
                 tmp = rb[0]-BD(d)[0];      tmp = rb[0]-BD(d)[0];
                 UTON(tmp,*b);      UTON(tmp,*b);
         } else    } else
                 *b = r;      *b = r;
 }  }
 void mulin(N n,unsigned int d,unsigned int *p)  void mulin(N n,unsigned int d,unsigned int *p)
 {  {
         unsigned int carry;    unsigned int carry;
         unsigned int *m,*me;    unsigned int *m,*me;
   
         bzero((char *)p,(int)((PL(n)+1)*sizeof(int)));    bzero((char *)p,(int)((PL(n)+1)*sizeof(int)));
         for ( carry = 0, m = BD(n), me = m+PL(n); m < me; p++, m++ ) {    for ( carry = 0, m = BD(n), me = m+PL(n); m < me; p++, m++ ) {
                 DMA2(*m,d,*p,carry,carry,*p)      DMA2(*m,d,*p,carry,carry,*p)
         }    }
         *p = carry;    *p = carry;
 }  }
   
 unsigned int divin(N n,unsigned int dvr,N *q)  unsigned int divin(N n,unsigned int dvr,N *q)
 {  {
         int d;    int d;
         unsigned int up;    unsigned int up;
         unsigned int *m,*mq,*md;    unsigned int *m,*mq,*md;
         N nq;    N nq;
   
         d = PL(n); m = BD(n);    d = PL(n); m = BD(n);
         if ( ( d == 1 ) && ( dvr > *m ) ) {    if ( ( d == 1 ) && ( dvr > *m ) ) {
                 *q = 0;      *q = 0;
                 return *m;      return *m;
         }    }
         *q = nq = NALLOC(d); INITRC(nq);    *q = nq = NALLOC(d); INITRC(nq);
         for ( md = m+d-1, mq = BD(nq)+d-1, up = 0; md >= m; md--, mq-- ) {    for ( md = m+d-1, mq = BD(nq)+d-1, up = 0; md >= m; md--, mq-- ) {
                 DSAB(dvr,up,*md,*mq,up)      DSAB(dvr,up,*md,*mq,up)
         }    }
         PL(nq) = (BD(nq)[d-1]?d:d-1);    PL(nq) = (BD(nq)[d-1]?d:d-1);
         if ( !PL(nq) )    if ( !PL(nq) )
                 FREEN(nq);      FREEN(nq);
         return ( up );    return ( up );
 }  }
   
 void bprintn(N n)  void bprintn(N n)
 {  {
         int l,i;    int l,i;
         unsigned int *b;    unsigned int *b;
   
         if ( !n )    if ( !n )
                 printf("0");      printf("0");
         else {    else {
                 l = PL(n); b = BD(n);      l = PL(n); b = BD(n);
                 for ( i = l-1; i >= 0; i-- )      for ( i = l-1; i >= 0; i-- )
                         printf("%010u|",b[i]);        printf("%010u|",b[i]);
         }    }
 }  }
   
 void bxprintn(N n)  void bxprintn(N n)
 {  {
         int l,i;    int l,i;
         unsigned int *b;    unsigned int *b;
   
         if ( !n )    if ( !n )
                 printf("0");      printf("0");
         else {    else {
                 l = PL(n); b = BD(n);      l = PL(n); b = BD(n);
                 for ( i = l-1; i >= 0; i-- )      for ( i = l-1; i >= 0; i-- )
                         printf("%08x|",b[i]);        printf("%08x|",b[i]);
         }    }
 }  }
   
 #if defined(VISUAL) || defined(i386)  #if (defined(_M_IX86) || defined(i386)) && !defined(__MINGW32__)
 void muln(N n1,N n2,N *nr)  void muln(N n1,N n2,N *nr)
 {  {
         unsigned int tmp,carry,mul;    unsigned int tmp,carry,mul;
         unsigned int *p1,*m1,*m2;    unsigned int *p1,*m1,*m2;
         int i,d1,d2,d;    int i,d1,d2,d;
         N r;    N r;
   
         if ( !n1 || !n2 )    if ( !n1 || !n2 )
                 *nr = 0;      *nr = 0;
         else if ( UNIN(n1) )    else if ( UNIN(n1) )
                 COPY(n2,*nr);      COPY(n2,*nr);
         else if ( UNIN(n2) )    else if ( UNIN(n2) )
                 COPY(n1,*nr);      COPY(n1,*nr);
         else if ( (PL(n1) == 1) && (PL(n2) == 1) ) {    else if ( (PL(n1) == 1) && (PL(n2) == 1) ) {
                 DM(*BD(n1),*BD(n2),carry,tmp)      DM(*BD(n1),*BD(n2),carry,tmp)
                 if ( carry ) {      if ( carry ) {
                         *nr = r = NALLOC(2); INITRC(r);        *nr = r = NALLOC(2); INITRC(r);
                         PL(r) = 2; p1 = BD(r); *p1++ = tmp; *p1 = carry;        PL(r) = 2; p1 = BD(r); *p1++ = tmp; *p1 = carry;
                 } else      } else
                         STON(tmp,*nr);        STON(tmp,*nr);
         } else {    } else {
                 d1 = PL(n1); d2 = PL(n2);      d1 = PL(n1); d2 = PL(n2);
                 d = d1+d2;      d = d1+d2;
                 *nr = r = NALLOC(d); INITRC(r);      *nr = r = NALLOC(d); INITRC(r);
                 bzero((char *)BD(r),(int)((d1+d2)*sizeof(int)));      bzero((char *)BD(r),(int)((d1+d2)*sizeof(int)));
                 for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )      for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
                         if ( mul = *m2 )        if ( mul = *m2 )
                                 muln_1(m1,d1,mul,BD(r)+i);          muln_1(m1,d1,mul,BD(r)+i);
                 PL(r) = (BD(r)[d-1]?d:d-1);      PL(r) = (BD(r)[d-1]?d:d-1);
         }    }
 }  }
   
 void _muln(N n1,N n2,N nr)  void _muln(N n1,N n2,N nr)
 {  {
         unsigned int mul;    unsigned int mul;
         unsigned int *m1,*m2;    unsigned int *m1,*m2;
         int i,d1,d2,d;    int i,d1,d2,d;
   
         if ( !n1 || !PL(n1) || !n2 || !PL(n2) )    if ( !n1 || !PL(n1) || !n2 || !PL(n2) )
                 PL(nr) = 0;      PL(nr) = 0;
         else if ( UNIN(n1) )    else if ( UNIN(n1) )
                 dupn(n2,nr);      dupn(n2,nr);
         else if ( UNIN(n2) )    else if ( UNIN(n2) )
                 dupn(n1,nr);      dupn(n1,nr);
         else {    else {
                 d1 = PL(n1); d2 = PL(n2);      d1 = PL(n1); d2 = PL(n2);
                 d = d1+d2;      d = d1+d2;
                 bzero((char *)BD(nr),(int)((d1+d2)*sizeof(int)));      bzero((char *)BD(nr),(int)((d1+d2)*sizeof(int)));
                 for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )      for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
                         if ( mul = *m2 )        if ( mul = *m2 )
                                 muln_1(m1,d1,mul,BD(nr)+i);          muln_1(m1,d1,mul,BD(nr)+i);
                 PL(nr) = (BD(nr)[d-1]?d:d-1);      PL(nr) = (BD(nr)[d-1]?d:d-1);
         }    }
 }  }
   
 void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)  void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
 {  {
         /* esi : p, edi : r, carry : ebx, s : ecx */    /* esi : p, edi : r, carry : ebx, s : ecx */
 #if defined(VISUAL)  #if defined(_M_IX86)
         __asm {    __asm {
         push esi    push esi
         push edi    push edi
         mov esi,p    mov esi,p
         mov edi,r    mov edi,r
         mov ecx,s    mov ecx,s
         xor ebx,ebx    xor ebx,ebx
         Lstart_muln:    Lstart_muln:
         mov eax,DWORD PTR [esi]    mov eax,DWORD PTR [esi]
         mul d    mul d
         add eax,DWORD PTR [edi]    add eax,DWORD PTR [edi]
         adc edx,0    adc edx,0
         add eax,ebx    add eax,ebx
         adc edx,0    adc edx,0
         mov DWORD PTR[edi],eax    mov DWORD PTR[edi],eax
         mov ebx,edx    mov ebx,edx
         lea     esi,DWORD PTR [esi+4]    lea  esi,DWORD PTR [esi+4]
         lea     edi,DWORD PTR [edi+4]    lea  edi,DWORD PTR [edi+4]
         dec ecx    dec ecx
         jnz Lstart_muln    jnz Lstart_muln
         mov DWORD PTR [edi],ebx    mov DWORD PTR [edi],ebx
         pop edi    pop edi
         pop esi    pop esi
         }    }
 #else  #else
         asm volatile("\    asm volatile("\
         movl    %0,%%esi;\    pushl  %%ebx;\
         movl    %1,%%edi;\    movl  %0,%%esi;\
         movl    $0,%%ebx;\    movl  %1,%%edi;\
         Lstart_muln:;\    movl  $0,%%ebx;\
         movl    (%%esi),%%eax;\    Lstart_muln:;\
         mull    %2;\    movl  (%%esi),%%eax;\
         addl    (%%edi),%%eax;\    mull  %2;\
         adcl    $0,%%edx;\    addl  (%%edi),%%eax;\
         addl    %%ebx,%%eax;\    adcl  $0,%%edx;\
         adcl    $0,%%edx;\    addl  %%ebx,%%eax;\
         movl    %%eax,(%%edi);\    adcl  $0,%%edx;\
         movl    %%edx,%%ebx;\    movl  %%eax,(%%edi);\
         leal    4(%%esi),%%esi;\    movl  %%edx,%%ebx;\
         leal    4(%%edi),%%edi;\    leal  4(%%esi),%%esi;\
         decl    %3;\    leal  4(%%edi),%%edi;\
         jnz             Lstart_muln;\    decl  %3;\
         movl    %%ebx,(%%edi)"\    jnz    Lstart_muln;\
         :\    movl  %%ebx,(%%edi);\
         :"m"(p),"m"(r),"m"(d),"m"(s)\    popl  %%ebx"\
         :"eax","ebx","edx","esi","edi");    :\
     :"m"(p),"m"(r),"m"(d),"m"(s)\
     :"eax","edx","esi","edi");
 #endif  #endif
 }  }
   
 void divnmain(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)  void divnmain(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)
 {  {
         int i,j;    int i,j;
         UL r,ltmp;    UL r,ltmp;
         unsigned int l,ur;    unsigned int l,ur;
         unsigned int *n1,*n2;    unsigned int *n1,*n2;
         unsigned int u,qhat;    unsigned int u,qhat;
         unsigned int v1,v2;    unsigned int v1,v2;
   
         v1 = m2[d2-1]; v2 = m2[d2-2];    v1 = m2[d2-1]; v2 = m2[d2-2];
 #if 0  #if 0
         if ( v1 == (1<<31) && !v2 ) {    if ( v1 == (1<<31) && !v2 ) {
                 divnmain_special(d1,d2,m1,m2,q);      divnmain_special(d1,d2,m1,m2,q);
                 return;      return;
         }    }
 #endif  #endif
         for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {    for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {
                 n1 = m1+d2; n2 = m1+d2-1;      n1 = m1+d2; n2 = m1+d2-1;
                 if ( *n1 == v1 ) {         if ( *n1 == v1 ) {
                         qhat = ULBASE - 1;              qhat = ULBASE - 1;
                         r = (UL)*n1 + (UL)*n2;        r = (UL)*n1 + (UL)*n2;
                 } else {         } else {
                         DSAB(v1,*n1,*n2,qhat,ur)        DSAB(v1,*n1,*n2,qhat,ur)
                         r = (UL)ur;        r = (UL)ur;
                 }      }
                 DM(v2,qhat,u,l)      DM(v2,qhat,u,l)
                 while ( 1 ) {         while ( 1 ) {
                         if ((r > (UL)u) || ((r == (UL)u) && (*(n1-2) >= l)))             if ((r > (UL)u) || ((r == (UL)u) && (*(n1-2) >= l)))
                                 break;          break;
                         if ( l >= v2 )        if ( l >= v2 )
                                 l -= v2;          l -= v2;
                         else {        else {
                                 l = (unsigned int)((UL)l+(ULBASE-(UL)v2)); u--;          l = (unsigned int)((UL)l+(ULBASE-(UL)v2)); u--;
                         }        }
                 r += v1; qhat--;            r += v1; qhat--;
                 }         }
                 if ( qhat ) {      if ( qhat ) {
                         u = divn_1(m2,d2,qhat,m1);        u = divn_1(m2,d2,qhat,m1);
                         if ( *(m1+d2) < u ) {        if ( *(m1+d2) < u ) {
                                 for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;          for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;
                                         i > 0; i--, n1++, n2++ ) {            i > 0; i--, n1++, n2++ ) {
                                         ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;            ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;
                                         *n1 = (unsigned int)(ltmp & BMASK);            *n1 = (unsigned int)(ltmp & BMASK);
                                         ur = (unsigned int)(ltmp >> BSH);            ur = (unsigned int)(ltmp >> BSH);
                                 }          }
                         }        }
                         *n1 = 0;        *n1 = 0;
                 }      }
                 *q = qhat;      *q = qhat;
         }    }
 }  }
   
 void divnmain_special(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)  void divnmain_special(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)
 {  {
         int i,j;    int i,j;
         UL ltmp;    UL ltmp;
         unsigned int ur;    unsigned int ur;
         unsigned int *n1,*n2;    unsigned int *n1,*n2;
         unsigned int u,qhat;    unsigned int u,qhat;
         unsigned int v1,v2;    unsigned int v1,v2;
   
         v1 = m2[d2-1]; v2 = 0;    v1 = m2[d2-1]; v2 = 0;
         for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {    for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {
                 n1 = m1+d2; n2 = m1+d2-1;      n1 = m1+d2; n2 = m1+d2-1;
                 qhat = ((*n1)<<1)|((*n2)>>31);      qhat = ((*n1)<<1)|((*n2)>>31);
                 if ( qhat ) {      if ( qhat ) {
                         u = divn_1(m2,d2,qhat,m1);        u = divn_1(m2,d2,qhat,m1);
                         if ( *(m1+d2) < u ) {        if ( *(m1+d2) < u ) {
                                 for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;          for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;
                                         i > 0; i--, n1++, n2++ ) {            i > 0; i--, n1++, n2++ ) {
                                         ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;            ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;
                                         *n1 = (unsigned int)(ltmp & BMASK);            *n1 = (unsigned int)(ltmp & BMASK);
                                         ur = (unsigned int)(ltmp >> BSH);            ur = (unsigned int)(ltmp >> BSH);
                                 }          }
                         }        }
                         *n1 = 0;        *n1 = 0;
                 }      }
                 *q = qhat;      *q = qhat;
         }    }
 }  }
   
 unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)  unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
 {  {
 /*  /*
         unsigned int borrow,l;    unsigned int borrow,l;
   
         for ( borrow = 0; s--; p++, r++ ) {    for ( borrow = 0; s--; p++, r++ ) {
                 DMA(*p,d,borrow,borrow,l)      DMA(*p,d,borrow,borrow,l)
                 if ( *r >= l )      if ( *r >= l )
                         *r -= l;        *r -= l;
                 else {      else {
                         *r = (unsigned int)((UL)*r+(ULBASE-(UL)l)); borrow++;        *r = (unsigned int)((UL)*r+(ULBASE-(UL)l)); borrow++;
                 }      }
         }    }
         return borrow;    return borrow;
 */  */
         /* esi : p, edi : r, borrow : ebx, s : ecx */    /* esi : p, edi : r, borrow : ebx, s : ecx */
 #if defined(VISUAL)  #if defined(_M_IX86)
         __asm {    __asm {
         push esi    push esi
         push edi    push edi
         mov esi,p    mov esi,p
         mov edi,r    mov edi,r
         mov ecx,s    mov ecx,s
         xor ebx,ebx    xor ebx,ebx
         Lstart_divn:    Lstart_divn:
         mov eax,DWORD PTR [esi]    mov eax,DWORD PTR [esi]
         mul d    mul d
         add eax,ebx    add eax,ebx
         adc edx,0    adc edx,0
         sub DWORD PTR [edi],eax    sub DWORD PTR [edi],eax
         adc edx,0    adc edx,0
         mov ebx,edx    mov ebx,edx
         lea     esi,DWORD PTR [esi+4]    lea  esi,DWORD PTR [esi+4]
         lea     edi,DWORD PTR [edi+4]    lea  edi,DWORD PTR [edi+4]
         dec ecx    dec ecx
         jnz Lstart_divn    jnz Lstart_divn
         mov eax,ebx    mov eax,ebx
         pop edi    pop edi
         pop esi    pop esi
         }    }
         /* return value is in eax. */    /* return value is in eax. */
 #else  #else
         unsigned int borrow;    unsigned int borrow;
   
         asm volatile("\    asm volatile("\
         movl    %1,%%esi;\    pushl  %%ebx;\
         movl    %2,%%edi;\    movl  %1,%%esi;\
         movl    $0,%%ebx;\    movl  %2,%%edi;\
         Lstart_divn:;\    movl  $0,%%ebx;\
         movl    (%%esi),%%eax;\    Lstart_divn:;\
         mull    %3;\    movl  (%%esi),%%eax;\
         addl    %%ebx,%%eax;\    mull  %3;\
         adcl    $0,%%edx;\    addl  %%ebx,%%eax;\
         subl    %%eax,(%%edi);\    adcl  $0,%%edx;\
         adcl    $0,%%edx;\    subl  %%eax,(%%edi);\
         movl    %%edx,%%ebx;\    adcl  $0,%%edx;\
         leal    4(%%esi),%%esi;\    movl  %%edx,%%ebx;\
         leal    4(%%edi),%%edi;\    leal  4(%%esi),%%esi;\
         decl    %4;\    leal  4(%%edi),%%edi;\
         jnz             Lstart_divn;\    decl  %4;\
         movl    %%ebx,%0"\    jnz    Lstart_divn;\
         :"=m"(borrow)\    movl  %%ebx,%0;\
         :"m"(p),"m"(r),"m"(d),"m"(s)\    popl  %%ebx"\
         :"eax","ebx","edx","esi","edi");    :"=m"(borrow)\
     :"m"(p),"m"(r),"m"(d),"m"(s)\
     :"eax","edx","esi","edi");
   
         return borrow;    return borrow;
 #endif  #endif
 }  }
   
Line 551  unsigned int divn_1(unsigned int *p,int s,unsigned int
Line 560  unsigned int divn_1(unsigned int *p,int s,unsigned int
   
 void muln(N n1,N n2,N *nr)  void muln(N n1,N n2,N *nr)
 {  {
         unsigned int tmp,carry,mul;    unsigned int tmp,carry,mul;
         unsigned int *p1,*pp,*m1,*m2;    unsigned int *p1,*pp,*m1,*m2;
         int i,j,d1,d2,d;    int i,j,d1,d2;
         N r;    N r;
   
         if ( !n1 || !n2 )    if ( !n1 || !n2 )
                 *nr = 0;      *nr = 0;
         else if ( UNIN(n1) )    else if ( UNIN(n1) )
                 COPY(n2,*nr);      COPY(n2,*nr);
         else if ( UNIN(n2) )    else if ( UNIN(n2) )
                 COPY(n1,*nr);      COPY(n1,*nr);
         else if ( (PL(n1) == 1) && (PL(n2) == 1) ) {    else if ( (PL(n1) == 1) && (PL(n2) == 1) ) {
                 DM(*BD(n1),*BD(n2),carry,tmp)      DM(*BD(n1),*BD(n2),carry,tmp)
                 if ( carry ) {      if ( carry ) {
                         *nr = r = NALLOC(2); INITRC(r);        *nr = r = NALLOC(2); INITRC(r);
                         PL(r) = 2; p1 = BD(r); *p1++ = tmp; *p1 = carry;        PL(r) = 2; p1 = BD(r); *p1++ = tmp; *p1 = carry;
                 } else      } else
                         STON(tmp,*nr);        STON(tmp,*nr);
         } else {    } else {
                 d1 = PL(n1); d2 = PL(n2); d = d1+d2;      d1 = PL(n1); d2 = PL(n2);
                 *nr = r = NALLOC(d); INITRC(r);      *nr = r = NALLOC(d1+d2); INITRC(r);
                 for ( i = 0, pp = BD(r); i < d; i++, pp++ ) *pp = 0;      bzero((char *)BD(r),(int)((d1+d2)*sizeof(int)));
                 for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )      for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
                         if ( mul = *m2 ) {        if ( mul = *m2 ) {
                                 for ( j = d1, carry = 0, p1 = m1, pp = BD(r)+i;          for ( j = d1, carry = 0, p1 = m1, pp = BD(r)+i;
                                         j--; pp++ ) {            j--; pp++ ) {
                                         DMA2(*p1++,mul,*pp,carry,carry,*pp)            DMA2(*p1++,mul,*pp,carry,carry,*pp)
                                 }          }
                                 *pp = carry;          *pp = carry;
                         }        }
                 PL(r) = (carry?d1+d2:d1+d2-1);      PL(r) = (carry?d1+d2:d1+d2-1);
         }    }
 }  }
   
 void _muln(N n1,N n2,N nr)  void _muln(N n1,N n2,N nr)
 {  {
         unsigned int tmp,carry,mul;    unsigned int carry=0,mul;
         unsigned int *p1,*pp,*m1,*m2;    unsigned int *p1,*pp,*m1,*m2;
         int i,j,d1,d2;    int i,j,d1,d2;
   
         if ( !n1 || !PL(n1) || !n2 || !PL(n2) )    if ( !n1 || !PL(n1) || !n2 || !PL(n2) )
                 PL(nr) = 0;      PL(nr) = 0;
         else if ( UNIN(n1) )    else if ( UNIN(n1) )
                 dupn(n2,nr);      dupn(n2,nr);
         else if ( UNIN(n2) )    else if ( UNIN(n2) )
                 dupn(n1,nr);      dupn(n1,nr);
         else {    else {
                 d1 = PL(n1); d2 = PL(n2);      d1 = PL(n1); d2 = PL(n2);
                 bzero((char *)BD(nr),(int)((d1+d2)*sizeof(int)));      bzero((char *)BD(nr),(int)((d1+d2)*sizeof(int)));
                 for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )      for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
                         if ( mul = *m2 ) {        if ( mul = *m2 ) {
                                 for ( j = d1, carry = 0, p1 = m1, pp = BD(nr)+i;          for ( j = d1, carry = 0, p1 = m1, pp = BD(nr)+i;
                                         j--; pp++ ) {            j--; pp++ ) {
                                         DMA2(*p1++,mul,*pp,carry,carry,*pp)            DMA2(*p1++,mul,*pp,carry,carry,*pp)
                                 }          }
                                 *pp = carry;          *pp = carry;
                         }        }
                 PL(nr) = (carry?d1+d2:d1+d2-1);      PL(nr) = (carry?d1+d2:d1+d2-1);
         }    }
 }  }
   
 /* r[0...s] = p[0...s-1]*d */  /* r[0...s] = p[0...s-1]*d */
   
 void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)  void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
 {  {
         unsigned int carry;    unsigned int carry;
   
         for ( carry = 0; s--; p++, r++ ) {    for ( carry = 0; s--; p++, r++ ) {
                 DMA2(*p,d,carry,*r,carry,*r)      DMA2(*p,d,carry,*r,carry,*r)
         }    }
         *r = carry;    *r = carry;
 }  }
   
 void divnmain(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)  void divnmain(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)
 {  {
         int i,j;    int i,j;
         UL r,ltmp;    UL r,ltmp;
         unsigned int l,ur,tmp;    unsigned int l,ur;
         unsigned int *n1,*n2;    unsigned int *n1,*n2;
         unsigned int u,qhat;    unsigned int u,qhat;
         unsigned int v1,v2;    unsigned int v1,v2;
   
         v1 = m2[d2-1]; v2 = m2[d2-2];    v1 = m2[d2-1]; v2 = m2[d2-2];
         for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {    for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {
                 n1 = m1+d2; n2 = m1+d2-1;      n1 = m1+d2; n2 = m1+d2-1;
                 if ( *n1 == v1 ) {         if ( *n1 == v1 ) {
                         qhat = ULBASE - 1;              qhat = ULBASE - 1;
                         r = (UL)*n1 + (UL)*n2;        r = (UL)*n1 + (UL)*n2;
                 } else {         } else {
                         DSAB(v1,*n1,*n2,qhat,ur)        DSAB(v1,*n1,*n2,qhat,ur)
                         r = (UL)ur;        r = (UL)ur;
                 }      }
                 DM(v2,qhat,u,l)      DM(v2,qhat,u,l)
                 while ( 1 ) {         while ( 1 ) {
                         if ((r > (UL)u) || ((r == (UL)u) && (*(n1-2) >= l)))             if ((r > (UL)u) || ((r == (UL)u) && (*(n1-2) >= l)))
                                 break;          break;
                         if ( l >= v2 )        if ( l >= v2 )
                                 l -= v2;          l -= v2;
                         else {        else {
                                 l = (unsigned int)((UL)l+(ULBASE-(UL)v2)); u--;          l = (unsigned int)((UL)l+(ULBASE-(UL)v2)); u--;
                         }        }
                 r += v1; qhat--;            r += v1; qhat--;
                 }         }
                 if ( qhat ) {      if ( qhat ) {
                         u = divn_1(m2,d2,qhat,m1);        u = divn_1(m2,d2,qhat,m1);
                         if ( *(m1+d2) < u ) {        if ( *(m1+d2) < u ) {
                                 for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;          for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;
                                         i > 0; i--, n1++, n2++ ) {            i > 0; i--, n1++, n2++ ) {
                                         ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;            ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;
                                         *n1 = (unsigned int)(ltmp & BMASK);            *n1 = (unsigned int)(ltmp & BMASK);
                                         ur = (unsigned int)(ltmp >> BSH);            ur = (unsigned int)(ltmp >> BSH);
                                 }          }
                         }        }
                         *n1 = 0;        *n1 = 0;
                 }      }
                 *q = qhat;      *q = qhat;
         }    }
 }  }
   
 unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)  unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
 {  {
         unsigned int borrow,l;    unsigned int borrow,l;
   
         for ( borrow = 0; s--; p++, r++ ) {    for ( borrow = 0; s--; p++, r++ ) {
                 DMA(*p,d,borrow,borrow,l)      DMA(*p,d,borrow,borrow,l)
                 if ( *r >= l )      if ( *r >= l )
                         *r -= l;        *r -= l;
                 else {      else {
                         *r = (unsigned int)((UL)*r+(ULBASE-(UL)l)); borrow++;        *r = (unsigned int)((UL)*r+(ULBASE-(UL)l)); borrow++;
                 }      }
         }    }
         return borrow;    return borrow;
 }  }
 #endif  #endif
Legend:
Removed from v.1.5  
changed lines
  Added in v.1.14
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>