![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to Z.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / engine|
Return to Z.c CVS log | Up to [local] / OpenXM_contrib2 / asir2000 / engine |
1.2 noro 1: #if 0
1.1 noro 2: #include "ca.h"
3: #include "inline.h"
1.2 noro 4: #endif
1.1 noro 5:
6: typedef struct oZ {
7: int p;
8: unsigned int b[1];
9: } *Z;
10:
1.3 ! noro 11: #define IMM_MAX 1073741823
! 12: #define IMM_MIN -1073741823
! 13:
! 14: /* immediate int -> Z */
! 15: #define IMMTOZ(c,n) ((c)>=IMM_MIN&&(c)<=IMM_MAX?(n)=(Z)(((c)<<1)|1):(n)=immtoz(c))
! 16: /* immediate Z ? */
! 17: #define IS_IMM(c) (((unsigned int)c)&1)
! 18: /* Z can be conver to immediate ? */
! 19: #define IS_IMMZ(n) (SL(n) == 1&&BD(n)[0]<=IMM_MAX)
! 20: /* Z -> immediate Z */
! 21: #define IMMZTOZ(n,z) (SL(n)<0?(z)=(Z)(((-BD(n)[0])<<1)|1):(Z)(((BD(n)[0])<<1)|1))
! 22: /* Z -> immediate int */
! 23: #define ZTOIMM(c) (((int)(c))>>1)
1.1 noro 24: #define ZALLOC(d) ((Z)MALLOC_ATOMIC(TRUESIZE(oZ,(d)-1,int)))
25: #define SL(n) ((n)->p)
26:
1.3 ! noro 27: Z immtoz(int c);
1.1 noro 28: Z qtoz(Q n);
29: Q ztoq(Z n);
30: Z chsgnz(Z n);
31: Z dupz(Z n);
32: Z addz(Z n1,Z n2);
33: Z subz(Z n1,Z n2);
34: Z mulz(Z n1,Z n2);
35: Z divsz(Z n1,Z n2);
36: Z divz(Z n1,Z n2,Z *rem);
37: Z gcdz(Z n1,Z n2);
38: Z gcdz_cofactor(Z n1,Z n2,Z *c1,Z *c2);
39: Z estimate_array_gcdz(Z *a,int n);
40: Z array_gcdz(Z *a,int n);
1.2 noro 41: void mkwcz(int k,int l,Z *t);
1.1 noro 42: int remzi(Z n,int m);
43: inline void _addz(Z n1,Z n2,Z nr);
44: inline void _subz(Z n1,Z n2,Z nr);
45: inline void _mulz(Z n1,Z n2,Z nr);
46: inline int _addz_main(unsigned int *m1,int d1,unsigned int *m2,int d2,unsigned int *mr);
47: inline int _subz_main(unsigned int *m1,int d1,unsigned int *m2,int d2,unsigned int *mr);
48:
1.3 ! noro 49: Z immtoz(int c)
! 50: {
! 51: Z z;
! 52:
! 53: z = ZALLOC(1);
! 54: if ( c < 0 ) {
! 55: SL(z) = -1; BD(z)[0] = -c;
! 56: } else {
! 57: SL(z) = 1; BD(z)[0] = c;
! 58: }
! 59: return z;
! 60: }
! 61:
! 62: int sgnz(Z n)
1.2 noro 63: {
64: if ( !n ) return 0;
1.3 ! noro 65: else if ( IS_IMM(n) ) return ZTOIMM(n)>0?1:1;
1.2 noro 66: else if ( SL(n) < 0 ) return -1;
67: else return 1;
68: }
69:
1.1 noro 70: z_mag(Z n)
71: {
1.3 ! noro 72: int n;
! 73:
! 74: if ( !n ) return 0;
! 75: else if ( IS_IMM(n) ) {
! 76: c = ZTOIMM(n);
! 77: if ( c < 0 ) c = -c;
! 78: for ( i = 0; c; c >>= 1, i++ );
! 79: return i;
! 80: }
! 81: else return n_bits((N)n);
1.1 noro 82: }
83:
84: Z qtoz(Q n)
85: {
1.3 ! noro 86: Z r,t;
1.1 noro 87:
88: if ( !n ) return 0;
89: else if ( !INT(n) )
90: error("qtoz : invalid input");
91: else {
1.3 ! noro 92: t = (Z)NM(n);
! 93: if ( IS_IMMZ(t) ) {
! 94: IMMZTOZ(t,r);
! 95: } else {
! 96: r = dupz((Z)t);
! 97: }
1.1 noro 98: if ( SGN(n) < 0 ) SL(r) = -SL(r);
99: return r;
100: }
101: }
102:
103: Q ztoq(Z n)
104: {
105: Q r;
106: Z nm;
1.3 ! noro 107: int sgn,c;
1.1 noro 108:
109: if ( !n ) return 0;
1.3 ! noro 110: else if ( IS_IMM(n) ) {
! 111: c = ZTOIMM(n);
! 112: STOQ(c,r);
! 113: return r;
! 114: } else {
1.1 noro 115: nm = dupz(n);
116: if ( SL(nm) < 0 ) {
117: sgn = -1;
118: SL(nm) = -SL(nm);
119: } else
120: sgn = 1;
121: NTOQ((N)nm,sgn,r);
122: return r;
123: }
124: }
125:
126: Z dupz(Z n)
127: {
128: Z nr;
129: int sd,i;
130:
131: if ( !n ) return 0;
1.3 ! noro 132: else if ( IS_IMM(n) ) return n;
! 133: else {
! 134: if ( (sd = SL(n)) < 0 ) sd = -sd;
! 135: nr = ZALLOC(sd);
! 136: SL(nr) = SL(n);
! 137: for ( i = 0; i < sd; i++ ) BD(nr)[i] = BD(n)[i];
! 138: return nr;
! 139: }
1.1 noro 140: }
141:
142: Z chsgnz(Z n)
143: {
144: Z r;
1.3 ! noro 145: int c;
1.1 noro 146:
147: if ( !n ) return 0;
1.3 ! noro 148: else if ( IS_IMM(n) ) {
! 149: c = -ZTOIMM(n);
! 150: IMMTOZ(c,r);
! 151: return r;
! 152: } else {
1.1 noro 153: r = dupz(n);
154: SL(r) = -SL(r);
155: return r;
156: }
157: }
158:
1.3 ! noro 159:
1.1 noro 160: Z addz(Z n1,Z n2)
161: {
1.3 ! noro 162: int d1,d2,d,c;
! 163: Z r,r1;
! 164: struct oZ t;
1.1 noro 165:
166: if ( !n1 ) return dupz(n2);
167: else if ( !n2 ) return dupz(n1);
1.3 ! noro 168: else if ( IS_IMM(n1) ) {
! 169: if ( IS_IMM(n2) ) {
! 170: c = ZTOIMM(n1)+ZTOIMM(n2);
! 171: IMMTOZ(c,r);
! 172: return r;
! 173: } else {
! 174: c = ZTOIMM(n1);
! 175: if ( c < 0 ) {
! 176: t.p = -1; t.b[0] = -c;
1.1 noro 177: } else {
1.3 ! noro 178: t.p = 1; t.b[0] = c;
1.1 noro 179: }
1.3 ! noro 180: if ( (d2 = SL(n2)) < 0 ) d2 = -d2;
! 181: r = ZALLOC(d2+1);
! 182: _addz(&t,n2,r);
! 183: return r;
! 184: }
! 185: } else if ( IS_IMM(n2) ) {
! 186: c = ZTOIMM(n2);
! 187: if ( c < 0 ) {
! 188: t.p = -1; t.b[0] = -c;
1.1 noro 189: } else {
1.3 ! noro 190: t.p = 1; t.b[0] = c;
! 191: }
! 192: if ( (d1 = SL(n1)) < 0 ) d1 = -d1;
! 193: r = ZALLOC(d1+1);
! 194: _addz(n1,&t,r);
! 195: return r;
! 196: } else {
! 197: if ( (d1 = SL(n1)) < 0 ) d1 = -d1;
! 198: if ( (d2 = SL(n2)) < 0 ) d2 = -d2;
! 199: d = MAX(d1,d2)+1;
! 200: r = ZALLOC(d);
! 201: _addz(n1,n2,r);
! 202: if ( !SL(r) ) r = 0;
! 203: else if ( IS_IMMZ(r) ) {
! 204: IMMZTOZ(r,r1); r = r1;
1.1 noro 205: }
1.3 ! noro 206: return r;
1.1 noro 207: }
208: }
209:
210: Z subz(Z n1,Z n2)
211: {
1.3 ! noro 212: int d1,d2,d,c;
! 213: Z r,r1;
! 214: struct oZ t;
1.1 noro 215:
1.3 ! noro 216: if ( !n1 ) {
! 217: r = dupz(n2);
! 218: SL(r) = -SL(r);
! 219: return r;
! 220: } else if ( !n2 ) return dupz(n1);
! 221: else if ( IS_IMM(n1) ) {
! 222: if ( IS_IMM(n2) ) {
! 223: c = ZTOIMM(n1)-ZTOIMM(n2);
! 224: IMMTOZ(c,r);
! 225: return r;
! 226: } else {
! 227: c = ZTOIMM(n1);
! 228: if ( c < 0 ) {
! 229: t.p = -1; t.b[0] = -c;
! 230: } else {
! 231: t.p = 1; t.b[0] = c;
! 232: }
! 233: if ( (d2 = SL(n2)) < 0 ) d2 = -d2;
! 234: r = ZALLOC(d2+1);
! 235: _subz(&t,n2,r);
! 236: return r;
! 237: }
! 238: } else if ( IS_IMM(n2) ) {
! 239: c = ZTOIMM(n2);
! 240: if ( c < 0 ) {
! 241: t.p = -1; t.b[0] = -c;
! 242: } else {
! 243: t.p = 1; t.b[0] = c;
! 244: }
! 245: if ( (d1 = SL(n1)) < 0 ) d1 = -d1;
! 246: r = ZALLOC(d1+1);
! 247: _subz(n1,&t,r);
! 248: return r;
! 249: } else {
! 250: if ( (d1 = SL(n1)) < 0 ) d1 = -d1;
! 251: if ( (d2 = SL(n2)) < 0 ) d2 = -d2;
! 252: d = MAX(d1,d2)+1;
! 253: r = ZALLOC(d);
! 254: _subz(n1,n2,r);
! 255: if ( !SL(r) ) r = 0;
! 256: else if ( IS_IMMZ(r) ) {
! 257: IMMZTOZ(r,r1); r = r1;
1.1 noro 258: }
1.3 ! noro 259: return r;
1.1 noro 260: }
261: }
262:
263: Z mulz(Z n1,Z n2)
264: {
1.3 ! noro 265: int d1,d2,sgn,i;
1.1 noro 266: unsigned int u1,u2,u,l;
1.3 ! noro 267: Z r;
1.1 noro 268:
269: if ( !n1 || !n2 ) return 0;
1.3 ! noro 270:
! 271: if ( IS_IMM(n1) ) {
! 272: sgn = 1;
! 273: u1 = ZTOIMM(n1); if ( u1 < 0 ) { sgn = -sgn; u1 = -u1; }
! 274: if ( IS_IMM(n2) ) {
! 275: u2 = ZTOIMM(n2); if ( u2 < 0 ) { sgn = -sgn; u2 = -u2; }
1.1 noro 276: DM(u1,u2,u,l);
1.3 ! noro 277: if ( !u ) {
! 278: IMMTOZ(l,r);
1.1 noro 279: } else {
1.3 ! noro 280: r = ZALLOC(2); SL(r) = 2; BD(r)[1] = u; BD(r)[0] = l;
1.1 noro 281: }
282: } else {
1.3 ! noro 283: if ( (d2 = SL(n2)) < 0 ) { sgn = -sgn; d2 = -d2; }
! 284: r = ZALLOC(d2+1);
! 285: for ( i = d2; i >= 0; i-- ) BD(r)[i] = 0;
! 286: muln_1(BD(n2),d2,u1,BD(r));
! 287: SL(r) = BD(r)[d2]?d2+1:d2;
1.1 noro 288: }
1.3 ! noro 289: } else if ( IS_IMM(n2) ) {
! 290: sgn = 1;
! 291: u2 = ZTOIMM(n2); if ( u2 < 0 ) { sgn = -sgn; u2 = -u2; }
! 292: if ( (d1 = SL(n1)) < 0 ) { sgn = -sgn; d1 = -d1; }
! 293: r = ZALLOC(d1+1);
! 294: for ( i = d1; i >= 0; i-- ) BD(r)[i] = 0;
! 295: muln_1(BD(n1),d1,u2,BD(r));
! 296: SL(r) = BD(r)[d1]?d1+1:d1;
! 297: } else {
! 298: sgn = 1;
! 299: if ( (d1 = SL(n1)) < 0 ) { sgn = -sgn; d1 = -d1; }
! 300: if ( (d2 = SL(n2)) < 0 ) { sgn = -sgn; d2 = -d2; }
! 301: r = ZALLOC(d1+d2);
! 302: _mulz(n1,n2,r);
1.1 noro 303: }
1.3 ! noro 304: if ( sgn < 0 ) SL(r) = -SL(r);
! 305: return r;
1.1 noro 306: }
307:
1.3 ! noro 308: /* kokokara */
1.1 noro 309:
310: Z divsz(Z n1,Z n2)
311: {
1.3 ! noro 312: int sgn,d1,d2;
! 313: Z q;
1.1 noro 314:
315: if ( !n2 ) error("divsz : division by 0");
1.3 ! noro 316: if ( !n1 ) return 0;
! 317:
! 318: if ( IS_IMM(n1) ) {
! 319: if ( !IS_IMM(n2) )
! 320: error("divsz : cannot happen");
! 321: c = ZTOIMM(n1)/ZTOIMM(n2);
! 322: IMMTOZ(c,q);
! 323: return q;
1.1 noro 324: }
1.3 ! noro 325: if ( IS_IMM(n2) ) {
! 326: sgn = 1;
! 327: u2 = ZTOIMM(n2); if ( u2 < 0 ) { sgn = -sgn; u2 = -u2; }
! 328: diviz(n1,u2,&q);
! 329: if ( sgn < 0 ) SL(q) = -SL(q);
! 330: return q;
! 331: }
! 332:
! 333: sgn = 1;
! 334: if ( (d2 = SL(n2)) < 0 ) { sgn = -sgn; d2 = -d2; }
! 335: if ( d2 == 1 ) {
! 336: diviz(n1,BD(u2)[0],&q);
! 337: if ( sgn < 0 ) SL(q) = -SL(q);
! 338: return q;
! 339: }
! 340: if ( (d1 = SL(n1)) < 0 ) { sgn = -sgn; d1 = -d1; }
! 341: if ( d1 < d2 ) error("divsz : cannot happen");
! 342: return q;
1.1 noro 343: }
344:
345: Z gcdz(Z n1,Z n2)
346: {
1.3 ! noro 347: int d1,d2;
1.1 noro 348: N gcd;
349:
350: if ( !n1 ) return n2;
351: else if ( !n2 ) return n1;
352:
353: n1 = dupz(n1);
354: if ( SL(n1) < 0 ) SL(n1) = -SL(n1);
355: n2 = dupz(n2);
356: if ( SL(n2) < 0 ) SL(n2) = -SL(n2);
357: gcdn((N)n1,(N)n2,&gcd);
358: return (Z)gcd;
359: }
360:
361: int remzi(Z n,int m)
362: {
363: unsigned int *x;
364: unsigned int t,r;
365: int i;
366:
367: if ( !n ) return 0;
368: i = SL(n);
369: if ( i < 0 ) i = -i;
370: for ( i--, x = BD(n)+i, r = 0; i >= 0; i--, x-- ) {
371: #if defined(sparc)
372: r = dsa(m,r,*x);
373: #else
374: DSAB(m,r,*x,t,r);
375: #endif
376: }
377: return r;
378: }
379:
380: Z gcdz_cofactor(Z n1,Z n2,Z *c1,Z *c2)
381: {
382: Z gcd;
383:
384: gcd = gcdz(n1,n2);
385: *c1 = divsz(n1,gcd);
386: *c2 = divsz(n2,gcd);
387: return gcd;
388: }
389:
390: Z estimate_array_gcdz(Z *b,int n)
391: {
392: int m,i,j,sd;
393: Z *a;
394: Z s,t;
395:
396: a = (Z *)ALLOCA(n*sizeof(Z));
397: for ( i = j = 0; i < n; i++ ) if ( b[i] ) a[j++] = b[i];
398: n = j;
399: if ( !n ) return 0;
400: if ( n == 1 ) return a[0];
401:
402: m = n/2;
403: for ( i = 0, s = 0; i < m; i++ ) {
404: if ( !a[i] ) continue;
405: else s = (SL(a[i])<0)?subz(s,a[i]):addz(s,a[i]);
406: }
407: for ( t = 0; i < n; i++ ) {
408: if ( !a[i] ) continue;
409: else t = (SL(a[i])<0)?subz(t,a[i]):addz(t,a[i]);
410: }
411: return gcdz(s,t);
412: }
413:
414: Z array_gcdz(Z *b,int n)
415: {
416: int m,i,j,sd;
417: Z *a;
418: Z gcd;
419:
420: a = (Z *)ALLOCA(n*sizeof(Z));
421: for ( i = j = 0; i < n; i++ ) if ( b[i] ) a[j++] = b[i];
422: n = j;
423: if ( !n ) return 0;
424: if ( n == 1 ) return a[0];
425: gcd = a[0];
426: for ( i = 1; i < n; i++ )
427: gcd = gcdz(gcd,a[i]);
428: return gcd;
429: }
430:
431: void _copyz(Z n1,Z n2)
432: {
433: int n,i;
434:
435: if ( !n1 || !SL(n1) ) SL(n2) = 0;
436: else {
437: n = SL(n2) = SL(n1);
438: if ( n < 0 ) n = -n;
439: for ( i = 0; i < n; i++ ) BD(n2)[i] = BD(n1)[i];
440: }
441: }
442:
443: void _addz(Z n1,Z n2,Z nr)
444: {
1.3 ! noro 445: int d1,d2;
1.1 noro 446:
447: if ( !n1 || !SL(n1) ) _copyz(n2,nr);
448: else if ( !n2 || !SL(n2) ) _copyz(n1,nr);
1.3 ! noro 449: else if ( (d1=SL(n1)) > 0 )
! 450: if ( (d2=SL(n2)) > 0 )
! 451: SL(nr) = _addz_main(BD(n1),d1,BD(n2),d2,BD(nr));
1.1 noro 452: else
1.3 ! noro 453: SL(nr) = _subz_main(BD(n1),d1,BD(n2),-d2,BD(nr));
! 454: else if ( (d2=SL(n2)) > 0 )
! 455: SL(nr) = _subz_main(BD(n2),d2,BD(n1),-d1,BD(nr));
1.1 noro 456: else
1.3 ! noro 457: SL(nr) = -_addz_main(BD(n1),-d1,BD(n2),-d2,BD(nr));
1.1 noro 458: }
459:
460: void _subz(Z n1,Z n2,Z nr)
461: {
1.3 ! noro 462: int d1,d2;
1.1 noro 463:
464: if ( !n1 || !SL(n1) ) _copyz(n2,nr);
465: else if ( !n2 || !SL(n2) ) {
466: _copyz(n1,nr);
467: SL(nr) = -SL(nr);
1.3 ! noro 468: } else if ( (d1=SL(n1)) > 0 )
! 469: if ( (d2=SL(n2)) > 0 )
! 470: SL(nr) = _subz_main(BD(n1),d1,BD(n2),d2,BD(nr));
1.1 noro 471: else
1.3 ! noro 472: SL(nr) = _addz_main(BD(n1),d1,BD(n2),-d2,BD(nr));
! 473: else if ( (d2=SL(n2)) > 0 )
! 474: SL(nr) = -_addz_main(BD(n1),-d1,BD(n2),d2,BD(nr));
1.1 noro 475: else
1.3 ! noro 476: SL(nr) = -_subz_main(BD(n1),-d1,BD(n2),-d2,BD(nr));
1.1 noro 477: }
478:
479: void _mulz(Z n1,Z n2,Z nr)
480: {
1.3 ! noro 481: int d1,d2;
1.1 noro 482: int i,d,sgn;
483: unsigned int mul;
484: unsigned int *m1,*m2;
485:
486: if ( !n1 || !SL(n1) || !n2 || !SL(n2) )
487: SL(nr) = 0;
488: else {
1.3 ! noro 489: d1 = SL(n1); d2 = SL(n2);
1.1 noro 490: sgn = 1;
1.3 ! noro 491: if ( d1 < 0 ) { sgn = -sgn; d1 = -d1; }
! 492: if ( d2 < 0 ) { sgn = -sgn; d2 = -d2; }
! 493: d = d1+d2;
1.1 noro 494: for ( i = d-1, m1 = BD(nr); i >= 0; i-- ) *m1++ = 0;
1.3 ! noro 495: for ( i = 0, m1 = BD(n1), m2 = BD(n2); i < d2; i++, m2++ )
! 496: if ( mul = *m2 ) muln_1(m1,d1,mul,BD(nr)+i);
1.1 noro 497: SL(nr) = sgn*(BD(nr)[d-1]?d:d-1);
498: }
499: }
500:
501: int _addz_main(unsigned int *m1,int d1,unsigned int *m2,int d2,unsigned int *mr)
502: {
503: int d,i;
504: unsigned int tmp,c;
505: unsigned int *t;
506:
507: if ( d2 > d1 ) {
508: t = m1; m1 = m2; m2 = t;
509: d = d1; d1 = d2; d2 = d;
510: }
511: #if defined(VISUAL)
512: __asm {
513: push esi
514: push edi
515: mov esi,m1
516: mov edi,m2
517: mov ebx,mr
518: mov ecx,d2
519: xor eax,eax
520: Lstart__addz:
521: mov eax,DWORD PTR [esi]
522: mov edx,DWORD PTR [edi]
523: adc eax,edx
524: mov DWORD PTR [ebx],eax
525: lea esi,DWORD PTR [esi+4]
526: lea edi,DWORD PTR [edi+4]
527: lea ebx,DWORD PTR [ebx+4]
528: dec ecx
529: jnz Lstart__addz
530: pop edi
531: pop esi
532: mov eax,0
533: adc eax,eax
534: mov c,eax
535: }
536: #elif defined(i386)
537: asm volatile("\
538: movl %1,%%esi;\
539: movl %2,%%edi;\
540: movl %3,%%ebx;\
541: movl %4,%%ecx;\
542: testl %%eax,%%eax;\
543: Lstart__addz:;\
544: movl (%%esi),%%eax;\
545: movl (%%edi),%%edx;\
546: adcl %%edx,%%eax;\
547: movl %%eax,(%%ebx);\
548: leal 4(%%esi),%%esi;\
549: leal 4(%%edi),%%edi;\
550: leal 4(%%ebx),%%ebx;\
551: decl %%ecx;\
552: jnz Lstart__addz;\
553: movl $0,%%eax;\
554: adcl %%eax,%%eax;\
555: movl %%eax,%0"\
556: :"=m"(c)\
557: :"m"(m1),"m"(m2),"m"(mr),"m"(d2)\
558: :"eax","ebx","ecx","edx","esi","edi");
559: #else
560: for ( i = 0, c = 0, mr = BD(nr); i < d2; i++, m1++, m2++, mr++ ) {
561: tmp = *m1 + *m2;
562: if ( tmp < *m1 ) {
563: tmp += c;
564: c = 1;
565: } else {
566: tmp += c;
567: c = tmp < c ? 1 : 0;
568: }
569: *mr = tmp;
570: }
571: #endif
572: for ( i = d2, m1 += d2, mr += d2; (i < d1) && c ; i++ ) {
573: tmp = *m1++ + c;
574: c = tmp < c ? 1 : 0;
575: *mr++ = tmp;
576: }
577: for ( ; i < d1; i++ )
578: *mr++ = *m1++;
579: *mr = c;
580: return (c?d1+1:d1);
581: }
582:
583: int _subz_main(unsigned int *m1,int d1,unsigned int *m2,int d2,unsigned int *mr)
584: {
585: int d,i,sgn;
586: unsigned int t,tmp,br;
587: unsigned int *m;
588:
589: if ( d1 > d2 ) sgn = 1;
590: else if ( d1 < d2 ) sgn = -1;
591: else {
592: for ( i = d1-1; i >= 0 && m1[i] == m2[i]; i-- );
593: if ( i < 0 ) return 0;
594: if ( m1[i] > m2[i] ) sgn = 1;
595: else if ( m1[i] < m2[i] ) sgn = -1;
596: }
597: if ( sgn < 0 ) {
598: m = m1; m1 = m2; m2 = m;
599: d = d1; d1 = d2; d2 = d;
600: }
601: #if defined(VISUAL)
602: __asm {
603: push esi
604: push edi
605: mov esi,m1
606: mov edi,m2
607: mov ebx,mr
608: mov ecx,d2
609: xor eax,eax
610: Lstart__subz:
611: mov eax,DWORD PTR [esi]
612: mov edx,DWORD PTR [edi]
613: sbb eax,edx
614: mov DWORD PTR [ebx],eax
615: lea esi,DWORD PTR [esi+4]
616: lea edi,DWORD PTR [edi+4]
617: lea ebx,DWORD PTR [ebx+4]
618: dec ecx
619: jnz Lstart__subz
620: pop edi
621: pop esi
622: mov eax,0
623: adc eax,eax
624: mov br,eax
625: }
626: #elif defined(i386)
627: asm volatile("\
628: movl %1,%%esi;\
629: movl %2,%%edi;\
630: movl %3,%%ebx;\
631: movl %4,%%ecx;\
632: testl %%eax,%%eax;\
633: Lstart__subz:;\
634: movl (%%esi),%%eax;\
635: movl (%%edi),%%edx;\
636: sbbl %%edx,%%eax;\
637: movl %%eax,(%%ebx);\
638: leal 4(%%esi),%%esi;\
639: leal 4(%%edi),%%edi;\
640: leal 4(%%ebx),%%ebx;\
641: decl %%ecx;\
642: jnz Lstart__subz;\
643: movl $0,%%eax;\
644: adcl %%eax,%%eax;\
645: movl %%eax,%0"\
646: :"=m"(br)\
647: :"m"(m1),"m"(m2),"m"(mr),"m"(d2)\
648: :"eax","ebx","ecx","edx","esi","edi");
649: #else
650: for ( i = 0, br = 0, mr = BD(nr); i < d2; i++, mr++ ) {
651: t = *m1++;
652: tmp = *m2++ + br;
653: if ( br > 0 && !tmp ) {
654: /* tmp = 2^32 => br = 1 */
655: }else {
656: tmp = t-tmp;
657: br = tmp > t ? 1 : 0;
658: *mr = tmp;
659: }
660: }
661: #endif
662: for ( i = d2, m1 += d2, mr += d2; (i < d1) && br; i++ ) {
663: t = *m1++;
664: tmp = t - br;
665: br = tmp > t ? 1 : 0;
666: *mr++ = tmp;
667: }
668: for ( ; i < d1; i++ )
669: *mr++ = *m1++;
670: for ( i = d1-1, mr--; i >= 0 && !*mr--; i-- );
671: return sgn*(i+1);
672: }
673:
674: /* XXX */
675:
676: void printz(Z n)
677: {
678: int sd;
679:
680: if ( !n )
681: fprintf(asir_out,"0");
682: else {
683: if ( (sd = SL(n)) < 0 ) { SL(n) = -SL(n); fprintf(asir_out,"-"); }
684: printn((N)n);
685: if ( sd < 0 ) SL(n) = -SL(n);
1.2 noro 686: }
687: }
688:
689: /*
690: * Dx^k*x^l = W(k,l,0)*x^l*Dx^k+W(k,l,1)*x^(l-1)*x^(k-1)*+...
691: *
692: * t = [W(k,l,0) W(k,l,1) ... W(k,l,min(k,l)]
693: * where W(k,l,i) = i! * kCi * lCi
694: */
695:
696: void mkwcz(int k,int l,Z *t)
697: {
698: int i,n,up,low;
699: N nm,d,c;
700:
701: n = MIN(k,l);
702: for ( t[0] = (Z)ONEN, i = 1; i <= n; i++ ) {
703: DM(k-i+1,l-i+1,up,low);
704: if ( up ) {
705: nm = NALLOC(2); PL(nm) = 2; BD(nm)[0] = low; BD(nm)[1] = up;
706: } else {
707: nm = NALLOC(1); PL(nm) = 1; BD(nm)[0] = low;
708: }
709: kmuln((N)t[i-1],nm,&d); divin(d,i,&c); t[i] = (Z)c;
1.1 noro 710: }
711: }