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

Annotation of OpenXM_contrib2/asir2000/engine/lmi.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/lmi.c,v 1.5 2017/01/08 03:05:39 noro Exp $
1.2       noro       49: */
1.1       noro       50: #include "ca.h"
                     51: #include "base.h"
                     52: #include "inline.h"
                     53:
                     54: typedef struct oGEN_LM {
1.6     ! noro       55:   int id;
        !            56:   N dense;
        !            57:   N sparse;
        !            58:   int bit;
        !            59:   unsigned int lower;
1.1       noro       60: } *GEN_LM;
                     61:
                     62: void pwrlm0(N,N,N *);
                     63: void gen_simpn(N,N*);
                     64: void gen_simpn_force(N,N*);
                     65:
                     66: int lm_lazy;
                     67:
                     68: static GEN_LM current_mod_lm;
                     69:
1.4       noro       70: void random_lm(LM *r)
1.1       noro       71: {
1.6     ! noro       72:   N n,t;
1.1       noro       73:
1.6     ! noro       74:   if ( !current_mod_lm )
        !            75:     error("random_lm : current_mod_lm is not set");
        !            76:   randomn(current_mod_lm->bit+1,&n);
        !            77:   gen_simpn_force(n,&t);
        !            78:   MKLM(t,*r);
1.1       noro       79: }
                     80:
1.4       noro       81: void ntosparsen(N p,N *bits)
1.1       noro       82: {
1.6     ! noro       83:   int l,i,j,nz;
        !            84:   N r;
        !            85:   unsigned int *pb,*rb;
        !            86:   if ( !p )
        !            87:     *bits = 0;
        !            88:   else {
        !            89:     l = n_bits(p);
        !            90:     pb = BD(p);
        !            91:     for ( i = l-1, nz = 0; i >= 0; i-- )
        !            92:       if ( pb[i>>5]&(1<<i&31) ) nz++;
        !            93:     *bits = r = NALLOC(nz); PL(r) = nz;
        !            94:     rb = BD(r);
        !            95:     for ( i = l-1, j = 0; i >= 0; i-- )
        !            96:       if ( pb[i>>5]&(1<<i&31) )
        !            97:         rb[j++] = i;
        !            98:   }
1.1       noro       99: }
                    100:
1.4       noro      101: void setmod_lm(N p)
1.1       noro      102: {
1.6     ! noro      103:   int i;
1.5       noro      104:     Q q;
1.1       noro      105:
1.6     ! noro      106:   if ( !current_mod_lm )
        !           107:     current_mod_lm = (GEN_LM)MALLOC(sizeof(struct oGEN_LM));
        !           108:   bzero((char *)current_mod_lm,sizeof(struct oGEN_LM));
        !           109:   current_mod_lm->dense = p;
        !           110:   ntosparsen(p,&current_mod_lm->sparse);
        !           111:   current_mod_lm->bit = n_bits(p);
        !           112:   current_mod_lm->id = UP2_DENSE; /* use ordinary rep. by default */
        !           113:   if ( !(current_mod_lm->bit%32) ) {
        !           114:     for ( i = PL(p)-1; i >= 1; i-- )
        !           115:       if ( BD(p)[i] != 0xffffffff )
        !           116:         break;
        !           117:     if ( !i ) {
        !           118:       current_mod_lm->lower = (unsigned int)((((UL)1)<<32)-((UL)BD(p)[0]));
        !           119:       current_mod_lm->id = UP2_SPARSE; /* XXX */
        !           120:     }
        !           121:   }
1.5       noro      122:     setmod_lf(p);
1.1       noro      123: }
                    124:
1.4       noro      125: void getmod_lm(N *p)
1.1       noro      126: {
1.6     ! noro      127:   if ( !current_mod_lm )
        !           128:     *p = 0;
        !           129:   else
        !           130:     *p = current_mod_lm->dense;
1.1       noro      131: }
                    132:
1.4       noro      133: void simplm(LM n,LM *r)
1.1       noro      134: {
1.6     ! noro      135:   N rem;
1.1       noro      136:
1.6     ! noro      137:   if ( !n )
        !           138:     *r = 0;
        !           139:   else if ( NID(n) != N_LM )
        !           140:     *r = n;
        !           141:   else {
        !           142:     gen_simpn(n->body,&rem);
        !           143:     MKLM(rem,*r);
        !           144:   }
1.1       noro      145: }
                    146:
