Annotation of OpenXM_contrib2/asir2000/asm/ddN.c, Revision 1.14
1.2 noro 1: /*
2: * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
3: * All rights reserved.
4: *
5: * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
6: * non-exclusive and royalty-free license to use, copy, modify and
7: * redistribute, solely for non-commercial and non-profit purposes, the
8: * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
9: * conditions of this Agreement. For the avoidance of doubt, you acquire
10: * only a limited right to use the SOFTWARE hereunder, and FLL or any
11: * third party developer retains all rights, including but not limited to
12: * copyrights, in and to the SOFTWARE.
13: *
14: * (1) FLL does not grant you a license in any way for commercial
15: * purposes. You may use the SOFTWARE only for non-commercial and
16: * non-profit purposes only, such as academic, research and internal
17: * business use.
18: * (2) The SOFTWARE is protected by the Copyright Law of Japan and
19: * international copyright treaties. If you make copies of the SOFTWARE,
20: * with or without modification, as permitted hereunder, you shall affix
21: * to all such copies of the SOFTWARE the above copyright notice.
22: * (3) An explicit reference to this SOFTWARE and its copyright owner
23: * shall be made on your publication or presentation in any form of the
24: * results obtained by use of the SOFTWARE.
25: * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
1.3 noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.2 noro 27: * for such modification or the source code of the modified part of the
28: * SOFTWARE.
29: *
30: * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
31: * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
32: * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
33: * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
34: * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
35: * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
36: * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
37: * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
38: * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
39: * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
40: * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
41: * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
42: * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
43: * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
44: * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
45: * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
46: * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
47: *
1.14 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/asm/ddN.c,v 1.13 2015/08/29 04:15:04 fujimoto Exp $
1.2 noro 49: */
1.1 noro 50: #ifndef FBASE
51: #define FBASE
52: #endif
53:
54: #include "ca.h"
55: #include "base.h"
56: #include "inline.h"
57:
1.12 ohara 58: #if defined(__GNUC__)
59: unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r) __attribute__ ((noinline));
60: void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r) __attribute__ ((noinline));
61: #endif
62:
1.4 noro 63: void divn(N n1,N n2,N *nq,N *nr)
1.1 noro 64: {
1.14 ! noro 65: int tmp,b;
! 66: int i,j,d1,d2,dd;
! 67: unsigned int *m1,*m2;
! 68: unsigned int uq,ur,msw;
! 69: N q,r,t1,t2;
! 70:
! 71: if ( !n2 )
! 72: error("divn: divsion by 0");
! 73: else if ( !n1 ) {
! 74: *nq = 0; *nr = 0;
! 75: } else if ( UNIN(n2) ) {
! 76: COPY(n1,*nq); *nr = 0;
! 77: } else if ( (d1 = PL(n1)) < (d2 = PL(n2)) ) {
! 78: COPY(n1,*nr); *nq = 0;
! 79: } else if ( d2 == 1 ) {
! 80: if ( d1 == 1 ) {
! 81: DQR(BD(n1)[0],BD(n2)[0],uq,ur)
! 82: STON(uq,*nq); STON(ur,*nr);
! 83: } else {
! 84: ur = divin(n1,BD(n2)[0],nq); STON(ur,*nr);
! 85: }
! 86: return;
! 87: } else if ( (d1 == d2) && ((tmp = cmpn(n1,n2)) <= 0) )
! 88: if ( !tmp ) {
! 89: *nr = 0; COPY(ONEN,*nq);
! 90: } else {
! 91: COPY(n1,*nr); *nq = 0;
! 92: }
! 93: else {
! 94: msw = BD(n2)[d2-1];
! 95: for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
! 96: b = BSH-i;
! 97: if ( b ) {
! 98: bshiftn(n1,-b,&t1); bshiftn(n2,-b,&t2);
! 99: } else {
! 100: COPY(n1,t1); COPY(n2,t2);
! 101: }
! 102: m1 = (unsigned int *)ALLOCA((PL(t1)+1)*sizeof(unsigned int));
! 103: bcopy(BD(t1),m1,PL(t1)*sizeof(unsigned int)); m1[PL(t1)] = 0;
! 104: m2 = BD(t2); dd = d1-d2;
! 105: *nq = q = NALLOC(dd+1); INITRC(q);
! 106: divnmain(d1,d2,m1,m2,BD(q)); FREEN(t1); FREEN(t2);
! 107: PL(q) = (BD(q)[dd]?dd+1:dd);
! 108: if ( !PL(q) ) {
! 109: FREEN(q);
! 110: }
! 111: for ( j = d2-1; (j >= 0) && (m1[j] == 0); j--);
! 112: if ( j == -1 )
! 113: *nr = 0;
! 114: else {
! 115: r = NALLOC(j+1); INITRC(r); PL(r) = j+1;
! 116: bcopy(m1,BD(r),(j+1)*sizeof(unsigned int));
! 117: if ( b ) {
! 118: bshiftn(r,b,nr); FREEN(r);
! 119: } else
! 120: *nr = r;
! 121: }
! 122: }
1.1 noro 123: }
124:
1.4 noro 125: void divsn(N n1,N n2,N *nq)
1.1 noro 126: {
1.14 ! noro 127: int d1,d2,dd,i,b;
! 128: unsigned int *m1,*m2;
! 129: unsigned int uq,msw;
! 130: N q,t1,t2;
! 131:
! 132: if ( !n1 )
! 133: *nq = 0;
! 134: else if ( UNIN(n2) ) {
! 135: COPY(n1,*nq);
! 136: } else if ( (d1 = PL(n1)) < (d2 = PL(n2)) )
! 137: error("divsn: cannot happen");
! 138: else if ( d2 == 1 ) {
! 139: if ( d1 == 1 ) {
! 140: uq = BD(n1)[0] / BD(n2)[0]; STON(uq,*nq);
! 141: } else
! 142: divin(n1,BD(n2)[0],nq);
! 143: } else {
! 144: msw = BD(n2)[d2-1];
! 145: for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
! 146: b = BSH-i;
! 147: if ( b ) {
! 148: bshiftn(n1,-b,&t1); bshiftn(n2,-b,&t2);
! 149: } else {
! 150: COPY(n1,t1); COPY(n2,t2);
! 151: }
! 152: m1 = (unsigned int *)ALLOCA((PL(t1)+1)*sizeof(unsigned int));
! 153: bcopy(BD(t1),m1,PL(t1)*sizeof(unsigned int)); m1[PL(t1)] = 0;
! 154: m2 = BD(t2); dd = d1-d2;
! 155: *nq = q = NALLOC(dd+1); INITRC(q);
! 156: divnmain(d1,d2,m1,m2,BD(q));
! 157: FREEN(t1); FREEN(t2);
! 158: PL(q) = (BD(q)[dd]?dd+1:dd);
! 159: if ( !PL(q) )
! 160: FREEN(q);
! 161: }
1.1 noro 162: }
163:
1.4 noro 164: void remn(N n1,N n2,N *nr)
1.1 noro 165: {
1.14 ! noro 166: int d1,d2,tmp;
! 167: unsigned int uq,ur;
! 168: N q;
! 169:
! 170: if ( !n2 )
! 171: error("remn: divsion by 0");
! 172: else if ( !n1 || UNIN(n2) )
! 173: *nr = 0;
! 174: else if ( (d1 = PL(n1)) < (d2 = PL(n2)) ) {
! 175: COPY(n1,*nr);
! 176: } else if ( d2 == 1 ) {
! 177: if ( d1 == 1 ) {
! 178: DQR(BD(n1)[0],BD(n2)[0],uq,ur)
! 179: STON(ur,*nr);
! 180: } else {
! 181: ur = rem(n1,BD(n2)[0]); UTON(ur,*nr);
! 182: }
! 183: return;
! 184: } else if ( (d1 == d2) && ((tmp = cmpn(n1,n2)) <= 0) )
! 185: if ( !tmp ) {
! 186: *nr = 0;
! 187: } else {
! 188: COPY(n1,*nr);
! 189: }
! 190: else
! 191: divn(n1,n2,&q,nr);
1.1 noro 192: }
193:
194: /* d = 2^(32*words)-lower */
195:
1.4 noro 196: void remn_special(N a,N d,int bits,unsigned int lower,N *b)
1.1 noro 197: {
1.14 ! noro 198: int words;
! 199: N r;
! 200: int rl,i;
! 201: unsigned int c,tmp;
! 202: unsigned int *rb,*src,*dst;
! 203:
! 204: if ( cmpn(a,d) < 0 ) {
! 205: *b = a;
! 206: return;
! 207: }
! 208: words = bits>>5;
! 209: r = NALLOC(PL(a)); dupn(a,r);
! 210: rl = PL(r); rb = BD(r);
! 211: while ( rl > words ) {
! 212: src = rb+words;
! 213: dst = rb;
! 214: i = rl-words;
! 215: for ( c = 0; i > 0; i--, src++, dst++ ) {
! 216: DMA2(*src,lower,c,*dst,c,*dst)
! 217: *src = 0;
! 218: }
! 219: for ( ; c; dst++ ) {
! 220: tmp = *dst + c;
! 221: c = tmp < c ? 1 : 0;
! 222: *dst = tmp;
! 223: }
! 224: rl = dst-rb;
! 225: }
! 226: for ( i = words-1; i >= 0 && !rb[i]; i-- );
! 227: PL(r) = i + 1;
! 228: if ( cmpn(r,d) >= 0 ) {
! 229: tmp = rb[0]-BD(d)[0];
! 230: UTON(tmp,*b);
! 231: } else
! 232: *b = r;
1.1 noro 233: }
1.4 noro 234: void mulin(N n,unsigned int d,unsigned int *p)
1.1 noro 235: {
1.14 ! noro 236: unsigned int carry;
! 237: unsigned int *m,*me;
1.1 noro 238:
1.14 ! noro 239: bzero((char *)p,(int)((PL(n)+1)*sizeof(int)));
! 240: for ( carry = 0, m = BD(n), me = m+PL(n); m < me; p++, m++ ) {
! 241: DMA2(*m,d,*p,carry,carry,*p)
! 242: }
! 243: *p = carry;
1.1 noro 244: }
245:
1.4 noro 246: unsigned int divin(N n,unsigned int dvr,N *q)
1.1 noro 247: {
1.14 ! noro 248: int d;
! 249: unsigned int up;
! 250: unsigned int *m,*mq,*md;
! 251: N nq;
! 252:
! 253: d = PL(n); m = BD(n);
! 254: if ( ( d == 1 ) && ( dvr > *m ) ) {
! 255: *q = 0;
! 256: return *m;
! 257: }
! 258: *q = nq = NALLOC(d); INITRC(nq);
! 259: for ( md = m+d-1, mq = BD(nq)+d-1, up = 0; md >= m; md--, mq-- ) {
! 260: DSAB(dvr,up,*md,*mq,up)
! 261: }
! 262: PL(nq) = (BD(nq)[d-1]?d:d-1);
! 263: if ( !PL(nq) )
! 264: FREEN(nq);
! 265: return ( up );
1.1 noro 266: }
267:
1.4 noro 268: void bprintn(N n)
1.1 noro 269: {
1.14 ! noro 270: int l,i;
! 271: unsigned int *b;
1.1 noro 272:
1.14 ! noro 273: if ( !n )
! 274: printf("0");
! 275: else {
! 276: l = PL(n); b = BD(n);
! 277: for ( i = l-1; i >= 0; i-- )
! 278: printf("%010u|",b[i]);
! 279: }
1.1 noro 280: }
281:
1.4 noro 282: void bxprintn(N n)
1.1 noro 283: {
1.14 ! noro 284: int l,i;
! 285: unsigned int *b;
1.1 noro 286:
1.14 ! noro 287: if ( !n )
! 288: printf("0");
! 289: else {
! 290: l = PL(n); b = BD(n);
! 291: for ( i = l-1; i >= 0; i-- )
! 292: printf("%08x|",b[i]);
! 293: }
1.1 noro 294: }
295:
1.13 fujimoto 296: #if (defined(_M_IX86) || defined(i386)) && !defined(__MINGW32__)
1.4 noro 297: void muln(N n1,N n2,N *nr)
1.1 noro 298: {
1.14 ! noro 299: unsigned int tmp,carry,mul;
! 300: unsigned int *p1,*m1,*m2;
! 301: int i,d1,d2,d;
! 302: N r;
! 303:
! 304: if ( !n1 || !n2 )
! 305: *nr = 0;
! 306: else if ( UNIN(n1) )
! 307: COPY(n2,*nr);
! 308: else if ( UNIN(n2) )
! 309: COPY(n1,*nr);
! 310: else if ( (PL(n1) == 1) && (PL(n2) == 1) ) {
! 311: DM(*BD(n1),*BD(n2),carry,tmp)
! 312: if ( carry ) {
! 313: *nr = r = NALLOC(2); INITRC(r);
! 314: PL(r) = 2; p1 = BD(r); *p1++ = tmp; *p1 = carry;
! 315: } else
! 316: STON(tmp,*nr);
! 317: } else {
! 318: d1 = PL(n1); d2 = PL(n2);
! 319: d = d1+d2;
! 320: *nr = r = NALLOC(d); INITRC(r);
! 321: bzero((char *)BD(r),(int)((d1+d2)*sizeof(int)));
! 322: for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
! 323: if ( mul = *m2 )
! 324: muln_1(m1,d1,mul,BD(r)+i);
! 325: PL(r) = (BD(r)[d-1]?d:d-1);
! 326: }
1.1 noro 327: }
328:
1.4 noro 329: void _muln(N n1,N n2,N nr)
1.1 noro 330: {
1.14 ! noro 331: unsigned int mul;
! 332: unsigned int *m1,*m2;
! 333: int i,d1,d2,d;
! 334:
! 335: if ( !n1 || !PL(n1) || !n2 || !PL(n2) )
! 336: PL(nr) = 0;
! 337: else if ( UNIN(n1) )
! 338: dupn(n2,nr);
! 339: else if ( UNIN(n2) )
! 340: dupn(n1,nr);
! 341: else {
! 342: d1 = PL(n1); d2 = PL(n2);
! 343: d = d1+d2;
! 344: bzero((char *)BD(nr),(int)((d1+d2)*sizeof(int)));
! 345: for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
! 346: if ( mul = *m2 )
! 347: muln_1(m1,d1,mul,BD(nr)+i);
! 348: PL(nr) = (BD(nr)[d-1]?d:d-1);
! 349: }
1.1 noro 350: }
351:
1.4 noro 352: void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
1.1 noro 353: {
1.14 ! noro 354: /* esi : p, edi : r, carry : ebx, s : ecx */
1.9 ohara 355: #if defined(_M_IX86)
1.14 ! noro 356: __asm {
! 357: push esi
! 358: push edi
! 359: mov esi,p
! 360: mov edi,r
! 361: mov ecx,s
! 362: xor ebx,ebx
! 363: Lstart_muln:
! 364: mov eax,DWORD PTR [esi]
! 365: mul d
! 366: add eax,DWORD PTR [edi]
! 367: adc edx,0
! 368: add eax,ebx
! 369: adc edx,0
! 370: mov DWORD PTR[edi],eax
! 371: mov ebx,edx
! 372: lea esi,DWORD PTR [esi+4]
! 373: lea edi,DWORD PTR [edi+4]
! 374: dec ecx
! 375: jnz Lstart_muln
! 376: mov DWORD PTR [edi],ebx
! 377: pop edi
! 378: pop esi
! 379: }
1.1 noro 380: #else
1.14 ! noro 381: asm volatile("\
! 382: pushl %%ebx;\
! 383: movl %0,%%esi;\
! 384: movl %1,%%edi;\
! 385: movl $0,%%ebx;\
! 386: Lstart_muln:;\
! 387: movl (%%esi),%%eax;\
! 388: mull %2;\
! 389: addl (%%edi),%%eax;\
! 390: adcl $0,%%edx;\
! 391: addl %%ebx,%%eax;\
! 392: adcl $0,%%edx;\
! 393: movl %%eax,(%%edi);\
! 394: movl %%edx,%%ebx;\
! 395: leal 4(%%esi),%%esi;\
! 396: leal 4(%%edi),%%edi;\
! 397: decl %3;\
! 398: jnz Lstart_muln;\
! 399: movl %%ebx,(%%edi);\
! 400: popl %%ebx"\
! 401: :\
! 402: :"m"(p),"m"(r),"m"(d),"m"(s)\
! 403: :"eax","edx","esi","edi");
1.1 noro 404: #endif
405: }
406:
1.4 noro 407: void divnmain(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)
1.1 noro 408: {
1.14 ! noro 409: int i,j;
! 410: UL r,ltmp;
! 411: unsigned int l,ur;
! 412: unsigned int *n1,*n2;
! 413: unsigned int u,qhat;
! 414: unsigned int v1,v2;
1.1 noro 415:
1.14 ! noro 416: v1 = m2[d2-1]; v2 = m2[d2-2];
1.1 noro 417: #if 0
1.14 ! noro 418: if ( v1 == (1<<31) && !v2 ) {
! 419: divnmain_special(d1,d2,m1,m2,q);
! 420: return;
! 421: }
1.1 noro 422: #endif
1.14 ! noro 423: for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {
! 424: n1 = m1+d2; n2 = m1+d2-1;
! 425: if ( *n1 == v1 ) {
! 426: qhat = ULBASE - 1;
! 427: r = (UL)*n1 + (UL)*n2;
! 428: } else {
! 429: DSAB(v1,*n1,*n2,qhat,ur)
! 430: r = (UL)ur;
! 431: }
! 432: DM(v2,qhat,u,l)
! 433: while ( 1 ) {
! 434: if ((r > (UL)u) || ((r == (UL)u) && (*(n1-2) >= l)))
! 435: break;
! 436: if ( l >= v2 )
! 437: l -= v2;
! 438: else {
! 439: l = (unsigned int)((UL)l+(ULBASE-(UL)v2)); u--;
! 440: }
! 441: r += v1; qhat--;
! 442: }
! 443: if ( qhat ) {
! 444: u = divn_1(m2,d2,qhat,m1);
! 445: if ( *(m1+d2) < u ) {
! 446: for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;
! 447: i > 0; i--, n1++, n2++ ) {
! 448: ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;
! 449: *n1 = (unsigned int)(ltmp & BMASK);
! 450: ur = (unsigned int)(ltmp >> BSH);
! 451: }
! 452: }
! 453: *n1 = 0;
! 454: }
! 455: *q = qhat;
! 456: }
1.1 noro 457: }
458:
1.4 noro 459: void divnmain_special(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)
1.1 noro 460: {
1.14 ! noro 461: int i,j;
! 462: UL ltmp;
! 463: unsigned int ur;
! 464: unsigned int *n1,*n2;
! 465: unsigned int u,qhat;
! 466: unsigned int v1,v2;
! 467:
! 468: v1 = m2[d2-1]; v2 = 0;
! 469: for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {
! 470: n1 = m1+d2; n2 = m1+d2-1;
! 471: qhat = ((*n1)<<1)|((*n2)>>31);
! 472: if ( qhat ) {
! 473: u = divn_1(m2,d2,qhat,m1);
! 474: if ( *(m1+d2) < u ) {
! 475: for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;
! 476: i > 0; i--, n1++, n2++ ) {
! 477: ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;
! 478: *n1 = (unsigned int)(ltmp & BMASK);
! 479: ur = (unsigned int)(ltmp >> BSH);
! 480: }
! 481: }
! 482: *n1 = 0;
! 483: }
! 484: *q = qhat;
! 485: }
1.1 noro 486: }
487:
1.4 noro 488: unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
1.1 noro 489: {
490: /*
1.14 ! noro 491: unsigned int borrow,l;
1.1 noro 492:
1.14 ! noro 493: for ( borrow = 0; s--; p++, r++ ) {
! 494: DMA(*p,d,borrow,borrow,l)
! 495: if ( *r >= l )
! 496: *r -= l;
! 497: else {
! 498: *r = (unsigned int)((UL)*r+(ULBASE-(UL)l)); borrow++;
! 499: }
! 500: }
! 501: return borrow;
1.1 noro 502: */
1.14 ! noro 503: /* esi : p, edi : r, borrow : ebx, s : ecx */
1.9 ohara 504: #if defined(_M_IX86)
1.14 ! noro 505: __asm {
! 506: push esi
! 507: push edi
! 508: mov esi,p
! 509: mov edi,r
! 510: mov ecx,s
! 511: xor ebx,ebx
! 512: Lstart_divn:
! 513: mov eax,DWORD PTR [esi]
! 514: mul d
! 515: add eax,ebx
! 516: adc edx,0
! 517: sub DWORD PTR [edi],eax
! 518: adc edx,0
! 519: mov ebx,edx
! 520: lea esi,DWORD PTR [esi+4]
! 521: lea edi,DWORD PTR [edi+4]
! 522: dec ecx
! 523: jnz Lstart_divn
! 524: mov eax,ebx
! 525: pop edi
! 526: pop esi
! 527: }
! 528: /* return value is in eax. */
1.1 noro 529: #else
1.14 ! noro 530: unsigned int borrow;
1.1 noro 531:
1.14 ! noro 532: asm volatile("\
! 533: pushl %%ebx;\
! 534: movl %1,%%esi;\
! 535: movl %2,%%edi;\
! 536: movl $0,%%ebx;\
! 537: Lstart_divn:;\
! 538: movl (%%esi),%%eax;\
! 539: mull %3;\
! 540: addl %%ebx,%%eax;\
! 541: adcl $0,%%edx;\
! 542: subl %%eax,(%%edi);\
! 543: adcl $0,%%edx;\
! 544: movl %%edx,%%ebx;\
! 545: leal 4(%%esi),%%esi;\
! 546: leal 4(%%edi),%%edi;\
! 547: decl %4;\
! 548: jnz Lstart_divn;\
! 549: movl %%ebx,%0;\
! 550: popl %%ebx"\
! 551: :"=m"(borrow)\
! 552: :"m"(p),"m"(r),"m"(d),"m"(s)\
! 553: :"eax","edx","esi","edi");
1.1 noro 554:
1.14 ! noro 555: return borrow;
1.1 noro 556: #endif
557: }
558:
559: #else
560:
1.4 noro 561: void muln(N n1,N n2,N *nr)
1.1 noro 562: {
1.14 ! noro 563: unsigned int tmp,carry,mul;
! 564: unsigned int *p1,*pp,*m1,*m2;
! 565: int i,j,d1,d2;
! 566: N r;
! 567:
! 568: if ( !n1 || !n2 )
! 569: *nr = 0;
! 570: else if ( UNIN(n1) )
! 571: COPY(n2,*nr);
! 572: else if ( UNIN(n2) )
! 573: COPY(n1,*nr);
! 574: else if ( (PL(n1) == 1) && (PL(n2) == 1) ) {
! 575: DM(*BD(n1),*BD(n2),carry,tmp)
! 576: if ( carry ) {
! 577: *nr = r = NALLOC(2); INITRC(r);
! 578: PL(r) = 2; p1 = BD(r); *p1++ = tmp; *p1 = carry;
! 579: } else
! 580: STON(tmp,*nr);
! 581: } else {
! 582: d1 = PL(n1); d2 = PL(n2);
! 583: *nr = r = NALLOC(d1+d2); INITRC(r);
! 584: bzero((char *)BD(r),(int)((d1+d2)*sizeof(int)));
! 585: for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
! 586: if ( mul = *m2 ) {
! 587: for ( j = d1, carry = 0, p1 = m1, pp = BD(r)+i;
! 588: j--; pp++ ) {
! 589: DMA2(*p1++,mul,*pp,carry,carry,*pp)
! 590: }
! 591: *pp = carry;
! 592: }
! 593: PL(r) = (carry?d1+d2:d1+d2-1);
! 594: }
1.1 noro 595: }
596:
1.4 noro 597: void _muln(N n1,N n2,N nr)
1.1 noro 598: {
1.14 ! noro 599: unsigned int carry=0,mul;
! 600: unsigned int *p1,*pp,*m1,*m2;
! 601: int i,j,d1,d2;
! 602:
! 603: if ( !n1 || !PL(n1) || !n2 || !PL(n2) )
! 604: PL(nr) = 0;
! 605: else if ( UNIN(n1) )
! 606: dupn(n2,nr);
! 607: else if ( UNIN(n2) )
! 608: dupn(n1,nr);
! 609: else {
! 610: d1 = PL(n1); d2 = PL(n2);
! 611: bzero((char *)BD(nr),(int)((d1+d2)*sizeof(int)));
! 612: for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
! 613: if ( mul = *m2 ) {
! 614: for ( j = d1, carry = 0, p1 = m1, pp = BD(nr)+i;
! 615: j--; pp++ ) {
! 616: DMA2(*p1++,mul,*pp,carry,carry,*pp)
! 617: }
! 618: *pp = carry;
! 619: }
! 620: PL(nr) = (carry?d1+d2:d1+d2-1);
! 621: }
1.1 noro 622: }
623:
624: /* r[0...s] = p[0...s-1]*d */
625:
1.4 noro 626: void muln_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
1.1 noro 627: {
1.14 ! noro 628: unsigned int carry;
1.1 noro 629:
1.14 ! noro 630: for ( carry = 0; s--; p++, r++ ) {
! 631: DMA2(*p,d,carry,*r,carry,*r)
! 632: }
! 633: *r = carry;
1.1 noro 634: }
635:
1.4 noro 636: void divnmain(int d1,int d2,unsigned int *m1,unsigned int *m2,unsigned int *q)
1.1 noro 637: {
1.14 ! noro 638: int i,j;
! 639: UL r,ltmp;
! 640: unsigned int l,ur;
! 641: unsigned int *n1,*n2;
! 642: unsigned int u,qhat;
! 643: unsigned int v1,v2;
! 644:
! 645: v1 = m2[d2-1]; v2 = m2[d2-2];
! 646: for ( j = d1-d2, m1 += j, q += j; j >= 0; j--, q--, m1-- ) {
! 647: n1 = m1+d2; n2 = m1+d2-1;
! 648: if ( *n1 == v1 ) {
! 649: qhat = ULBASE - 1;
! 650: r = (UL)*n1 + (UL)*n2;
! 651: } else {
! 652: DSAB(v1,*n1,*n2,qhat,ur)
! 653: r = (UL)ur;
! 654: }
! 655: DM(v2,qhat,u,l)
! 656: while ( 1 ) {
! 657: if ((r > (UL)u) || ((r == (UL)u) && (*(n1-2) >= l)))
! 658: break;
! 659: if ( l >= v2 )
! 660: l -= v2;
! 661: else {
! 662: l = (unsigned int)((UL)l+(ULBASE-(UL)v2)); u--;
! 663: }
! 664: r += v1; qhat--;
! 665: }
! 666: if ( qhat ) {
! 667: u = divn_1(m2,d2,qhat,m1);
! 668: if ( *(m1+d2) < u ) {
! 669: for ( i = d2, qhat--, ur = 0, n1 = m1, n2 = m2;
! 670: i > 0; i--, n1++, n2++ ) {
! 671: ltmp = (UL)*n1 + (UL)*n2 + (UL)ur;
! 672: *n1 = (unsigned int)(ltmp & BMASK);
! 673: ur = (unsigned int)(ltmp >> BSH);
! 674: }
! 675: }
! 676: *n1 = 0;
! 677: }
! 678: *q = qhat;
! 679: }
1.1 noro 680: }
681:
1.4 noro 682: unsigned int divn_1(unsigned int *p,int s,unsigned int d,unsigned int *r)
1.1 noro 683: {
1.14 ! noro 684: unsigned int borrow,l;
1.1 noro 685:
1.14 ! noro 686: for ( borrow = 0; s--; p++, r++ ) {
! 687: DMA(*p,d,borrow,borrow,l)
! 688: if ( *r >= l )
! 689: *r -= l;
! 690: else {
! 691: *r = (unsigned int)((UL)*r+(ULBASE-(UL)l)); borrow++;
! 692: }
! 693: }
! 694: return borrow;
1.1 noro 695: }
696: #endif
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>