Annotation of OpenXM_contrib2/asir2000/engine/igcdhack.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/engine/igcdhack.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
! 2: #if 0
! 3: #include "base.h"
! 4:
! 5: #define SWAP(a,b,Type) { Type temp=a; a = b; b = temp; }
! 6: #define SIGNED_VAL(a,s) ((s)>0?(a): -(a))
! 7: #endif
! 8:
! 9: #undef FULLSET
! 10:
! 11:
! 12: static int negate_nbd( b, l )
! 13: int *b, l;
! 14: /* It is assumed that l > 0 and b[0] != 0 */
! 15: {
! 16: int i, nz;
! 17:
! 18: *b = BASE - *b;
! 19: for ( nz = i = 1 ; i < l; ) {
! 20: b++, i++;
! 21: if ( (*b ^= BMASK) != 0 ) nz = i;
! 22: }
! 23: return nz;
! 24: }
! 25:
! 26: /****************/
! 27: /****************/
! 28:
! 29: static int abs_U_V_maxrshift( u, lu, v, lv, a )
! 30: int *u, lu, *v, lv, *a;
! 31: /*
! 32: * Compute | U - V | >> (as much as possible) in a[], where U = u[0:lu-1]
! 33: * and V = v[0:lv-1]. Value returned is the length of the result
! 34: * multiplied by the sign of (U-V).
! 35: */
! 36: {
! 37: int sign, rsh, t, b, i, i_1, l, *p = a;
! 38:
! 39: /* first, determine U <, ==, > V */
! 40: sign = 1;
! 41: if ( lu == lv ) {
! 42: int i;
! 43:
! 44: for ( i = lu -1; i >= 0; i-- ) {
! 45: if ( u[i] == v[i] ) continue;
! 46: lu = lv = i + 1;
! 47: if ( u[i] > v[i] ) goto U_V;
! 48: else goto swap;
! 49: }
! 50: return 0;
! 51: } else if ( lu > lv ) goto U_V;
! 52: swap:
! 53: SWAP( lu, lv, int ); SWAP( u, v, int * );
! 54: sign = -1;
! 55: U_V: /* U > V, i.e., (lu > lv) || (lu == lv && u[lu-1] > v[lv-1]) */
! 56: for ( i = 0; i < lv; i++ ) if ( u[i] != v[i] ) break;
! 57: if ( i != 0 ) u += i, lu -= i, v += i, lv -= i;
! 58: if ( lv == 0 ) { /* no more pv[], and lu > 0 */
! 59: i = trailingzeros_nbd( u, &rsh );
! 60: i = rshift_nbd( u, lu, rsh, i, a );
! 61: return SIGNED_VAL( i, sign );
! 62: }
! 63: /**/
! 64: if ( (t = u[0] - v[0]) >= 0 ) b = 0;
! 65: else b = 1, t += BASE;
! 66: TRAILINGZEROS( t, rsh );
! 67: /**/
! 68: if ( rsh != 0 ) {
! 69: int w, lsh = BSH - rsh;
! 70:
! 71: /* subtract v[*] from u[*] */
! 72: for ( i = 1, l = 0; i < lv; i++, t = w >> rsh ) {
! 73: if ( (w = u[i] - v[i] - b) >= 0 ) b = 0;
! 74: else b = 1, w += BASE;
! 75: if ( (*p++ = t | ((w << lsh) & BMASK)) != 0 ) l = i;
! 76: }
! 77: /* subtraction of V completed. Note there may exist non-0 borrow b */
! 78: if ( i >= lu ) /* lu == lv, and b==0 is implied */ goto chk_nz;
! 79: /* we still have (lu-i) elms in U */
! 80: if ( b != 0 ) { /* non-zero borrow */
! 81: if ( (w = u[i++] -1) != 0 ) t |= ((w << lsh) & BMASK);
! 82: if ( i >= lu ) { /* lu == lv+1. The M.S. word w may beome 0 */
! 83: if ( w != 0 ) {
! 84: *p++ = t;
! 85: l = lv;
! 86: t = w >> rsh;
! 87: } else lu--;
! 88: goto chk_nz;
! 89: }
! 90: *p++ = t;
! 91: if ( w < 0 ) {
! 92: while ( (w = u[i++]) == 0 ) *p++ = BMASK;
! 93: w -= 1;
! 94: *p++ = (BMASK >> rsh) | ((w << lsh) & BMASK);
! 95: }
! 96: t = w >> rsh;
! 97: }
! 98: for ( ; i < lu; i++, t = w >> rsh ) *p++ = t | (((w = u[i]) << lsh) & BMASK);
! 99: l = i -1;
! 100: chk_nz:
! 101: if ( t != 0 ) *p = t, l = /*1 + (p - a)*/ lu;
! 102: goto ret_nz;
! 103: }
! 104: /* rsh == 0 : need not be right-shifted */
! 105: *p++ = t, l = 1;
! 106: for ( i = 1; i < lv; ) {
! 107: if ( (t = u[i] - v[i] - b) >= 0 ) b = 0;
! 108: else b = 1, t += BASE;
! 109: i++;
! 110: if ( (*p++ = t) != 0 ) l = i;
! 111: }
! 112: if ( i >= lu ) /* lu == lv, and b==0 is implied */ goto ret_nz;
! 113: if ( b != 0 ) { /* non-zero borrow */
! 114: t = u[i++] -1;
! 115: if ( i >= lu ) /* lu == lv+1. The M.S. word w may beome 0 */ goto chk_nz;
! 116: if ( t < 0 ) {
! 117: do *p++ = BMASK; while ( (t = u[i++]) == 0 );
! 118: t -= 1;
! 119: if ( i >= lu ) {
! 120: l = lu;
! 121: if ( t == 0 ) l--;
! 122: else *p = t;
! 123: goto ret_nz;
! 124: }
! 125: }
! 126: *p++ = t;
! 127: }
! 128: for ( ; i < lu; i++ ) *p++ = u[i];
! 129: l = lu;
! 130: /**/
! 131: ret_nz:
! 132: return SIGNED_VAL( l, sign );
! 133: }
! 134:
! 135: /****************/
! 136: /****************/
! 137:
! 138: #ifdef CALL
! 139: static int aV_plus_c(), aV_plus_c_sh();
! 140: #endif
! 141: static int abs_U_aV_c(), abs_U_aV_c_sh();
! 142: static int abs_aVplusc_U(), abs_aVplusc_U_sh();
! 143:
! 144: static int abs_U_aV_maxrshift( u, lu, a, v, lv, ans )
! 145: int *u, lu, a, *v, lv, *ans;
! 146: {
! 147: int i, t, c, rsh, l, *p = ans, w, lsh;
! 148: int hi;
! 149:
! 150: if ( a == 1 ) return( (l = abs_U_V_maxrshift( u, lu, v, lv, ans )) >= 0 ? l : -l );
! 151: /**/
! 152: for ( l = lu <= lv ? lu : lv, i = c = 0; i < l; ) {
! 153: DMA( a, v[i], c, hi, t, not_used );
! 154: c = hi;
! 155: if ( (t -= u[i++]) != 0 ) goto non0w;
! 156: }
! 157: /**/
! 158: if ( i < lv ) { /* no more u[] elm.'s, and we still have non-zero v[] elm's */
! 159: do {
! 160: DMA( a, v[i], c, hi, t, not_used );
! 161: i++, c = hi;
! 162: } while ( t == 0 );
! 163: TRAILINGZEROS( t, rsh );
! 164: v += i, lv -= i;
! 165: #ifdef CALL
! 166: if ( rsh == 0 ) {
! 167: *p++ = t;
! 168: return( aV_plus_c( a, v, lv, c, p ) + 1 );
! 169: } else return aV_plus_c_sh( a, v, lv, c, rsh, t, ans );
! 170: #else
! 171: i = 0;
! 172: if ( rsh == 0 ) {
! 173: *p++ = t;
! 174: for ( ; i < lv; i++, c = hi, *p++ = t ) { DMA( a, v[i], c, hi, t, not_used ); }
! 175: if ( c == 0 ) return( i + 1 );
! 176: *p = c;
! 177: return( i + 2 );
! 178: } else {
! 179: for ( lsh = BSH - rsh; i < lv; i++, c = hi, t = w >> rsh ) {
! 180: DMA( a, v[i], c, hi, w, not_used );
! 181: *p++ = t | ((w << lsh) & BMASK);
! 182: }
! 183: if ( c != 0 ) {
! 184: *p++ = t | ((c << lsh) & BMASK);
! 185: i++;
! 186: t = c >> rsh;
! 187: }
! 188: if ( t == 0 ) return i;
! 189: *p = t;
! 190: return( i + 1 );
! 191: }
! 192: #endif
! 193: }
! 194: /* no more v[] elm.'s */
! 195: if ( i >= lu ) {
! 196: if ( c == 0 ) return 0;
! 197: c_sh:
! 198: while ( (c&1) == 0 ) c >>= 1;
! 199: *p = c;
! 200: return 1;
! 201: }
! 202: t = u[i++] - c;
! 203: if ( i >= lu ) {
! 204: if ( t == 0 ) return 0;
! 205: c = t < 0 ? -t : t;
! 206: goto c_sh;
! 207: } else if ( t == 0 ) {
! 208: while ( (t = u[i++]) == 0 ) ;
! 209: c = 0;
! 210: } else if ( t < 0 ) t += BASE, c = 1;
! 211: else c = 0;
! 212: /* diff turned out to be > 0 */
! 213: u += i, lu -= i;
! 214: TRAILINGZEROS( t, rsh );
! 215: i = 0;
! 216: if ( rsh == 0 ) {
! 217: *p++ = t;
! 218: if ( c != 0 ) {
! 219: for ( ; i < lu; *p++ = BMASK ) if ( (t = u[i++]) != 0 ) break;
! 220: t--;
! 221: if ( i >= lu && t == 0 ) return i;
! 222: *p++ = t;
! 223: }
! 224: for ( ; i < lu; i++ ) *p++ = u[i];
! 225: } else { /* rsh != 0 */
! 226: lsh = BSH - rsh;
! 227: if ( c != 0 ) { /* there still exists u[] elms.'s because diff is known to be > 0 */
! 228: if ( (w = u[i++]) == 0 ) {
! 229: *p++ = t | ((BMASK << lsh) & BMASK);
! 230: for ( ; (w = u[i++]) == 0; *p++ = BMASK ) ;
! 231: t = BMASK >> rsh;
! 232: }
! 233: w--;
! 234: *p++ = t | ((w << lsh) & BMASK);
! 235: t = w >> rsh;
! 236: }
! 237: for ( ; i < lu; i++, t = w >> rsh ) *p++ = t | (((w = u[i]) << lsh) & BMASK);
! 238: if ( t == 0 ) return( i );
! 239: *p = t;
! 240: }
! 241: return( i + 1 );
! 242: /**/
! 243: non0w:
! 244: u += i, lu -= i, v += i, lv -= i;
! 245: /* | t + BASE*(c + a*v[0:lv-1] - u[0:lu-1] | */
! 246: if ( lv < lu ) { /* a*V < U if lv<lu-1, and very likely if lv==lu-1 */
! 247: t = -t; /* We shall compute U - a*V */
! 248: if ( t < 0 ) t += BASE, c++;
! 249: } else /* a*V > U if lv>lu, and very likely even if lv==lu */
! 250: if ( t < 0 ) t += BASE, c--;
! 251: TRAILINGZEROS( t, rsh );
! 252: return( rsh == 0 ? (lv < lu ? abs_U_aV_c( t, u, lu, a, v, lv, c, ans ) :
! 253: abs_aVplusc_U( t, a, v, lv, c, u, lu, ans )) :
! 254: lv < lu ? abs_U_aV_c_sh( t, u, lu, a, v, lv, c, rsh, ans ) :
! 255: abs_aVplusc_U_sh( t, a, v, lv, c, u, lu, rsh, ans ) );
! 256: }
! 257:
! 258: /*********/
! 259:
! 260: static int aV_plus_c( a, v, lv, c, p )
! 261: int a, *v, lv, c, *p;
! 262: /*
! 263: * Compute (a*V + c) in p[], where V = v[0:lv-1].
! 264: */
! 265: {
! 266: int i, t;
! 267: int hi;
! 268:
! 269: for ( i = 0; i < lv; i++, *p++ = t, c = hi ) { DMA( a, v[i], c, hi, t, not_used); }
! 270: if ( c == 0 ) return i;
! 271: *p = c;
! 272: return( i + 1 );
! 273: }
! 274:
! 275: static int aV_plus_c_sh( a, v, lv, c, rsh, t, p )
! 276: int a, *v, lv, c, rsh, *p;
! 277: /*
! 278: * Compute (t + BASE*((a*V + c) >> rsh)) in p[], where V = v[0:lv-1].
! 279: */
! 280: {
! 281: int i, w, lsh = BSH - rsh;
! 282: int hi;
! 283:
! 284: for ( i = 0; i < lv; i++, c = hi, t = w >> rsh ) {
! 285: DMA( a, v[i], c, hi, w, not_used );
! 286: *p++ = t | ((w << lsh) & BMASK);
! 287: }
! 288: if ( c != 0 ) {
! 289: *p++ = t | ((c << lsh) & BMASK);
! 290: i++;
! 291: t = c >> rsh;
! 292: }
! 293: if ( t == 0 ) return i;
! 294: *p = t;
! 295: return( i + 1 );
! 296: }
! 297:
! 298: /*********/
! 299:
! 300: static int abs_aVplusc_U( t, a, v, lv, c, u, lu, ans )
! 301: int t, a, *v, lv, c, *u, lu, *ans;
! 302: /* compute | t + BASE*(a*V+c - U) | in ans[], where V=v[0:lv-1], U=u[0:lu-1],
! 303: * and a > 1, -1 <= c < BASE, lv >= lu >=0 && t != 0 are assumed.
! 304: */
! 305: {
! 306: int i, l, b, *p = ans;
! 307: int hi;
! 308:
! 309: *p++ = t, l = 1;
! 310: if ( c < 0 ) c = 0, b = 1;
! 311: else b = 0;
! 312: for ( i = 0; i < lu; *p++ = t ) {
! 313: DMA( a, v[i], c, hi, t, not_used );
! 314: t -= u[i++], c = hi;
! 315: if ( b != 0 ) t--;
! 316: b = 0;
! 317: if ( t > 0 ) l = i + 1;
! 318: else if ( t < 0 ) {
! 319: t += BASE, l = i + 1;
! 320: if ( c != 0 ) c--;
! 321: else b = 1;
! 322: }
! 323: }
! 324: /* at this point, i == lu and we have written (i+1) entries in ans[] */
! 325: if ( i >= lv ) { /* no more elms in v[] */
! 326: if ( b != 0 ) {
! 327: if ( i == 0 ) { ans[i] = BASE - t; return 1; }
! 328: l = negate_nbd( ans, i ); /* <================ */
! 329: /* if ( (t = BMASK ^ ans[i]) == 0 ) return l;*/
! 330: if ( (t ^= BMASK) == 0 ) return l;
! 331: ans[i] = t;
! 332: return( i + 1 );
! 333: }
! 334: if ( c == 0 ) return l;
! 335: *p = c;
! 336: return( i + 2 );
! 337: }
! 338: /* diff turned out to be > 0 */
! 339: if ( b != 0 ) { /* c == 0 */
! 340: while ( (b = v[i++]) == 0 ) *p++ = BMASK;
! 341: DM( a, b, hi, t );
! 342: if ( t != 0 ) t--, c = hi;
! 343: else t = BMASK, c = hi - 1;
! 344: *p++ - t;
! 345: }
! 346: #ifdef CALL
! 347: return( aV_plus_c( a, &v[i], lv - i, c, p ) + i +1 );
! 348: #else
! 349: for ( ; i < lv; i++, *p++ = t, c = hi ) { DMA( a, v[i], c, hi, t, not_used );}
! 350: if ( c == 0 ) return( i + 1 );
! 351: *p = c;
! 352: return( i + 2 );
! 353: #endif
! 354: }
! 355:
! 356: static int abs_aVplusc_U_sh( t, a, v, lv, c, u, lu, rsh, ans )
! 357: int t, a, *v, lv, c, *u, lu, rsh, *ans;
! 358: {
! 359: int i, l, w, b, lsh = BSH - rsh, *p = ans;
! 360: int hi;
! 361:
! 362: if ( c < 0 ) c = 0, b = 1;
! 363: else b = 0;
! 364: for ( i = 0; i < lu; t = w >> rsh ) {
! 365: DMA( a, v[i], c, hi, w, not_used );
! 366: w -= u[i++], c = hi;
! 367: if ( b != 0 ) w--;
! 368: b = 0;
! 369: if ( w < 0 ) {
! 370: w += BASE;
! 371: if ( c != 0 ) c--;
! 372: else b = 1;
! 373: }
! 374: if ( (*p++ = t | ((w << lsh) & BMASK)) != 0 ) l = i;
! 375: }
! 376: /* at this point, i == lu and we have written i entries in ans[] */
! 377: if ( i >= lv ) {
! 378: if ( b != 0 ) { /* diff turned out to be < 0 */
! 379: if ( i == 0 ) {
! 380: *p = (BASE >> rsh) - t;
! 381: return 1;
! 382: }
! 383: l = negate_nbd( ans, i ); /* <================ */
! 384: if ( (t ^= (BMASK >> rsh)) == 0 ) return l;
! 385: *p = t;
! 386: return( i + 1 );
! 387: }
! 388: if ( c == 0 ) {
! 389: if ( t == 0 ) return l;
! 390: *p = t;
! 391: return( i + 1 );
! 392: }
! 393: *p++ = t | ((c << lsh) & BMASK);
! 394: if ( (t = c >> rsh) == 0 ) return( i + 1 );
! 395: *p = t;
! 396: return( i + 2 );
! 397: }
! 398: /* lv >= lu+1 = i+1, and diff is > 0 */
! 399: if ( b != 0 ) { /* c == 0 */
! 400: if ( (w = v[i++]) == 0 ) {
! 401: *p++ = t | ((BMASK << lsh) & BMASK);
! 402: while ( (w = v[i++]) == 0 ) *p++ = BMASK;
! 403: t = BMASK >> rsh;
! 404: } else {
! 405: DM( a, w, c, b );
! 406: if ( b != 0 ) b--;
! 407: else b = BMASK, c--;
! 408: *p++ = t | ((b << lsh) & BMASK);
! 409: t = b >> rsh;
! 410: }
! 411: }
! 412: /**/
! 413: #ifdef CALL
! 414: return( aV_plus_c_sh( a, &v[i], lv - i, c, rsh, t, p ) + i );
! 415: #else
! 416: for ( ; i < lv; i++, t = w >> rsh, c = hi ) {
! 417: DMA( a, v[i], c, hi, w, not_used );
! 418: *p++ = t | ((w << lsh) & BMASK);
! 419: }
! 420: if ( c != 0 ) {
! 421: *p++ = t | ((c << lsh) & BMASK);
! 422: i++, t = c >> rsh;
! 423: }
! 424: if ( t == 0 ) return i;
! 425: *p = t;
! 426: return( i + 1 );
! 427: #endif
! 428: }
! 429:
! 430: /*********/
! 431:
! 432: static int abs_U_aV_c( t, u, lu, a, v, lv, c, ans )
! 433: int t, *u, lu, a, *v, lv, c, *ans;
! 434: /* compute | t + BASE*(U - a*V - c) | in ans[], where U = u[0:lu-1], V = v[0:lv-1],
! 435: * c is <= BASE-1, and lu >= lv+1 && t != 0 assumed.
! 436: */
! 437: {
! 438: int i, w, l, *p = ans;
! 439: int hi;
! 440:
! 441: *p++ = t, l = 1;
! 442: for ( i = 0; i < lv; *p++ = t ) {
! 443: DMA( a, v[i], c, hi, w, not_used );
! 444: t = u[i] - w;
! 445: i++, c = hi;
! 446: if ( t != 0 ) l = i + 1;
! 447: if ( t < 0 ) t += BASE, c++;
! 448: }
! 449: /* at this point, i == lv */
! 450: if ( lu == lv+1 ) {
! 451: if ( (t = u[i]) == c ) return l;
! 452: else if ( t > c ) { *p = t - c; return( lu + 1 ); }
! 453: l = negate_nbd( ans, lu ); /* <================ */
! 454: if ( (c -= (t + 1)) == 0 ) return l;
! 455: *p = c;
! 456: return( lu + 1 );
! 457: }
! 458: /* lu >= lv+2 and the diff turned out to be >= 0 */
! 459: if ( c != 0 ) {
! 460: /* note c <= BASE. If c == BASE, we must care about the length of the result
! 461: of the case when lu == lv+2 && u[lu-2:lu-1]=u[i:i+1] == BASE */
! 462: if ( c >= BASE && i+2 == lu && u[i+1] == 1 ) {
! 463: /* m.s. word of the diff becomes 0 */
! 464: if ( (t = u[i] - (c & BMASK)) == 0 ) return l;
! 465: *p = t;
! 466: return lu;
! 467: }
! 468: for ( ; 1; c = 1, *p++ = t + BASE ) if ( (t = u[i++] - c) >= 0 ) break;
! 469: if ( t == 0 && i >= lu ) return lu;
! 470: *p++ = t;
! 471: }
! 472: for ( ; i < lu; i++ ) *p++ = u[i];
! 473: return( lu + 1 );
! 474: }
! 475:
! 476: static int abs_U_aV_c_sh( t, u, lu, a, v, lv, c, rsh, ans )
! 477: int t, *u, lu, a, *v, lv, c, rsh, *ans;
! 478: /* r-shift version of the above */
! 479: {
! 480: int i, w, l, lsh = BSH - rsh, *p = ans;
! 481: int hi;
! 482:
! 483: for ( l = i = 0; i < lv; t = w >> rsh ) {
! 484: DMA( a, v[i], c, hi, w, not_used );
! 485: w = u[i] - w, c = hi;
! 486: i++;
! 487: if ( w < 0 ) w += BASE, c++;
! 488: if ( (*p++ = t | ((w << lsh) & BMASK)) != 0 ) l = i;
! 489: }
! 490: /* at this point, i == lv, and we have written lv elm.s in ans[] */
! 491: if ( lu == lv+1 ) {
! 492: if ( (w = u[i] - c) == 0 ) {
! 493: if ( t == 0 ) return l;
! 494: *p = t;
! 495: return lu; /* == lv+1 */
! 496: } else if ( w > 0 ) {
! 497: l0: *p++ = t | ((w << lsh) & BMASK);
! 498: if ( (w >>= rsh) == 0 ) return lu; /* == lv+1 */
! 499: *p = w;
! 500: return( lu + 1 );
! 501: }
! 502: /* diff turned out to be < 0 */
! 503: w = - w - 1;
! 504: if ( lv == 0 ) { /* lu == 1 */
! 505: t = (BASE >> rsh) - t;
! 506: if ( w != 0 ) goto l0;
! 507: *p = t;
! 508: return 1;
! 509: }
! 510: l = negate_nbd( ans, lv ); /* <================ */
! 511: t ^= (BMASK >> rsh);
! 512: if ( w != 0 ) goto l0;
! 513: if ( t == 0 ) return l;
! 514: *p = t;
! 515: return lu; /* == lv + 1 */
! 516: }
! 517: /* now, lu >= lv+2 == i+2 */
! 518: if ( c != 0 ) {
! 519: /* note c <= BASE. If c == BASE, we must care about the length of the result
! 520: of the case when lu == lv+2 && u[lu-2:lu-1]=u[i:i+1] == BASE */
! 521: if ( c >= BASE && i+2 == lu && t == 0 && u[i] == 0 && u[i+1] == 1 )
! 522: return l; /* m.s. word of the diff has become 0 */
! 523: for ( ; 1; t = w >> rsh, c = 1 ) {
! 524: if ( (w = u[i++] - c) >= 0 ) break;
! 525: w += BASE;
! 526: *p++ = t | ((w << lsh) & BMASK);
! 527: }
! 528: t |= ((w << lsh) & BMASK);
! 529: w >>= rsh;
! 530: if ( i >= lu ) {
! 531: if ( w != 0 ) {
! 532: *p++ = t;
! 533: *p = w;
! 534: return( lu + 1 );
! 535: } else if ( t == 0 ) return( i - 1 );
! 536: else { *p = t; return i; }
! 537: }
! 538: *p++ = t;
! 539: t = w;
! 540: }
! 541: for ( ; i < lu; i++, t = w >> rsh ) *p++ = t | (((w = u[i]) << lsh) & BMASK);
! 542: if ( t == 0 ) return i;
! 543: *p = t;
! 544: return( i + 1 );
! 545: }
! 546:
! 547: /****************/
! 548: /****************/
! 549: #ifdef FULLSET
! 550: static int abs_axU_V_maxrshift(), abs_axU_bxV_maxrshift();
! 551:
! 552: static int abs_aU_bV_maxrshift( a, u, lu, b, v, lv, ans )
! 553: int a, *u, lu, b, *v, lv, *ans;
! 554: /*
! 555: * Compute | a*U - b*V | >> (as much as possible) in ans[], where U=u[0:lu-1]
! 556: * and V=v[0:lv-1].
! 557: */
! 558: {
! 559: if ( a == 1 ) return abs_U_aV_maxrshift( u, lu, b, v, lv, ans );
! 560: else if ( b == 1 ) return abs_U_aV_maxrshift( v, lv, a, u, lu, ans );
! 561: else if ( a == b ) return abs_axU_V_maxrshift( a, u, lu, v, lv, ans );
! 562: return abs_axU_bxV_maxrshift( a, u, lu, b, v, lv, ans );
! 563: }
! 564:
! 565: static int abs_axU_V_maxrshift( a, u, lu, v, lv, ans )
! 566: int a, *u, lu, *v, lv, *ans;
! 567: /*
! 568: * Compute a*| U - V | >> (as much as possible) in ans[], where U=u[0:lu-1]
! 569: * and V=v[0:lv-1], 0 < a < 2^BASE.
! 570: */
! 571: {
! 572: int i, l, b, c, t, rsh, *p = ans, lsh, w, x;
! 573: int hi;
! 574:
! 575: if ( lu > lv ) goto U_V;
! 576: else if ( lu < lv ) goto swap;
! 577: else {
! 578: for ( i = lu; i-- > 0; ) {
! 579: if ( u[i] == v[i] ) continue;
! 580: lu = lv = i + 1;
! 581: if ( u[i] > v[i] ) goto U_V;
! 582: else goto swap;
! 583: }
! 584: return 0;
! 585: }
! 586: swap:
! 587: SWAP( lu, lv, int ); SWAP( u, v, int * );
! 588: U_V:
! 589: for ( b = c = i = 0; i < lv; ) {
! 590: if ( (w = u[i] - v[i] - b) < 0 ) w += BASE, b = 1;
! 591: else b = 0;
! 592: DMA( a, w, c, hi, t, not_used );
! 593: i++, c = hi;
! 594: if ( t != 0 ) goto non0w;
! 595: }
! 596: while ( i < lu ) {
! 597: if ( (w = u[i] - b) < 0 ) w += BASE, b = 1;
! 598: else b = 0;
! 599: DMA( a, w, c, hi, t, not_used );
! 600: i++, c = hi;
! 601: if ( t != 0 ) goto non0w;
! 602: }
! 603: if ( b != 0 ) c -= a;
! 604: else if ( c == 0 ) return 0;
! 605: while ( (c&1) == 0 ) c >>= 1;
! 606: *p = c;
! 607: return 1;
! 608: /**/
! 609: non0w:
! 610: u += i, lu -= i, v += i, lv -= i;
! 611: TRAILINGZEROS( t, rsh );
! 612: i = 0;
! 613: if ( rsh == 0 ) {
! 614: *p++ = t, l = 0;
! 615: for ( ; i < lv; c = hi ) {
! 616: if ( (w = u[i] - v[i] - b) < 0 ) w += BASE, b = 1;
! 617: else b = 0;
! 618: DMA( a, w, c, hi, t, not_used );
! 619: i++;
! 620: if ( (*p++ = t) != 0 ) l = i;
! 621: }
! 622: for ( ; i < lu; c = hi ) {
! 623: if ( (w = u[i] - b) < 0 ) w += BASE, b = 1;
! 624: else b = 1;
! 625: DMA( a, w, c, hi, t, not_used );
! 626: i++;
! 627: if ( (*p++ = t) != 0 ) l = i;
! 628: }
! 629: if ( b != 0 ) c -= a;
! 630: if ( c == 0 ) return( l + 1 );
! 631: *p = c;
! 632: return( i + 2 );
! 633: } else {
! 634: for ( lsh = BSH - rsh, l = 0; i < lv; c = hi, t = x >> rsh ) {
! 635: if ( (w = u[i] - v[i] - b) < 0 ) w += BASE, b = 1;
! 636: else b = 0;
! 637: DMA( a, w, c, hi, x, not_used );
! 638: i++;
! 639: if ( (*p++ = t | ((x << lsh) & BMASK)) != 0 ) l = i;
! 640: }
! 641: for ( ; i < lu; i++, c = hi, t = x >> rsh ) {
! 642: if ( (w = u[i] - b) < 0 ) w += BASE, b = 1;
! 643: else b = 0;
! 644: DMA( a, w, c, hi, x, not_used );
! 645: i++;
! 646: if ( (*p++ = t | ((x << lsh) & BMASK)) != 0 ) l = i;
! 647: }
! 648: if ( b != 0 ) c -= a;
! 649: if ( c != 0 ) {
! 650: *p++ = t | ((c << lsh) & BMASK);
! 651: i++, t = c >> rsh;
! 652: } else if ( t == 0 ) return l;
! 653: if ( t == 0 ) return( i );
! 654: *p = t;
! 655: return( i + 1 );
! 656: }
! 657: }
! 658: #endif /* FULLSET */
! 659:
! 660: static int abs_axU_bxV_maxrshift( a, u, lu, b, v, lv, ans )
! 661: int a, *u, lu, b, *v, lv, *ans;
! 662: /*
! 663: * Compute | a*U - b*V | >> (as much as possible) in ans[], where U=u[0:lu-1]
! 664: * and V=v[0:lv-1], 0 < a, b < 2^BASE.
! 665: */
! 666: {
! 667: int i, l, c, d, s, t, w, rsh, *p = ans, lsh;
! 668: int hi;
! 669: static int bw_int32();
! 670:
! 671: BitWidth( a, c ); BitWidth( u[lu -1], s );
! 672: BitWidth( b, d ); BitWidth( v[lv -1], t );
! 673: if ( lu < lv || lu == lv && (c + s) < (d + t) ) {
! 674: SWAP( a, b, int ); SWAP( lu, lv, int ); SWAP( u, v, int * );
! 675: }
! 676: for ( i = c = d = 0; i < lv; ) {
! 677: DMA( a, u[i], c, hi, t, not_used ); /* 0 <= c <= 2^BSH -2 */
! 678: c = hi;
! 679: DMA( b, v[i], d, hi, s, not_used ); /* 0 <= d <= 2^BSH -1 */
! 680: i++, d = hi;
! 681: if ( (t -= s) == 0 ) continue;
! 682: else if ( t < 0 ) {
! 683: t += BASE;
! 684: if ( c != 0 ) c--;
! 685: else d++;
! 686: }
! 687: goto non0w;
! 688: }
! 689: if ( i >= lu ) { /* no more elm.'s in u[] and v[] */
! 690: only_c:
! 691: if ( (c -= d) == 0 ) return 0;
! 692: else if ( c < 0 ) c = -c;
! 693: sh_onlyc:
! 694: /* TRAILINGZEROS( c, rsh ); */
! 695: while ( (c&1) == 0 ) c >>= 1;
! 696: *p = c;
! 697: return 1;
! 698: }
! 699: /* there is at least one more elm. in u[] */
! 700: DMA( a, u[i], c, hi, t, not_used );
! 701: i++, c = hi;
! 702: if ( i >= lu ) { /* in this case, diff still can be <= 0 */
! 703: if ( hi == 0 ) { c = t; goto only_c; }
! 704: if ( (t -= d) == 0 ) { c = hi; goto sh_onlyc; }
! 705: else if ( t < 0 ) t += BASE, hi--;
! 706: TRAILINGZEROS( t, rsh );
! 707: if ( rsh == 0 ) *p++ = t;
! 708: else {
! 709: *p++ = t | ((hi << (BSH - rsh)) & BMASK);
! 710: hi >>= rsh;
! 711: }
! 712: if ( hi == 0 ) return 1;
! 713: *p = hi;
! 714: return 2;
! 715: } /* diff turned out to be > 0 */ else if ( (t -= d) > 0 ) d = 0;
! 716: else if ( t < 0 ) t += BASE, d = 1;
! 717: else {
! 718: while ( i < lu ) {
! 719: DMA( a, u[i], c, hi, t, not_used );
! 720: i++, c = hi;
! 721: if ( t != 0 ) break;
! 722: }
! 723: }
! 724: u += i, lu -= i;
! 725: TRAILINGZEROS( t, rsh );
! 726: l = i = 0;
! 727: if ( rsh == 0 ) {
! 728: *p++ = t;
! 729: goto no_more_v;
! 730: } else goto no_more_v_sh;
! 731: /**/
! 732: non0w:
! 733: u += i, lu -= i, v += i, lv -= i;
! 734: TRAILINGZEROS( t, rsh );
! 735: i = 0;
! 736: if ( rsh == 0 ) {
! 737: *p++ = t, l = 0;
! 738: for ( ; i < lv; *p++ = t ) {
! 739: DMA( a, u[i], c, hi, t, not_used );
! 740: c = hi;
! 741: DMA( b, v[i], d, hi, s, not_used );
! 742: i++, d = hi;
! 743: if ( (t -= s) == 0 ) continue;
! 744: else if ( t < 0 ) {
! 745: t += BASE;
! 746: if ( c != 0 ) c--;
! 747: else d++;
! 748: }
! 749: l = i;
! 750: }
! 751: no_more_v:
! 752: if ( i >= lu ) {
! 753: final_c_d:
! 754: if ( (c -= d) == 0 ) return( l + 1 );
! 755: else if ( c < 0 ) {
! 756: l = negate_nbd( ans, i + 1 ); /* <================ */
! 757: if ( (c = -c - 1) == 0 ) return l;
! 758: }
! 759: *p = c;
! 760: return( i + 2 );
! 761: }
! 762: /* a*u[i:lu-1] + c - d */
! 763: if ( c >= d ) c -= d;
! 764: else {
! 765: DMA( a, u[i], c, hi, t, not_used );
! 766: if ( i+1 >= lu ) {
! 767: if ( hi == 0 ) { c = t; goto final_c_d; }
! 768: if ( (t -= d) < 0 ) t += BASE, c = hi - 1;
! 769: else c = hi;
! 770: i++, *p++ = t;
! 771: goto final_c;
! 772: }
! 773: i++, c = hi;
! 774: if ( (t -= d) >= 0 ) *p++ = t;
! 775: else {
! 776: *p++ = t + BASE;
! 777: if ( c != 0 ) c--;
! 778: else {
! 779: for ( ; i < lu; *p++ = BMASK ) {
! 780: DMA( a, u[i], c, hi, t, not_used );
! 781: i++, c = hi;
! 782: if ( t == 0 ) continue;
! 783: *p++ = t - 1;
! 784: goto aU;
! 785: }
! 786: c--;
! 787: goto final_c;
! 788: }
! 789: }
! 790: }
! 791: /*** return( aV_plus_c( a, &u[i], lu - i, c, p ) + i + 1 ); ***/
! 792: aU:
! 793: for ( ; i < lu; i++, c = hi, *p++ = t ) { DMA( a, u[i], c, hi, t, not_used ); }
! 794: final_c:
! 795: if ( c == 0 ) return( i + 1 );
! 796: *p = c;
! 797: return( i + 2 );
! 798: } else { /* rsh != 0 */
! 799: for ( lsh = BSH - rsh; i < lv; t = w >> rsh ) {
! 800: DMA( a, u[i], c, hi, w, not_used );
! 801: c = hi;
! 802: DMA( b, v[i], d, hi, s, not_used );
! 803: i++, d = hi;
! 804: if ( (w -= s) < 0 ) {
! 805: w += BASE;
! 806: if ( c != 0 ) c--;
! 807: else d++;
! 808: }
! 809: if ( (*p++ = t | ((w << lsh) & BMASK)) != 0 ) l = i;
! 810: }
! 811: no_more_v_sh:
! 812: if ( i >= lu ) {
! 813: final_c_d_sh:
! 814: if ( (c -= d) == 0 )
! 815: if ( t == 0 ) return l;
! 816: else { *p = t; return( i + 1 ); }
! 817: else if ( c < 0 ) {
! 818: if ( i == 0 ) t = (BASE >> rsh) - t;
! 819: else {
! 820: l = negate_nbd( ans, i ); /* <================ */
! 821: /* t = (BASE >> rsh) - t;*/
! 822: t ^= ((BMASK >> rsh));
! 823: }
! 824: c = -c -1;
! 825: }
! 826: *p++ = t | ((c << lsh) & BMASK);
! 827: if ( (c >>= rsh) == 0 ) return( i + 1 );
! 828: *p = c;
! 829: return( i + 2 );
! 830: } else if ( c >= d ) c -= d;
! 831: else {
! 832: DMA( a, u[i], c, hi, w, not_used );
! 833: if ( i+1 >= lu ) {
! 834: if ( hi == 0 ) { c = w; goto final_c_d_sh; }
! 835: if ( (w -= d) < 0 ) w += BASE, c = hi - 1;
! 836: else c = hi;
! 837: i++, *p++ = t | ((w << lsh) & BMASK);
! 838: t = w >> rsh;
! 839: goto final_c_sh;
! 840: }
! 841: i++, c = hi;
! 842: if ( (w -= d) >= 0 ) { *p++ = t | ((w << lsh) & BMASK); t = w >> rsh; }
! 843: else {
! 844: w += BASE;
! 845: *p++ = t | ((w << lsh) & BMASK);
! 846: t = w >> rsh;
! 847: if ( c != 0 ) c--;
! 848: else {
! 849: while ( i < lu ) {
! 850: DMA( a, u[i], c, hi, w, not_used );
! 851: i++, c = hi;
! 852: if ( w == 0 ) {
! 853: *p++ = t | ((BMASK << lsh) & BMASK);
! 854: t = BMASK >> rsh;
! 855: continue;
! 856: } else {
! 857: w--;
! 858: *p++ = t | ((w << lsh) & BMASK);
! 859: t = w >> rsh;
! 860: goto aU_sh;
! 861: }
! 862: }
! 863: c--;
! 864: goto final_c_sh;
! 865: }
! 866: }
! 867: }
! 868: /*** return( aV_plus_c_sh( a, &u[i], lu - i, c, rsh, t ) + i ); ***/
! 869: aU_sh:
! 870: for ( ; i < lu; i++, c = hi, t = w >> rsh ) {
! 871: DMA( a, u[i], c, hi, w, not_used );
! 872: *p++ = t | ((w << lsh) & BMASK);
! 873: }
! 874: final_c_sh:
! 875: if ( c != 0 ) {
! 876: *p++ = t | ((c << lsh) & BMASK);
! 877: i++, t = c >> rsh;
! 878: }
! 879: if ( t == 0 ) return i;
! 880: *p = t;
! 881: return( i + 1 );
! 882: }
! 883: }
! 884: /****************/
! 885: /****************/
! 886: #ifdef FULLSET
! 887: static int UplusV_maxrshift(), aUplusV_maxrshift(), axUplusV_maxrshift(), aUplusbV_maxrshift();
! 888:
! 889: static int aU_plus_bV_maxrshift( a, u, lu, b, v, lv, ans )
! 890: int a, *u, lu, b, *v, lv, *ans;
! 891: {
! 892: if ( a == 1 )
! 893: return( b == 1 ? UplusV_maxrshift( u, lu, v, lv, ans ) :
! 894: aUplusV_maxrshift( b, v, lv, u, lu, ans ) );
! 895: else if ( b == 1 ) return aUplusV_maxrshift( a, u, lu, v, lv, ans );
! 896: else if ( a == b ) return axUplusV_maxrshift( a, u, lu, v, lv, ans );
! 897: return aUplusbV_maxrshift( a, u, lu, b, v, lv, ans );
! 898: }
! 899:
! 900: static int UplusV_maxrshift( u, lu, v, lv, ans )
! 901: int *u, lu, *v, lv, *ans;
! 902: /*
! 903: * Compute ( (U + V) >> (as much as possible) ) in ans[], where U=u[0:lu-1],
! 904: * V=v[0:lv-1], and 0 < a < BASE.
! 905: */
! 906: {
! 907: int i, t, c, rsh, *p = ans;
! 908: int hi;
! 909:
! 910: if ( lu < lv ) { SWAP( lu, lv, int ); SWAP( u, v, int * ); }
! 911: for ( c = i = 0; i < lv; ) {
! 912: t = (c += (u[i] + v[i])) & BMASK;
! 913: /* 0 <= c <= 2*(2^BSH -1) + 1 = 2^BSH + (2^BSH -1) */
! 914: i++, c >>= BSH;
! 915: if ( t != 0 ) goto non0w;
! 916: }
! 917: while ( i < lu ) {
! 918: t = (c += u[i]) & BMASK;
! 919: i++, c >>= BSH;
! 920: if ( t != 0 ) goto non0w;
! 921: }
! 922: *p = c;
! 923: return 1;
! 924: /**/
! 925: non0w:
! 926: u += i, lu -= i, v += i, lv -= i;
! 927: TRAILINGZEROS( t, rsh );
! 928: i = 0;
! 929: if ( rsh == 0 ) {
! 930: *p++ = t;
! 931: for ( ; i < lv; i++, c >>= BSH ) *p++ = (c += (u[i] + v[i])) & BMASK;
! 932: if ( c != 0 ) {
! 933: while ( i < lu )
! 934: if ( (c = u[i++] + 1) >= BASE ) *p++ = c & BMASK;
! 935: else {
! 936: *p++ = c;
! 937: goto cpu;
! 938: }
! 939: *p = 1; /* == c */
! 940: return( lu + 2 );
! 941: }
! 942: cpu:
! 943: for ( ; i < lu; i++ ) *p++ = u[i];
! 944: return( lu + 1 );
! 945: } else {
! 946: int lsh = BSH - rsh;
! 947:
! 948: for ( ; i < lv; i++, c >>= BSH ) {
! 949: c += (u[i] + v[i]);
! 950: *p++ = t | ((c << lsh) & BMASK);
! 951: t = (c & BMASK) >> rsh;
! 952: }
! 953: if ( c != 0 ) {
! 954: while ( i < lu ) {
! 955: c = u[i++] + 1;
! 956: *p++ = t | ((c << lsh) & BMASK);
! 957: t = (c & BMASK) >> rsh;
! 958: if ( (c >>= BSH) == 0 ) goto cpu_sh;
! 959: }
! 960: *p = t | (1 << lsh);
! 961: return( lu + 1 );
! 962: }
! 963: cpu_sh:
! 964: for ( ; i < lu; t = c >> rsh ) *p++ = t | (((c = u[i++]) << lsh) & BMASK);
! 965: if ( t == 0 ) return( lu );
! 966: *p = t;
! 967: return( lu + 1 );
! 968: }
! 969: }
! 970:
! 971: static int aUplusV_maxrshift( a, u, lu, v, lv, ans )
! 972: int a, *u, lu, *v, lv, *ans;
! 973: /*
! 974: * Compute ( (a*U + V) >> (as much as possible) ) in ans[], where U=u[0:lu-1],
! 975: * V=v[0:lv-1], and 1 < a < BASE.
! 976: */
! 977: {
! 978: int i, l, t, c, rsh, *p = ans, w, lsh;
! 979: int hi;
! 980:
! 981: for ( l = lu < lv ? lu : lv, c = i = 0; i < l; ) {
! 982: /* Note that (carry in a*U[*]) <= 2^BSH-2 because {(2^j-1)}^2 + (2^j-2)
! 983: * = 2^j*(2^j -2)+(2^j -1), and the following c is <= 2^j -1 because
! 984: * {(2^j -1)}^2 + (2^j -1) + (2^j -1) = 2^j*(2^j -1) + (2^j -1).
! 985: * ^^^^^^carry^^^^^^ pv[?].
! 986: */
! 987: c += v[i];
! 988: DMA( a, u[i], c, hi, t, not_used );
! 989: i++, c = hi;
! 990: if ( t != 0 ) goto non0sum;
! 991: }
! 992: if ( i >= lu ) { /* compute (c + v[i:lv-1]) */
! 993: while ( i < lv ) {
! 994: t = (c += v[i++]) & BMASK;
! 995: c >>= BSH;
! 996: if ( t != 0 ) goto non0w_v;
! 997: }
! 998: *p = 1;
! 999: return 1;
! 1000: non0w_v:
! 1001: v += i, lv -= i;
! 1002: i = 0;
! 1003: TRAILINGZEROS( t, rsh );
! 1004: if ( rsh == 0 ) {
! 1005: *p++ = t;
! 1006: L_addv:
! 1007: if ( c != 0 ) {
! 1008: while ( i < lv )
! 1009: if ( (c = u[i++] + 1) >= BASE ) *p++ = c & BMASK;
! 1010: else {
! 1011: *p++ = c;
! 1012: goto cpv;
! 1013: }
! 1014: *p = 1; /* == c */
! 1015: return( lv + 2 );
! 1016: }
! 1017: cpv:
! 1018: for ( ; i < lv; i++ ) *p++ = v[i];
! 1019: return( lv + 1 );
! 1020: } else {
! 1021: lsh = BSH - rsh;
! 1022: L_addv_sh:
! 1023: if ( c != 0 ) {
! 1024: while ( i < lv ) {
! 1025: c += v[i++];
! 1026: *p++ = t | ((c << lsh) & BMASK);
! 1027: t = (c & BMASK) >> rsh;
! 1028: if ( (c >>= BSH) == 0 ) goto cpv_sh;
! 1029: }
! 1030: *p = t | (1 << lsh);
! 1031: return( lv + 1 );
! 1032: }
! 1033: cpv_sh:
! 1034: for ( ; i < lv; t = c >> rsh ) *p++ = t | (((c = v[i++]) << lsh) & BMASK);
! 1035: if ( t == 0 ) return( lv );
! 1036: *p = t;
! 1037: return( lv + 1 );
! 1038: }
! 1039: } else { /* i >= lv, and compute (c + a*u[i:lu-1]) */
! 1040: while ( i < lu ) {
! 1041: DMA( a, u[i], c, hi, t, not_used );
! 1042: i++, c = hi;
! 1043: if ( t != 0 ) goto non0w_u;
! 1044: }
! 1045: TRAILINGZEROS( c, rsh );
! 1046: *p = c;
! 1047: return 1;
! 1048: /**/
! 1049: non0w_u:
! 1050: u += i, lu -= i;
! 1051: TRAILINGZEROS( t, rsh );
! 1052: i = 0;
! 1053: if ( rsh == 0 ) {
! 1054: *p++ = t;
! 1055: L_addu:
! 1056: #ifdef CALL
! 1057: return( aV_plus_c( a, u, lv, c, p ) + 1 );
! 1058: #else
! 1059: for ( ; i < lu; i++, *p++ = t, c = hi ) { DMA( a, u[i], c, hi, t, not_used ); }
! 1060: if ( c == 0 ) return( lu + 1 );
! 1061: *p = c;
! 1062: return( lu + 2 );
! 1063: #endif
! 1064: } else {
! 1065: lsh = BSH - rsh;
! 1066: L_addu_sh:
! 1067: #ifdef CALL
! 1068: return aV_plus_c_sh( a, u, lu, c, rsh, t, p );
! 1069: #else
! 1070: for ( ; i < lu; i++, c = hi, t = w >> rsh ) {
! 1071: DMA( a, u[i], c, hi, w, not_used );
! 1072: *p++ = t | ((w << lsh) & BMASK);
! 1073: }
! 1074: if ( c != 0 ) {
! 1075: *p++ = t | ((c << lsh) & BMASK);
! 1076: i++, t = c >> rsh;
! 1077: }
! 1078: if ( t == 0 ) return i;
! 1079: *p = t;
! 1080: return( i + 1 );
! 1081: #endif
! 1082: }
! 1083: }
! 1084: /**/
! 1085: non0sum:
! 1086: u += i, lu -= i, v += i, lv -= i, l -= i;
! 1087: TRAILINGZEROS( t, rsh );
! 1088: i = 0;
! 1089: if ( rsh == 0 ) {
! 1090: *p++ = t;
! 1091: for ( ; i < l; i++, *p++ = t, c = hi ) {
! 1092: c += v[i];
! 1093: DMA( a, u[i], c, hi, t, not_used );
! 1094: }
! 1095: if ( i >= lu ) /* compute (c + v[i:lv-1]). note i may be >= lv */ goto L_addv;
! 1096: else /* i >= lv, and compute (c + u[i:lu-1]) */ goto L_addu;
! 1097: } else {
! 1098: lsh = BSH - rsh;
! 1099: for ( ; i < l; i++, c = hi, t = w >> rsh ) {
! 1100: c += v[i];
! 1101: DMA( a, u[i], c, hi, w, not_used );
! 1102: *p++ = t | ((w << lsh) & BMASK);
! 1103: }
! 1104: if ( i >= lu ) /* compute (c + v[i:lv-1]) >> rsh. note i may be >= lv */ goto L_addv_sh;
! 1105: else /* i >= lv, and compute (c + u[i:lu-1]) >> rsh */ goto L_addu_sh;
! 1106: }
! 1107: }
! 1108:
! 1109:
! 1110: static int axUplusV_maxrshift( a, u, lu, v, lv, ans )
! 1111: int a, *u, lu, *v, lv, *ans;
! 1112: /*
! 1113: * Compute ( a*(U + V) >> (as much as possible) ) in ans[], where U=u[0:lu-1],
! 1114: * V=v[0:lv-1], and 1 < a < BASE.
! 1115: */
! 1116: {
! 1117: int i, t, c, d, rsh, w, *p = ans, lsh, x;
! 1118: int hi;
! 1119:
! 1120: if ( lu < lv ) { SWAP( lu, lv, int ); SWAP( u, v, int * ); }
! 1121: for ( i = c = d = 0; i < lv; ) {
! 1122: w = (d += (u[i] + v[i])) & BMASK;
! 1123: DMA( a, w, c, hi, t, not_used );
! 1124: i++, d >>= BSH, c = hi;
! 1125: if ( t != 0 ) goto non0w;
! 1126: }
! 1127: while ( i < lu ) {
! 1128: w = (d += u[i++]) & BMASK;
! 1129: DMA( a, w, c, hi, t, not_used );
! 1130: d >>= BSH, c = hi;
! 1131: if ( t != 0 ) goto non0w;
! 1132: }
! 1133: if ( d != 0 ) c += a;
! 1134: TRAILINGZEROS( c, rsh );
! 1135: if ( rsh != 0 || c <= BMASK ) { *p = c; return 1; }
! 1136: *p++ = c & BMASK;
! 1137: *p = c >> BSH;
! 1138: return 2;
! 1139: /**/
! 1140: non0w:
! 1141: u += i, lu -= i, v += i, lv -= i;
! 1142: TRAILINGZEROS( t, rsh );
! 1143: i = 0;
! 1144: if ( rsh == 0 ) {
! 1145: *p++ = t;
! 1146: for ( ; i < lv; i++, d >>= BSH, *p++ = t, c = hi ) {
! 1147: w = (d += (u[i] + v[i])) & BMASK;
! 1148: DMA( a, w, c, hi, t, not_used );
! 1149: }
! 1150: for ( ; i < lu; d >>= BSH, *p++ = t, c = hi ) {
! 1151: w = (d += u[i++]) & BMASK;
! 1152: DMA( a, w, c, hi, t, not_used );
! 1153: }
! 1154: if ( d != 0 ) c += a;
! 1155: if ( c == 0 ) return( lu + 1 );
! 1156: *p++ = c & BMASK;
! 1157: if ( (c >>= BSH) == 0 ) return( lu + 2 );
! 1158: *p = c;
! 1159: return( lu + 3 );
! 1160: } else {
! 1161: for ( lsh = BSH - rsh; i < lv; i++, d >>= BSH, t = x >> rsh, c = hi ) {
! 1162: w = (d += (u[i] + v[i])) & BMASK;
! 1163: DMA( a, w, c, hi, x, not_used );
! 1164: *p++ = t | ((x << lsh) & BMASK);
! 1165: }
! 1166: for ( ; i < lu; d >>= BSH, t = x >> rsh, c = hi ) {
! 1167: w = (d += u[i++]) & BMASK;
! 1168: DMA( a, w, c, hi, x, not_used );
! 1169: *p++ = t | ((x << lsh) & BMASK);
! 1170: }
! 1171: if ( d != 0 ) c += a;
! 1172: t |= ((c << lsh) & BMASK);
! 1173: c >>= rsh;
! 1174: if ( c != 0 ) {
! 1175: *p++ = t;
! 1176: *p = c;
! 1177: return( lu + 2 );
! 1178: } else if ( t == 0 ) return lu;
! 1179: *p = t;
! 1180: return( lu + 1 );
! 1181: }
! 1182: }
! 1183: #endif /* FULLSET */
! 1184:
! 1185: static int aUplusbV_maxrshift( a, u, lu, b, v, lv, ans )
! 1186: int a, *u, lu, b, *v, lv, *ans;
! 1187: /*
! 1188: * Compute ( (a*U + b*V) >> (as much as possible) ) in ans[], where U=u[0:lu-1],
! 1189: * V=v[0:lv-1], and 0 < a,b < BASE.
! 1190: */
! 1191: {
! 1192: int i, c, d, s, t, w, rsh, *p = ans, lsh;
! 1193: int hi;
! 1194:
! 1195: if ( lu < lv ) { SWAP( a, b, int ); SWAP( lu, lv, int ); SWAP( u, v, int * ); }
! 1196: for ( c = d = i = 0; i < lv; ) {
! 1197: DMA( a, u[i], c, hi, t, not_used ); /* 0 <= c <= 2^BSH -1 */
! 1198: c = hi;
! 1199: DMA( b, v[i], d, hi, s, not_used ); /* 0 <= d <= 2^BSH -2 */
! 1200: i++, d = hi, t += s;
! 1201: if ( t > BMASK ) c++, t &= BMASK;
! 1202: if ( t != 0 ) goto non0w;
! 1203: }
! 1204: c += d;
! 1205: for ( d = 0; i < lu; ) {
! 1206: DMA( a, u[i], c, hi, t, not_used );
! 1207: i++, c = hi;
! 1208: if ( t != 0 ) goto non0w;
! 1209: }
! 1210: TRAILINGZEROS( c, rsh );
! 1211: if ( rsh != 0 || c <= BMASK ) { *p = c; return 1; }
! 1212: *p++ = c & BMASK;
! 1213: *p = 1;
! 1214: return 2;
! 1215: /**/
! 1216: non0w:
! 1217: u += i, lu -= i, v += i, lv -= i;
! 1218: TRAILINGZEROS( t, rsh );
! 1219: i = 0;
! 1220: if ( rsh == 0 ) {
! 1221: *p++ = t;
! 1222: for ( ; i < lv; i++, *p++ = t ) {
! 1223: DMA( a, u[i], c, hi, t, not_used );
! 1224: c = hi;
! 1225: DMA( b, v[i], d, hi, s, not_used );
! 1226: d = hi, t += s;
! 1227: if ( t > BMASK ) c++, t &= BMASK;
! 1228: }
! 1229: c += d;
! 1230: for ( ; i < lu; i++, *p++ = t, c = hi ) {
! 1231: DMA( a, u[i], c, hi, t, not_used );
! 1232: }
! 1233: if ( c == 0 ) return( lu + 1 );
! 1234: else if ( c <= BMASK ) { *p = c; return( lu + 2 ); }
! 1235: *p++ = c & BMASK;
! 1236: *p = 1;
! 1237: return( lu + 3 );
! 1238: } else {
! 1239: for ( lsh = BSH - rsh; i < lv; i++, t = w >> rsh ) {
! 1240: DMA( a, u[i], c, hi, w, not_used );
! 1241: c = hi;
! 1242: DMA( b, v[i], d, hi, s, not_used );
! 1243: d = hi, w += s;
! 1244: if ( w > BMASK ) c++, w &= BMASK;
! 1245: *p++ = t | ((w << lsh) & BMASK);
! 1246: }
! 1247: c += d; /* <= 2^BSH + (2^BSH - 3) */
! 1248: for ( ; i < lu; i++, c = hi, t = w >> rsh ) {
! 1249: DMA( a, u[i], c, hi, w, not_used );
! 1250: *p++ = t | ((w << lsh) & BMASK);
! 1251: }
! 1252: if ( c == 0 ) {
! 1253: if ( t == 0 ) return lu;
! 1254: *p = t;
! 1255: return( lu + 1 );
! 1256: }
! 1257: *p++ = t | ((c << lsh) & BMASK);
! 1258: if ( (c >>= rsh) == 0 ) return( lu + 1 );
! 1259: *p = c;
! 1260: return( lu + 2 );
! 1261: }
! 1262: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to igcdhack.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / engine