[BACK]Return to int.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2018 / builtin

Annotation of OpenXM_contrib2/asir2018/builtin/int.c, Revision 1.4

1.1       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@sec.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:  *
1.4     ! noro       48:  * $OpenXM: OpenXM_contrib2/asir2018/builtin/int.c,v 1.3 2018/09/25 07:36:01 noro Exp $
1.1       noro       49: */
                     50: #include "ca.h"
                     51: #include "parse.h"
                     52: #include "base.h"
                     53:
                     54: void Pidiv(), Pirem(), Pigcd(), Pilcm(), Pfac(), Prandom(), Pinv();
                     55: void Pup2_inv(),Pgf2nton(), Pntogf2n();
                     56: void Pup2_init_eg(), Pup2_show_eg();
                     57: void Piqr(), Pprime(), Plprime(), Pinttorat();
                     58: void Piand(), Pior(), Pixor(), Pishift();
                     59: void Pisqrt();
                     60: void Plrandom();
                     61: void Pset_upkara(), Pset_uptkara(), Pset_up2kara(), Pset_upfft();
                     62: void Pmt_save(), Pmt_load();
                     63: void Psmall_jacobi();
                     64: void Pdp_set_mpi();
                     65: void Pntoint32(),Pint32ton();
                     66:
                     67: void Pigcdbin(), Pigcdbmod(), PigcdEuc(), Pigcdacc(), Pigcdcntl();
                     68:
                     69: void Pihex();
                     70: void Pimaxrsh(), Pilen();
                     71: void Ptype_t_NB();
                     72:
                     73: struct ftab int_tab[] = {
                     74:   {"dp_set_mpi",Pdp_set_mpi,-1},
                     75:   {"isqrt",Pisqrt,1},
                     76:   {"idiv",Pidiv,2},
                     77:   {"irem",Pirem,2},
                     78:   {"iqr",Piqr,2},
                     79:   {"igcd",Pigcd,-2},
                     80:   {"ilcm",Pilcm,2},
                     81:   {"up2_inv",Pup2_inv,2},
                     82:   {"up2_init_eg",Pup2_init_eg,0},
                     83:   {"up2_show_eg",Pup2_show_eg,0},
                     84:   {"type_t_NB",Ptype_t_NB,2},
                     85:   {"gf2nton",Pgf2nton,1},
                     86:   {"ntogf2n",Pntogf2n,1},
                     87:   {"set_upkara",Pset_upkara,-1},
                     88:   {"set_uptkara",Pset_uptkara,-1},
                     89:   {"set_up2kara",Pset_up2kara,-1},
                     90:   {"set_upfft",Pset_upfft,-1},
                     91:   {"inv",Pinv,2},
                     92:   {"inttorat",Pinttorat,3},
                     93:   {"fac",Pfac,1},
                     94:   {"prime",Pprime,1},
                     95:   {"lprime",Plprime,1},
                     96:   {"random",Prandom,-1},
                     97:   {"lrandom",Plrandom,1},
                     98:   {"iand",Piand,2},
                     99:   {"ior",Pior,2},
                    100:   {"ixor",Pixor,2},
                    101:   {"ishift",Pishift,2},
                    102:   {"small_jacobi",Psmall_jacobi,2},
                    103:
                    104:   {"igcdbin",Pigcdbin,2},    /* HM@CCUT extension */
                    105:   {"igcdbmod",Pigcdbmod,2},  /* HM@CCUT extension */
                    106:   {"igcdeuc",PigcdEuc,2},    /* HM@CCUT extension */
                    107:   {"igcdacc",Pigcdacc,2},    /* HM@CCUT extension */
                    108:   {"igcdcntl",Pigcdcntl,-1},  /* HM@CCUT extension */
                    109:
                    110:   {"mt_save",Pmt_save,1},
                    111:   {"mt_load",Pmt_load,1},
                    112:   {"ntoint32",Pntoint32,1},
                    113:   {"int32ton",Pint32ton,1},
                    114:   {0,0,0},
                    115: };
                    116:
                    117: static int is_prime_small(unsigned int);
                    118: static unsigned int gcd_small(unsigned int,unsigned int);
                    119: int TypeT_NB_check(unsigned int, unsigned int);
                    120: int mpi_mag;
                    121:
                    122: void Pntoint32(NODE arg,USINT *rp)
                    123: {
                    124:   Z q,z;
                    125:   unsigned int t;
                    126:
                    127:   asir_assert(ARG0(arg),O_N,"ntoint32");
                    128:   q = (Z)ARG0(arg);
                    129:   if ( !q ) {
                    130:     MKUSINT(*rp,0);
                    131:     return;
                    132:   }
                    133:   if ( !INT(q) || !smallz(q) )
                    134:     error("ntoint32 : invalid argument");
                    135:   absz(q,&z);
1.4     ! noro      136:   t = ZTOS(z);
1.1       noro      137:   if ( sgnz(q) < 0 )
                    138:     t = -(int)t;
                    139:   MKUSINT(*rp,t);
                    140: }
                    141:
                    142: void Pint32ton(NODE arg,Z *rp)
                    143: {
                    144:   int t;
                    145:
                    146:   asir_assert(ARG0(arg),O_USINT,"int32ton");
                    147:   t = (int)BDY((USINT)ARG0(arg));
1.4     ! noro      148:   STOZ(t,*rp);
1.1       noro      149: }
                    150:
                    151: void Pdp_set_mpi(NODE arg,Z *rp)
                    152: {
                    153:   if ( arg ) {
                    154:     asir_assert(ARG0(arg),O_N,"dp_set_mpi");
1.4     ! noro      155:     mpi_mag = ZTOS((Q)ARG0(arg));
1.1       noro      156:   }
1.4     ! noro      157:   STOZ(mpi_mag,*rp);
1.1       noro      158: }
                    159:
                    160: void Psmall_jacobi(NODE arg,Z *rp)
                    161: {
                    162:   Z a,m;
                    163:   int a0,m0,s;
                    164:
                    165:   a = (Z)ARG0(arg);
                    166:   m = (Z)ARG1(arg);
                    167:   asir_assert(a,O_N,"small_jacobi");
                    168:   asir_assert(m,O_N,"small_jacobi");
                    169:   if ( !a )
                    170:      *rp = ONE;
                    171:   else if ( !m || !INT(m) || !INT(a)
                    172:     || !smallz(m) || !smallz(a) || sgnz(m) < 0 || evenz(m) )
                    173:     error("small_jacobi : invalid input");
                    174:   else {
1.4     ! noro      175:     a0 = ZTOS(a); m0 = ZTOS(m);
1.1       noro      176:     s = small_jacobi(a0,m0);
1.4     ! noro      177:     STOZ(s,*rp);
1.1       noro      178:   }
                    179: }
                    180:
                    181: int small_jacobi(int a,int m)
                    182: {
                    183:   int m4,m8,a4,j1,i,s;
                    184:
                    185:   a %= m;
                    186:   if ( a == 0 || a == 1 )
                    187:     return 1;
                    188:   else if ( a < 0 ) {
                    189:     j1 = small_jacobi(-a,m);
                    190:     m4 = m%4;
                    191:     return m4==1?j1:-j1;
                    192:   } else {
                    193:     for ( i = 0; a && !(a&1); i++, a >>= 1 );
                    194:     if ( i&1 ) {
                    195:       m8 = m%8;
                    196:       s = (m8==1||m8==7) ? 1 : -1;
                    197:     } else
                    198:       s = 1;
                    199:     /* a, m are odd */
                    200:     j1 = small_jacobi(m%a,a);
                    201:     m4 = m%4; a4 = a%4;
                    202:     s *= (m4==1||a4==1) ? 1 : -1;
                    203:     return j1*s;
                    204:   }
                    205: }
                    206:
                    207: void Ptype_t_NB(NODE arg,Z *rp)
                    208: {
1.4     ! noro      209:   if ( TypeT_NB_check(ZTOS((Q)ARG0(arg)),ZTOS((Q)ARG1(arg))))
1.1       noro      210:     *rp = ONE;
                    211:   else
                    212:     *rp = 0;
                    213: }
                    214:
                    215: int TypeT_NB_check(unsigned int m, unsigned int t)
                    216: {
                    217:   unsigned int p,k,u,h,d;
                    218:
                    219:   if ( !(m%8) )
                    220:     return 0;
                    221:   p = t*m+1;
                    222:   if ( !is_prime_small(p) )
                    223:     return 0;
                    224:   for ( k = 1, u = 2%p; ; k++ )
                    225:     if ( u == 1 )
                    226:       break;
                    227:     else
                    228:       u = (2*u)%p;
                    229:   h = t*m/k;
                    230:   d = gcd_small(h,m);
                    231:   return d == 1 ? 1 : 0;
                    232: }
                    233:
                    234: /*
                    235:  * a simple prime checker
                    236:  * return value: 1  ---  prime number
                    237:  *               0  ---  composite number
                    238:  */
                    239:
                    240: static int is_prime_small(unsigned int a)
                    241: {
                    242:   unsigned int b,t,i;
                    243:
                    244:   if ( !(a%2) ) return 0;
                    245:   for ( t = a, i = 0; t; i++, t >>= 1 );
                    246:   /* b >= sqrt(a) */
                    247:   b = 1<<((i+1)/2);
                    248:
                    249:   /* divisibility test by all odd numbers <= b */
                    250:   for ( i = 3; i <= b; i += 2 )
                    251:     if ( !(a%i) )
                    252:       return 0;
                    253:   return 1;
                    254: }
                    255:
                    256: /*
                    257:  * gcd for unsigned int as integers
                    258:  * return value: GCD(a,b)
                    259:  *
                    260:  */
                    261:
                    262:
                    263: static unsigned int gcd_small(unsigned int a,unsigned int b)
                    264: {
                    265:   unsigned int t;
                    266:
                    267:   if ( b > a ) {
                    268:     t = a; a = b; b = t;
                    269:   }
                    270:   /* Euclid's algorithm */
                    271:   while ( 1 )
                    272:     if ( !(t = a%b) ) return b;
                    273:     else {
                    274:       a = b; b = t;
                    275:     }
                    276: }
                    277:
                    278: void Pmt_save(NODE arg,Z *rp)
                    279: {
                    280:   int ret;
                    281:
                    282:   ret = mt_save(BDY((STRING)ARG0(arg)));
1.4     ! noro      283:   STOZ(ret,*rp);
1.1       noro      284: }
                    285:
                    286: void Pmt_load(NODE arg,Z *rp)
                    287: {
                    288:   int ret;
                    289:
                    290:   ret = mt_load(BDY((STRING)ARG0(arg)));
1.4     ! noro      291:   STOZ(ret,*rp);
1.1       noro      292: }
                    293:
                    294: void isqrt(Z a,Z *r);
                    295:
                    296: void Pisqrt(NODE arg,Z *rp)
                    297: {
                    298:   Z a;
                    299:   Z r;
                    300:
                    301:   a = (Z)ARG0(arg);
                    302:   asir_assert(a,O_N,"isqrt");
                    303:   if ( !a )
                    304:     *rp = 0;
                    305:   else if ( sgnz(a) < 0 )
                    306:     error("isqrt : negative argument");
                    307:   else {
                    308:     isqrt(a,rp);
                    309:   }
                    310: }
                    311:
                    312: void Pidiv(NODE arg,Z *rp)
                    313: {
                    314:   Z r;
                    315:   Z dnd,dvr;
                    316:
                    317:   dnd = (Z)ARG0(arg); dvr = (Z)ARG1(arg);
                    318:   asir_assert(dnd,O_N,"idiv");
                    319:   asir_assert(dvr,O_N,"idiv");
                    320:   if ( !dvr )
                    321:     error("idiv: division by 0");
                    322:   else if ( !dnd )
                    323:     *rp = 0;
                    324:   else
                    325:     divqrz(dnd,dvr,rp,&r);
                    326: }
                    327:
                    328: void Pirem(NODE arg,Z *rp)
                    329: {
                    330:   Z q,dnd,dvr;
                    331:
                    332:   dnd = (Z)ARG0(arg); dvr = (Z)ARG1(arg);
                    333:   asir_assert(dnd,O_N,"irem");
                    334:   asir_assert(dvr,O_N,"irem");
                    335:   if ( !dvr )
                    336:     error("irem: division by 0");
                    337:   else if ( !dnd )
                    338:     *rp = 0;
                    339:   else
                    340:     divqrz(dnd,dvr,&q,rp);
                    341: }
                    342:
                    343: void iqrv(VECT a,Z dvr,LIST *rp);
                    344:
                    345: void Piqr(NODE arg,LIST *rp)
                    346: {
                    347:   Z dnd,dvr,a,b;
                    348:   NODE node;
                    349:
                    350:   dnd = (Z)ARG0(arg); dvr = (Z)ARG1(arg);
                    351:   if ( !dvr )
                    352:     error("iqr: division by 0");
                    353:   else if ( !dnd )
                    354:     a = b = 0;
                    355:   else if ( OID(dnd) == O_VECT ) {
                    356:     iqrv((VECT)dnd,dvr,rp); return;
                    357:   } else {
                    358:     asir_assert(dnd,O_N,"iqr");
                    359:     asir_assert(dvr,O_N,"iqr");
                    360:     divqrz(dnd,dvr,&a,&b);
                    361:   }
                    362:   node = mknode(2,a,b); MKLIST(*rp,node);
                    363: }
                    364:
                    365: void Pinttorat(NODE arg,LIST *rp)
                    366: {
1.2       noro      367:   int ret;
1.3       noro      368:   Z c,m,t,b,nm,dn;
1.1       noro      369:   NODE node;
                    370:
                    371:   asir_assert(ARG0(arg),O_N,"inttorat");
                    372:   asir_assert(ARG1(arg),O_N,"inttorat");
                    373:   asir_assert(ARG2(arg),O_N,"inttorat");
                    374:   c = (Z)ARG0(arg); m = (Z)ARG1(arg); b = (Z)ARG2(arg);
1.3       noro      375:   remz(c,m,&t); c = t;
1.2       noro      376:   ret = inttorat(c,m,b,&nm,&dn);
                    377:   if ( !ret )
                    378:     *rp = 0;
                    379:   else {
                    380:     node = mknode(2,nm,dn); MKLIST(*rp,node);
                    381:   }
1.1       noro      382: }
                    383:
                    384: void Pigcd(NODE arg,Z *rp)
                    385: {
                    386:   Z n1,n2;
                    387:
                    388:   if ( argc(arg) == 1 ) {
                    389:     gcdvz((VECT)ARG0(arg),rp);
                    390:     return;
                    391:   }
                    392:   n1 = (Z)ARG0(arg); n2 = (Z)ARG1(arg);
                    393:   asir_assert(n1,O_N,"igcd");
                    394:   asir_assert(n2,O_N,"igcd");
                    395:   gcdz(n1,n2,rp);
                    396: }
                    397:
                    398: void iqrv(VECT a,Z dvr,LIST *rp)
                    399: {
                    400:   int i,n;
                    401:   VECT q,r;
                    402:   Z *b;
                    403:   NODE n0;
                    404:
                    405:   if ( !dvr )
                    406:     error("iqrv: division by 0");
                    407:   n = a->len; b = (Z *)BDY(a);
                    408:   MKVECT(q,n); MKVECT(r,n);
                    409:   for ( i = 0; i < n; i++ )
                    410:     divqrz(b[i],dvr,(Z *)&BDY(q)[i],(Z *)&BDY(r)[i]);
                    411:   n0 = mknode(2,q,r); MKLIST(*rp,n0);
                    412: }
                    413:
                    414: /*
                    415:  * gcd = GCD(a,b), ca = a/g, cb = b/g
                    416:  */
                    417:
                    418: void igcd_cofactor(Z a,Z b,Z *gcd,Z *ca,Z *cb)
                    419: {
                    420:   Z g;
                    421:
                    422:   if ( !a ) {
                    423:     if ( !b )
                    424:       error("igcd_cofactor : invalid input");
                    425:     else {
                    426:       *ca = 0;
                    427:       *cb = ONE;
                    428:       *gcd = b;
                    429:     }
                    430:   } else if ( !b ) {
                    431:     *ca = ONE;
                    432:     *cb = 0;
                    433:     *gcd = a;
                    434:   } else {
                    435:     gcdz(a,b,&g);
                    436:     divsz(a,g,ca);
                    437:     divsz(b,g,cb);
                    438:     *gcd = g;
                    439:   }
                    440: }
                    441:
                    442: void Pilcm(NODE arg,Z *rp)
                    443: {
                    444:   Z n1,n2,g,q,l;
                    445:
                    446:   n1 = (Z)ARG0(arg); n2 = (Z)ARG1(arg);
                    447:   asir_assert(n1,O_N,"ilcm");
                    448:   asir_assert(n2,O_N,"ilcm");
                    449:   if ( !n1 || !n2 )
                    450:     *rp = 0;
                    451:   else {
                    452:     gcdz(n1,n2,&g); divsz(n1,g,&q);
                    453:     mulz(q,n2,&l); absz(l,rp);
                    454:   }
                    455: }
                    456:
                    457: void Piand(NODE arg,Z *rp)
                    458: {
                    459:   mpz_t t;
                    460:   Z n1,n2;
                    461:
                    462:   n1 = (Z)ARG0(arg); n2 = (Z)ARG1(arg);
                    463:   asir_assert(n1,O_N,"iand");
                    464:   asir_assert(n2,O_N,"iand");
                    465:   if ( !n1 || !n2 ) *rp = 0;
                    466:   else {
                    467:     mpz_init(t);
                    468:     mpz_and(t,BDY(n1),BDY(n2));
                    469:     MPZTOZ(t,*rp);
                    470:   }
                    471: }
                    472:
                    473: void Pior(NODE arg,Z *rp)
                    474: {
                    475:   Z n1,n2;
                    476:   mpz_t t;
                    477:
                    478:   n1 = (Z)ARG0(arg); n2 = (Z)ARG1(arg);
                    479:   asir_assert(n1,O_N,"ior");
                    480:   asir_assert(n2,O_N,"ior");
                    481:   if ( !n1 ) *rp = n2;
                    482:   else if ( !n2 ) *rp = n1;
                    483:   else {
                    484:     mpz_init(t);
                    485:     mpz_ior(t,BDY(n1),BDY(n2));
                    486:     MPZTOZ(t,*rp);
                    487:   }
                    488: }
                    489:
                    490: void Pixor(NODE arg,Z *rp)
                    491: {
                    492:   Z n1,n2;
                    493:   mpz_t t;
                    494:
                    495:   n1 = (Z)ARG0(arg); n2 = (Z)ARG1(arg);
                    496:   asir_assert(n1,O_N,"ixor");
                    497:   asir_assert(n2,O_N,"ixor");
                    498:   if ( !n1 ) *rp = n2;
                    499:   else if ( !n2 ) *rp = n1;
                    500:   else {
                    501:     mpz_init(t);
                    502:     mpz_xor(t,BDY(n1),BDY(n2));
                    503:     MPZTOZ(t,*rp);
                    504:   }
                    505: }
                    506:
                    507: void Pishift(NODE arg,Z *rp)
                    508: {
                    509:   int i;
                    510:   Z n1,s;
                    511:   mpz_t t;
                    512:
                    513:   n1 = (Z)ARG0(arg);
                    514:   s = (Z)ARG1(arg);
                    515:   asir_assert(n1,O_N,"ixor");
                    516:   asir_assert(s,O_N,"ixor");
1.4     ! noro      517:   bshiftz(n1,ZTOS(s),rp);
1.1       noro      518: }
                    519:
                    520: void isqrt(Z a,Z *r)
                    521: {
                    522:   int k;
                    523:   Z x,t,x2,xh,quo,rem;
                    524:
                    525:   if ( !a )
                    526:     *r = 0;
                    527:   else if ( UNIZ(a) )
                    528:     *r = ONE;
                    529:   else {
                    530:     k = z_bits((Q)a); /* a <= 2^k-1 */
                    531:     bshiftz(ONE,-((k>>1)+(k&1)),&x); /* a <= x^2 */
                    532:     while ( 1 ) {
                    533:       mulz(x,x,&t);
                    534:       if ( cmpz(t,a) <= 0 ) {
                    535:         *r = x; return;
                    536:       } else {
                    537:         if ( tstbitz(x,0) )
                    538:           addz(x,a,&t);
                    539:         else
                    540:           t = a;
                    541:         bshiftz(x,-1,&x2); divqrz(t,x2,&quo,&rem);
                    542:         bshiftz(x,1,&xh); addz(quo,xh,&x);
                    543:       }
                    544:     }
                    545:   }
                    546: }
                    547:
                    548: void Pup2_init_eg(Obj *rp)
                    549: {
                    550:   up2_init_eg();
                    551:   *rp = 0;
                    552: }
                    553:
                    554: void Pup2_show_eg(Obj *rp)
                    555: {
                    556:   up2_show_eg();
                    557:   *rp = 0;
                    558: }
                    559:
                    560: void Pgf2nton(NODE arg,Z *rp)
                    561: {
                    562:   if ( !ARG0(arg) )
                    563:     *rp = 0;
                    564:   else
                    565:     up2toz(((GF2N)ARG0(arg))->body,rp);
                    566: }
                    567:
                    568: void Pntogf2n(NODE arg,GF2N *rp)
                    569: {
                    570:   UP2 t;
                    571:
                    572:   ztoup2((Z)ARG0(arg),&t);
                    573:   MKGF2N(t,*rp);
                    574: }
                    575:
                    576: void Pup2_inv(NODE arg,P *rp)
                    577: {
                    578:   UP2 a,b,t;
                    579:
                    580:   ptoup2(ARG0(arg),&a);
                    581:   ptoup2(ARG1(arg),&b);
                    582:   invup2(a,b,&t);
                    583:   up2top(t,rp);
                    584: }
                    585:
                    586: void Pinv(NODE arg,Num *rp)
                    587: {
                    588:   Num n;
                    589:   Z mod;
                    590:   MQ r;
                    591:   int inv;
                    592:
                    593:   n = (Num)ARG0(arg); mod = (Z)ARG1(arg);
                    594:   asir_assert(n,O_N,"inv");
                    595:   asir_assert(mod,O_N,"inv");
                    596:   if ( !n || !mod )
                    597:     error("inv: invalid input");
                    598:   else
                    599:     switch ( NID(n) ) {
                    600:       case N_Q:
                    601:         invz((Z)n,mod,(Z *)rp);
                    602:         break;
                    603:       case N_M:
1.4     ! noro      604:         inv = invm(CONT((MQ)n),ZTOS(mod));
1.1       noro      605:         STOMQ(inv,r);
                    606:         *rp = (Num)r;
                    607:         break;
                    608:       default:
                    609:         error("inv: invalid input");
                    610:     }
                    611: }
                    612:
                    613: void Pfac(NODE arg,Z *rp)
                    614: {
                    615:   asir_assert(ARG0(arg),O_N,"fac");
1.4     ! noro      616:   factorialz(ZTOS((Q)ARG0(arg)),rp);
1.1       noro      617: }
                    618:
                    619: void Plrandom(NODE arg,Z *rp)
                    620: {
                    621:   asir_assert(ARG0(arg),O_N,"lrandom");
1.4     ! noro      622:   randomz(ZTOS((Q)ARG0(arg)),rp);
1.1       noro      623: }
                    624:
                    625: void Prandom(NODE arg,Z *rp)
                    626: {
                    627:   unsigned int r;
                    628:
                    629: #if 0
                    630: #if defined(_PA_RISC1_1)
                    631:   r = mrand48()&BMASK;
                    632: #else
                    633:   if ( arg )
1.4     ! noro      634:     srandom(ZTOS((Q)ARG0(arg)));
1.1       noro      635:   r = random()&BMASK;
                    636: #endif
                    637: #endif
                    638:   if ( arg )
1.4     ! noro      639:     mt_sgenrand(ZTOS((Q)ARG0(arg)));
1.1       noro      640:   r = mt_genrand();
1.4     ! noro      641:   UTOZ(r,*rp);
1.1       noro      642: }
                    643:
                    644: #if defined(VISUAL) || defined(__MINGW32__)
                    645: void srandom(unsigned int);
                    646:
                    647: static unsigned int R_Next;
                    648:
                    649: unsigned int random() {
                    650:         if ( !R_Next )
                    651:                 R_Next = 1;
                    652:         return R_Next = (R_Next * 1103515245 + 12345);
                    653: }
                    654:
                    655: void srandom(unsigned int s)
                    656: {
                    657:     if ( s )
                    658:       R_Next = s;
                    659:         else if ( !R_Next )
                    660:             R_Next = 1;
                    661: }
                    662: #endif
                    663:
                    664: void Pprime(NODE arg,Z *rp)
                    665: {
                    666:   int i;
                    667:
                    668:   asir_assert(ARG0(arg),O_N,"prime");
1.4     ! noro      669:   i = ZTOS((Q)ARG0(arg));
1.1       noro      670:   if ( i < 0 || i >= 1900 )
                    671:     *rp = 0;
                    672:   else
1.4     ! noro      673:     STOZ(sprime[i],*rp);
1.1       noro      674: }
                    675:
                    676: void Plprime(NODE arg,Z *rp)
                    677: {
                    678:   int i,p;
                    679:
                    680:   asir_assert(ARG0(arg),O_N,"lprime");
1.4     ! noro      681:   i = ZTOS((Q)ARG0(arg));
1.1       noro      682:   if ( i < 0 )
                    683:     *rp = 0;
                    684:   else
                    685:     p = get_lprime(i);
1.4     ! noro      686:   STOZ(p,*rp);
1.1       noro      687: }
                    688:
                    689: extern int up_kara_mag, up_tkara_mag, up_fft_mag;
                    690:
                    691: void Pset_upfft(NODE arg,Z *rp)
                    692: {
                    693:   if ( arg ) {
                    694:     asir_assert(ARG0(arg),O_N,"set_upfft");
1.4     ! noro      695:     up_fft_mag = ZTOS((Q)ARG0(arg));
1.1       noro      696:   }
1.4     ! noro      697:   STOZ(up_fft_mag,*rp);
1.1       noro      698: }
                    699:
                    700: void Pset_upkara(NODE arg,Z *rp)
                    701: {
                    702:   if ( arg ) {
                    703:     asir_assert(ARG0(arg),O_N,"set_upkara");
1.4     ! noro      704:     up_kara_mag = ZTOS((Q)ARG0(arg));
1.1       noro      705:   }
1.4     ! noro      706:   STOZ(up_kara_mag,*rp);
1.1       noro      707: }
                    708:
                    709: void Pset_uptkara(NODE arg,Z *rp)
                    710: {
                    711:   if ( arg ) {
                    712:     asir_assert(ARG0(arg),O_N,"set_uptkara");
1.4     ! noro      713:     up_tkara_mag = ZTOS((Q)ARG0(arg));
1.1       noro      714:   }
1.4     ! noro      715:   STOZ(up_tkara_mag,*rp);
1.1       noro      716: }
                    717:
                    718: extern int up2_kara_mag;
                    719:
                    720: void Pset_up2kara(NODE arg,Z *rp)
                    721: {
                    722:   if ( arg ) {
                    723:     asir_assert(ARG0(arg),O_N,"set_up2kara");
1.4     ! noro      724:     up2_kara_mag = ZTOS((Q)ARG0(arg));
1.1       noro      725:   }
1.4     ! noro      726:   STOZ(up2_kara_mag,*rp);
1.1       noro      727: }
                    728:
                    729: void Pigcdbin(NODE arg,Z *rp)
                    730: {
                    731:   Pigcd(arg,rp);
                    732: }
                    733:
                    734: void Pigcdbmod(NODE arg,Z *rp)
                    735: {
                    736:   Pigcd(arg,rp);
                    737: }
                    738:
                    739: void Pigcdacc(NODE arg,Z *rp)
                    740: {
                    741:   Pigcd(arg,rp);
                    742: }
                    743:
                    744: void PigcdEuc(NODE arg,Z *rp)
                    745: {
                    746:   Pigcd(arg,rp);
                    747: }
                    748:
                    749: void Pigcdcntl(NODE arg,Z *rp)
                    750: {
                    751:   *rp = ONE;
                    752: }

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