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