1.5       noro      147: void simplm_force(LM n,LM *r)
                    148: {
1.6     ! noro      149:   N rem;
1.5       noro      150:
1.6     ! noro      151:   if ( !n )
        !           152:     *r = 0;
        !           153:   else if ( NID(n) != N_LM )
        !           154:     *r = n;
        !           155:   else {
        !           156:     gen_simpn_force(n->body,&rem);
        !           157:     MKLM(rem,*r);
        !           158:   }
1.5       noro      159: }
                    160:
                    161:
1.4       noro      162: void qtolm(Q q,LM *l)
1.1       noro      163: {
1.6     ! noro      164:   N rn;
        !           165:   LM a,b,c;
1.1       noro      166:
1.6     ! noro      167:   if ( !q || (OID(q)==O_N && ((NID(q) == N_LM) || (NID(q) == N_GFPN))) ) { /* XXX */
        !           168:     *l = (LM)q;
        !           169:   } else if ( OID(q) == O_N && NID(q) == N_Q ) {
        !           170:     gen_simpn(NM(q),&rn); MKLM(rn,a);
        !           171:     if ( SGN(q) < 0 ) {
        !           172:       chsgnlm(a,&c); a = c;
        !           173:     }
        !           174:     if ( DN(q) ) {
        !           175:       gen_simpn_force(DN(q),&rn); MKLM(rn,b);
        !           176:       if ( !b )
        !           177:         error("qtolm : invalid denominator");
        !           178:       divlm(a,b,l);
        !           179:     } else
        !           180:       *l = a;
        !           181:   } else
        !           182:     error("qtolm : invalid argument");
1.1       noro      183: }
                    184:
                    185: #define NZLM(a) ((a)&&(NID(a)==N_LM))
                    186:
1.4       noro      187: void addlm(LM a,LM b,LM *c)
1.1       noro      188: {
1.6     ! noro      189:   N t,t1;
        !           190:   LM z;
1.1       noro      191:
1.6     ! noro      192:   qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z;
        !           193:   if ( !a )
        !           194:     *c = b;
        !           195:   else if ( !b )
        !           196:     *c = a;
        !           197:   else {
        !           198:     addn(a->body,b->body,&t);
        !           199:     gen_simpn(t,&t1);
        !           200:     MKLM(t1,*c);
        !           201:   }
1.1       noro      202: }
                    203:
1.4       noro      204: void sublm(LM a,LM b,LM *c)
1.1       noro      205: {
1.6     ! noro      206:   int s;
        !           207:   N t;
        !           208:   LM z;
        !           209:
        !           210:   qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z;
        !           211:   if ( !b )
        !           212:     *c = a;
        !           213:   else if ( !a )
        !           214:     chsgnlm(b,c);
        !           215:   else {
        !           216:     s = subn(a->body,b->body,&t);
        !           217:     if ( !s )
        !           218:       *c = 0;
        !           219:     else {
        !           220:       MKLM(t,z);
        !           221:       if ( s < 0 )
        !           222:         chsgnlm(z,c);
        !           223:       else
        !           224:         *c = z;
        !           225:     }
        !           226:   }
1.1       noro      227: }
                    228:
1.4       noro      229: void mullm(LM a,LM b,LM *c)
1.1       noro      230: {
1.6     ! noro      231:   N t,r;
        !           232:   LM z;
1.1       noro      233:
1.6     ! noro      234:   qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z;
        !           235:   if ( !a || !b )
        !           236:     *c = 0;
        !           237:   else {
        !           238:     kmuln(a->body,b->body,&t);
        !           239:     gen_simpn(t,&r);
        !           240:     MKLM(r,*c);
        !           241:   }
1.1       noro      242: }
                    243:
