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

Annotation of OpenXM_contrib2/asir2000/engine/igcdhack.c, Revision 1.6

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

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>