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

Diff for /OpenXM_contrib2/asir2000/engine/gfs.c between version 1.18 and 1.19

version 1.18, 2017/02/27 05:14:54 version 1.19, 2018/03/29 01:32:52
Line 45 
Line 45 
  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,   * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.   * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
  *   *
  * $OpenXM: OpenXM_contrib2/asir2000/engine/gfs.c,v 1.17 2003/05/07 06:26:49 noro Exp $   * $OpenXM: OpenXM_contrib2/asir2000/engine/gfs.c,v 1.18 2017/02/27 05:14:54 noro Exp $
 */  */
 #include "ca.h"  #include "ca.h"
 #include "inline.h"  #include "inline.h"
Line 61  int *current_gfs_ntoi;
Line 61  int *current_gfs_ntoi;
 int *current_gfs_iton;  int *current_gfs_iton;
   
 struct prim_root_info {  struct prim_root_info {
         int p;    int p;
         int extdeg;    int extdeg;
         int defpoly;    int defpoly;
         int prim_root;    int prim_root;
 };  };
   
 /* if p >= SF_THRESHOLD, usual representation is used */  /* if p >= SF_THRESHOLD, usual representation is used */
Line 477  struct prim_root_info prim_root_info_tab[] = {
Line 477  struct prim_root_info prim_root_info_tab[] = {
   
 void dec_um(int p,int a,UM u)  void dec_um(int p,int a,UM u)
 {  {
         int i;    int i;
   
         for ( i = 0; a; i++, a /= p )    for ( i = 0; a; i++, a /= p )
                 COEF(u)[i] = a%p;      COEF(u)[i] = a%p;
         DEG(u) = i-1;    DEG(u) = i-1;
 }  }
   
 /*  /*
Line 493  void dec_um(int p,int a,UM u)
Line 493  void dec_um(int p,int a,UM u)
   
 void setmod_sf(int p,int n)  void setmod_sf(int p,int n)
 {  {
         int r,i,q1,q,t,t1;    int r,i,q1,q,t,t1;
         UM dp;    UM dp;
   
         if ( p >= SF_THRESHOLD ) {    if ( p >= SF_THRESHOLD ) {
                 if ( n > 1 )      if ( n > 1 )
                         error("setmod_ff : p^n is too large");        error("setmod_ff : p^n is too large");
                 current_gfs_p = p;      current_gfs_p = p;
                 current_gfs_q = p;      current_gfs_q = p;
                 current_gfs_q1 = p-1;      current_gfs_q1 = p-1;
                 current_gfs_ext = 0;      current_gfs_ext = 0;
                 current_gfs_ntoi = 0;      current_gfs_ntoi = 0;
                 current_gfs_iton = 0;      current_gfs_iton = 0;
                 current_gfs_plus1 = 0;      current_gfs_plus1 = 0;
                 return;      return;
         }    }
   
         for ( i = 0, q = 1; i < n; i++ )    for ( i = 0, q = 1; i < n; i++ )
                 q *= p;      q *= p;
     if ( p == current_gfs_p && q == current_gfs_q ) return;      if ( p == current_gfs_p && q == current_gfs_q ) return;
         dp = UMALLOC(n);    dp = UMALLOC(n);
         r = search_defpoly_and_primitive_root(p,n,dp);    r = search_defpoly_and_primitive_root(p,n,dp);
         if ( !r ) {    if ( !r ) {
                 generate_defpoly_um(p,n,dp);      generate_defpoly_um(p,n,dp);
                 r = generate_primitive_root_enc(p,n,dp);      r = generate_primitive_root_enc(p,n,dp);
                 if ( !r )      if ( !r )
                         error("setmod_sf : primitive root not found");        error("setmod_sf : primitive root not found");
         }    }
         current_gfs_p = p;    current_gfs_p = p;
         current_gfs_q = q;    current_gfs_q = q;
         current_gfs_q1 = q1 = q-1;    current_gfs_q1 = q1 = q-1;
         if ( n > 1 )    if ( n > 1 )
                 umtop(CO->v,dp,&current_gfs_ext);      umtop(CO->v,dp,&current_gfs_ext);
         else    else
                 current_gfs_ext = 0;      current_gfs_ext = 0;
         current_gfs_iton = (int *)MALLOC(q1*sizeof(int));    current_gfs_iton = (int *)MALLOC(q1*sizeof(int));
         current_gfs_iton[0] = 1;    current_gfs_iton[0] = 1;
         for ( i = 1; i < q1; i++ )    for ( i = 1; i < q1; i++ )
                 current_gfs_iton[i] = mulremum_enc(p,n,dp,current_gfs_iton[i-1],r);      current_gfs_iton[i] = mulremum_enc(p,n,dp,current_gfs_iton[i-1],r);
   
         current_gfs_ntoi = (int *)MALLOC(q*sizeof(int));    current_gfs_ntoi = (int *)MALLOC(q*sizeof(int));
         current_gfs_ntoi[0] = -1;    current_gfs_ntoi[0] = -1;
         for ( i = 0; i < q1; i++ )    for ( i = 0; i < q1; i++ )
                 current_gfs_ntoi[current_gfs_iton[i]] = i;      current_gfs_ntoi[current_gfs_iton[i]] = i;
   
         current_gfs_plus1 = (int *)MALLOC(q*sizeof(int));    current_gfs_plus1 = (int *)MALLOC(q*sizeof(int));
         for ( i = 0; i < q1; i++ ) {    for ( i = 0; i < q1; i++ ) {
                 t = current_gfs_iton[i];      t = current_gfs_iton[i];
                 /* add 1 to the constant part */      /* add 1 to the constant part */
                 t1 = (t/p)*p+((t+1)%p);      t1 = (t/p)*p+((t+1)%p);
                 current_gfs_plus1[i] = current_gfs_ntoi[t1];      current_gfs_plus1[i] = current_gfs_ntoi[t1];
         }    }
 }  }
   
 int search_defpoly_and_primitive_root(int p,int n,UM dp)  int search_defpoly_and_primitive_root(int p,int n,UM dp)
 {  {
         int l,min,max,mid,p1,i,ind,t;    int l,min,max,mid,p1,i,ind,t;
   
         l = sizeof(prim_root_info_tab)/sizeof(struct prim_root_info);    l = sizeof(prim_root_info_tab)/sizeof(struct prim_root_info);
         min = 0; max = l-1;    min = 0; max = l-1;
         ind = -1;    ind = -1;
         while ( max - min > 1 ) {    while ( max - min > 1 ) {
                 mid = (max+min)/2;      mid = (max+min)/2;
                 p1 = prim_root_info_tab[mid].p;      p1 = prim_root_info_tab[mid].p;
                 if ( p1 == p ) {      if ( p1 == p ) {
                         ind = mid; break;        ind = mid; break;
                 } else if ( p1 > p )      } else if ( p1 > p )
                         max = mid;        max = mid;
                 else      else
                         min = mid;        min = mid;
         }    }
         if ( ind < 0 ) {    if ( ind < 0 ) {
                 if ( prim_root_info_tab[min].p == p )      if ( prim_root_info_tab[min].p == p )
                         ind = min;        ind = min;
                 else if ( prim_root_info_tab[max].p == p )      else if ( prim_root_info_tab[max].p == p )
                         ind = max;        ind = max;
                 else      else
                         return 0; /* XXX */        return 0; /* XXX */
         }    }
         /* now prim_root_info_tab[ind].p = p */    /* now prim_root_info_tab[ind].p = p */
         t = ind - (prim_root_info_tab[ind].extdeg-1);    t = ind - (prim_root_info_tab[ind].extdeg-1);
         /* now prim_root_info_tab[t].extdeg = 1 */    /* now prim_root_info_tab[t].extdeg = 1 */
         for ( i = t; prim_root_info_tab[i].p == p; i++ )    for ( i = t; prim_root_info_tab[i].p == p; i++ )
                 if ( prim_root_info_tab[i].extdeg == n )      if ( prim_root_info_tab[i].extdeg == n )
                         break;        break;
         if ( prim_root_info_tab[i].p != p )    if ( prim_root_info_tab[i].p != p )
                 return 0;      return 0;
         dec_um(p,prim_root_info_tab[i].defpoly,dp);    dec_um(p,prim_root_info_tab[i].defpoly,dp);
         return prim_root_info_tab[i].prim_root;    return prim_root_info_tab[i].prim_root;
 }  }
   
 void generate_defpoly_um(int p,int n,UM dp)  void generate_defpoly_um(int p,int n,UM dp)
 {  {
         int i,j,a,q;    int i,j,a,q;
         UM wf,wdf,wgcd;    UM wf,wdf,wgcd;
   
         wf = W_UMALLOC(n);    wf = W_UMALLOC(n);
         wdf = W_UMALLOC(n);    wdf = W_UMALLOC(n);
         wgcd = W_UMALLOC(n);    wgcd = W_UMALLOC(n);
         COEF(dp)[n] = 1;    COEF(dp)[n] = 1;
         DEG(dp) = n;    DEG(dp) = n;
         for ( i = 0, q = 1; i < n; i++ )    for ( i = 0, q = 1; i < n; i++ )
                 q *= p;      q *= p;
         for ( i = 0; i < q; i++ ) {    for ( i = 0; i < q; i++ ) {
                 for ( j = 0, a = i; a; j++, a /= p )      for ( j = 0, a = i; a; j++, a /= p )
                         COEF(dp)[j] = a%p;        COEF(dp)[j] = a%p;
                 for ( ; j < n; j++ )      for ( ; j < n; j++ )
                         COEF(dp)[j] = 0;        COEF(dp)[j] = 0;
                 cpyum(dp,wf);      cpyum(dp,wf);
                 diffum(p,dp,wdf);      diffum(p,dp,wdf);
                 gcdum(p,wf,wdf,wgcd);      gcdum(p,wf,wdf,wgcd);
                 if ( DEG(wgcd) >= 1 )      if ( DEG(wgcd) >= 1 )
                         continue;        continue;
                 mini(p,dp,wf);      mini(p,dp,wf);
                 if ( DEG(wf) <= 0 )      if ( DEG(wf) <= 0 )
                         return;        return;
         }    }
 }  }
   
 int generate_primitive_root_enc(int p,int n,UM dp)  int generate_primitive_root_enc(int p,int n,UM dp)
 {  {
         int i,r,rj,j,q;    int i,r,rj,j,q;
   
         if ( p == 2 && n == 1 )    if ( p == 2 && n == 1 )
                 return 1;      return 1;
   
         for ( i = 0, q = 1; i < n; i++ )    for ( i = 0, q = 1; i < n; i++ )
                  q *= p;       q *= p;
         for ( r = n==1?2:p; r < q; r++ ) {    for ( r = n==1?2:p; r < q; r++ ) {
                 rj = r;      rj = r;
                 for ( j = 1; j < q-1 && rj != 1; j++ )      for ( j = 1; j < q-1 && rj != 1; j++ )
                         rj = mulremum_enc(p,n,dp,rj,r);        rj = mulremum_enc(p,n,dp,rj,r);
                 if ( j == q-1 )      if ( j == q-1 )
                         return r;        return r;
         }    }
         /* not found */    /* not found */
         return 0;    return 0;
 }  }
   
 /* [a(p)]*[b(p)] in GF(p^n) -> [a(x)*b(x) mod dp(x)]_{x->p} */  /* [a(p)]*[b(p)] in GF(p^n) -> [a(x)*b(x) mod dp(x)]_{x->p} */
   
 int mulremum_enc(int p,int n,UM dp,int a,int b)  int mulremum_enc(int p,int n,UM dp,int a,int b)
 {  {
         int i,dr,r;    int i,dr,r;
         UM wa,wb,wc,wq;    UM wa,wb,wc,wq;
   
         if ( n == 1 )    if ( n == 1 )
                 return (a*b)%p;      return (a*b)%p;
         if ( !a || !b )    if ( !a || !b )
                 return 0;      return 0;
   
         wa = W_UMALLOC(n);    wa = W_UMALLOC(n);
         dec_um(p,a,wa);    dec_um(p,a,wa);
   
         wb = W_UMALLOC(n);    wb = W_UMALLOC(n);
         dec_um(p,b,wb);    dec_um(p,b,wb);
   
         wc = W_UMALLOC(2*n);    wc = W_UMALLOC(2*n);
         wq = W_UMALLOC(2*n);    wq = W_UMALLOC(2*n);
         mulum(p,wa,wb,wc);    mulum(p,wa,wb,wc);
         dr = divum(p,wc,dp,wq);    dr = divum(p,wc,dp,wq);
         for ( i = dr, r = 0; i >= 0; i-- )    for ( i = dr, r = 0; i >= 0; i-- )
                 r = r*p+COEF(wc)[i];      r = r*p+COEF(wc)[i];
         return r;    return r;
 }  }
   
 /* sigma : alpha -> alpha^q */  /* sigma : alpha -> alpha^q */
   
 void gfs_galois_action(GFS a,Q e,GFS *c)  void gfs_galois_action(GFS a,Q e,GFS *c)
 {  {
         Q p;    Q p;
         int i,k;    int i,k;
         GFS t,s;    GFS t,s;
   
         t = a;    t = a;
         k = QTOS(e);    k = QTOS(e);
         STOQ(current_gfs_p,p);    STOQ(current_gfs_p,p);
         for ( i = 0; i < k; i++ ) {    for ( i = 0; i < k; i++ ) {
                 pwrgfs(t,p,&s); t = s;      pwrgfs(t,p,&s); t = s;
         }    }
         *c = t;    *c = t;
 }  }
   
 /* GF(pn)={0,1,a,a^2,...} -> GF(pm)={0,1,b,b^2,...}; a->b^k */  /* GF(pn)={0,1,a,a^2,...} -> GF(pm)={0,1,b,b^2,...}; a->b^k */
   
 void gfs_embed(GFS z,int k,int pm,GFS *c)  void gfs_embed(GFS z,int k,int pm,GFS *c)
 {  {
         int t;    int t;
   
         if ( !z )    if ( !z )
                 *c = 0;      *c = 0;
         else {    else {
                 t = dmar(k,CONT(z),0,pm-1);      t = dmar(k,CONT(z),0,pm-1);
                 MKGFS(t,*c);      MKGFS(t,*c);
         }    }
 }  }
   
 /* 0 <= index <= q-1 */  /* 0 <= index <= q-1 */
   
 void indextogfs(int index, GFS *c)  void indextogfs(int index, GFS *c)
 {  {
         if ( index == 0 )    if ( index == 0 )
                 *c = 0;      *c = 0;
         else if ( index >= current_gfs_q )    else if ( index >= current_gfs_q )
                 error("indextogfs : exhausted");      error("indextogfs : exhausted");
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 MKGFS(index,*c);      MKGFS(index,*c);
         } else {    } else {
                 MKGFS(index-1,*c);      MKGFS(index-1,*c);
         }    }
 }  }
   
 void itogfs(int n, GFS *c)  void itogfs(int n, GFS *c)
 {  {
         n = _itosf(n);    n = _itosf(n);
         if ( !n )    if ( !n )
                 *c = 0;      *c = 0;
         else {    else {
                 n = IFTOF(n);      n = IFTOF(n);
                 MKGFS(n,*c);      MKGFS(n,*c);
         }    }
 }  }
   
 void iftogfs(int n, GFS *c)  void iftogfs(int n, GFS *c)
 {  {
         if ( !n )    if ( !n )
                 *c = 0;      *c = 0;
         else {    else {
                 MKGFS(IFTOF(n),*c);      MKGFS(IFTOF(n),*c);
         }    }
 }  }
   
 void qtogfs(Q a,GFS *c)  void qtogfs(Q a,GFS *c)
 {  {
         int s;    int s;
   
         if ( a && (SGN(a) < 1) )    if ( a && (SGN(a) < 1) )
                 error("qtogfs : invalid argument");      error("qtogfs : invalid argument");
         s = QTOS(a)%current_gfs_q;    s = QTOS(a)%current_gfs_q;
         itogfs(s,c);    itogfs(s,c);
 }  }
   
 void mqtogfs(MQ a,GFS *c)  void mqtogfs(MQ a,GFS *c)
 {  {
         if ( !a )    if ( !a )
                 *c = 0;      *c = 0;
         else    else
                 itogfs(CONT(a),c);      itogfs(CONT(a),c);
 }  }
   
 void gfstomq(GFS a,MQ *c)  void gfstomq(GFS a,MQ *c)
 {  {
         if ( !a )    if ( !a )
                 *c = 0;      *c = 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 UTOMQ(CONT(a),*c);      UTOMQ(CONT(a),*c);
         } else {    } else {
                 UTOMQ(current_gfs_iton[CONT(a)],*c);      UTOMQ(current_gfs_iton[CONT(a)],*c);
         }    }
 }  }
   
 void gfstopgfs(GFS a,V v,P *c)  void gfstopgfs(GFS a,V v,P *c)
 {  {
         MQ t;    MQ t;
         Q q;    Q q;
   
         if ( !a )    if ( !a )
                 *c = 0;      *c = 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 UTOMQ(CONT(a),t);      UTOMQ(CONT(a),t);
                 STOQ(CONT(t),q);      STOQ(CONT(t),q);
                 *c = (P)q;      *c = (P)q;
         } else    } else
                 enc_to_p(current_gfs_p,current_gfs_iton[CONT(a)],v,c);      enc_to_p(current_gfs_p,current_gfs_iton[CONT(a)],v,c);
 }  }
   
 void ntogfs(Obj a,GFS *b)  void ntogfs(Obj a,GFS *b)
 {  {
         P t;    P t;
   
         if ( !current_gfs_q1 )    if ( !current_gfs_q1 )
                 error("addgfs : current_gfs_q is not set");      error("addgfs : current_gfs_q is not set");
         if ( !a || (OID(a)==O_N && NID(a) == N_GFS) )    if ( !a || (OID(a)==O_N && NID(a) == N_GFS) )
                 *b = (GFS)a;      *b = (GFS)a;
         else if ( OID(a) == O_N && NID(a) == N_M )    else if ( OID(a) == O_N && NID(a) == N_M )
                 mqtogfs((MQ)a,b);      mqtogfs((MQ)a,b);
         else if ( OID(a) == O_N && NID(a) == N_Q ) {    else if ( OID(a) == O_N && NID(a) == N_Q ) {
                 ptomp(current_gfs_p,(P)a,&t); mqtogfs((MQ)t,b);      ptomp(current_gfs_p,(P)a,&t); mqtogfs((MQ)t,b);
         } else    } else
                 error("ntogfs : invalid argument");      error("ntogfs : invalid argument");
 }  }
   
 void addgfs(GFS a,GFS b,GFS *c)  void addgfs(GFS a,GFS b,GFS *c)
 {  {
         int ai,bi,ci;    int ai,bi,ci;
         GFS z;    GFS z;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         ntogfs((Obj)b,&z); b = z;    ntogfs((Obj)b,&z); b = z;
         if ( !a )    if ( !a )
                 *c = b;      *c = b;
         else if ( !b )    else if ( !b )
                 *c = a;      *c = a;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 ai = CONT(a); bi = CONT(b);      ai = CONT(a); bi = CONT(b);
                 ci = ai+bi-current_gfs_q;      ci = ai+bi-current_gfs_q;
                 if ( ci == 0 )      if ( ci == 0 )
                         *c = 0;        *c = 0;
                 else {      else {
                         if ( ci < 0 )        if ( ci < 0 )
                                 ci += current_gfs_q;          ci += current_gfs_q;
                         MKGFS(ci,*c);        MKGFS(ci,*c);
                 }      }
         } else {    } else {
                 ai = CONT(a); bi = CONT(b);      ai = CONT(a); bi = CONT(b);
                 if ( ai > bi ) {      if ( ai > bi ) {
                         /* tab[ai]+tab[bi] = tab[bi](tab[ai-bi]+1) */        /* tab[ai]+tab[bi] = tab[bi](tab[ai-bi]+1) */
                         ci = current_gfs_plus1[ai-bi];        ci = current_gfs_plus1[ai-bi];
                         if ( ci < 0 )        if ( ci < 0 )
                                 *c = 0;          *c = 0;
                         else {        else {
                                 ci += bi;          ci += bi;
                                 if ( ci >= current_gfs_q1 )          if ( ci >= current_gfs_q1 )
                                         ci -= current_gfs_q1;            ci -= current_gfs_q1;
                                 MKGFS(ci,*c);          MKGFS(ci,*c);
                         }        }
                 } else {      } else {
                         /* tab[ai]+tab[bi] = tab[ai](tab[bi-ai]+1) */        /* tab[ai]+tab[bi] = tab[ai](tab[bi-ai]+1) */
                         ci = current_gfs_plus1[bi-ai];        ci = current_gfs_plus1[bi-ai];
                         if ( ci < 0 )        if ( ci < 0 )
                                 *c = 0;          *c = 0;
                         else {        else {
                                 ci += ai;          ci += ai;
                                 if ( ci >= current_gfs_q1 )          if ( ci >= current_gfs_q1 )
                                         ci -= current_gfs_q1;            ci -= current_gfs_q1;
                                 MKGFS(ci,*c);          MKGFS(ci,*c);
                         }        }
                 }      }
         }    }
 }  }
   
 void subgfs(GFS a,GFS b,GFS *c)  void subgfs(GFS a,GFS b,GFS *c)
 {  {
         GFS t,z;    GFS t,z;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         ntogfs((Obj)b,&z); b = z;    ntogfs((Obj)b,&z); b = z;
         if ( !b )    if ( !b )
                 *c = a;      *c = a;
         else {    else {
                 chsgngfs(b,&t);      chsgngfs(b,&t);
                 addgfs(a,t,c);      addgfs(a,t,c);
         }    }
 }  }
   
 void mulgfs(GFS a,GFS b,GFS *c)  void mulgfs(GFS a,GFS b,GFS *c)
 {  {
         int ai,bi,ci;    int ai,bi,ci;
         GFS z;    GFS z;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         ntogfs((Obj)b,&z); b = z;    ntogfs((Obj)b,&z); b = z;
         if ( !a || !b )    if ( !a || !b )
                 *c = 0;      *c = 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 ai = CONT(a); bi = CONT(b);      ai = CONT(a); bi = CONT(b);
                 DMAR(ai,bi,0,current_gfs_q,ci);      DMAR(ai,bi,0,current_gfs_q,ci);
                 MKGFS(ci,*c);      MKGFS(ci,*c);
         } else {    } else {
                 ai = CONT(a) + CONT(b);      ai = CONT(a) + CONT(b);
                 if ( ai >= current_gfs_q1 )      if ( ai >= current_gfs_q1 )
                         ai -= current_gfs_q1;        ai -= current_gfs_q1;
                 MKGFS(ai,*c);      MKGFS(ai,*c);
         }    }
 }  }
   
 void divgfs(GFS a,GFS b,GFS *c)  void divgfs(GFS a,GFS b,GFS *c)
 {  {
         int ai,bi,ci;    int ai,bi,ci;
         GFS z;    GFS z;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         ntogfs((Obj)b,&z); b = z;    ntogfs((Obj)b,&z); b = z;
         if ( !b )    if ( !b )
                 error("divgfs : division by 0");      error("divgfs : division by 0");
         else if ( !a )    else if ( !a )
                 *c = 0;      *c = 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 ai = CONT(a); bi = invm(CONT(b),current_gfs_q);      ai = CONT(a); bi = invm(CONT(b),current_gfs_q);
                 DMAR(ai,bi,0,current_gfs_q,ci);      DMAR(ai,bi,0,current_gfs_q,ci);
                 MKGFS(ci,*c);      MKGFS(ci,*c);
         } else {    } else {
                 ai = CONT(a) - CONT(b);      ai = CONT(a) - CONT(b);
                 if ( ai < 0 )      if ( ai < 0 )
                         ai += current_gfs_q1;        ai += current_gfs_q1;
                 MKGFS(ai,*c);      MKGFS(ai,*c);
         }    }
 }  }
   
 void chsgngfs(GFS a,GFS *c)  void chsgngfs(GFS a,GFS *c)
 {  {
         int ai;    int ai;
         GFS z;    GFS z;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         if ( !a )    if ( !a )
                 *c = 0;      *c = 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 ai = current_gfs_q - CONT(a);      ai = current_gfs_q - CONT(a);
                 MKGFS(ai,*c);      MKGFS(ai,*c);
         } else if ( current_gfs_q1&1 )    } else if ( current_gfs_q1&1 )
                 *c = a;      *c = a;
         else {    else {
                 /* r^((q-1)/2) = -1 */      /* r^((q-1)/2) = -1 */
                 ai = CONT(a)+(current_gfs_q1>>1);      ai = CONT(a)+(current_gfs_q1>>1);
                 if ( ai >= current_gfs_q1 )      if ( ai >= current_gfs_q1 )
                         ai -= current_gfs_q1;        ai -= current_gfs_q1;
                 MKGFS(ai,*c);      MKGFS(ai,*c);
         }    }
 }  }
   
 void pwrgfs(GFS a,Q b,GFS *c)  void pwrgfs(GFS a,Q b,GFS *c)
 {  {
         N an,tn,rn;    N an,tn,rn;
         GFS t,s,z;    GFS t,s,z;
         int ai;    int ai;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         if ( !b )    if ( !b )
                 itogfs(1,c);      itogfs(1,c);
         else if ( !a )    else if ( !a )
                 *c = 0;      *c = 0;
         else if ( !current_gfs_ntoi) {    else if ( !current_gfs_ntoi) {
                 ai = pwrm(current_gfs_q,CONT(a),QTOS(b));      ai = pwrm(current_gfs_q,CONT(a),QTOS(b));
                 MKGFS(ai,*c);      MKGFS(ai,*c);
         } else {    } else {
                 STON(CONT(a),an); muln(an,NM(b),&tn);      STON(CONT(a),an); muln(an,NM(b),&tn);
                 STON(current_gfs_q1,an); remn(tn,an,&rn);      STON(current_gfs_q1,an); remn(tn,an,&rn);
                 if ( !rn )      if ( !rn )
                         itogfs(1,c);        itogfs(1,c);
                 else if ( SGN(b) > 0 )      else if ( SGN(b) > 0 )
                         MKGFS(BD(rn)[0],*c);        MKGFS(BD(rn)[0],*c);
                 else {      else {
                         itogfs(1,&t);        itogfs(1,&t);
                         MKGFS(BD(rn)[0],s);        MKGFS(BD(rn)[0],s);
                         divgfs(t,s,c);        divgfs(t,s,c);
                 }      }
         }    }
 }  }
   
 int cmpgfs(GFS a,GFS b)  int cmpgfs(GFS a,GFS b)
 {  {
         GFS z;    GFS z;
   
         ntogfs((Obj)a,&z); a = z;    ntogfs((Obj)a,&z); a = z;
         if ( !a )    if ( !a )
                 return !b ? 0 : -1;      return !b ? 0 : -1;
         else    else
                 if ( !b )      if ( !b )
                         return 1;        return 1;
                 else {      else {
                         if ( CONT(a) > CONT(b) )        if ( CONT(a) > CONT(b) )
                                 return 1;          return 1;
                         else if ( CONT(a) < CONT(b) )        else if ( CONT(a) < CONT(b) )
                                 return -1;          return -1;
                         else        else
                                 return 0;          return 0;
                 }      }
 }  }
   
 void pthrootgfs(GFS a,GFS *b)  void pthrootgfs(GFS a,GFS *b)
 {  {
         Q p;    Q p;
         int e,i;    int e,i;
         GFS t,s;    GFS t,s;
   
         STOQ(characteristic_sf(),p);    STOQ(characteristic_sf(),p);
         e = extdeg_sf()-1;    e = extdeg_sf()-1;
         t = a;    t = a;
         for ( i = 0; i < e; i++ ) {    for ( i = 0; i < e; i++ ) {
                 pwrgfs(t,p,&s); t = s;      pwrgfs(t,p,&s); t = s;
         }    }
         *b = t;    *b = t;
 }  }
   
 void randomgfs(GFS *r)  void randomgfs(GFS *r)
 {  {
         unsigned int t;    unsigned int t;
   
         if ( !current_gfs_q1 )    if ( !current_gfs_q1 )
                 error("addgfs : current_gfs_q is not set");      error("addgfs : current_gfs_q is not set");
         t = mt_genrand()%current_gfs_q;    t = mt_genrand()%current_gfs_q;
         indextogfs(t,r);    indextogfs(t,r);
 }  }
   
 /* arithmetic operations for 'immediate values of GFS */  /* arithmetic operations for 'immediate values of GFS */
   
 int _addsf(int a,int b)  int _addsf(int a,int b)
 {  {
         if ( !a )    if ( !a )
                 return b;      return b;
         else if ( !b )    else if ( !b )
                 return a;      return a;
   
         a = IFTOF(a); b = IFTOF(b);    a = IFTOF(a); b = IFTOF(b);
   
         if ( !current_gfs_ntoi ) {    if ( !current_gfs_ntoi ) {
                 a = a+b-current_gfs_q;      a = a+b-current_gfs_q;
                 if ( a == 0 )      if ( a == 0 )
                         return 0;        return 0;
                 else {      else {
                         if ( a < 0 )        if ( a < 0 )
                                 a += current_gfs_q;          a += current_gfs_q;
                         return FTOIF(a);        return FTOIF(a);
                 }      }
         }    }
   
         if ( a > b ) {    if ( a > b ) {
                 /* tab[a]+tab[b] = tab[b](tab[a-b]+1) */      /* tab[a]+tab[b] = tab[b](tab[a-b]+1) */
                 a = current_gfs_plus1[a-b];      a = current_gfs_plus1[a-b];
                 if ( a < 0 )      if ( a < 0 )
                         return 0;        return 0;
                 else {      else {
                         a += b;        a += b;
                         if ( a >= current_gfs_q1 )        if ( a >= current_gfs_q1 )
                                 a -= current_gfs_q1;          a -= current_gfs_q1;
                         return FTOIF(a);        return FTOIF(a);
                 }      }
         } else {    } else {
                 /* tab[a]+tab[b] = tab[a](tab[b-a]+1) */      /* tab[a]+tab[b] = tab[a](tab[b-a]+1) */
                 b = current_gfs_plus1[b-a];      b = current_gfs_plus1[b-a];
                 if ( b < 0 )      if ( b < 0 )
                         return 0;        return 0;
                 else {      else {
                         b += a;        b += a;
                         if ( b >= current_gfs_q1 )        if ( b >= current_gfs_q1 )
                                 b -= current_gfs_q1;          b -= current_gfs_q1;
                         return FTOIF(b);        return FTOIF(b);
                 }      }
         }    }
 }  }
   
 int _chsgnsf(int a)  int _chsgnsf(int a)
 {  {
         if ( !a )    if ( !a )
                 return 0;      return 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 a = current_gfs_q-IFTOF(a);      a = current_gfs_q-IFTOF(a);
                 return FTOIF(a);      return FTOIF(a);
         } else if ( current_gfs_q1&1 )    } else if ( current_gfs_q1&1 )
                 return a;      return a;
         else {    else {
                 /* r^((q-1)/2) = -1 */      /* r^((q-1)/2) = -1 */
                 a = IFTOF(a);      a = IFTOF(a);
                 a += (current_gfs_q1>>1);      a += (current_gfs_q1>>1);
                 if ( a >= current_gfs_q1 )      if ( a >= current_gfs_q1 )
                         a -= current_gfs_q1;        a -= current_gfs_q1;
                 return FTOIF(a);      return FTOIF(a);
         }    }
 }  }
   
 int _subsf(int a,int b)  int _subsf(int a,int b)
 {  {
         if ( !a )    if ( !a )
                 return _chsgnsf(b);      return _chsgnsf(b);
         else if ( !b )    else if ( !b )
                 return a;      return a;
         else    else
                 return _addsf(a,_chsgnsf(b));      return _addsf(a,_chsgnsf(b));
 }  }
   
 int _mulsf(int a,int b)  int _mulsf(int a,int b)
 {  {
         int c;    int c;
   
         if ( !a || !b )    if ( !a || !b )
                 return 0;      return 0;
         else if ( !current_gfs_ntoi ) {    else if ( !current_gfs_ntoi ) {
                 a = IFTOF(a); b = IFTOF(b);      a = IFTOF(a); b = IFTOF(b);
                 DMAR(a,b,0,current_gfs_q,c);      DMAR(a,b,0,current_gfs_q,c);
                 return FTOIF(c);      return FTOIF(c);
         } else {    } else {
                 a = IFTOF(a) + IFTOF(b);      a = IFTOF(a) + IFTOF(b);
                 if ( a >= current_gfs_q1 )      if ( a >= current_gfs_q1 )
                         a -= current_gfs_q1;        a -= current_gfs_q1;
                 return FTOIF(a);      return FTOIF(a);
         }    }
 }  }
   
 int _invsf(int a)  int _invsf(int a)
 {  {
         if ( !a ) {    if ( !a ) {
                 error("_invsf : division by 0");      error("_invsf : division by 0");
                 /* NOTREACHED */      /* NOTREACHED */
                 return -1;      return -1;
         } else if ( !current_gfs_ntoi ) {    } else if ( !current_gfs_ntoi ) {
                 a = invm(IFTOF(a),current_gfs_q);      a = invm(IFTOF(a),current_gfs_q);
                 return FTOIF(a);      return FTOIF(a);
         } else {    } else {
                 a = current_gfs_q1 - IFTOF(a);      a = current_gfs_q1 - IFTOF(a);
                 return FTOIF(a);      return FTOIF(a);
         }    }
 }  }
   
 int _divsf(int a,int b)  int _divsf(int a,int b)
 {  {
         int c;    int c;
   
         if ( !b ) {    if ( !b ) {
                 error("_divsf : division by 0");      error("_divsf : division by 0");
                 /* NOTREACHED */      /* NOTREACHED */
                 return -1;      return -1;
         } else if ( !current_gfs_ntoi ) {    } else if ( !current_gfs_ntoi ) {
                 b = invm(IFTOF(b),current_gfs_q);      b = invm(IFTOF(b),current_gfs_q);
                 a = IFTOF(a);      a = IFTOF(a);
                 DMAR(a,b,0,current_gfs_q,c);      DMAR(a,b,0,current_gfs_q,c);
                 return FTOIF(c);      return FTOIF(c);
         } else if ( !a )    } else if ( !a )
                 return 0;      return 0;
         else {    else {
                 a = IFTOF(a) - IFTOF(b);      a = IFTOF(a) - IFTOF(b);
                 if ( a < 0 )      if ( a < 0 )
                         a += current_gfs_q1;        a += current_gfs_q1;
                 return FTOIF(a);      return FTOIF(a);
         }    }
 }  }
   
 int _pwrsf(int a,int b)  int _pwrsf(int a,int b)
 {  {
         GFS at,ct;    GFS at,ct;
         Q bt;    Q bt;
         int c;    int c;
   
         if ( !b )    if ( !b )
                 return _onesf();      return _onesf();
         else if ( !a )    else if ( !a )
                 return 0;      return 0;
         if ( !current_gfs_ntoi ) {    if ( !current_gfs_ntoi ) {
                 a = pwrm(current_gfs_q,IFTOF(a),b);      a = pwrm(current_gfs_q,IFTOF(a),b);
                 return FTOIF(a);      return FTOIF(a);
         } else {    } else {
                 iftogfs(a,&at);      iftogfs(a,&at);
                 STOQ(b,bt);      STOQ(b,bt);
                 pwrgfs(at,bt,&ct);      pwrgfs(at,bt,&ct);
                 c = CONT(ct);      c = CONT(ct);
                 return FTOIF(c);      return FTOIF(c);
         }    }
 }  }
   
 int _onesf()  int _onesf()
 {  {
         return !current_gfs_ntoi ? FTOIF(1) : FTOIF(0);    return !current_gfs_ntoi ? FTOIF(1) : FTOIF(0);
 }  }
   
 int _itosf(int n)  int _itosf(int n)
 {  {
         int i;    int i;
   
         /* XXX */    /* XXX */
 #if 0  #if 0
         n %= current_gfs_p;    n %= current_gfs_p;
 #else  #else
         n %= current_gfs_q;    n %= current_gfs_q;
 #endif  #endif
         if ( !n )    if ( !n )
                 return 0;      return 0;
         i = !current_gfs_ntoi ? n : current_gfs_ntoi[n];    i = !current_gfs_ntoi ? n : current_gfs_ntoi[n];
         i = FTOIF(i);    i = FTOIF(i);
         if ( n < 0 )    if ( n < 0 )
                 i = _chsgnsf(i);      i = _chsgnsf(i);
         return i;    return i;
 }  }
   
 int _isonesf(int a)  int _isonesf(int a)
 {  {
         return a == _onesf();    return a == _onesf();
 }  }
   
 int _randomsf()  int _randomsf()
 {  {
         int t;    int t;
   
         t = (int) (mt_genrand() % current_gfs_q);    t = (int) (mt_genrand() % current_gfs_q);
         if ( !current_gfs_ntoi )    if ( !current_gfs_ntoi )
                 return t ? FTOIF(t) : 0;      return t ? FTOIF(t) : 0;
         else    else
                 return t!=current_gfs_q1 ? FTOIF(t) : 0;      return t!=current_gfs_q1 ? FTOIF(t) : 0;
 }  }
   
 int field_order_sf()  int field_order_sf()
 {  {
         return current_gfs_q;    return current_gfs_q;
 }  }
   
 int characteristic_sf()  int characteristic_sf()
 {  {
         return current_gfs_p;    return current_gfs_p;
 }  }
   
 int extdeg_sf()  int extdeg_sf()
 {  {
         if ( !current_gfs_ext )    if ( !current_gfs_ext )
                 return 1;      return 1;
         else    else
                 return UDEG(current_gfs_ext);      return UDEG(current_gfs_ext);
 }  }
Legend:
Removed from v.1.18  
changed lines
  Added in v.1.19
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>