1.4       noro      244: void divlm(LM a,LM b,LM *c)
1.1       noro      245: {
1.6     ! noro      246:   LM r;
        !           247:   Q t,m,i;
        !           248:   LM z;
        !           249:
        !           250:   qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z;
        !           251:   if ( !b )
        !           252:     error("divlm: division by 0");
        !           253:   else if ( !a )
        !           254:     *c = 0;
        !           255:   else {
        !           256:     NTOQ(b->body,1,t); NTOQ(current_mod_lm->dense,1,m);
        !           257:     invl(t,m,&i);
        !           258:     MKLM(NM(i),r);
        !           259:     mullm(a,r,c);
        !           260:   }
1.1       noro      261: }
                    262:
1.4       noro      263: void chsgnlm(LM a,LM *c)
1.1       noro      264: {
1.6     ! noro      265:   LM t;
        !           266:   N s,u;
1.1       noro      267:
1.6     ! noro      268:   if ( !a )
        !           269:     *c = a;
        !           270:   else {
        !           271:     qtolm((Q)a,&t); a = t;
        !           272:     gen_simpn_force(a->body,&s);
        !           273:     subn(current_mod_lm->dense,s,&u);
        !           274:     MKLM(u,*c);
        !           275:   }
1.1       noro      276: }
                    277:
1.4       noro      278: void pwrlm(LM a,Q b,LM *c)
1.1       noro      279: {
1.6     ! noro      280:   LM t;
        !           281:   N s;
1.1       noro      282:
1.6     ! noro      283:   if ( !b ) {
        !           284:     MKLM(ONEN,*c);
        !           285:   } else if ( !a )
        !           286:     *c = 0;
        !           287:   else {
        !           288:     qtolm((Q)a,&t); a = t;
        !           289:     pwrlm0(a->body,NM(b),&s);
        !           290:     MKLM(s,*c);
        !           291:   }
1.1       noro      292: }
                    293:
1.4       noro      294: void pwrlm0(N a,N n,N *c)
1.1       noro      295: {
1.6     ! noro      296:   N n1,t,t1,t2,r1;
        !           297:   int r;
1.1       noro      298:
1.6     ! noro      299:   if ( !n )
        !           300:     *c = ONEN;
        !           301:   else if ( UNIN(n) )
        !           302:     *c = a;
        !           303:   else {
        !           304:     r = divin(n,2,&n1);
        !           305:     pwrlm0(a,n1,&t);
        !           306:     kmuln(t,t,&t1); gen_simpn(t1,&r1);
        !           307:     if ( r ) {
        !           308:       kmuln(r1,a,&t2); gen_simpn(t2,&r1);
        !           309:     }
        !           310:     *c = r1;
        !           311:   }
1.1       noro      312: }
                    313:
1.4       noro      314: int cmplm(LM a,LM b)
1.1       noro      315: {
1.6     ! noro      316:   LM z;
1.1       noro      317:
1.6     ! noro      318:   qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z;
        !           319:   if ( !a )
        !           320:     if ( !b )
        !           321:       return 0;
        !           322:     else
        !           323:       return -1;
        !           324:   else if ( !b )
        !           325:       return 1;
        !           326:   else
        !           327:     return cmpn(a->body,b->body);
1.1       noro      328: }
                    329:
                    330: void remn_special(N,N,int,unsigned int ,N *);
                    331:
1.4       noro      332: void gen_simpn(N a,N *b)
1.1       noro      333: {
1.6     ! noro      334:   if ( !current_mod_lm )
        !           335:     error("gen_simpn: current_mod_lm is not set");
        !           336:   if ( lm_lazy )
        !           337:     *b = a;
        !           338:   else if ( current_mod_lm->id == UP2_SPARSE )
        !           339:     remn_special(a,current_mod_lm->dense,current_mod_lm->bit,current_mod_lm->lower,b);
        !           340:   else
        !           341:     remn(a,current_mod_lm->dense,b);
1.1       noro      342: }
                    343:
1.4       noro      344: void gen_simpn_force(N a,N *b)
1.1       noro      345: {
1.6     ! noro      346:   if ( !current_mod_lm )
        !           347:     error("gen_simpn_force: current_mod_lm is not set");
        !           348:   else if ( current_mod_lm->id == UP2_SPARSE )
        !           349:     remn_special(a,current_mod_lm->dense,current_mod_lm->bit,current_mod_lm->lower,b);
        !           350:   else
        !           351:     remn(a,current_mod_lm->dense,b);
1.1       noro      352: }
                    353:

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